(Demeter main page) (Cookbook main page)
| /01DataHandling /02Fitting /03AtomsAndFeff /04PlottingAndConfiguration /05HgDNA /06Uranyl |
Fitting
This section covers the parts of Demeter that expand upon the functionality of Artemis.
Fitting basics
How do I perform a simple fit?
Simple fits using Demeter resemble feffit input files in many ways. It is hard to appreciate the full capabilities of Demeter with a very simple fit, but still it demonstrates many of Demeter's capabilities.
1 #!/usr/bin/perl ### simple fit to Cu foil data
2 ## -----------------------------------------------------------
3 ## boilerplate required in every script:
4 use Demeter; ## strict and warnings automatically imported with Demeter
5 ## -----------------------------------------------------------
6 my $dobject = Demeter::Data -> new();
7 @dobject->set_mode(screen => 0, ifeffit => 1, file => ">cufit.iff");
8
9 ## set up your data
10 $dobject ->set(file => "cu10k.chi",
11 name => '10K copper data',
12 fft_kmin => 3, fft_kmax => 14,
13 fit_space => 'r',
14 fit_k1 => 1, fit_k3 => 1,
15 bft_rmin => 1.6, bft_rmax => 4.3,
16 );
17
18 ## set up several parameters
19 my @gdsobjects = (Demeter::GDS -> new(type => 'guess', name => 'alpha', mathexp => 0),
20 Demeter::GDS -> new(type => 'guess', name => 'amp', mathexp => 1),
21 Demeter::GDS -> new(type => 'guess', name => 'enot', mathexp => 0),
22 Demeter::GDS -> new(type => 'guess', name => 'theta', mathexp => 500),
23 Demeter::GDS -> new(type => 'set', name => 'temp', mathexp => 300),
24 Demeter::GDS -> new(type => 'set', name => 'sigmm', mathexp => 0.00052),
25 );
26
27 ## set up the first 5 scattering paths
28 my @pobjects = ();
29 foreach my $i (0 .. 4) {
30 my $j = $i+1;
31 $pobjects[$i] = Demeter::Path -> new();
32 $pobjects[$i]->set(data => $dobject,
33 folder => './',
34 file => "feff000$j.dat",
35 s02 => 'amp',
36 e0 => 'enot',
37 delr => 'alpha*reff',
38 sigma2 => 'debye(temp, theta) + sigmm',
39 );
40 };
41
42 ## make a Fit object, which is just a collection of GDS, Data, and Path objects
43 my $fitobject = Demeter::Fit -> new(gds => \@gdsobjects,
44 data => [$dobject],
45 paths => \@pobjects,
46 );
47
48 ## do the fit (or the sum of paths)
49 $fitobject -> fit;
50 ##$fitobject -> sum;
51
52 ## Make a plot
53 $fitobject->po->set(plot_fit => 1, plot_win => 0, plot_res => 0,
54 kweight => 2, r_pl =>'m','q_pl' =>'r');
55 my $s = 0; # stack the plot interestingly...
56 foreach my $obj ($dobject, @pobjects,) {
57 $obj -> plot('R');
58 $s -= 0.8;
59 $dobject -> set('y_offset'=>$s);
60 };
61
62 ## save the results of the fit
63 $dobject->save("fit", "cufit.fit");
64
65 ## write a log file
66 $fitobject -> logfile("cufit.log", "Fit to copper data", q{});
At line 7, we set up Demeter's disposal channels. In this case, commands will be sent to Ifeffit for processing and written to a file for saving, but the commands will not be echoed to the screen as you run the script.
The data set is defined at lines 10-16 in a manner we are already familiar with.
At lines 19-25, the various parameters of the fit are defined. This syntax is, admittedly, a bit bulky. The alternative syntax is
1 $gdsobject = $any_object->simpleGDS("guess alpha = 0");
The paths and path parameters are set up in lines 28-40. Again, this bears a strong resemblance to feffit. Note that perl is used to loop over 5 very similar paths, using a counter $j to change the name of the feffNNNN.dat file.
Finally, the fit is defined in lines 43-46. A Fit is simply a collection of data sets, parameters, and paths. There is enough information in the definitions of the parts of the fit to completely define how the fit should be translated into ifeffit commands. Each argument of the fit definition is a list reference.
Until line 40, Demeter is not significantly less wordy than a feffit input file or an ifeffit script, although you get to use high-level language features to automate some parts of the problem. There is not much that can be done about the wordiness -- you have to define all the data processing and path parameters in some manner. Demeter begins to shine after line 40. Once everything is set up, other chores can be performed with considerable economy of language.
Once armed with a Fit object, many things are easy. The fit is performed at line 49. The fitting results are saved as chi(k) to a file at line 63. A log file is written at line 66. The results of the fit are plotted in an interesting way at lines 54-61.
How do I perform a multiple data set fit?
Multiple data set fits are very natural in Demeter. You simply define more data, paths, and parameters then group them all together as a Fit object. Here is an example for measurements of a copper foil at two temperatures.
1 #!/usr/bin/perl ### multiple data set fit to Cu foil data at two temperatures
2 ## -----------------------------------------------------------
3 ## boilerplate required in every script:
4 use Demeter; ## strict and warnings automatically imported with Demeter
5 ## -----------------------------------------------------------
6 my $demeter = Demeter->new; # a utility object
7 $demeter->set_mode(screen => 0, ifeffit => 1, file => ">mdsfit.iff",);
8
9 my @common = (fft_kmin => 3, fft_kmax => 14,
10 bft_rmax => 1.0, bft_rmax => 4.3,
11 fit_k1 => 1, fit_k3 => 1,);
12
13 my $data_010k = Demeter::Data -> new();
14 $data_010k -> set(@common);
15 $data_010k -> set(file => "cu10k.chi",
16 name => '10 K copper data',
17 );
18 my $data_150k = Demeter::Data -> new();
19 $data_150k -> set(@common);
20 $data_150k -> set(file => "cu150k.chi",
21 name => '150 K copper data',
22 );
23
24 ###--- make GDS objects for an isotropic expansion, correlated Debye model
25 my @gdsobjects =
26 (Demeter::GDS -> new(type => 'guess', name => 'alpha010', mathexp => 0,),
27 Demeter::GDS -> new(type => 'guess', name => 'alpha150', mathexp => 0,),
28 ## here is some syntactic sugar
29 $demeter -> simpleGDS("guess amp = 1"),
30 Demeter::GDS -> new(type => 'guess', name => 'enot', mathexp => 0,),
31 Demeter::GDS -> new(type => 'guess', name => 'theta', mathexp => 500,),
32 Demeter::GDS -> new(type => 'set', name => 'sigmm', mathexp => 0.0005,),
33 );
34
35 ## make Path objects for the first 5 paths in copper (3 shell fit)
36 my @paths_010k = ();
37 foreach my $i (0 .. 4) {
38 my $j = $i+1;
39 $paths_010k[$i] = Demeter::Path -> new();
40 $paths_010k[$i]->set(data => $data_010k,
41 folder => './',
42 file => "feff000$j.dat",
43 name => "10K, path $j",
44 s02 => 'amp',
45 e0 => 'enot',
46 delr => 'alpha010*reff',
47 sigma2 => 'debye(10, theta) + sigmm',
48 );
49 };
50
51 ## clone these 5 paths for the second data set
52 my @paths_150k = ();
53 foreach my $i (0 .. 4) {
54 my $j = $i+1;
55 $paths_150k[$i] = $paths_010k[$i] -> clone(data => $data_150k,
56 name => "150K, path $j",
57 delr => 'alpha150*reff',
58 sigma2 => 'debye(150, theta) + sigmm',
59 );
60 };
61
62 ## make a Fit object
63 my $fitobject = Demeter::Fit -> new(gds => \@gdsobjects,
64 data => [$data_010k, $data_150k],
65 paths => [@paths_010k, @paths_150k],
66 );
67 ## do the fit
68 $fitobject -> fit;
69
70 ## plot the results of the fit
71 $data_010k -> po -> set(plot_fit => 1);
72 $data_010k -> set(y_offset => 0.5) -> plot('r');
73 $data_150k -> plot('r');
At lines 18-22 a second data set is defined. Then at line 27 a new parameter is defined for the second data set. Finally, the paths used for the first data set are cloned and altered for the second data set at lines 52-60. From there, everything proceeds the same as for any other fit. The Fit is a collection of data, paths, and parameters. The paths are associated with particular data sets at lines 40 and 55.
Note that two of the tricks mentioned in earlier recipes are used to streamline the definition of all these objects. An array of common data parameters is defined at lines 9-11 and used at lines 14 and 19. Then the object cloning trick is used at lines 57-61 to define the paths for the second data set.
Model building features
What is the characteristic value of a data set and how do I use it in a fit?
One awkward aspect of the previous example is the need to define the sigma2 parameter differently for the two data sets. There is an obvious relationship between the what the data is an how the sigma2 parameter is defined. Demeter allows you to define a characteristic value for every data set. This value can be anything of interest to you as you create your fiting model. In this case, the characteristic value is the temperature at which the data were measured.
First, add a line to the data definition by editing lines 15-24 in the recipe above to read:
1 my $data_010k = Demeter::Data -> new();
2 $data_010k -> set(@common);
3 $data_010k -> set(file => "cu10k.chi",
4 name => '10 K copper data',
5 cv => 10,
6 );
7 my $data_150k = Demeter::Data -> new();
8 $data_150k -> set(@common);
9 $data_150k -> set(file => "cu150k.chi",
10 name => '150 K copper data',
11 cv => 150,
12 );
Then modify the definitions of the paths for the first data set at lines 38-51:
1 my @paths_010k = ();
2 foreach my $i (0 .. 4) {
3 my $j = $i+1;
4 $paths_010k[$i] = Demeter::Path -> new();
5 $paths_010k[$i]->set(data => $data_010k,
6 folder => './',
7 file => "feff000$j.dat",
8 name => "10K, path $j",
9 s02 => 'amp',
10 e0 => 'enot',
11 delr => 'alpha010*reff',
12 sigma2 => 'debye([cv], theta) + sigmm',
13 );
14 };
Finally, you can exclude line 60 from the definition of the paths for the second data set:
1 my @paths_150k = ();
2 foreach my $i (0 .. 4) {
3 my $j = $i+1;
4 $paths_150k[$i] = $paths_010k[$i] -> clone(data => $data_150k,
5 name => "150K, path $j",
6 delr => 'alpha150*reff',
7 );
8 };
As Demeter sets up the fit and begins generating Ifeffit commands, the string [cv] is replaced by the appropriate characteristic value. A value of 10 will be used for the paths associated with the first data set and of 150 for the second data set.
For those who know the old feffit program, the cv is much like feffit's local parameter.
What is a local guess and how do I use it in a fit?
The second awkward aspect of the multiple data set is the need to define two very similar guess parameters at lines 28 and 29. The alpha_010 and alpha_150 parameters serve the same purpose in the math expressions for the delr path parameter, while allowing that parameter to float freely for each data set.
To address this, Demeter offers the lguess parameter type. First, you should define characteristic values for each data set as in the previous recipe. Then, replace lines 28 and 29 with this one line:
1 Demeter::GDS -> new(type => 'lguess', name => 'alpha', mathexp => 0,),
then alter line 48 like so:
1 delr => 'alpha*reff',
and simply remove line 59 altogether so that the lines defining the paths for the second data set read like this:
1 $paths_150k[$i] = $paths_010k[$i] -> clone(data => $data_150k,
2 name => "150K, path $j",
3 sigma2 => 'debye([cv], theta) + sigmm',
4 );
When Demeter sets up the fit, an actual guess parameter will be automatically generated for each data set. In this case, these automatic guess parameters will be called alpha_10 and alpha_150. (If you do not define a characteristic value, a unique four-character string will be used as the suffix of the automatic guess parameter.) The Path parameters using the lguess parameter will be rewritten to use the proper automatic guess parameter and the fit will proceed accordingly.
You can use lguess parameters in any math expression for any parameter or path parameter. Any def parameters which depend on an lguess parameter will be expanded into automatic def parameters in the same manner.
Bringing together lguess parameters and characteristic values
Bringing together lguess values and characteristic values makes Demeter very powerful and expressive. Here is the two-temperature copper script from the Data Handling recipe page rewritten to use characteristic values and an lguess parameter:
1 #!/usr/bin/perl ### multiple data set fit to Cu foil data at two temperatures
2 ## -----------------------------------------------------------
3 ## boilerplate required in every script:
4 use Demeter; ## strict and warnings automatically imported with Demeter
5 ## -----------------------------------------------------------
6 my $demeter = Demeter -> new();
7 $demeter->set_mode({screen => 0, ifeffit => 1, file => ">mdsfit.iff",});
8
9 my @common = (fft_kmin => 3, fft_kmax => 14,
10 bft_rmax => 1.0, bft_rmax => 4.3,
11 fit_k1 => 1, fit_k3 => 1,);
12
13 my $data_010k = Demeter::Data -> new();
14 $data_010k -> set(@common);
15 $data_010k -> set(file => "cu10k.chi",
16 name => '10 K copper data',
17 cv => 10,
18 );
19 my $data_150k = Demeter::Data -> new();
20 $data_150k -> set(@common);
21 $data_150k -> set(file => "cu150k.chi",
22 name => '150 K copper data',
23 cv => 150,
24 );
25
26 ###--- make GDS objects for an isotropic expansion, correlated Debye model
27 my @gdsobjects =
28 (Demeter::GDS -> new(type => 'lguess', name => 'alpha', mathexp => 0,),
29 $demeter -> simpleGDS("guess amp = 1"), ## here is some syntactic sugar provided by the Tools module
30 Demeter::GDS -> new(type => 'guess', name => 'enot', mathexp => 0,),
31 Demeter::GDS -> new(type => 'guess', name => 'theta', mathexp => 500,),
32 Demeter::GDS -> new(type => 'set', name => 'sigmm', mathexp => 0.0005,),
33 );
34
35 ## make Path objects for the first 5 paths in copper (3 shell fit)
36 my @paths_010k = ();
37 foreach my $i (0 .. 4) {
38 my $j = $i+1;
39 $paths_010k[$i] = Demeter::Path -> new();
40 $paths_010k[$i]->set(data => $data_010k,
41 folder => './',
42 file => "feff000$j.dat",
43 name => "10K, path $j",
44 s02 => 'amp',
45 e0 => 'enot',
46 delr => 'alpha*reff',
47 sigma2 => 'debye([cv], theta) + sigmm',
48 );
49 };
50
51 ## clone these 5 paths for the second data set
52 my @paths_150k = ();
53 foreach my $i (0 .. 4) {
54 my $j = $i+1;
55 $paths_150k[$i] = $paths_010k[$i] -> clone(data => $data_150k,
56 name => "150K, path $j",
57 );
58 };
59
60 ## make a Fit object
61 my $fitobject = Demeter::Fit -> new(gds => \@gdsobjects,
62 data => [$data_010k, $data_150k],
63 paths => [@paths_010k, @paths_150k],
64 );
65 ## do the fit
66 $fitobject -> fit;
67
68 ## plot the fit results
69 $data_010k -> po -> set(plot_fit => 1);
70 $data_010k -> y_offset(0.5);
71 $_ -> plot('r') foreach ($data_010k, $data_150k);
How do I use more than one Feff calculation in a fit?
This is quite simple. In every example we have seen so far, each feffNNNN.dat file came from a single Feff calculation. But there is no explicit definition of a Feff calculation in any of the examples we have seen so far. (But see later recipes for how to use Demeter to control your Feff caluclations.) Each Path object has folder and file attributes. The folder attribute points to a location on disk, which the file attribute specifies which feffNNNN.dat file to use. If you have different Feff calculations sitting in different folders on disk, simple point your Data objects at the calculation of your choice. From there, everything is the same as on previous examples.
See the Ag-Au example in examples/ag-au in the Demeter distribution.
Here is a script to fit the AgBrCl data that is distributed with Horae as an Artemis example:
![/!\](/moin_static185/ifeffit/img/alert.png "/!\"){kind=small}
Example not written!
1 #!/usr/bin/perl ### AgBrCl fit with two feff calculations
2 use Demeter; ## strict and warnings imported automatically with Demeter
3
4 ## -----------------------------------------------------------
5
6 ## ...
Working with Fit objects
How do I save a fit?
After running a fit, you can save it like so:
1 $fit -> fit;
2 $fit -> freeze("myfit.dpj");
This write a normal zip file (but with the .dpj extension) containing [http://www.yaml.org YAML] serializations of the Fit object and the results of the fit. You can import a fit result like so:
1 my $fit = Demeter::Fit->thaw("myfit.dpj");
2 $fit -> po -> set(plot_fit=>1);
3 my @data = @{ $fit -> data };
4 foreach my $d (@data) {
5 $d -> plot('R');
6 };
These lines import the fit and plot the data and fit for each data set in the fit.
How do I compare fits?
Once you have thawed a fit, you can access the results of the fit or generate a log file in the normal manner. To prepare a report on multiple fits, simply thaw each of them and use the reporting methods of the Fit object. For instance:
1 my @fits = (Demeter::Fit->thaw("fit1.dpj"),
2 Demeter::Fit->thaw("fit2.dpj"),
3 Demeter::Fit->thaw("fit3.dpj"),
4 );
5 print " fit R-factor\n";
6 print "----------------------\n";
7 print " 1 ", $fits[0]->statistic("rfactor"), "\n";
8 print " 2 ", $fits[1]->statistic("rfactor"), "\n";
9 print " 3 ", $fits[2]->statistic("rfactor"), "\n";