NLLSF example

From MikePeel.net

Jump to: navigation, search

Below is example code. See GSL#Nonlinear Least Squares Fitting for the explanation.

int fit_f(const gsl_vector * parameters, void * data, gsl_vector * f)
{
	unsigned int i;
	
	/*
	 Get pointers to the actual data.
	*/
	unsigned int n = ((struct data_def *)data)->n;
	double *x = ((struct data_def *)data)->x;
	double *delta_T = ((struct data_def *)data)->dT;
	double *sigma_T = ((struct data_def *)data)->sigma_T;
	
	/*
	 Get values for the parameters
	*/
	double y = gsl_vector_get(parameters, 0);
	double ps_t1 = gsl_vector_get(parameters, 1);
	double ps_t2 = gsl_vector_get(parameters, 2);
	double calc_delta_T;
	for (i = 0; i < n; i++)
	{
		/*
		 Calculate the expected delta_T for these parameters
		 */
		calc_delta_T = get_Delta_T(x[i], y, ps_t1, ps_t2);
		
		/*
		 Compare the calculated value for delta_T with the measured delta_T, and output the difference to vector f.
		 Enforce keeping the variables positive.
		*/
		if (y < 0)
		{
			gsl_vector_set(parameters, 0, 0.0);
			y = 0;
		}
		if (ps_t1 < 0)
		{
			gsl_vector_set(parameters, 1, 0.0);
			ps_t1 = 0;
		}
		if (ps_t2 < 0)
		{
			gsl_vector_set(parameters, 2, 0.0);
			ps_t2 = 0;
		}

		gsl_vector_set(f, i, (calc_delta_T - delta_T[i]) / sigma_T[i]);
	};
	return GSL_SUCCESS;
};