function [ve, Npeak, width, pos, contrast, width_err, pos_err, contrast_err] = esrfit_Npeak(f, data, Npeak, smth) if find(isnan(data)) ve = [0,0,0,0,0,0,0]; return; end %smooting the curve n = length(data)/smth; %averging over n neighbors 35 datasmth = smooth(data,n)'; %find peaks (contrast,position) if length(findpeaks(-datasmth)) >= Npeak [peakMag,peakLoc] = findpeaks(-datasmth,'NPeaks',Npeak,'SortStr','descend'); [peakLoc,Ind] = sort(peakLoc); peakMag = peakMag(Ind); else Npeak = length(findpeaks(-datasmth)); [peakMag,peakLoc] = findpeaks(-datasmth,'NPeaks',Npeak,'SortStr','descend'); [peakLoc,Ind] = sort(peakLoc); peakMag = peakMag(Ind); end %defining init parameters for fit bkg = mean(maxk(data,10)); %background if normalized close to one pos = f(peakLoc); %peakLoc*(f(end)-f(1))/length(f)+f(1); amp = bkg+peakMag; width = 0.0001.*ones(1,Npeak); %compile init parameters for Lorentz fitting if Npeak>1 initPar = vertcat(amp,pos,width); %concatenate initPar = initPar(1:end); %make line vector initPar = [bkg,initPar]; %add bkg else initPar = vertcat(amp,pos,width); %concatenate initPar = initPar(1:end); %make line vector initPar = [bkg,initPar']; %add bkg end % fit vs = initPar; % fct = @(v,x) Lorentz_xN_Func(Npeak, v, x); fct = @(v,x) Gaussian_xN_Func(Npeak, v, x); [ve,r,J,COVB,mse] = nlinfit(f,data,fct,vs); %extract error bars errfitres_abs = nlparci(ve,r,'jacobian',J); %absolute value boundaries errfitres_diff = abs((errfitres_abs(:,1) - errfitres_abs(:,2))/2); %difference values --> error bar bkg = ve(1); bkg_err = errfitres_diff(1); for i=1:Npeak width(i) = 1000*ve(4+(i-1)*3); %width in MHz width_err(i) = 1000*errfitres_diff(4+(i-1)*3); %error bar in MHz pos(i) = ve(3+(i-1)*3); %peak position in GHz pos_err(i) = errfitres_diff(3+(i-1)*3); %error bar in MHz contrast(i) = 100*ve(2+(i-1)*3)/bkg; %contrast in % contrast_err(i) = 100*(errfitres_diff(2+(i-1)*3)/bkg + ve(2+(i-1)*3)/bkg^2*bkg_err); %error bar in MHz end end