%-------------------------------------------------------------------------------
% Project 3
% Group 8
% Isaac Gerg (idg101)
% CSE 486
%-------------------------------------------------------------------------------

clc;
close all;

%-------------------------------------------------------------------------------
% Read in images
	imgwheelnoise = imread('wheelnoise.gif', 'gif');

%-------------------------------------------------------------------------------
% Part a - Sobel Edge Detection
    imgwheelnoise_SobelEdgeDetection = imdetectedges( imgwheelnoise, 'sobel', 0.29);
	imwrite(imgwheelnoise_SobelEdgeDetection, 'wheelnoise_sobel.jpg', 'jpg');
    
%-------------------------------------------------------------------------------
% Part b - Hough Transform
    disp('Running voting algorithm...');

    E = imgwheelnoise_SobelEdgeDetection;
    d_rho = 3;          % 0 <= rho <= 363
    d_theta = 1;        % -pi/2 < theta < pi
    Rmin = -256;
    Rmax = 365;
    R = Rmax - Rmin;    % Adding 0 instead of 1 so array A has integer height.
    Tmin = -89;
    Tmax = 179;
    T= Tmax - Tmin + 1;     % Add 1 to account for zero.
    % Conversion factors
    toRadians = pi/180;
    toDegrees = 180/pi;
    
   % Initialize A(R,T), P
   % A is accumulator, P is the list of pixels that voted for a particular point.
    A = zeros(T / d_theta , R / d_rho);
    P = cell(T / d_theta , R / d_rho);
    
    for i = 1:258,
        for j = 1:258,
            for h = Tmin:Tmax,
                if E(i,j) == 255
                    rho = i*cos(h*d_theta*toRadians) + j*sin(h*d_theta*toRadians) - Rmin;
                    k = round((rho)/d_rho) + 1;
                    A(h + 90, k) = A(h + 90, k) + 1;
                    mxTmp = P{h + 90, k};
                    P{h + 90, k} = [mxTmp; i j];
                end
            end
        end
    end
    
   disp('Done.');
   figure('Name', 'Hough Space', 'NumberTitle', 'off', 'MenuBar', 'none');
   imagesc(A);
   title('Hough Space');
   xlabel('((rho + 256) / 3) (pixels)');
   ylabel('theta + 90 (degrees)');
 
%-------------------------------------------------------------------------------
% Part c - Display Hough Transform
    disp('Creating image from hough space...');    

    iMin = min(min(A));
    iMax = max(max(A));
    
    % ---- f = 0.0 ----
    f = 0.00;
    intT = iMin + f * (iMax - iMin);
    
    imgHough0_0 = zeros(258, 258);

    for theta = 1:(T / d_theta),
        for rho = 1:(R / d_rho),
            if A(theta, rho) > intT
                mxTmp = P{theta, rho};
                for i = 1:A(theta, rho),
                    imgHough0_0(mxTmp(i,1), mxTmp(i,2)) = 255;
                end
            end
        end
    end

    % ---- f = 0.25 ----
    f = 0.25;
    intT = iMin + f * (iMax - iMin);
    
    imgHough0_25 = zeros(258, 258);
    
    for theta = 1:(T / d_theta),
        for rho = 1:(R / d_rho),
            if A(theta,rho) > intT
                mxTmp = P{theta, rho};
                for i = 1:A(theta, rho),
                    imgHough0_25(mxTmp(i,1), mxTmp(i,2)) = 255;
                end
            end
        end
    end
    
    % ---- f = 0.50 ----
    f = 0.50;
    intT = iMin + f * (iMax - iMin);
    
    imgHough0_50 = zeros(258, 258);
    
    for theta = 1:(T / d_theta),
        for rho = 1:(R / d_rho),
            if A(theta,rho) > intT
                mxTmp = P{theta, rho};
                for i = 1:A(theta, rho),
                    imgHough0_50(mxTmp(i,1), mxTmp(i,2)) = 255;
                end
            end
        end
    end

    % ---- f = 0.75 ----
    f = 0.75;
    intT = iMin + f * (iMax - iMin);
    
    imgHough0_75 = zeros(258, 258);
    
    for theta = 1:(T / d_theta),
        for rho = 1:(R / d_rho),
            if A(theta,rho) > intT
                mxTmp = P{theta, rho};
                for i = 1:A(theta, rho),
                    imgHough0_75(mxTmp(i,1), mxTmp(i,2)) = 255;
                end
            end
        end
    end
    
   disp('Done.');
   figure('Name', 'Images After Hough Transform', 'NumberTitle', 'off', 'MenuBar', 'none');
   colormap('gray');
   subplot(2,2,1);
   imagesc(imgHough0_0);
   title('f = 0.00');
   subplot(2,2,2);
   imagesc(imgHough0_25);
   title('f = 0.25');
   subplot(2,2,3);
   imagesc(imgHough0_50);
   title('f = 0.50');
   subplot(2,2,4);
   imagesc(imgHough0_75);
   title('f = 0.75');

   
   imwrite(imgHough0_0, 'wheelnoise_hough_transform0_0.jpg', 'jpg');
   imwrite(imgHough0_25, 'wheelnoise_hough_transform0_25.jpg', 'jpg');
   imwrite(imgHough0_50, 'wheelnoise_hough_transform0_50.jpg', 'jpg');
   imwrite(imgHough0_75, 'wheelnoise_hough_transform0_75.jpg', 'jpg');