Skip to content
Documentation
Reference
Curves & Fitting

Curves & Fitting

A curve file holds a fit model and its results. For the full list of c* commands see c — Fitting. It is created with cc (attached to a data file) or cca (standalone, no data reference). The curve is defined by an expression in the dummy variable t.

Introductory example 

 > fm 101 3 testdata
0 > =ox (i-50)/10
0 > =oy exp(-(x-j)^2/2/(1.4^j)^2)*ran(.8,1.2+0*i)
0 > cc p0*exp(-(t-p1)^2/2/p2^2)
1 > op0 1
1 > op1 1
1 > op2 1
1 > cf
1 > 0,1 p 1
0,1 > cp

Line by line:

  1. Create a test data set (101 points, 3 spectra, named testdata).
  2. Set x-scale from −5 to 5 in steps of 0.1.
  3. Set y to a noisy Gaussian whose width grows with spectrum index j.
  4. Create a Gaussian fit curve with three parameters. 5–7. Set initial parameter values (the default is 1 anyway).
  5. Fit — Levenberg-Marquardt least squares.
  6. Plot data (file 0) and fit (file 1) together.
  7. Print the fitted parameters.

The dummy variable t

Inside the curve expression, t is the running variable (the independent variable). It plays the role that x plays for data. All other expression syntax — file constants, z coordinates, parameters — works exactly as described in Expressions.

Convolution

Two special expression forms enable fitting with instrumental resolution:

conv(theory)

Convolutes theory with the resolution function set by cv. The resolution file must be set beforehand:

2 > cv 0      # use file 0 as the resolution
2 > cf        # fit with convolution

An optional shift argument shifts the convolution centre: conv(theory, shift).

resol(shift)

Evaluates to a copy of the resolution function (a Dirac delta peak in the limit of perfect resolution). Used to model an elastic line:

cc p0 + p1*resol(p2) + p3*conv(kwwc(t,p4,p5), p2)

Here:

  • p0 is a flat background.
  • p1*resol(p2) is a delta peak (elastic line) shifted by p2.
  • p3*conv(...) is the quasi-elastic line, shifted by p2.

Parameter management

Command Effect
op0 1.5 Set parameter p0 to 1.5
op0 pi/z0 Set p0 using an expression
cx 2 Fix p2 (exclude from fit)
cx 2+ Fix p2 and all parameters with higher indices
cu 2 Unfix p2
cg 1 Make p1 global (same value for all spectra in the file)
cp Print all parameters

Enforcing positive parameters

If a fit converges to a physically meaningless negative value, reset that parameter to its absolute value and refit:

11-13 > cf
         ...  (p1 is negative)
11-13 > op1 abs(p1)
11-13 > cf

Fit settings

Command Effect
cd <k> Set data file to fit (expression for k allowed)
cv <k> Set resolution/convolution file
cv- Disable convolution
cr <expr> Restrict fit to points where <expr> is true
cwc Constant weights
cwv Weights = reciprocal variance (default)
cf Run fit
cfs Run fit, allow slow (N²) convolution
cst 30 Set fit timeout to 30 seconds
cse 1e-13 Set fit accuracy (sqrt of differentiation step)

Extracting parameters

After fitting, use ci to extract a parameter into a new workspace:

13 > ci 1     # extract p1 vs z0 into a new workspace
13 > ci *     # extract all parameters

The result is a data file where y = parameter value and x = the corresponding z0 (e.g. scattering angle). Alternatively, oi p0 does the same for p0.

Slow convolution error

If you see:

not equidistant and slow convolution not allowed
exception in fit function

the energy grid is not equidistant (typically caused by mpa with a factor that produces a non-integer result). Three solutions:

  1. Use mpaf <n> instead of mpa — bins by a fixed integer factor, always produces equidistant channels, and enables fast convolution.
  2. Use cfs instead of cf — allows slow (O(N²)) convolution at the cost of speed.
  3. Reduce data to the elastic region first (mr abs(x)<6), then the number of points is small enough that cfs is fast.