%**************************************************************************
%Initialize variables
t=0:5e-6:1e-3;              % Time
A=2;                        % Amplitude
f=1000;                     % Frequency
w=(2*pi)*f;                 % Omega
Vdc=0;                      % DC component

timepoints=length(t);       % Number of data points to plot each harmonic
lastharmonic=171;           % Last harmonic number
endslow=11;                 % Last harmonic in slow part of video
HA=zeros(lastharmonic,timepoints);     % Harmonic array

%**************************************************************************
%Variables associated with the getframe and movie functions
num=1;                      % Frame counter for video
nframes = 1;
Frames = moviein(nframes);  % Movie frames
mov = avifile('D:\Data\Research\StudentProjects\FourierConceptIllustration\Superposition.avi','compression','Cinepak')   %movie frame to avi conversion
%}
%**************************************************************************
% Calculate harmonics 1..endslow and make slow part of movie
cla                                     % Clear plot window
for (x=1:2:lastharmonic)                % Only calculate odd harmonics
    HA(x,:)=((2*A)/(pi*x))*sin(x*w*t);  % Calculate harmonic x
    HAtotal=sum(HA)+Vdc;                % Add harmonic x to the total
    % Put harmonic sum in upper plot window
    subplot(2,1,1);                     % Use upper plot window
    plot(t,HAtotal);                    % Plot new sum of harmonic series
    xlabel('time (s)');  ylabel('Voltage (V)'); % Label axes
    xlim([-1e-4 1.1e-3]);               % Keep axes constant and large
    ylim([-2,2]);                       %   enough to show entire plot
    
    % Show this harmonic in lower plot window
    subplot(2,1,2);                     % Use lower plot window
    plot(t,HA(x,:));                    % Plot this harmonic
    xlabel('time (s)');  ylabel('Voltage (V)'); % Label axes
    xlim([-1e-4 1.1e-3]);               % Keep axes constand and large
    ylim([-2,2]);                       %   enough to show largest harmonic
%{    
    if x==1      % Save graphic of fundamental to file
        print -r300 -djpeg 'D:\Data\Research\StudentProjects\FourierConceptIllustration\fundamental.jpg';
    end;
%}   
    %Capture extra frames to make first part of video slower
    if x <= endslow                          % Last harmonic to show in slow-motion
        for (y=1:1:9)
            Frames(num)=getframe(gcf);               %getframe function
            num=num+1;
        end
        pause(1);                           % Pause for 1 second so viewer can see
%        if x=endslow     % Save plot of fundamental plus harmonics 3..endslow
%            print -r300 -djpeg 'D:\Data\Research\StudentProjects\FourierConceptIllustration\serieswithfewharmonics.jpg';
%        end
    end
    Frames(num)=getframe(gcf);               %getframe function
    num = num +1;
end
%**************************************************************************
% Plot final square wave and superimposed component harmonics
xsquare=[0, 0, 0.5e-3, 0.5e-3, 1e-3, 1e-3]; % x-axis data points for square wave
ysquare=[0, 1,   1,      -1,    -1,    0];  % y-axis data points for square wave
% Square wave in upper plot window
subplot(2,1,1);             
plot(xsquare,ysquare);
xlabel('time (s)');  ylabel('Voltage (V)');
xlim([-1e-4 1.1e-3]);
ylim([-2 2]);
% Superimposed harmonics in lower plot window
subplot(2,1,2);             
plot(t,HA);
xlabel('time (s)');  ylabel('Voltage (V)');
xlim([-1e-4 1.1e-3]);
ylim([-2,2]);
% Save extra frames to hold final image
for x=1:20    
    Frames(num)=getframe(gcf);
    num = num +1;
end
% Save graphic of square wave and superimposed harmonics
%print -r300 -djpeg 'D:\Data\Research\StudentProjects\FourierConceptIllustration\squarewave.jpg';
mov = addframe(mov,Frames);
mov = close(mov);
%**************************************************************************