File gaussian_lorentzian_peak_generator.sce: creates a data set and saves it in the noisy_gaussian_lorentzian_peak.txt file.
chdir("mypath\");
exec("gaussian_lorentzian_functions.sci")
// *************
// * Constants *
// *************
paramlorentz(1) = 5; // height of the lorentzian curve
paramlorentz(2) = 2; // "width" of the lorentzian curve
paramgauss(1) = 10; // height of the gaussian curve
paramgauss(2) = 3; // "width" of the gaussian curve
var=0.5; // variance of the noise normal law
nbpts = 100 // number of points
halfwidth = max(3*paramgauss(2), 3*paramlorentz(2)) // for x range
step = 2*halfwidth/nbpts;
// ******************
// * Initialisation *
// ******************
X = zeros(nbpts,1);
Y = zeros(nbpts,1);
// ****************
// * Main program *
// ****************
// Data generation
i=(1:nbpts)'
X = step*i - halfwidth;
Y = profile([paramlorentz;paramgauss], X) + var*rand(X, "normal");
end
// Saving
csvWrite([X, Y], "noisy_gaussian_lorentzian_peak.txt")
File gauss_lorentz_desummation.sce: process the data to extract the lorentzian and the gaussian components.
// ******************// * Initialisation *// ******************clf;chdir("mypath\");exec("gaussian_lorentzian_functions.sci")// *************// * Constants *// *************Ainit=[1;1;1;1];// initial parameters of the model function// *************// * Functions *// *************// Residual functionfunction[e]=res(A, x, y)e=profil(A,x)-y;endfunction// ****************// * Main program *// ****************// Data readingdata=csvRead("noisy_gaussian_lorentzian_peak.txt")X=data(:,1);Y=data(:,2);// Data processing[S,Aopt]=leastsq(list(res,X,Y),Ainit)Yopt=profil(Aopt,X);YLopt=lorentz(Aopt(1:2),X);YGopt=gauss(Aopt(3:4),X);// Displayplot(X,Yopt,"-r")plot(X,YLopt,"-c")plot(X,YGopt,"-g")plot(X,Y,"+b")height=max(Y);xmin=min(X)xstring(xmin,height*0.95,..."lorentzian: Al = "+string(0.01*round(100*Aopt(1)))...+"; Hl = "+string(0.01*round(100*Aopt(2))))xstring(xmin,3*height/4,..."gaussian: Ag = "+string(0.01*round(100*Aopt(3)))...+"; Hg = "+string(0.01*round(100*Aopt(4))))xstring(xmin,height/2,..."variance: S = "+string(0.01*round(100*S)))
רישיון
אני, בעל זכויות היוצרים על עבודה זו, מפרסם בזאת את העבודה תחת הרישיון הבא:
ייחוס – יש לתת ייחוס הולם, לתת קישור לרישיון, ולציין אם נעשו שינויים. אפשר לעשות את זה בכל צורה סבירה, אבל לא בשום צורה שמשתמע ממנה שמעניק הרישיון תומך בך או בשימוש שלך.
שיתוף זהה – אם תיצרו רמיקס, תשנו, או תבנו על החומר, חובה עליכם להפיץ את התרומות שלך לפי תנאי רישיון זהה או תואם למקור.