function [] = SPPplus_figures_v1(resp_wave,Rec_wav,Gp_t,Gpp_t,G_speed,...
    Gp_t_dot,Gpp_t_dot,time_wave,f_domain,ft_amp_sign_fft_resp_s,M,W,...
    print_figs,analtype,filename,figsvect)
% This function produces figures using the SPP analysis; if specified by
% the user, the figures can also be saved to .jpg files
    %Inputs: resp_wave = Lx3 matrix of the strain, rate and stress data
                %input to the analysis, with each row representing a
                %measuring point
            %Rec_Wav = Lx3 matrix of the strain, rate and stress data
                %output by the analysis, with each row representing a
                %measuring point (if there is no difference between this
                %and resp_wave, enter 0)
            %Gp_t = Lx1 vector contaning G'_t values output by the analysis
            %Gpp_t = Lx1 vector contaning G''_t values output by the analysis
            %G_speed = Lx1 vector contaning modulus speeds output by the analysis
            %Gp_t_dot = Lx1 vector contaning dG'_t/dt values output by the analysis
            %Gpp_t_dot = Lx1 vector contaning dG''_t/dt values output by the analysis
            %time_wave = Lx1 vector of time at each measurement point
            %f_domain = Lx1 vector listing the range of fourier harmonics,
                %(enter 0 if not defined)
            %ft_amp_sign_fft_resp_s =  Lx1 vector of normalized fourier
                %magnitudes(enter 0 if not defined)
            %M = number of harmonics for stress (enter 0 if not defined)
            %W = number of obtainable harmonics (enter 0 if not defined)
            %print_figs = binary value determining whether to save figures
                %as .jpg files
            %analtype = analysis method used (as a string)
            %filename = name of file from which the data originated
            %figsvect = vector that specifies which figures to create

%Produce figure title
fig_name=sprintf('SPPplus_%s_v1 by Rogers group@UIUC',analtype);

%Standard SPP metrics figure
if figsvect(1) == 1
    figure('Name',fig_name)
    %Stress-Strain plot
    subplot(2,3,1)
    plot(resp_wave(:,1),resp_wave(:,3),'k--')
    if Rec_wav == 0
    else
        hold on
        plot(Rec_wav(:,1),Rec_wav(:,3),'r--o','MarkerSize',1)
    end
    hold on
    xlim1 = get(gca,'xlim');  %Get x range 
    ylim2 = get(gca,'ylim');  %Get y range 
    maxX=max(abs(xlim1));
    maxY=max(abs(ylim2));
    axis([-1.2*maxX 1.2*maxX -1.2*maxY 1.2*maxY]);
    hold on
    plot([-1.2*maxX 1.2*maxX],[0 0],'k:','LineWidth',0.1)
    hold on
    plot([0 0],[-1.2*maxY 1.2*maxY],'k:','LineWidth',0.1)
    xlabel('Strain [-]')
    ylabel('Stress [Pa]')
    title('Stress-strain')
    if Rec_wav == 0
        legend('Raw stress output','Location','northeast')
    else
        legend('Raw stress output', 'Reconstructed stress','Location','northeast')
    end
    hold off
    
    %Stress-Rate plot
    subplot(2,3,4)
    plot(resp_wave(:,2),resp_wave(:,3),'k--')
    if Rec_wav == 0
    else
        hold on
        plot(Rec_wav(:,2),Rec_wav(:,3),'r--o','MarkerSize',1)
    end
    xlim1 = get(gca,'xlim');  %Get x range 
    ylim2 = get(gca,'ylim');  %Get y range 
    maxX=max(abs(xlim1));
    maxY=max(abs(ylim2));
    axis([-1.2*maxX 1.2*maxX -1.2*maxY 1.2*maxY]);
    hold on
    plot([-1.2*maxX 1.2*maxX],[0 0],'k:','LineWidth',0.1)
    hold on
    plot([0 0],[-1.2*maxY 1.2*maxY],'k:','LineWidth',0.1)
    xlabel('Rate [1/s]')
    ylabel('Stress [Pa]')
    title('Stress-rate')
    hold off
    
    %Cole-Cole plot
    subplot(2,3,3)
    plot(Gp_t,Gpp_t,'r','LineWidth',1.5)
    xlim1 = get(gca,'xlim');  %Get x range 
    ylim2 = get(gca,'ylim');  %Get y range 
    maxX=max(Gp_t);
    minX=min(Gp_t);
    maxY=max(Gpp_t);
    minY=min(Gpp_t);
    axis([1.2*minX 1.2*maxX 1.2*minY 1.2*maxY]);
    hold on
    plot([-10^10 10^10],[0 0],'k:','LineWidth',0.1)
    hold on
    plot([0 0],[-10^10 10^10],'k:','LineWidth',0.1)
    hold on
    plot([-10^10 10^10],[-10^10 10^10],'k:','LineWidth',0.1)
    xlabel('G''_{t} [Pa]')
    ylabel('G"_{t} [Pa]')
    title('Cole-Cole plot')
    hold off
    
    %dG"_{t}/dt-dG'_{t}/dt plot
    subplot(2,3,5)
    plot(Gp_t_dot,Gpp_t_dot,'r','LineWidth',1.5)
    xlim1 = get(gca,'xlim');  %Get x range 
    ylim2 = get(gca,'ylim');  %Get y range 
    maxX=max(Gp_t_dot);
    minX=min(Gp_t_dot);
    maxY=max(Gpp_t_dot);
    minY=min(Gpp_t_dot);
    axis([1.2*minX 1.2*maxX 1.2*minY 1.2*maxY]);
    hold on
    plot([-10^10 10^10],[0 0],'k:','LineWidth',0.1)
    hold on
    plot([0 0],[-10^10 10^10],'k:','LineWidth',0.1)
    xlabel('dG''_{t}/dt [Pa/s]')
    ylabel('dG"_{t}/dt [Pa/s]')
    title('dG"_{t}/dt-dG''_{t}/dt')
    
    %Speed-G'_{t} plot
    subplot(2,3,6)
    plot(Gp_t,G_speed,'r','LineWidth',1.5)
    xlim1 = get(gca,'xlim');  %Get x range 
    ylim2 = get(gca,'ylim');  %Get y range 
    maxX=max(Gp_t);
    minX=min(Gp_t);
    maxY=max(G_speed);
    minY=min(G_speed);
    axis([1.2*minX 1.2*maxX 1.2*minY 1.2*maxY]);
    hold on
    plot([-10^10 10^10],[0 0],'k:','LineWidth',0.1)
    hold on
    plot([0 0],[-10^10 10^10],'k:','LineWidth',0.1)
    xlabel('G''_{t} [pa]')
    ylabel('Speed [Pa/s]')
    title('Speed-G''_{t}')
    
    %Speed-G"_{t} plot
    subplot(2,3,2)
    plot(G_speed,Gpp_t,'r','LineWidth',1.5)
    xlim1 = get(gca,'xlim');  %Get x range 
    ylim2 = get(gca,'ylim');  %Get y range 
    maxX=max(G_speed);
    minX=min(G_speed);
    maxY=max(Gpp_t);
    minY=min(Gpp_t);
    axis([1.2*minX 1.2*maxX 1.2*minY 1.2*maxY]);
    hold on
    plot([-10^10 10^10],[0 0],'k:','LineWidth',0.1)
    hold on
    plot([0 0],[-10^10 10^10],'k:','LineWidth',0.1)
    ylabel('G"_{t} [Pa]')
    xlabel('Speed [Pa/s]')
    title('Speed-G"_{t}')
    
    %Save figure if requested
    if print_figs == 1
        fignum = gcf;
        fignum.Position = [0,0,1920,1080];
        image1_name = sprintf('%s_SPP_LAOS_%s_PLOT.jpg',filename,analtype);
        print(image1_name,'-djpeg')
    end
end

%FT-spectra for stress
if figsvect(2) == 1
    figure('Name',fig_name)
    %Stress-Time plot
    subplot(2,1,1)
    plot(time_wave,resp_wave(:,3),'r','LineWidth',2)
    xlabel('Time [s]')
    ylabel('Strain [-]')
    title('Stress-time')
    
    %Fourier Harmonic Magnitudes plot
    subplot(2,1,2)
    stem(f_domain,ft_amp_sign_fft_resp_s,'b.','LineWidth',0.5)
    hold on
    plot([M M],[10^-10 0.7],'r-','LineWidth',2)
    set(gca,'yscale','log');
    xlim([0 W+2]);
    ylabel('I_{n}/I_{1} [-]')
    xlabel('Number of harmonics [-]')
    title('Fourier spectrum','fontweight','bold')
    
    %Save figure if requested
    if print_figs == 1
        fignum = gcf;
        fignum.Position = [0,0,1920,1080];
        image2_name = sprintf('%s_SPP_LAOS_%s_HARMONICS.jpg',filename,analtype);
        print(image2_name,'-djpeg')
    end
end

end

