% Motion Estimation by BlockMatching

% aviread is remplaced by VideoReader in Matlab R20015a
vidObj = VideoReader('Vectra21Frames.mpg');

% sequence loading

vidHeight  = vidObj.Height;
vidWidth   = vidObj.Width;
seq = struct('cdata', zeros(vidHeight, vidWidth, 3, 'uint8'), ...
    'colormap', []); 

k = 1;
while hasFrame(vidObj)
    seq(k).cdata = readFrame(vidObj);
    k=k+1;
end

% Video sequence display
figure (1); 
movie(seq);

%  We consider 2 successive images
im1 = seq(51).cdata;
im2 = ; %%% TO COMPLETE

% Conversion : RGB -> YUV, we use the Y component
YCBCR1 = rgb2ycbcr(im1);
YCBCR2 = ; %%% TO COMPLETE

Y1 = double( YCBCR1(:, :, 1) ); % the reference : Y(t-1)
Y2 = ; %%% A COMPLETER % the target : Y(t)

% We consider only a sub-image (a part of the image) => less computation
% Rq : The image dimension has to be modulo the block size (namely 8)
Y1 = Y1(1:120, 81:160);
Y2 = Y2(1:120, 81:160);

[nblignes, nbcolonnes] = size( Y1 );

% MSE Definition
MSE = @(x1, x2) mean2((x1 - x2).^2);

% Loop on the blocks "b" of the target Y(t), namely here Y2
%   => taking into account of the border effect 
%       (such as the research window is always included into the reference  
%       Y(t-1), namely here Y1)           
%       ("we remove " 8 rows and 8 colums around Y2)

% Display of Y2
 
figure(2);
imshow( uint8( Y2 ) ); 

    
    for i = 9 : 8 : nblignes-9
        for j = 9 : 8 : nbcolonnes-9
            
            % taking of b
            b = Y2(i:i+7, j:j+7);
            
            % taking of the research window "B" into the reference 
            % Y(t-1), namely here Y1 
            B = Y1(i-8:i+7, j-8:j+7);
            
            % MSE computation
            % Rq : nlfilter computes for all the blocks (8x8) of B
            % without overlapping
            EQM = ; %%% TO COMPLETE with "nlfilter"
            
            % selecr the minimum MSE
            EQM_min ; %%% TO COMPLETER with "min"
            
            % position of the block with the minimal MSE
            % Rq : if they are numerous, 
            %   take the 1st occurency
            [I, J] = ; %%% TO COMPLETE with "find"
            
            % Display of the motion vectors on Y2
            hold on;
            X = [ i+4 i+I(1)-4 ];
            Y = [ j+4 j+J(1)-4 ];
            line(X, Y);                  
                       
        end
    end