-
Notifications
You must be signed in to change notification settings - Fork 5
fit_curve_model
This module is intended for the fitting of spectral bell-shaped curves that one often encounters in experimental measurements.
Inside the Utils module are available different bell-shaped curves such as Gaussian, Lorentzian or Voigt profile to modelling different phenomena.
Here, we provide a small method to use those models and fit them in order to obtain useful information about the underlying processes. It uses the mean(abs2.(ydata .- ymodel)./σ) as the objective function, where ydata and ymodel are the experimental and modelled spectra, and σ is the standard deviation of each element of ydata.
Examples are posted in the following pages, and the full code for them can be obtained from here.
In order to use this module we can invoke as follow:
solFC = fit_curve_model(
model,
x_data,
y_data,
seed;
options=Optim.Options(),
alg=SAMIN(),
lb=-2.0.*seed,
ub=2.0.*seed,
σ::Array{T1,1}=ones(length(ydata)),
)model::Symbol tells the function which type of model to use for the fitting. There are three models so far that you can check in the Utils module: :gaussian, :lorentzian and :voigtian. Depending on the model you need to pass the right seed.
x_data::Array{Float64,1} is the x-axis data of the spectrum. Normally, here you have the wavelength, the frequency or energy, depending on the measurement.
y_data::Array{Float64,1} is the y-axis data of the spectrum. It has the same length as x_data.
Be sure that length(xdata) == length(ydata).
seed::Array{Array{Float64,1},1} are the initial guess of the parameters to optimise. The way to construct this input parameter depends on the number of peaks to fit. For instance, to fit one Gaussian peak we can have:
seed = [[0.0], # offset of the curve
[15.0, 0.0, sqrt(0.2)]] # parameters of the Gaussian peakNow, say you what to fit two Voigtian peaks, you can write:
seed = [[1.0], # offset of the curve
[10.0, -1.6, 0.0, 1.53], # parameters of the first Voigtian peak
[24.0, 5.0, 0.0, 1.53]] # parameters of the second Voigtian peakThis way of input allows determining how many peaks you want to fit. Remember that the first element of seed is the offset.
The parameters allowed in each model can be found in the Utils module.
options is an Optim.options() structure from the Optim.jl package. You can set here the options to run the algorithm, for instance, the number of iterations.
alg = SAMIN(), is the optimisation algorithm to use. By default, and it is suggested not to change it, it uses the simulated annealing with bounds. The reason not to change this is that this algorithm allows input bounds. It is an optional parameter in case you want to pass specific arguments to it. Check Optim.jl.
lb and ub are the lower and upper bounds that define the region where the algorithm will search the optimum. You need to be careful to pass bounds that include the seed, otherwise, you will be prompted an error. By default, these parameters take about lb=-2.0.*seed and ub=2.0.*seed. Both have the same type as the seed parameter, Array{Array{Float64,1},1}.
σ::Array{Float64,1} sets the standard deviation of the experimental ydata to use in the optimisation process.
solFC::NamedTuple is a structure returned from the fitting process with the following fields:
-
solution: The status info from the Optim solver -
model: Symbol with the model used in the fitting procedure -
xdata: x-axis experimental data input -
ydata: y-axis experimental data input -
seed: Initial guess input in the optimisation process -
lb: The lower bounds input -
ub: The upper bounds input -
ymodel: The final reproduction of the total line shape curve (ymodel) -
optParams: The optimal parameters found by the optimiser -
sigma: The standard deviation used in the objective function