4.2.2.13 Fitting with Integral using NAG Library

1 Summary
2 What you will learn
3 Example and Steps

Summary

Origin allows user to define an Origin C fitting function which involves an integral. You can call NAG functions to perform the integration while defining the fitting function. There are built-in functions in Origin C which perform integration. For the current example, the NAG solution is recommended. It has a better performance compared to the built-in integration algorithm. Note that an infinite NAG integrator is used here.

Minimum Origin Version Required: Origin 8.0 SR6

What you will learn

This tutorial will show you how to:

Create a fitting function using the Fitting Function Organizer
Create a fitting function with a Definite Integral using a NAG integration routine
Set up the Initial Code for the fitting function

Example and Steps

We will fit the following model:

$y=y_0+\int_{-\infty}^{x} \frac{A}{w\sqrt{\frac{\pi}{2}}} e^{-2\frac{(t-x_c)^2}{w^2}} dt$

Here $y_0$ $A$ , $xc$ and $w$ are the model parameters we want to obtain from the data fitting. The fitting procedure can be outlined into the following steps:

Define the Function

Press F9 to open the Fitting Function Organizer and then create a new Category named FittingWithIntegral. Define a new fitting function nag_integration_fitting in the new category as follow:

Function Name:	nag_integration_fitting
Function Type:	User-Defined
Independent Variables:	x
Dependent Variables:	y
Parameter Names:	y0, A, xc, w
Function Form:	Origin C
Function:

Click the button beside the Function box to open the code builder then follow the steps below:

Scroll up to go to line

//add the header file for the NAG functions here.

and paste the following code below this line.

#include <OC_nag.h>

Go to line

// and access in your fitting function.

and paste the following code below.

struct user   // parameters in the integrand
{
	double amp, center, width;
 
};
// Function supplied by user, return the value of the integrand at a given x.
static double NAG_CALL f_callback(double x, Nag_User *comm) 
{
        struct user *sp = (struct user *)(comm->p);
 
        double amp, center, width;    // temp variable to accept the parameters in the Nag_User communication struct
        amp = sp->amp;
        center = sp->center;
        width = sp->width;
 
	return amp * exp( -2*(x - center)*(x - center)/width/width ) / (width*sqrt(PI/2));
}

Go to line

// Beginning of editable part.

and paste the following code below.

// Through the absolute accuracy epsabs, relative accuracy epsrel and max_num_subint you can  
	// control the precision of the integration you need 
	// if epsrel is set negative, the absolute accuracy will be used. 
	// Similarly, you can control only relative accuracy by set the epsabs negative
	double epsabs = 0.0, epsrel = 0.0001;
	// The max number of sub-intervals needed to evaluate the function in the integral
	// The more difficult the integrand the larger max_num_subint should be
	// For most problems 200 to 500 is adequate and recommended
	Integer max_num_subint = 200;
	 
	// Result keeps the approximate integral value returned by the algorithm
	// abserr is an estimate of the error which should be an upper bound for the |I - result| 
	// where I is the integral value
	double result, abserr;
	 
	// The structure of type Nag_QuadProgress, 
	// it contains pointers allocated memory internally with max_num_subint elements
	Nag_QuadProgress qp;
	 
	// The NAG error parameter (structure)
	static NagError fail;
	 
	// Parameters passed to integrand by Nag_User communication struct
	        Nag_User comm;	
	struct user s;
	        s.amp = A;
	        s.center = xc;
	        s.width = w;
	        comm.p = (Pointer)&s;
	 
	// Perform integration
	// There are 3 kinds of infinite boundary types you can use in Nag infinite integrator
	// Nag_LowerSemiInfinite, Nag_UpperSemiInfinite, Nag_Infinite
	d01smc(f_callback, Nag_LowerSemiInfinite, x, epsabs, epsrel, max_num_subint, &result, &abserr, &qp, &comm, &fail);
	 
	        // you may want to exam the error by printing out error message, just uncomment the following lines
	// if (fail.code != NE_NOERROR)
	        // printf("%s\n", fail.message);
	 
	// For the error other than the following three errors which are due to bad input parameters 
	// or allocation failure  NE_INT_ARG_LT  NE_BAD_PARAM   NE_ALLOC_FAIL
	// You will need to free the memory allocation before calling the integration routine again to avoid memory leakage
	if (fail.code != NE_INT_ARG_LT && fail.code != NE_BAD_PARAM && fail.code != NE_ALLOC_FAIL)
	{
		NAG_FREE(qp.sub_int_beg_pts);
		NAG_FREE(qp.sub_int_end_pts);
		NAG_FREE(qp.sub_int_result);
		NAG_FREE(qp.sub_int_error);
	}
	 
	// Calculate the fitted value
	y = y0 + result;

Click Compile to compile the file.

In the above code, we firstly define the integrand as a callback function f_callback just outside the fitting function body _nlsfnag_integration_fitting. Note that we parametrize the integrand function with the variables amp, center and width, and pass them into the callback funtion through the Nag_User struct. Inside the fitting function, we perform the integration using NAG integrator d01smc.

Calling NAG functions should be more efficient than writing your own routines. Using an analogous method, you can perform finite, infinite, one-dimension and multi-dimension quadrature in your fitting function. Please read the NAG Quadrature page and select a proper routine.

Set the Initial Values for the Parameters or Setup the Initial Code

As it is a user-defined fitting function, you have to supply the initial guess values for the parameters. You may do it later by setting them manually in the Parameter tab in Nonlinear Curve Fit dialog.

Simulate the Function

After entering the function body codes, you can click the Compile button in Code Builder to check syntax errors. And then click Return to Dialog button to go back Fitting Function Organizer dialog box. Now click the Save button to generate the .FDF file (Function definition file).

Once you have a .FDF file, you can click the Simulate button to simulate a curve, this will be very helpful to evaluate the initial values. In the simcurve dialog, enter some proper parameter values and X range, and see what the curve looks like in the Preview panel.

Fit the Curve

Before you start to fit the curve, it is very helpful to simulate the function first. Performing integration may take some time, if there is any mistake, you may see Origin "freeze" after you click the Fit button. So in the Fitting Function Organizer dialog, select the function we defined and click the Simulate button. This will bring up the simcurve X-Function. Enter some "guess" values and click the Apply button. If the simulated curve looks like your source data, you can go further to fit.

To test the fitting function:

Import \Samples\Curve Fitting\Replicate Response Data.dat to Origin.
Next we are going to use log scale of the data in Col(A). To do this, in the F(x) = column label row of column A, type in the formula Col(A) = log(Col(A)) and hit Enter once to change the data.
Highlight column A and B and create a scatter plot, the shape is a sigmoid curve.
Bring up the NLFit dialog from Analysis: Fitting: Nonlinear Curve Fit menu item.
Select the fitting function we defined in the previous sections and go to the Parameters tab, initialize all parameters by 1 and fit (click Fit button).
You are supposed to see these results: