
X = imread('moon.png');

% Roberts partial derivatives
dx =[-1 0; 0 1];
dy =[0 1; -1 0];

X = double(X);
DX = imfilter(X, dx);
% rescaled in [0, 255]
DX2 = (DX-min(DX(:)))/(max(DX(:))-min(DX(:)))*255;

DY = imfilter(X, dy);
% rescaled in [0, 255]
DY2 = (DY-min(DY(:)))/(max(DY(:))-min(DY(:)))*255;

% magnitude of the gradient 
DXY = sqrt(DX.^2+DY.^2);

fprintf('DX (x derivative) range: [%d, %d]\n', [min(DX(:)) max(DX(:))]);
fprintf('DY (y derivative) range: [%d, %d]\n', [min(DY(:)) max(DY(:))]);
fprintf('DXY (gradient) range: [%f, %f]\n', [min(DXY(:)) max(DXY(:))]);

DXY = uint8(DXY);
DXY2 = histeq(DXY);

X = uint8(X);
DX = uint8(DX);
DX2 = uint8(DX2);
DY = uint8(DY);
DY2 = uint8(DY2);
DXY = uint8(DXY);

figure; imshow(X); title('Original image');
figure; imshow(DX); title('X derivative (Roberts)');
figure; imshow(DX2); title('X derivative rescaled (Roberts)');
figure; imshow(DY); title('Y derivative (Roberts)');
figure; imshow(DY2); title('Y derivative rescaled (Roberts)');
figure; imshow(DXY); title('Gradient (Roberts)');
figure; imshow(DXY2); title('Equalized gradient (Roberts)');
