% This is the program for vignetting for ee362
% Title: fstop.m 
% written by Stephen Kaye
% Started on 3/04/06
%
% This program takes a cut of the image data in a radial
% direction.  The middle of the picture is where the vector
% starts and goes toward the edge by a slope of 1 x-pixel
% to 1 y-pixel.  This could be changed to take a different
% cut along the data.
% 
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%

%clear all can't clear all, need the vectors to stay around 

filename = 'sd500_01_7.1f.jpg';
image = imread(filename);
%imshow(image);

db_image = double(image);
Y = 0.30.*db_image(:,:,1) + 0.59.*db_image(:,:,2) + 0.11.*db_image(:,:,3);

% Now we have the luminance of the image.  Let's take a cut in the radial direction
% find center of picture

gray_im = Y;  % this is the intensity, use it 
gray_doub = double(gray_im);

clear 'Y';

% let's resize the image so that Matlab can compute the contour
im_size = size(gray_doub);
aspect_ratio = im_size(2)/im_size(1);
rebinned = imresize(gray_doub, [750 (750*aspect_ratio)]);

re_size = size(rebinned);
x_center = re_size(1)/2;
y_center = re_size(2)/2;

slope = aspect_ratio; 

% Now filter the image
% We're interested in low frequency light falloff, not pixel to pixel noise
filter_2d = fspecial('gaussian', [3 3], 0.3);
filtered_gray_im = imfilter(rebinned, filter_2d);

max_filt_gray_im = max(max(filtered_gray_im));
norm_filt_gray = filtered_gray_im./max_filt_gray_im; 

for i = 1:x_center 
	inten_dist(i) = filtered_gray_im(x_center+i,y_center+i);
end

% normalize the intensity distribution
max_inten = max(max(inten_dist));
norm_inten_71 = inten_dist./max_inten;
cos_index = [1:numel(norm_inten_71)]*2.65e-2*pi/180;
index = [1:numel(norm_inten_71)]*5.22;

plot(index,norm_inten_71, index, cos(cos_index)); 
print -djpeg test.jpg;


% All this plotting needs to take place after the script
% has been run on images with various fstops.
% After that, then we can plot the data all on one graph.

%plot(norm_inten);
%title('Intensity from middle of picture');
%xlabel('Distance');
%ylabel('Intensity'); 
%print -djpeg intensity60.jpg;
 
