Wrapper of the cMPFIT library for Nim.
For the (mostly autogenerated) documentation, see:
https://vindaar.github.io/nim-mpfit/
The following example shows the basic usage of the library. Note that
the actual code related to mpfit-nim is only the definition of proc
expH(...) and proc fitHalfLife(...).
In addition an echoResult proc is defined, which can be used to
pretty print the result of the fit. Alternatively, you can call
pretty taking the final parameters and the mp_result (i.e. the
tuple returned by fit as separate arguments) to get a string
representation of it.
Note: we use the mpfit/plotting submodule here to generate a plot in
one line. It depends on ggplotnim.
import std / [strutils, sequtils, strformat]
import pkg / [zero_functional, seqmath]
import mpfit
const
filename = "data/half_life_muon.txt"
func expH(p: seq[float], x: float): float =
## the function we'd like to fit. Any user defined function needs to
## be of the signature
## `proc[T](p: seq[T], x: T): T`
## i.e. conform to the `FuncProto[T]` type defined in mpfit_nim
result = p[0] * exp(-p[1] * x)
proc parseHalfLifeData(filename: string): (seq[float], seq[float]) =
## Parse the input file. First create seq of tuple of floats
## then convert that to tuple of seq[float]
let s = readFile(filename).splitLines --> filter('#' notin it and it.len > 0).
map(it.splitWhitespace).
map((it[0].parseFloat,
it[1].parseFloat))
result[0] = s --> map(it[0])
result[1] = s --> map(it[1])
proc fitHalfLife(bins, counts, countsErr: seq[float]): (seq[float], mp_result) =
## the actual code which performs the fitting. Call the `fit` proc
## with the user defined function to be fitted as the first argument,
## the initial parameter guess as the second and finally x, y and y_err
# start parameters
let p = [1400.0, 1.0]
# now just call fit
let (pRes, res) = fit(expH, p, bins, counts, countsErr)
echoResult(pRes, res = res)
result = (pRes, res)
echo &"The lifetime of the muon is ~ {1.0 / pRes[1]:.2f} µs"
when isMainModule:
# first parse the data from the file
let (bins, counts) = parseHalfLifeData(filename)
# calculates errors: poisson errors on the counts
let countsErr = counts.mapIt(sqrt(it))
# perform the fit and echo results
let (pRes, res) = fitHalfLife(bins, counts, countsErr)
# plot the data and the fit
import mpfit / plotting # import plotting convenience function
plot(
expH, pRes, # the function we fit and resulting fit params
bins, counts, countsErr, # the input data & errors
res, # the `mp_result` returned from the `fit` call
xMin = 0.0, xMax = 10.0, # customize range of plot
xlabel = "time / μs", ylabel = "# counts", title = "Muon half life measurement", # and labels
outfile = "../media/muon_lifetime_measurement.png", # save as png
verbose = false # we set `verbose` to false, as we already `echoResult` manually
)which outputs the following:
χ² = 74.0799 (8 DOF)
χ²/dof = 9.25999
NPAR = 2
NFREE = 2
NPEGGED = 0
NITER = 11
NFEV = 33
P[0] = 1937.1 /- 45.6327
P[1] = 0.515508 /- 0.00876106
The lifetime of the muon is ~ 1.94 µsand creates this plot:
The library depends on the cMPFIT library as a shared object. Either get the source code from here or use the code in c_src/. Compile the C library as follows:
gcc -c -Wall -Werror -fpic mpfit.c mpfit.h
gcc -shared -o libmpfit.so mpfit.owhich should create a libmpfit.so file in the same directory. The
Nim library will link against it. Either copy the shared library to
the location of your Nim code, in which you use mpfit-nim, or install
it system wide, depending on your system it may look like the
following (Ubuntu x64):
cp libmpfit.so /usr/lib/x86_64-linux-gnu
Once the shared library is available, there shouldn’t be anything else to do to use the library. Note: the example given in examples/fit_half_life.nim requires a few Nim libraries, which are not dependencies, since they are only used in the example, notably:
seqmath(for linspace)plotly(chroma) (to plot the data and fit)zero_functional(to parse the data)
The library consists of a single exported fit procedure, which has
the following signature:
proc fit*[T](f: FuncProto[T], pS: openArray[T], x, y, ey: openArray[T]): (seq[T], mp_result) =the first argument is a user defined function (see below), the following arguments are:
pS: the first guess for the parametersx: data for xy: data for yey: errors for y
Note: currently the ey may not be an empty sequence, nor 0, since we
use it as a weight. (TODO: change that!)
The mp_result object contains the chi^2 values for the fit, the
errors on the parameters and additional information about the internal
fitting process (e.g. number of times the user defined function was
called).
The type is defined in src/wrapper/mpfit_wrapper.nim.
The FuncProto[T] type is the following:
proc [T](p: seq[T], x: T): Tdefined in src/mpfit_nim.nim. The user defined function needs to conform to that (see the example above).
The C code is governed by the licence as shown in c_src/DISCLAIMER. The Nim code is published under the MIT license.