function J = mydemosaic(I)
    J = I;
    
    M = size(I, 1);
    N = size(I, 2);
    %channel = size(size(I), 2);    % == 2 => one channel
                                    % == 3 => three channel


    red_mask = uint8(repmat([1 0; 0 0], M/2, N/2));
    green_mask = uint8(repmat([0 1; 1 0], M/2, N/2));
    blue_mask = uint8(repmat([0 0; 0 1], M/2, N/2));

    R=I.*red_mask;
    G=I.*green_mask;
    B=I.*blue_mask;

    % Interpolation for the green at the missing points
    G= G + imfilter(G, [0 1 0; 1 0 1; 0 1 0]/4);
    
    % Interpolation for the blue at red points
    B1 = imfilter(B,[1 0 1; 0 0 0; 1 0 1]/4);
    
    % Interpolation for blue at green points
    B2 = imfilter(B+B1,[0 1 0; 1 0 1; 0 1 0]/4);
    B = B + B1 + B2;
    
    % Interpolation for the red at blue points
    R1 = imfilter(R,[1 0 1; 0 0 0; 1 0 1]/4);
    
    % Interpolation for red at green points 
    R2 = imfilter(R+R1,[0 1 0; 1 0 1; 0 1 0]/4);
    R = R + R1 + R2;
    
    J(:,:,1)=R; 
    J(:,:,2)=G; 
    J(:,:,3)=B;
    
end