fc = 195.943e12; #Center frequency (center frequency of ONA and TWLM must be the same).
fs = 2e+13; #Sample rate, which is inverse of time step (sample rate in the INTERCONNECT project file must be the same).
maxcoef = 7; #Gain fitting parameters.
tol = 1e-3;
maxiter = 100;
rolloff = 0.05;

fit_out=mczfit("sampleGainData.txt","sampleGainData",fc,fs,maxcoef,tol,maxiter,true,rolloff);
fit_out.fc = fc;
fit_out.fs = fs;

#Plot fitted and original gain
ncPlotInd = [1:12];
numcurves=length(ncPlotInd);
leg=cell(2*numcurves);
for (i=1:numcurves) {
    leg{i} = "Data Nc= "+num2str(fit_out.operatingState(ncPlotInd(i)));
    leg{i+numcurves} = "Fit Nc= "+num2str(fit_out.operatingState(ncPlotInd(i)));
}

plot(fit_out.frequency/1e12, fit_out.input(:,ncPlotInd),"frequency [THz]","Gain [1/m]");
holdon;
plot(fit_out.frequency/1e12, fit_out.output(:,ncPlotInd),"frequency [THz]", " ","Gain [1/m]","pen=--");
legend(leg);
holdoff;

#expand x-axis to accomodate legend
xx = fit_out.frequency/1e12;
xmin = xx(1);
xmax = xx(1) + 1.7*(xx(length(xx))-xx(1));
setplot("x min",xmin);
setplot("x max",xmax);
setplot("legend position",6);