clc;
clear all;

fid1 = fopen('nissl-rat.vol');
status1 = fseek(fid1, 0 , 0);
x_size1 = fread (fid1, 1, 'int32');
y_size1 = fread (fid1, 1, 'int32');
z_size1 = fread (fid1 , 1, 'int32');

Total_size = x_size1 * y_size1 * z_size1;
% fileInfo = dir('nissl_rat.vol');
% fileSizeInBytes = fileInfo.bytes;

x_size = double(x_size1);
y_size = double(y_size1);
z_size = double(z_size1);

position1 = ftell(fid1);
M = (fread(fid1,[x_size1 y_size1*z_size1], 'uint8'));
A = reshape(M, [x_size1 y_size1 z_size1]);
b = A(:,:,2);
Img = mat2gray(b);
figure;
imshow(Img);

position2 = ftell(fid1);
fclose(fid1);

Im = mat2gray(A(101:228, 101:228, 101:228));
%% 3D gradient
sx = [-1 0 1; -2 0 2; -1 0 1];
sy = [1 2 1; 0 0 0; -1 -2 -1];
szx = cat(3, sx, sx, sx);
szy = cat(3, sy, sy, sy);
szz = cat(3, ones(3,3), zeros(3,3), -ones(3,3));
Gx = convn(Im, szx, 'same');
Gy = convn(Im, szy, 'same');
Gz = convn(Im, szz, 'same');
grdmg = sqrt(Gx.^2 + Gy.^2 + Gz.^2);