4.4 NLFit and Peak Analyzer


Python functions can be used for performing nonlinear curve fitting. Peak functions defined with Python can also be used in Peak Analyzer. A fitting function file (FDF file) will need to be created which includes the Python function and script commands to install any Python packages that are needed for your Python function.

The following example outlines how to create an FDF with Python function.

Fitting Function for Mass Diffusion

This section outlines the procedure to create a fitting function for the Mass Diffusion equation:
MassDiffusionFittingFunction.png
using a Python function.

  1. First we need to create a new FDF file. Open the Fitting Function Builder tool using the Origin menu Tools: Fitting Function Builder...
  2. On the first page, select Create a New Function.
  3. On the second page, set Function Name to MassDiffusion, then select Python Function(Vector) for function type.
  4. On the third page, set independent variables to x, dependent variables to y, and set parameters to y0, D, L
  5. On the fourth page, enter
    from mpmath import nsum, exp, inf
    import numpy as np
    
    def MassDiffuse(x, y0, D, L):
    	sm = [float((nsum(lambda ii: 1/(2*ii+1)**2*exp(-D*(2*ii+1)**2*np.pi**2*t/(4*L**2)),[0, inf]))) for t in x]
    	return [y0*(1-8/np.pi**2*t) for t in sm]
    

    in the Python Function(Vector) edit box. Enter

    y = MassDiffuse(x, y0, D, L)
    

    in the Function Body edit box.

  6. On the fifth page, enter
    if(Python.chk("mpmath numpy") > 1)
    	return 1;
    return 0;
    

    in the Python Package Check(Labtalk Script) edit box. This function requires the mpmath and numpy packages. These need to be installed prior to using the function. To install the packages, you can use the Python Packages tool accessible from the Connectivity menu in Origin. When you share the function with another colleague, the packages can be automatically installed when that user adds the function to their Origin.

  7. Click Next button then click Finish button.
  8. You can test the fitting function with the following data for X and Y:
    0.00	-0.078
    0.20	 0.604
    0.40	 0.842
    0.60	 1.101
    0.80	 1.029
    1.00	 1.083
    1.20	 0.828
    1.40	 0.884
    1.60	 0.991
    1.80	 1.005
    2.00	 0.915
    

    Then use the initial parameter values:

    y0 = 1, D = 0.5, L = 0.5
    

    and perform the fitting. The final parameter results will be:

    Fitted values: y0 = 0.97846, D = 0.37563, L = 0.43781
    

Evaluate FDF in Python code

This shows how to define a fitting function with evaluate FDF in Python code.

  1. First we need to create a new FDF file. Open the Fitting Function Builder tool using the Origin menu Tools: Fitting Function Builder...
  2. On the first page, select Create a New Function.
  3. On the second page, set Function Name to HgPyEx, then select Python Function(Vector) for Function Type.
  4. On the third page, set Independent Variables to x, Dependent Variables to y. Set Parameters to:
    dx196,dx198,dx199,dx200,dx201,dx202,Vbl,Vamp,a1V,a2V,wLor,dFwm,fshift,P,T,L,nA

    set Derived Parameters to:

    ConcRef,wLorPval,fshiftPval
  5. On the fourth page in Parameters tab, set the values to coressponding parameters:
    Names = Vbl,Vamp,a1V,a2V,wLor,dFwm,fshift,P,T,L,nA
    Initial Values = 0(F),0.99629(V),0.000231(V),0(F),8.54(F),1(V),-2.54(F),1(F),1(F),89.87071(F),1.15e-05(V)
    Number Of Significant Digits = 0,0,0,0,6,7,6,0,0,0,15
    Unit = ,,,,MHz/Torr,,MHz/Torr,,,,
    

    In the Python Function(Vector) edit box, enter

    import numpy as np
    import originpro as op
    
    def Voigt(x, xc, A, alpha, gamma):
        """
        Return the Voigt line shape at x with Lorentzian component FWHM gamma
        and Gaussian component FWHM alpha.
        """
        #vp = ['y0','xc','A','wG','wL']
        va = [0, xc, A, alpha, gamma]
        x1=x.tolist()
        return np.array( op.evaluate_FDF( 'Voigt', x1, va ) )
    
    def HgOFunc( x, dx196, dx198, dx199, dx200, dx201, dx202, Vbl, Vamp, a1V, a2V, wLor, dFwm, fshift, P, T, L, nA):
    	# define transition frequencies
    	f1 = 1181.55994207E3
    	f2 = 1181.55583492E3
    	f3 = 1181.54042558E3
    	f4 = 1181.56256226E3
    	f5 = 1181.55102924E3
    	f6 = 1181.54117507E3
    	f7 = 1181.55515739E3
    	f8 = 1181.56270616E3
    	f9 = 1181.54573191E3
    	f10 = 1181.54052152E3
    	
    	f1_9 = 0
    	f2_9 = 0
    	f3_9 = 0
    	f4_9 = 0
    	f5_9 = 0
    	f6_9 = 0
    	f7_9 = 0
    	f8_9 = 0
    	f9_9 = 0
    	f10_9 = 0
    	
    	Pval=P
    	wLorP=wLor*P/1000 # units GHz = MHz/Torr*Torr/1000
    	fshiftP=fshift*P/1000 # units GHz = MHz/Torr*Torr/1000
    	
    	#x=x+f9+dFwm
    	vx = np.array( x )
    	vx = vx+f9+dFwm
    	
    	m196 = 195.96581  
    	m198 = 197.96674
    	m199 = 198.96825
    	m200 = 199.96825
    	m201 = 200.97028
    	m202 = 201.97062
    	m204 = 203.97347
    	
    	'''
    	w1 = 7.162326E-07*f1*(T/m196)^0.5   
    	w2 = 7.162326E-07*f2*(T/m198)^0.5
    	w3 = 7.162326E-07*f3*(T/m199)^0.5
    	w4 = 7.162326E-07*f4*(T/m199)^0.5
    	w5 = 7.162326E-07*f5*(T/m200)^0.5
    	w6 = 7.162326E-07*f6*(T/m201)^0.5
    	w7 = 7.162326E-07*f7*(T/m201)^0.5
    	w8 = 7.162326E-07*f8*(T/m201)^0.5
    	w9 = 7.162326E-07*f9*(T/m202)^0.5
    	w10= 7.162326E-07*f10*(T/m204)^0.5
    	'''
    	
    	w1 = 7.162326E-07*f1*(T/m196)**0.5   
    	w2 = 7.162326E-07*f2*(T/m198)**0.5
    	w3 = 7.162326E-07*f3*(T/m199)**0.5
    	w4 = 7.162326E-07*f4*(T/m199)**0.5
    	w5 = 7.162326E-07*f5*(T/m200)**0.5
    	w6 = 7.162326E-07*f6*(T/m201)**0.5
    	w7 = 7.162326E-07*f7*(T/m201)**0.5
    	w8 = 7.162326E-07*f8*(T/m201)**0.5
    	w9 = 7.162326E-07*f9*(T/m202)**0.5
    	w10= 7.162326E-07*f10*(T/m204)**0.5
    	
    	
    	g1 = 3 
    	g2 = 3 
    	g3 = 1 
    	g4 = 2 
    	g5 = 3 
    	g6 = 1.5 
    	g7 = 1 
    	g8 = 0.5
    	g9 = 3 
    	g10 = 3   
    	
    	x196 = 0.0015 + dx196  
    	x198 = 0.1004 + dx198
    	x199 = 0.1694 + dx199
    	x200 = 0.2314 + dx200
    	x201 = 0.1317 + dx201
    	x202 = 0.2974 + dx202
    	x204 = (1 - x196-x198-x199-x200-x201-x202)
    		
    	
    	'''dF1=x-f1-fshiftP;
    	dF2=x-f2-fshiftP;
    	dF3=x-f3-fshiftP;
    	dF4=x-f4-fshiftP;
    	dF5=x-f5-fshiftP;
    	dF6=x-f6-fshiftP;
    	dF7=x-f7-fshiftP;
    	dF8=x-f8-fshiftP;
    	dF9=x-f9-fshiftP;
    	dF10=x-f10-fshiftP;'''
    	
    	dF1=vx-f1-fshiftP
    	dF2=vx-f2-fshiftP
    	dF3=vx-f3-fshiftP
    	dF4=vx-f4-fshiftP
    	dF5=vx-f5-fshiftP
    	dF6=vx-f6-fshiftP
    	dF7=vx-f7-fshiftP
    	dF8=vx-f8-fshiftP
    	dF9=vx-f9-fshiftP
    	dF10=vx-f10-fshiftP
    		
    	'''A1 = (g1*x196/f1^2)*nlf_Voigt(dF1,0,f1_9,1,w1,wLorP) + (g2*x198/f2^2)*nlf_Voigt(dF2,0,f2_9,1,w2,wLorP) + (g3*x199/f3^2)*nlf_Voigt(dF3,0,f3_9,1,w3,wLorP);
    	A2 = (g4*x199/f4^2)*nlf_Voigt(dF4,0,f4_9,1,w4,wLorP) + (g5*x200/f5^2)*nlf_Voigt(dF5,0,f5_9,1,w5,wLorP) + (g6*x201/f6^2)*nlf_Voigt(dF6,0,f6_9,1,w6,wLorP);
    	A3 = (g7*x201/f7^2)*nlf_Voigt(dF7,0,f7_9,1,w7,wLorP) + (g8*x201/f8^2)*nlf_Voigt(dF8,0,f8_9,1,w8,wLorP) + (g9*x202/f9^2)*nlf_Voigt(dF9,0,f9_9,1,w9,wLorP);
    	A4 = (g10*x204/f10^2)*nlf_Voigt(dF10,0,f10_9,1,w10,wLorP);'''
    	
    	A1 = (g1*x196/f1**2)*Voigt(dF1,f1_9,1,w1,wLorP) + (g2*x198/f2**2)*Voigt(dF2,f2_9,1,w2,wLorP) + (g3*x199/f3**2)*Voigt(dF3,f3_9,1,w3,wLorP)
    	A2 = (g4*x199/f4**2)*Voigt(dF4,f4_9,1,w4,wLorP) + (g5*x200/f5**2)*Voigt(dF5,f5_9,1,w5,wLorP) + (g6*x201/f6**2)*Voigt(dF6,f6_9,1,w6,wLorP)
    	A3 = (g7*x201/f7**2)*Voigt(dF7,f7_9,1,w7,wLorP) + (g8*x201/f8**2)*Voigt(dF8,f8_9,1,w8,wLorP) + (g9*x202/f9**2)*Voigt(dF9,f9_9,1,w9,wLorP)
    	A4 = (g10*x204/f10**2)*Voigt(dF10,f10_9,1,w10,wLorP)
    	
    	# y value
    	#y = (Vbl + (Vamp + a1V*dF9 + a2V*dF9^2))*exp(-(nA*L*(2.99792458E8)^2*(1/(8*3.14159265358979))*(A1+A2+A3+A4)));
    	y = (Vbl + (Vamp + a1V*dF9 + a2V*dF9**2))*np.exp(-(nA*L*(2.99792458E8)**2*(1/(8*3.14159265358979))*(A1+A2+A3+A4)))
    	
    	return y.tolist()
    

    In the Function Body edit box, enter

    y=HgOFunc( x, dx196, dx198, dx199, dx200, dx201, dx202, Vbl, Vamp, a1V, a2V, wLor, dFwm, fshift, P, T, L, nA)
  6. On the fifth page, in the Python Package Check(Labtalk Script) edit box, enter
    if(Python.chk("numpy") > 1)
            return 1;
    return 0;
    if(Python.chk("numpy") > 1)
            return 1;
    return 0;
    
  7. On the sixth page, check Use Custom Code, and in Initialization Code box, enter
    //Code to be executed to initialize parameters
    //Get the current worksheet page.
    range rpage=ry.getpage()$;  
    //Get the data worksheet
    range rlayer=ry.getlayer()$; 
    //Get the data worksheet index
    int inext=rlayer.index;
    //Make sure the data workbook is active
    win -a %(ry.getpage()$);
    //Make sure the data worksheet is active
    page.active=inext;
     
    double Temp = mean(col(12)); 
    T=Temp; 
    double Press = mean(col(13));
    P=Press;
    Vamp=0.99629;
    a1V=2.31e-4;
    wLor=8.95;
    fshift=-2.54;
    dFwm=1;
    L=89.87071;
    nA=1.95e-7*P/100*100;
    
  8. On the eighth page in Derived Parameters tab, set the values to coressponding parameters:
    Unit = ug/m3
    Names = ConcRef,wLorPval,fshiftPval
    Meanings = Ref conc at T
    

    In Derived Parameters Equations box, enter

    ConcRef=(760/P)*((T)/293.15)*((nA/((10/1000)*8411821.73)*(1E+27/1000000))/(6.022140857E+23))*200.59*(100)^3*1000000
    wLorPval=wLor*P/1000
    fshiftPval=fshift*P/1000
    
  9. Click Finish button.