Construct a 256 x 256 image of the function gaus(r) = exp(-pi*r^2) over the domain -1.5 < (x,y) < 1.5. Then construct an image of the Laplacian of gaus(r) over the same domain. Find the rms difference between this image and those constructed by applying a Laplacian convolution filter to the original image. Plot the rms difference as a function of alpha. Find the value of alpha that minimizes this difference.
global m1 m2 z lz range steps; range = 3; steps = 256; x = linspace(-range/2,range/2,steps); y = -x; x = ones(steps,1)*x; y = y' * ones(1,steps); r2 = x.*x + y.*y; z = exp(-pi*r2); view = imnorm(z); imshow(view); imwrite(view,'gaus.jpg');
The image constructed is the negative of the Laplacian.
lz = -4*pi*(pi*r2-1).*exp(-pi*r2);
view = imnorm(lz);
imshow(view);
imwrite(view,'lgaus.jpg');
m1 = [0 1 0; 1 -4 1; 0 1 0];
m2 = [1 0 1; 0 -4 0; 1 0 1];
clear rms_error
lmin = -.5; lmax = 0.5; lsteps = 31;
for i = 0:(lsteps-1)
a = lmin+ (lmax-lmin)*i/(lsteps-1);
rms_error(i+1) = func(a);
end
%Plot the Error Graph
plot(linspace(lmin,lmax,lsteps),rms_error);
title('RMS Error vs. Alpha');
xlabel('Alpha');
ylabel('RMS Error');
%Determine the Best Alpha based upon RMS Error amin = fminsearch(@func,0) %Show the Best Laplacian Matrix lap = ((1-amin).*m1+amin.*m2)/(1+amin)
The results are
amin =
-0.2174
lap =
-0.2779 1.5557 -0.2779
1.5557 -5.1114 1.5557
-0.2779 1.5557 -0.2779
function error = func(v) global m1 m2 z lz range steps a = v(1); lap = ((1-a)*m1+a*m2)/(1+a); zlap = -filter2(lap,z,'valid')/(range/(steps-1))^2; diff = zlap-lz; error = sqrt(mean(diff(:).^2));
function out = imnorm(inp) % IMNORM % out = imnorm(inp) % normalizes input image between (0..1) % for unipolar image divide by the maximum value % bipolar image mapped such that zero -> 0.5 and image scales between (0..1) maxv = max(inp(:)); minv = min(inp(:)); if (minv<0.0) vn = max(maxv,-minv); out = (inp/vn+1)/2; else out = inp/maxv; end
Maintained by John Loomis, last updated 14 March 2000