function D3M_2Dfun(RR,videoname, foldername, T, Tshift, step)
% D3M-2Dfun(RR) - creates 2D contour D-3M movie frames
%
% input --> RR: two-column matrix containing
%           in the first column the time instants of the R peaks (in s)
%           in the second column the RR intervals (in s)
%           videoname: string containing the name for the 
%           video
%           foldername (optional): string containing the name of the
%           folder where the frames will be saved
%           T (default=10): time interval to display in one frame (in min)
%           Tshift (default=1): time shift between frames (in min)
%           step (default=0.15): step between plot contour lines
%
% Developed by A Burykin (burykin@gmail.com), L Citi, 
% MD Costa and AL Goldberger
% and modified by S Mariani (sara.mariani@wyss.harvard.edu)
% Wyss Institute at Harvard
% when using this function, please reference: Burykin, Anton et al.
% "Dynamical density delay maps: simple, new method for visualising 
% the behaviour of complex systems." BMC medical informatics and 
% decision making 14.1 (2014): 6.

narginchk(2, 6);

if (nargin>=3 && ~isempty (foldername))
    % save frames in folder specified by user
    savepic=1;
    mkdir(foldername);       
    else savepic=0;
end

if (nargin<4 || isempty (T))
    % use default interval
    T=10;
end
    
if (nargin<5 || isempty (Tshift))
    % use default shift
    Tshift=1;
end

if (nargin<6 || isempty (step))
    % use default step
    step=0.15;
end

vidObj = VideoWriter(videoname);
vidObj.FrameRate = 10;
open(vidObj);


t=RR(:,1);
t=t-t(1);
t=t/60.0; % sec2min
rr=RR(:,2);

scrsz = get(0,'ScreenSize');

k1=1;
k2=find(t>=T,1);

lolim=min(RR(:,2));
hilim=max(RR(:,2));

k=1;
while t(k2)<=t(end)
    
    t1=t(k1);
    t2=t(k2);
    
    t_disp=t(k1:k2);
    rr_disp=rr(k1:k2);

    y1=rr_disp(1:end-1);
    y2=rr_disp(2:end);

    [ctrs1, ctrs2, F, hAxes] = dscatter2(y1,y2,'PLOTTYPE','contour');
    close
    F = F./max(F(:));
 
    fig1=figure('Position',...
    [0.05*scrsz(3) 0.05*scrsz(4) 0.3*scrsz(3) 0.8*scrsz(4)],...
    'Color',[1 1 1]);

    subplot(4,1,[1 3])
    contour(ctrs1, ctrs2, F,'LineWidth',2,'LevelStep',step)
    colorbar;
    xlim([0.3 1.5])
    ylim([0.3 1.5])
    grid on
    axis square
    title('Contour Map');
    xlabel('RR(n), sec');
    ylabel('RR(n+1), sec');   
    
    subplot(4,1,4)
    plot(t_disp,rr_disp,'r','LineWidth',1.5)
    grid on
    xlim([min(t_disp) max(t_disp)])
    ylim([lolim hilim])
    title('Cardiac Interbeat Interval Time Series');
    ylabel('RR, sec');
    xlabel('Time, min');
    
    % capture frame here  
    Fr1=getframe(fig1);
    Im1=frame2im(Fr1);
    
    if savepic
    im_file=[foldername '/' int2str(k) '.bmp'];
    imwrite(Im1,im_file);
    end

    writeVideo(vidObj,Fr1);
    clear Fr1 Im1 Im1_map im_file;
    close(fig1);
    
    k1=find(t>=t1+Tshift,1);
    k2=find(t>=t2+Tshift,1);
    k=k+1;
end
close(vidObj);