Spectroscopy -- First Steps

 

                                   

 

7. Spectral Fitting in SPEX.

We must now go into the spectral fitting section of SPEX. Do this with the following command:

SPEX> /fit
 

You don't actually have to type this command to execute the following commands and do the spectral fits but it has the advantage of allowing you to type "list" to see what values the various parameters have and what commands are available in this new section of SPEX.

The first thing to do is define the type of function you want to use for the fit. The default is a single thermal plus a broken power-law function, f_vth_bpow. This can be considered as the workhorse model because you can control the 6 free parameters in such a way that you can make it work as a single thermal, a single or a double power-law, or any combination of these functions. For this demonstration we will use it as a thermal plus a single power law by fixing the parameters of the higher energy power law so that it does not contribute significantly to the total flux in the energy range of interest below 100 keV.

The full list of available functions can be seen by looking at the procedure names in the Solarsoft distribution under packages/xray. You can read the documentation for each function by typing xdoc,' followed by the name of the function, e.g.

SPEX> xdoc,'f_vth_bpow
 

You can get more information than you probably need at this point on the available functions by typing

SPEX> xdoc,'model_components.pro
 

The workhorse function f_vth_bpow returns the differential photon spectrum seen at the Earth for a two component model comprised of an optically-thin thermal-bremsstrahlung function normalized to the Sun-Earth distance plus a non-thermal broken power-law function normalized to the incident flux on the detector.

F_VTH_BPOW Parameters
a(0) emission measure in units of 1049 cm-3
a(1) kT plasma temperature in keV
For the "non-thermal" function it uses the broken power law function, F_BPOW.PRO:
a(2) normalization of broken power-law at Epivot (usually 50 keV)
a(3) negative power law index below break
a(4) break energy
a(5) negative power law index above break
a(6) low energy cutoff for power-law
a(7) spectral index below cutoff, between -1.5 and -2.0
default parameters - a(6) and a(7) default to 10 keV and -1.5 can be overwritten by values in common function_com

To see the spectra in count space before the fitting begins, enter

SPEX> count
 

and the neglected window comes to life.


 


 

The yellow curve is the flare count spectrum and the orange curve is the background. The light yellow line represents the input photon power-law spectrum, and the green line is the thermal component. The blue line is the predicted count rate response of the detectors to that input photon spectrum. When the fit is complete, the blue curves should match the yellow, more or less. But for now, we have arbitrary values for the function parameters.

It is important to choose good starting parameters for the chosen function. For the interval we have used in this demonstration, good starting parameters are as follows:

SPEX> apar, 3, 1.3, 0.05, 4, 400, 2
 

Note: Do not try to change these parameters with this command after doing the least-squares fit with the "fit" command. Instead, use the "photon" command as described below.

How would you choose these parameters yourself? One way is to start by choosing the first power-law parameters so that the function equals the measured photon spectrum at 50 keV and has about the right slope. Choose the other parameters such that the other components don't contribute significantly to the total function in the energy range of interest. Fix these other parameters with the "free" command introduced below and then do the fit. Once you have good power-law parameters, try estimating the thermal parameters you need to improve the fit at low energies. If you need a different power-law slope at higher energies, try guessing the required parameters for the second power-law function and freeing them prior to repeating the "fit" command. Finally, let all of the parameters of all the required components go free and obtain your final best fit.

The array "free" allows parameters to be fixed if the corresponding number is set to 0 or to vary if set to 1.We want to free the first four parameters and fix the last two as follows;

SPEX> free, 1, 1, 1, 1, 0, 0
 

There are three ways now to select the energy bins to be included in any spectral fits in SPEX.

  1. Use the two element parameter ERANGE. They are in the units of the energy bins (nominally keV).
  2. Use the ENRANGE command to make a graphic selection of 1 or multiple ranges of energy bins.
  3. Use the ENRANGE command with arguments, e.g. enrange, 10., 40., 60, 200. would select channels from 10-40 keV and 60-200 keV

For this demonstration we will keep it simple and use the ERANGE command to set the energy interval used in the fit to be 10 - 100 keV. This is done by entering the following:

SPEX> erange,10.,100.
 

You can set this manually with the "enrange" command. This brings up a spectral window and allows you to click on the start and end of each energy range that you want to use. But to keep it simple, we will stick with 10 - 100 keV.

The f_vth_bpow function has a cutoff of the power laws such that below 10 keV, the power-law index is set to -1.5. This can be changed to 1 keV and -1.5 with the following command:

SPEX> a_cutoff = [1.0,1.5]
 

where the first number is the cutoff energy and the second number is the power-law spectral index below that energy. Note the equals sign that allows SPEX to interpret this as an IDL command.

The SPEX command "fit" obtains best-fit spectra for all intervals between "ifirst" and "ilast." In our case since we have only selected one interval, both of these parameters are zero.

In preparation for showing the fitted spectrum and residuals, enlarge the Spectral Window to the full height of your screen. Before typing "fit" to obtain the best-fit spectra, it's a good idea to first type "list" and review carefully all of the settings and parameters. Otherwise, SPEX has a tendency to go into an unending loop if you have some parameter wrong, and the only recourse is to kill IDL and start SPEX all over again.

When you are ready, type "fit" to carry out the least-squares fit and determine the best-fit parameters that give the smallest value of chi-squared. After some calculations the following window appears:


 


 

You can see that we are getting a better fit between the blue curve (the predicted count spectrum) and the yellow curve (the measured count spectrum).

Then type "photon" to see the results of the fit. The following window will appear:


 


 

Click on the CHANGE PARAMETERS button in the resulting pop-up window to see the parameters used in the latest fit.


 


 

You can also change them for the plot and as starting parameters for a new fit. Clicking on the DONE button brings up a plot of the photon spectrum showing the data points and the different components of the fitting function using the newly changed parameters.

The best fit spectrum that I was able to obtain for our 1-minute time interval at the peak of the flare is shown below:


 


 

The plot shows the data points in yellow, the components of the fitting function in green and the summed fitting function in red. The best fit parameters are shown at the bottom. The thermal component has an emission measure of 2.6 x 1049 cm2 and a temperature of 1.3 keV; the power-law component has a slope of -5 and a flux at 50 keV of 0.04 photons s-1 cm-2 keV-1. The higher energy power-law component has fixed parameters chosen to ensure that it does not contribute significantly to the total function in the energy range of interest.

Note the points in the photon spectrum plot at energies above ~70 keV that lie above the best-fit function. The cause of this discrepancy is currently unknown.

Another way to display the spectrum is to type "count" to see the figure shown below.


 


 

In this plot, the measured count-rate spectrum is shown as counts s-1 cm-2 keV-1 with the background added back in. Also shown are the fitted function and its components converted to count rates with background included. The difference spectrum between the measured and calculated rates is also shown.

If more than one time interval is analyzed, the power-law index can be plotted vs. time as follows:

SPEX> idl, base=getutbase()
SPEX> idl,wset,0
SPEX> idl,utplot,avg(xselect,0),apar_arr[3,*], base
 

Note that "xselect" has the start and end times of the intervals that have been analyzed. For different parameters of the fits, choose the appropriate numbers for the first dimension of apar_arr.

 

Previous | Start | Next

 

 

Responsible NASA Official:
Brian Dennis
Web Design:
Merrick Berg

Solar Physics Laboratory, Goddard Space Flight Center

Space Science Laboratory, University of California Berkeley
 
Responsible Berkeley Official:
Hugh Hudson
Systems Admin:
Jon Loran

This page last updated: June 27, 2011