clear;
load("BlockModeTimeWindowEffects.icp");
mxFact = [1.0; 5.0; 25.0];
numTimeWindows = length(mxFact);

pulseAmplitude=cell(numTimeWindows);
pulseTime = cell(numTimeWindows);
tw = matrix(numTimeWindows,1);
for (i=1:numTimeWindows){
    extFact = mxFact(i);
    feval("calculateParameters.lsf");
    run;
    waituntildone;
    
    #extract results
    myr=getresult("OOSC_2","mode 1/signal");
    x=myr.getparameter(myr.getparameter);
    y=myr.getattribute(myr.getattribute);
    pulseTime{i} = x;
    pulseAmplitude{i}=y;
    tw(i) = getnamed("","time window");
    
    
    #extract input pulse
    if(i==1){
         myr=getresult("OOSC_1","mode 1/signal");
         inputPulseTime=myr.getparameter(myr.getparameter);
         inputPulseAmplitude=myr.getattribute(myr.getattribute);
    }
        
}
#plot results
I=sqrt(-1);
width=getnamed("GAUSS_1","width");
inputPower = getnamed("GAUSS_1","power");
fc=getnamed("GAUSS_1","frequency");
chirp=getnamed("GAUSS_1","chirp");


Tfwhm=width*timeWindow;
T0g=Tfwhm/(2*sqrt(log(2)));

#plot input pulse
myLegend = cell(numTimeWindows+1);
tc= ( min(inputPulseTime) + max(inputPulseTime) ) /2.0;
plot((inputPulseTime-tc)/T0g,abs(inputPulseAmplitude)^2/inputPower);
holdon;
myLegend{1} = "z=0 (i.e, input pulse)";

#plot output pulses
format short;
for(i=1:numTimeWindows){
    tc= ( min(pulseTime{i}) + max(pulseTime{i}) ) /2.0;
    plot((pulseTime{i}-tc)/T0g,abs(pulseAmplitude{i})^2/inputPower);
    holdon;
    myLegend{i+1} = "z=L, time window=" + num2str( tw(i)/1e-12 ) +" ps";
       
}
legend(myLegend);
setplot("x min", -tw(numTimeWindows)/2/T0g);
setplot("x max",tw(numTimeWindows)/2/T0g);
setplot("x label","(t-L/v_g)/T_0");
setplot("y label", "(Pulse Power)/(Input Pulse Power)");
format long;
holdoff;