%% Walk through the image of when pixels were first seen by shadow edge and
%% interpolate where the A and B points were on a scanline of the
%% groundplane. A & B are found by interpolating the A,B coordinates of the
%% real frames before and after this intermediate frame.

str = sprintf('output_%05d.jpg',startimg);
img = imread(str);
imshow(img);

%% pick the scanlines
[upperrow dummy] = ginput(1);
[lowerrow dummy] = ginput(1);

Idist = double(zeros(sizey,sizex));

%% reconstruct things inside the scanlines
for i = upperrow+1:lowwerrow-1
    for j = 1:sizex
        t = Iedgetime(i,j);
        if (t ~= endimg+1)
            ceilt = ceil(t);
            floort = floor(t);
        
            alpha = mod(t,1);
                    
            A = getLinePoint(upperrow,ceilt,Ishadow,Icontrast,sizex);
            B = getLinePoint(lowerrow,ceilt,Ishadow,Icontrast,sizex);
            Ap = getLinePoint(upperrow,floort,Ishadow,Icontrast,sizex);
            Bp = getLinePoint(lowerrow,floort,Ishadow,Icontrast,sizex);
            if (A ~= [-1 -1]) 
                if (B ~= [-1 -1]) 
                    if (Ap ~= [-1 -1]) 
                        if (Bp ~= [-1 -1])
                            
                            Amix = (1-alpha)*Ap+alpha*A;
                            Bmix = (1-alpha)*Bp+alpha*B;
                            
                            %% Cast rays to Amix, Bmix to CGroundPlane to
                            %% get camera space coordinates
                            AmixRay = normalize(Amix',fc,cc,kc,alpha_c);
                            AmixRay = [AmixRay(1) AmixRay(2) 1];
                            t = rayplane([0,0,0],AmixRay,CXaxis,CYaxis,COrigin);
                            CAmix = AmixRay*t;
                            BmixRay = normalize(Bmix',fc,cc,kc,alpha_c);
                            BmixRay = [BmixRay(1) BmixRay(2) 1];
                            t = rayplane([0,0,0],BmixRay,CXaxis,CYaxis,COrigin);
                            CBmix = BmixRay*t;

                            %% Make Plane from CAmix CBmix LightSource
                            
                            %% Find intersection between this pixels ray,
                            %% and the shadow plane
                            PixelRay = normalize([i,j]',fc,cc,kc,alpha_c);
                            PixelRay = [PixelRay(1) PixelRay(2) 1];
                            
                            t = rayplane([0,0,0],PixelRay,CAmix,CBmix,LightSource);
                            
                            %% The 3d point we just located, for plotting
                            %% or anything else...
                            Point = t*PixelRay;
                            
                            %We just want the projection z
                            %distance
                            Idist(i,j) = t;
                        end 
                    end 
                end 
            end
        end
    end
end
