Peak Fitting in Demeter
I have recently added an experimental new feature to Demeter that uses Fityk as a backend for peak fitting.
Fitting arbitrary lineshapes to data is something that can be done with Ifeffit's minimize() function. Indeed, it is how the existing peak fitting utility in Athena works. That, however, remains the weakest part of Athena. It has always been somewhat unstable and often fails to find decent fits even given sensible starting positions. Is it Ifeffit's fault? Is it Athena's fault? I have looked into it many times and never come to a clear conclusion, so I thought I would pursue another avenue.
Fityk is a robust and powerful program specifically geared toward the problem of fitting arbitrary lineshapses to arbitrary data. It has a nice GUI as well as a command line interface for the power-user. The author also provides SWIG support which can be used to make a perl wrapper.
Fityk has a large set of built-in line shapes, all of which are available to the Demeter script. Demeter adds two "user defined functions" -- step-like functions based on an arctangent and an error function -- specifically for use with XANES data.
Example fit
I have incorporated this wrapper in Demeter using the PeakFit and LineShape objects. Here I show an example of a fit to some Ti K-edge data using this new interface. This picture is a screenshot taken after running a simple fit:
Example script
The script that made this fit is short and simple:
1 #!/usr/bin/perl
2 use Demeter qw(:ui=screen :plotwith=gnuplot);
3
4 ### import some XANES data and process it
5 my $data = Demeter::Data->new(file => 'tipb.305',
6 energy => '$1',
7 numerator => '$2',
8 denominator => '$3',
9 ln => 1,
10 );
11
12 ### define a peak fitting object and assign the XANES data to it
13 my $peak = Demeter::PeakFit->new(xmin=>-15, xmax=>5, screen => 1);
14 $peak -> data($data);
15
16 ### add some lineshapes to the fit
17 $peak -> add('atan', center=>4976.5, name=>'arctangent');
18 $peak -> add('gaussian', center=>4969.8, name=>'Peak1');
19 $peak -> add('lorentzian', center=>4966.2, name=>'Peak2');
20
21 ### do the fit and write the results to the screen
22 $peak -> fit;
23 print $peak -> report;
24
25 ### make a pretty plot
26 $data -> po -> set(e_norm=>1, emin=>-20, emax=>30, plot_res=>0);
27 $_ -> plot('e') foreach ($data, $peak, @{$peak->lineshapes});
28 $peak -> pause;