**RHESSI Spectroscopy \-- Second Steps**

 

 

Last modified on 25 February 2003

**List of topics\*:**

1.  [How to analyze a flare with different attenuator states](#topic_1)

```{=html}
<!-- -->
```
1.  [How to select multiple analysis intervals](#topic_2)

2.  [How to change fitting models  ](#topic_3)

3.  [How to save and restore the fitting parameters and spectral
    data](#topic_4)

4.  [How to plot the fitted spectra outside SPEX](#topic_5)

5.  [How to make a movie of the spectra](#topic_6)

6.  [How to choose and evaluate a background model](#topic_7)

7.  [How to create and customize postscript plots within
    SPEX](#Method_1)

\*: Throughout this tutorial, the X4.8 flare on July 23 2003 (1.-6. what
does this mean?) and the M1.5 flare on June 1 2002 (7.-8. what does this
mean?) are used as examples to explain the topics. Words in italics
represent SPEX commands

**1. **[]{#topic_1 .anchor}**How to analyze a flare with different
attenuator states?**

> During most strong flares, the attenuator state will be changed. The
> following plot is an uncorrected quicklook light curve of the X4.8
> flare on July 23 2003 plotted in the GUI.
>
>  
>
>  
>
> As indicated by the purple horizontal bars, before 00:26:00 UT, the
> attenuator state was A1 (thin in). From 00:26:00 UT to 00:59:24 UT,
> the attenuator state was mostly A3 (both thick and thin in) with four
> short periods in the A1 state. After 00:59:24 UT, the attenuator state
> returned to A1 state. 
>
> To analyze count flux data in SPEX, a count rate spectrum file and a
> spectrometer response matrix (SRM) file must be read into the active
> session of SPEX.  The response matrix is stored in the array DRM (for
> detector response matrix) within the SPEX common blocks.  When using
> the RHESSI spectrum object or GUI to create the count rate spectrum
> file and response matrix file, the name of the response matrix file is
> written in the RESPFILE parameter in the header of the spectrum FITS
> file.  [A different srm file]{.underline} [is needed to analyze the
> count rate data]{.underline} [obtained]{.underline} [in]{.underline}
> [each attenuator state]{.underline} and detector
> combination[.]{.underline} The srm file and count spectra files are
> written out from the Spectrum widget of the GUI as explained in
> Spectroscopy -- First Steps. Alternatively, the srm files can be
> produced [from the spectrum object (sp_obj-\>filewrite,
> /fits,/build,\... ). They are appropriate for the attenuator state in
> the first interval after the start time.]{.underline} [If the
> attenator state changes during the flare, then a new srm file must be
> generated for each state by choosing the first time interval with the
> new attenuator state.]{.underline}
>
> The srm file is specified in SPEX with *dfile*, and then use
> *read_drm* to create the DRM array used with SPEX.
>
> For this 23 July 2002 X4.8 event, one spectrum file and two srmfiles
> were saved with the GUI. They are:
>
>                hsi_spectrum_20020723_001000.fits                    
> \-- spectrum file
>
> hsi_srm_20020723_001000_atten1.fits                  \-- srm file for
> A1 state
>
> hsi_srm_20020723_001000_atten3.fits                  \-- srm file for
> A3 state
>
> If we want to analyze the spectra before 00:26:00 UT, we have to first
> read the spectrum file, and then read the srm file appropriate for the
> A1 state. 
>
> ***SPEX\>  \_1file,
> hsi_spectrum_20020723_001000.fits***                   ***  ***
>
> ***SPEX\>  preview***
>
> ***SPEX\> graph***
>
> ***SPEX\> dfile,hsi_srm_20020723_001000_atten1.fits***
>
> ***SPEX\> read_drm***
>
> This should give the following response if everything is OK:
>
> Command:
> read_drm                                                              
>
> MRDFITS: Null image, NAXIS=0
>
> MRDFITS: Binary table.  6 columns by  92 rows.
>
> MRDFITS: Binary table.  3 columns by  1 rows.
>
> minmax(drm)     0.000000   
> 0.0738023                                          
>
> Change parameters or enter a command

              

> If we want to analyze the spectra between 00:26:00 UT and 00:59:24 UT
> (for now, because of the pileup problem we will not try to analyze the
> time intervals during the four short periods in the A1 state), we have
> to read the srm file appropriate for the A3 state:
>
> ***SPEX\>  dfile,hsi_srm_20020723_001000_atten3.fits***
>
> ***SPEX\> read_drm***
>
> It will give a similar output as above. Again, for the analysis of
> spectra after 00:59:24 UT, we need to switch back to the A1 srm file:
>
> ***SPEX\>  dfile,hsi_srm_20020723_001000_atten1.fits***
>
> ***SPEX\> read_drm***
>
> ***Once the correct spectruma and srm files have been entered into
> SPEX, the analysis can proceed as described in Spectroscopy --First
> Steps.***

**2. **[]{#topic_2 .anchor}**     How to select multiple analysis
intervals?**

**        **

> SPEX is capable to doing spectral fitting for multiple time intervals.
> In order to do this, we need to first select multiple analysis
> intervals. Since the July 23 event lasts for a very long time (about 1
> hour), it is very difficult to select time intervals in the crowded
> time history window, even though the rates are binned every 4 s.
>  First, we can enlarge the time history widow by clicking on one edge
> and dragging it to produce the desired size. The lightcurve with
> background subtracted then appears as follows:
>
>  
>
>  
>
>  
>
>  
>
>  An enlarged version of the time of interest can be obtained with the
> *zoom* command. First type *zoom*, then select the "graphic input",
> "continue by selecting from the window", then select the start time
> and end time of the zoom-in region. A zoomed in time history will
> appear as follows: 
>
>  
>
>  
>
> It is much more convenient to select time intervals from the window
> above to analyze the main impulsive peak of the flare (from 00:26:20
> UT to 00:39:40 UT).
>
> After you type *select,* the following window will pop up:
>
>  
>
>  
>
> This window lists three ways to select time intervals.
>
> []{#Method_1 .anchor}**Method 1**: Boundaries for contiguous
> intervals. This means the end of time one interval is the start time
> of the following interval. To use this method, click on the top button
> in the above window, click "continue by selecting from the window" in
> the next widget, then select the start of the first interval by moving
> the curser to the desitred time in the Time History Window and
> clicking, and then select the end of the first , second, third ...
> intervals until you have selected all the intervals.  To show what you
> have selected, type the following in SPEX:
>
> ***SPEX\>graph***
>
> ***SPEX\>display***
>
> The following selected time history window will appear:
>
>  
>
>  
>
> We selected four contiguous intervals. We can check the variable
> "xselect" to see what the start and end time are for each time
> interval in unit of seconds:
>
> ***SPEX\>
> help,xselect                                                              ***
>
> XSELECT (SPEX_PROC_COM)
>
>                 DOUBLE    = Array\[2, 4\]
>
> ***SPEX\>print,  xselect***
>
>        1647.9990       1711.9990
>
>        1711.9990       1747.9990
>
>        1747.9990       1784.0000
>
>        1784.0000       1819.9990
>
>  
>
> If we want to know the xselect in units of UT, then we need to get UT
> base time (in unit of seconds) by issuing an IDL command:
>
> ***SPEX\> idl, ut_base = getutbase()***
>
>  
>
> We check the ut_base in UT time, we type:
>
> ***SPEX\> idl, ptim, ut_base***
>
>             23-Jul-2002 00:00:00.000
>
>  
>
> Then we print xselect in UT time format:
>
> ***SPEX\> idl, ptim, xselect + base***
>
> 2002/07/23 00:27:27.999 2002/07/23 00:28:31.999
>
> 2002/07/23 00:28:31.999 2002/07/23 00:29:07.999
>
> 2002/07/23 00:29:07.999 2002/07/23 00:29:44.000
>
> 2002/07/23 00:29:44.000 2002/07/23 00:30:19.999
>
> Variable *ifirst* and *ilast* are the start and end time interval
> number to fit. The default of *ifirst* is 0, *ilast* is the last
> interval, which is 3 for this case.
>
> ***SPEX\> print,ifirst,
> ilast                                                        ***
>
> ***       0           3***
>
> If we type *count* to see the count flux spectrum, the spectrum is for
> the first time interval (00:27:27.999 \-- 00:28:31.999 UT):
>
>  
>
>  
>
> Using what we have learned from [Brian's Spectroscopy - First
> Steps](http://hesperia.gsfc.nasa.gov/~dennis/spectroscopy/first_steps.htm),
> we set the function we want to fit (f_vth_bpow model is the default),
> the energy range to fit, initial estimates of the fitting parameters,
> and which ones should be free to vary. We then type *fit* to fit the
> spectra for all the intervals sequentially.. As the fitting proceeds,
> the spectral window will show the fitted spectrum, one by one, and
> stop at the final time interval.  All the fitting parameters are saved
> in variable array, apar_arr:
>
> ***SPEX\> print, apar_arr***
>
>       1.24413      3.35509      14.1803      2.93961      81.1886     
> 3.32381
>
>       2.68748      3.31781      14.0138      2.70112      78.9045     
> 3.14351
>
>       3.57966      3.21793      13.6295      2.60562      81.8491     
> 3.09796
>
>       3.91184      3.17164      8.94853      2.71290      76.7049     
> 2.99297
>
> The reduced chi-square values are stored in the array chi:
>
> ***SPEX\> print, chi***
>
>      0.642303      1.05372     0.780159      1.05493
>
> If we only want to fit the spectrum for interval 0, we need change the
> value of *ilast* to 0
>
> ***SPEX\> ilast, 0***
>
> ***SPEX\>print, ifirst, ilast***

       0          0

> ***SPEX\>count***
>
> ***SPEX\>fit***
>
> The following plot of the best fit f_vth_bpow photon flux spectrum
> appears along with the residuals expressed as the number of sigma
> above or below the model values:
>
>  
>
>  
>
> If we want to fit the spectrum at the second time interval, we need to
> reset *ifirst* and * ilast:*
>
> ***SPEX\>ifirst, 1***
>
> ***SPEX\>ilast, 1***
>
> Then we can do the fit for the second interval. Similarly, we can do
> the third and forth interval and so on.
>
> **Method 2**: Boundaries for discrete intervals. Since the intervals
> are discrete, we need to select the start and end time of each
> interval.  Click on the middle button, click "continue by selecting
> from the window" in the next widget, then select the start and the end
> of each time interval by moving the curser to the appropriate times in
> the Time History Window and clicking with any mouse button. As many
> intervals can be selected as needed but the end of one interval will
> not necessarily be the start of the next.  To show what we have
> selected, we should type the following in SPEX:
>
> ***SPEX\>graph***
>
> ***SPEX\>display***
>
> The selected discrete intervals are shown in the following:
>
>  
>
>  
>
> The vertical orange dotted/broken lines mark the start/end of each
> interval, respectively. To get the time in seconds for each time
> interval:
>
> ***SPEX\> print, xselect***
>
>        1663.9990       1692.0000
>
>        1727.9990       1743.9990
>
>        1760.0000       1779.9990
>
>        1824.0000       1851.9990
>
> Type *fit* to fit all intervals one by one. If we want to fit a
> specific interval, set the *ifirst* and *ilast* to that interval.
> These steps are similar to those in [Mmethod 1](#Method_1) above.
>
> **Method 3**: As limits on sample groups of N (to be chosen). This
> method allows us to choose one long interval and divide it up into
> equal intervals of length N times the basic binning time of the data.
> In our example, we have used 4-s binning times Thus, if we want to fit
> spectra in every 4-s interval, we would set N to 1. If wanted to use
> 8-s intervals we would set N to 2. Here, we set N to 5, so we are
> fitting the spectrum contiguously with 20-s time intervals.
>
> Click on the button marked AS LIMITS ON SAMPLE GROUPS OF N (TO BE
> CHOSEN). The following window pops up: 
>
>  
>
>  
>
>  
>
> Select 5 in the righthand column, click the "ready" button, then click
> "continue by selecting from the window." Finally in the Time History
> Window, click the start and the end of time of the period in which you
> want to make your 20-s intervals. In order to see the intervals we
> selected, type the following:
>
> ***SPEX\>graph***
>
> ***SPEX\>display***
>
> And this window appears:
>
>  
>
>  
>
> We can see evenly selected contiguous 20-s intervals. By checking the
> size of the "xselect" array, we know how many intervals we have
> selected:
>
> ***SPEX\> help, xselect***
>
> XSELECT (SPEX_PROC_COM)
>
> DOUBLE = Array\[2, 37\]
>
> We selected 37 time intervals. We can fit the 37 spectra sequentially
> by assuring that ifirst = 0 and ilast =- 36 and typing fit. After it
> finishes, we can check the best-fit parameters:
>
> ***SPEX\> fit***
>
> ***SPEX\> help,apar_arr***
>
> APAR_ARR (SPEX_PROC_COM)
>
> FLOAT = Array\[6, 37\]
>
> ***SPEX\>print, apar_arr\[\*, 0:3\]***
>
>      0.756227    3.20582      9.46661      3.31592      75.3729     
> 3.83593
>
>     1.09641      3.35124      15.3057      2.90970      86.1484     
> 3.37075
>
>     1.60972      3.38625      16.0934      2.82412      85.3917     
> 3.15549
>
>     2.13309      3.36383      13.8587      2.75644      67.5149     
> 3.19990
>
> ***SPEX\>print, chi\[0:3\]***

                                               0.928469    1.00548    
0.788376     0.934216   

> The window with the fitted spectra and residuals cycles through very
> rapidly in the Spectral Window. However, we can use *create_ps* to
> save all the spectral and residual plots as postscript files
>
> :
>
> ***SPEX\>psstyle ='fullportrait'***
>
> ***SPEX\>create_ps***
>
> Then postscript files named like "fitting_100212_110946.ps" will be
> saved in your IDL working directory. It has 37 pages of spectral
> plots. Each one is one interval we select. This *create_ps* can also
> be used with the former two methods. For more detailed explanation see
> [Topic 8](#Method_1).
>
> If the postscript files we saved are not appropriate for publication
> purposes, or we need plots of the time history of some fitting
> parameters, then we have to save all the parameters and use them
> outside SPEX to make publication-quality plots. 

**3.     ** []{#topic_3 .anchor}**How to change fitting models?**

As explained in [Spectroscopy- First
Steps](http://hesperia.gsfc.nasa.gov/~dennis/spectroscopy/first_steps.htm),
there are many spectrum models available in SPEX. The default is the
photon flux model named \"f_vth_bpow\". It is a combination of an
isothermal function plus a double power-law. Let us try another model
named \"f_vth_thin\". 

 

> Model \"f_vth_thin\" assumes a thin-target bremsstrahlung spectrum. By
> fitting with this model, we are able to get several parameters of the
> electron ensemble that produced the observed X-ray photon flux
> assuming thin-target interactions. This is incontrast to the
> parameters of the measured photon spectrum obtained using the
> \"f_vth_bpow\" model. 
>
>  
>
> First, let us change the model to \"f_vth_thin\" to fit one interval
> (#8) of the above 37 intervals we selected.
>
> S ***PEX\> ifirst, 8***
>
> ***SPEX\> ilast,*** ***8***
>
> set to one interval at 00:30:00 \-- 00:30:19.999 UT
>
> ***SPEX\>  f_model, f_vth_thin***
>
> ***SPEX\> photon***
>
>  
>
> Then we will need to input the initial gauss of the 8 parameters. Here
> are the brief explanation of those 8 parameters:
>
> a(0)= emission measure in units of 10\^49 cm-3
>
> a(1)= KT plasma temperature in keV
>
> a(2) - normalization factor (in 1.0d55), i.e. plasma density \*
> nonthermal electron density \* volume of source
>
> a(3) - Power-law index of the electron distribution function below
> break energy.
>
> a(4) - Break energy in the electron distribution function (in keV)
>
> a(5) - Power-law index of the electron distribution function above
> break energy
>
> a(6) - Low energy cutoff in the electron distribution function (in
> keV)
>
> a(7) - High energy cutoff in the electron distribution function (in
> keV)
>
>  
>
>  Since the computation of this model involves an integration, the
> fitting speed is much slow than that of \"f_vth_bpow\" model. One
> thing we can do it to decrease the high energy cutoff, i.e. a(7). The
> default is 32000 keV, for most spectrum (fitting energy range below
> 300 keV), we can cut down to 5000 keV without affecting the high
> energy flux. 
>
>  
>
> As a rule of thumb, we should always obtain good starting parameters
> before we do any fitting. After changing the parameters, you can check
> the spectral window to compare the model spectrum and the data points.
> If there is not good agreement, the model parameters should be changed
> by typing *the* *photon* command*.* When the starting parameters give
> good agreement with the data points, then you can expect the fitting
> program to converge to good best-fit values. This is achieved with the
> *fit* command.
>
> ***SPEX\>fit***
>
>  
>
> The following plot is produced for interval 0 in the spectrum window:
>
>  
>
>  
>
> We can check the fitting parameters as follows
>
> :
>
> SPEX\>*print, apar*
>
> 4.25845 3.19432 17.1155 1.58129 131.000 2.47378 10.0000 5000.00
>
> SPEX\>*print, ch*i
>
> 1.09850
>
> There is another model called \"f_vth_thick\"that fits a thick-target
> bremsstrahlung spectrum. It also has 8 parameters with the same
> physical meaning as those of \"f_vth_thin\" except for the
> normalization factor, i.e. a(2). For the thick-target model, a(2) is
> the nonthermal electron density \* source area (in units of 1.0e25
> electrons cm^-1^). Since \"f_vth_thick\".involves two integrations,
> the fitting speed is much slower than for the other functions. It
> takes a few minutes to fit one interval. The fitting time depends on
> how many free parameter you set and the value of the high energy
> cut-off -- a(7). This should be at least a factor of ??? greater than
> the highest photon energy to be fit.

**4.     ** []{#topic_4 .anchor}**How to save and restore the fitting
parameters and spectral data?**

> If you want to save your actual SPEX session, simply type
>
> SPEX\>*save_event*
>
> This command creates the files \"summary.dat\" and \"spex_dflts.dat\"
> in your IDL working directory. If you call SPEX again, you can
> continue where you left off by restoring the data with
>
> SPEX\>*restore_event, /verb*
>
> SPEX allows you to export all the resulting parameters to an IDL save
> file, so that you can do any further analysis you are interested in
> outside of SPEX but within the IDL environment. An easy way to
> accomplish this is to call
>
> *** SPEX\> idl, spex_save_data , /all,
> \'spex07202*\_bpowfit_spectra.*sav'***
>
> Whenever we type *idl*, we use the command out of SPEX. The detail of
> SPEX commands can be found at [Analyzing Flare X-ray Spectra Using
> SPEX](http://hesperia.gsfc.nasa.gov/ssw/packages/spex/doc/spex_help.html). 
>
> Using the ***spex_save_data*** command, a file named
> 'spex0723_bpowfit_spectra.sav' will be saved in our IDL working
> directory. It contains all the spectral data and fitting parameters we
> need for later use. We can get the detailed explanation of all saved
> variables from Richard's
> [help_spex.txt](ftp://sohoftp.nascom.nasa.gov/solarsoft/packages/spex/doc/help_spex.txt).
> Here I explain some of them for later use:
>
>   APAR_ARR        FLOAT     = Array(6, 37)
>
> Model best-fit parameters for each interval.
>
> APAR_SIGMA      FLOAT     = Array(6, 37)
>
>         Curvefit uncertainties on the best fit parameters in APAR_ARR
>
> CHI             FLOAT     = Array(1, 37)
>
> Reduced chi square for each accumulation interval   

XSELECT         DOUBLE    = Array(2, 37)

> Time boundaries on accumulation intervals referenced to the same base
> UT.

UTBASE          DOUBLE    =   7.4338560e+008

The UT base

EDGES           FLOAT     = Array(2, 82)

Detector energy channel edges

OBSI            FLOAT     = Array(82, 37)

> Count flux in accumulation interval I.

EOBSI           FLOAT     = Array(82, 37)

Uncertainties on OBSI.

> BACKI           FLOAT     = Array(82, 37)
>
> Background for accumulation interval.
>
> CONVI          FLOAT     = Array(82, 37)
>
> Factor for conversion from count spectrum to photon spectrum.
>
> When we want to restore the data in an IDL session, type

***IDL\> restore,  'spex0723_bpowfit_spctra.sav'***

> Use the help command to see which variables are defined

***IDL\> help***

**5.     ** []{#topic_5 .anchor}**How to plot the fitted spectra outside
SPEX?**

                  

After restoring the data file, in order to plot the spectra, we first
need to compute the photon flux, error bars, background flux, model flux
and so on.

To compute the photon flux (=(obsi-backi)/convi) :

***IDL\>phflux = (obsi-backi)/convi***

To compute the error of the photon flux (=eobsi/convi):

***IDL\>errflux = eobsi/convi***

To compute the background flux (=backi/convi):

***IDL\>backflux = backi/convi***

***IDL\>help, phflux, errflux, backflux***

> PHFLUX          FLOAT     = Array\[82, 37\]
>
> ERRFLUX         FLOAT     = Array\[82, 37\]
>
> BACKFLUX        FLOAT     = Array\[82, 37\]

Since the "edges" is the edges of energy channel, we need to compute the
mean value of each channel:

***Idl\>print, edges\[\*, 0:2\]***

>       3.00000      4.00000
>
>       4.00000      5.00000

5.  6.00000

6.  

*** Idl\>edge_products, edges, mean = eph***

***Idl\>print, eph\[0:2\]***

3.50000      4.50000      5.50000

> So far, we have all the spectral data for plots. We also need the
> f_vth_bpow model flux so that we can overplot on the data like SPEX
> does. By using *xdoc* introduced in [Spectroscopy - First
> Steps](http://hesperia.gsfc.nasa.gov/~dennis/spectroscopy/first_steps.htm),
> we know f_vth_pow model uses f_vth.pro and f_bpow.pro to compute the
> thermal and nonthermal photon flux. The total photon flux is the sum
> of the two fluxes. Let us compute the model flux for the first
> interval using the saved fitting parameters from apar_arr.
>
> First, we get the fitting parameters for the first spectral fit:
>
> ***IDL\> fitpara = apar_arr\[\*, 0\]***
>
> ***IDL\> print, fitpara***
>
>                 0.756227      3.20582      9.46661      3.31592     
> 75.3729      3.83593
>
> Secondly, we can check the UT time of that interval to see whether
> they are correct by comparing the time history window:
>
> ***IDL\>uttime = anytim(xselect\[\*, 0\] + utbase, /ecs)***
>
> ***IDL\>print, uttime***
>
> 2002/07/23 00:27:24.000 2002/07/23 00:27:43.999
>
> Thirdly, compute the thermal photon flux:

***IDL\>thermflux = f_vth(edges, fitpara\[0:1\])***

> ***IDL\> help, thermflux***
>
> THERMFLUX       FLOAT     = Array\[82\]
>
> Finally, we comput the nonthermal photon flux. But we should do
> something beforehand. There is always a low energy cut-off for the
> power-law photon distribution. From [Spectroscopy - First
> Steps](http://hesperia.gsfc.nasa.gov/~dennis/spectroscopy/first_steps.htm),
> we know the default is 1 keV cut-off, and spectral index below 1 keV
> is -1.5. The nonthermal f\_ bpow model uses these two parameters, but
> in the saved fitting parameters, they are not saved. So we have to add
> them at the end of the variable array "fitpara" array:
>
> ***IDL\> a_cutoff = \[1, 1.5\]***
>
> ***IDL\> fitpara = \[fitpara, a_cutoff\]***
>
> ***IDL\> print, fitpara***
>
> ***             ***  0.756227      3.20582      9.46661     
> 3.31592      75.3729      3.83593      1.00000     1.50000
>
> ***IDL\> nonthermflux = f_bpow(eph, fitpara\[2:7\])***
>
> ***IDL\> help, nonthermflux***
>
> NONTHERMFLUX    FLOAT     = Array\[82\]
>
> Now we have all the data we want to plot the spectrum of the first
> interval. The spectrum can be plotted as follows:
>
> ***IDL\> loadct, 4***
>
> ***IDL\> window, 10, xsize = 600, ysize = 500***
>
> ***IDL\> plot, eph, thermflux, /xlog, /ylog, xtitle = \'Photon Energy
> (keV)\', ytitle = \'Flux (photons s!u-1!n cm!u-2!n keV!u-1!n)\',
> xrange = \[10.0,*** ***300.\], yrange=\[1.0e-3, 1.0e6\], xstyle = 1,
> ystyle = 1, /nodata***
>
> ***IDL\> oplot, eph, thermflux, lineStyle = 1, color = 47, thick =
> 4*** ***,***plot thermal flux
>
> ***IDL\> oplot, eph, nonthermFlux, lineStyle = 2, color = 110, thick =
> 4***;plot nonthermal flux***    ***
>
> ***IDL\> oplot, eph, thermflux + nonthermflux, linestyle = 0, color
> =255 , thick = 4***;plot total flux
>
> ***IDL\> oplot, eph, phflux\[\*, 0\], psym = 1, color = 144, symsize =
> 1, thick = 2                                                
>   ***;plot spectrum data
>
> ***IDL\> oplot_err, eph, spexFlux\[\*, 0\], yerr = errFlux\[\*,
> indplot\], psym = 1, symsize = 1, color =*** ***144*** ;plot spectrum
> data error bar
>
> ***IDL\> oplot, eph, backflux\[\*, 0\], psym = 5, symsize = 0.5, color
> =*** ***208*** ;plot the background flux
>
> The following is the spectrum we just plotted:
>
>  
>
>  
>
>  
>
>  

**6.    ** []{#topic_6 .anchor}**How to make a** **movie of the**
**spectra?**

> We already know how to plot the spectra in IDL (see [Topic
> #5](#topic_5)). We can use the IDL routine *TVRD* to read the current
> widow containing the spectrum for each interval. Then we use
> *WRITE_JPEG* to save it as a jpeg file. People who want to know the
> details of this should refer to Fanning's book [IDL Programming
> Techniques](http://www.dfanning.com/documents/books.html).
>
> The following is a sample code:
>
>  
>
> ***IDL\>image24 = tvrd(true = 1)*** ;read current IDL window
>
> ***IDL\>write_jpeg,  jpegName, image24, quality = 75, true =***
> ***1*** ;write to a jpeg file, jpegname is the name of the file
>
> Each spectrum is thus saved as a jpeg file. There are many ways to
> make a movie from multiple jpeg files. Here are two of them:
>
> \(1\) Using movie maker [VideoMach](http://www.gromada.com/)
> (shareware, noncommercial use is free !). It is very easy to make MPEG
> or AVI movies.  Here is the [sample
> movie](../movies/spex_movie_0723_peak.mpeg). As you can see, we add
> the flare light curve at the top of the spectrum window and mark the
> time of each spectrum in lightcurve.
>
> \(2\) using JSMOVIE in SSW. It is good for web publishing. Here is the
> [sample movie](../movies/spex_movie_0723_peak.html).

7.[]{#topic_7 .anchor} **How to choose and evaluate a background model**

> There are lots of useful hints and explanations about background
> subtraction in [RHESSI Spectroscopy - First
> Steps](http://hesperia.gsfc.nasa.gov/~dennis/spectroscopy/first_steps.htm).
> Since a reasonable background subtraction is a fundamental
> pre-requisite for RHESSI hard X-ray spectroscopy, here is something
> more about it.\
> \
> When you analyze the flare of interest, be sure to first check it's
> the long-term variation of the background from the quicklook observing
> summary plots, not only during the time of the flare but also for at
> least one RHESSI orbit before and after the flare. Here we will
> concentrate on the event that occurred on June 1, 2002 that peaked at
> 03:56 UT. The following window shows the observing summary from 01:00
> to 07:00 UT, revealing a sinusoidal variation of the background
> emission due to changes in the magnetic latitude of the RHESSI
> spacecraft along its orbit (for details on the RHESSI background
> components see David Smith\'s [Problems, Pecularities and Phenomena in
> RHESSI
> Spectroscopy](http://hesperia.gsfc.nasa.gov/rhessidatacenter/spectroscopy/ppp1.pdf)
> document). The flare occurs around the maximum of the background
> emission when the spacecraft was at high geomagnetic latitudes. This
> makes the background subtraction a rather subtle issue. 
>
> For the analysis, set the observation time interval in the GUI to
> 1-Jun-03:20:00 - 1-Jun-04:15:00:\
>   
>
> Then go to \"Retrieve/Process Data\" -\> \"Spectrum\...\". Choose all
> front detectors except 2 and 7 (the peculiarities of individual
> detectors are described in the [Problems, Pecularities and Phenomena
> in RHESSI
> Spectroscopy](http://hesperia.gsfc.nasa.gov/rhessidatacenter/spectroscopy/ppp1.pdf)
> document by David Smith). Set the time resolution to 4.0 s and choose
> bin code 7 (197 bins, 1-keV bands covering 3 to 100 keV, then 5-keV
> bands to 600 keV - check the button \"Show Binning Codes\" to see what
> choices you have). This is the best choice for this flare since there
> is emission detectable above 100 keV.
>
> Click the "Plot Spectrum" button to accumulate the spectra and display
> them for all the intervals. Once the spectra have been obtained, click
> the button \"Write output file\...\" to save the spectral data. Select
> the full response matrix (\"0-Full calculation of diagonal,
> off-diagonal terms\") -\> \"Write FITS file\". Your window should look
> like:
>
> Then call SPEX from an IDL session with the command \"spex_proc\".
>  Read the spectral data and the spectrometer response (srm) matrix
> into SPEX and plot the count rates as follows:\
> \
> ***SPEX\> data,hessi,front\
> SPEX\> \_1file,hsi_spectrum_20020601_032000.fits\
> SPEX\> dfile,hsi_srm_20020601_032000.fits\
> SPEX\> graph\
> SPEX\> energy_bands,6,12,12,25,25,50,50,300\
> SPEX\> scale_bands,1,2,6,8\
> SPEX\> th_ytype, 1\
> SPEX\> graph***
>
>                             
>
> The sinusoidal background behavior tends to show up more prominently
> at high energies. So, you might try to use a higher-order polynomial
> fit at these energies. Let us apply a linear fit for the global energy
> range as well as for the energy band 6-12 keV, and a 3rd order
> polynomial fit for the higher energy bands (12-25, 25-50, 50-300
> keV).\
> \
> ***SPEX\> back_order, 1,1, 3,3,3***\
> \
> Now, starting with the first band we have to choose reasonable
> intervals before and after the flare during which the background is
> calculated.\
> \
> ***SPEX\> use_band= -1\
> SPEX\> background ***    ; select graphically the interval to be used
> for background determination\
> \
> The same for the other energy bands. And you probably want to change
> the y-range.\
> \
> ***SPEX\> th_yrange,1e-3,10\
> SPEX\> use_band=0\
> SPEX\> back                  *** ; select background intervals*\
> **SPEX\> use_band=1\
> SPEX\> back  **                 ; * select background intervals***\
> SPEX\> use_band=2\
> SPEX\> back                  *** ;  select background intervals***\
> SPEX\> use_band=3\
> SPEX\> back                 **  ;  select background intervals*\
> \
> The two images illustrate possible background selections for
> use_band=0 (3-12 keV band, chosen background model: order 1) and
> use_band=2 (25-50 keV, chosen background model: order 3) .\
> \
> \
> \
> \
> \
> There is no easy way within SPEX to plot the original count rates
> together with the modeled background. However, as always in SPEX you
> can execute every IDL-command from the SPEX prompt by typing \"idl,
> command\". Let us use this method to get a basic idea whether the
> derived background model is reasonable. The count rates before
> background subtraction are stored in the two-dimensional variable
> \"flux\" and the background data in \"back\":\
> \
> ***SPEX\> help, flux, back\
> ***\
> The first index refers to energy bins, the second to time bins. See
> the
> [help_spex.txt](ftp://sohoftp.nascom.nasa.gov/solarsoft/packages/spex/doc/help_spex.txt)
> document by Richard Schwartz for a description of parameters used
> within SPEX.\
> \
> Let us first check the energy band 6-12 keV where the chosen model is
> a 1st order polynomial (and remember that our data start at 3 keV).\
> \
> ***SPEX\> idl, plot, total ( flux \[ 6-3:12-3-1, \* \], 1 ), /ylog,
> /xst\
> SPEX\> idl, oplot, total ( back \[ 6-3:12-3-1, \* \], 1 )\
> ***\
> \
> \
> Probably you\'d prefer to plot the count rates as function of time, so
> let us do this for the 25-50 keV band:\
> \
> ***SPEX\> idl, utplot, ut \[ 0, \* \], total( flux \[ 25-3:50-3-1, \*
> \] , 1 ), /ylog, /xst, col=170\
> SPEX\> idl, outplot, ut \[0, \*\], total( back\[ 25-3:50-3-1, \* \]
> ,1)***\
> \
> \
> \
> In this case, the 3rd order plynomial fit yields a very good result.
> After evaluating the background model for the other energy bands, you
> may apply a 3rd order fit to all of them. To do this, set
> \"back_order, 3, 3, 3, 3, 3\" , and redo the background subtraction
> for use_band = -1 and use_band = 0.\
> If your background estimates turn out to be rather bad, then you may
> wish to start from scratch. To undefine any background you have
> already subtracted, use the command \"background, clear\" (sets the
> background to zero everywhere). Afterwards, don\'t forget to start (as
> always) with \"use_band= -1\".\
> \
> Sometimes you may find that you don\'t have a long enough continuous
> observation interval before and/or after the flare (due to RHESSI
> night times or the spacecraft passing through the South Atlantic
> Anomaly) to choose a background model of an order as high as 3 (by the
> way: 3 is the highest order SPEX allows). In these cases, it\'s
> probably the most safe solution to choose for each energy band
> individually a short interval before and after the flare (1 min or so,
> depending on the actual context) and apply a linear fit, i.e. use
> \"backorder, 1, 1, 1, 1, 1\". In particular, in short events 1st order
> models give reasonable results also when the flare is located on top
> of a high-latitude variation. When evaluating your chosen background
> model, keep in mind that inaccuracies in the background subtraction
> have a relatively stronger effect at high energies where the flare
> count rates are low. For alternative approaches to the background
> determination (e.g., using the data of one RHESSI orbit before/after
> the flare) see the suggestion in David Smith\'s [Behind, Beneath and
> Before HESSI Spectroscopy or
> \"BBB\"](http://hesperia.gsfc.nasa.gov/Data_Analysis/hedaws/goddard/smith/index.htm).
> You can also use data from 15 orbits before and/or after the flare to
> try and duplicate the same range of geomagnetic parameters as for the
> flare observation itself.\
> \
> Finally, select a time interval around the maximum of the 25-50 keV
> emission and plot the count spectrum.\
> \
> ***SPEX\> select   ;*** select graphically a time range***\
> SPEX\> count***\
> \
>  
>
> Comparing the background subtracted count spectrum (yellow) with the
> background spectrum (orange) you see that at the chosen time interval
> your flare spectrum is reliable up to about 200 keV. Above,  the
> measured count spectrum becomes dominated by the background spectrum.

**8.** **How to create and customize postscript plots within SPEX**

> SPEX allows you to create a postscript plot of the last plot(s) that
> you have done (i.e. in the time history or the spectral window) . The
> command is ***create_ps \[, filename\]***. So, for instance, if your
> last plot has been the flare count spectrum, then create a ps-file
> with\
> \
> ***SPEX\> create_ps, countspec_01June02.ps***\
> \
> If you have selected multiple time ranges to make spectral fits at all
> these time steps (as discussed in [Topic 2](#topic_2).), then the
> ***create_ps*** command makes a plot that contains all these spectral
> fits (and residuals), each on one page. This feature may be
> interesting to you when you want to have a quicklook to see if the
> fits you obtained in all the intervals are reasonable (check the fits
> as well as the residuals). You can enlarge the plot range to the full
> page size with calling:\
> \
> ***SPEX\> psstyle=\"FULLPORTRAIT\"***\
> \
> before the ***create_ps*** command.\
> \
> The standard setting of the psstyle parameter is \"PORTRAIT\", another
> option is  \"LANDSCAPE\".\
> \
> If you want to switch to color mode in the ps-plot, call first\
> \
> ***SPEX\> idl, sps, /color***\
> \
> The ***create_ps*** command provides you with an easy tool within SPEX
> to make hardcopy outputs of your results. However, if you want to make
> plots with specific settings then you have to export your data to IDL
> (see [Topic 3](#topic_3).) and use your own routines.