20 January 2005 Example using OSPEX
Natalia describes her problems in the following email:
********************************************
Dear colleagues!
Excuse me for this question. But it`s really important for us and I
hope you will be able to help us.
We need to build the spectrum for event 20 Jan 2005
06:46:00-06:47:00.
I`ve found the ready spectrum ( .png file) for this time on Hessi site.
To my regret I couldn`t attach figures for you. Your server doesn`t take such
messages.
For making spectrum I used
time 06:30-06:50
Time for background 06:30-06:32
Energy 2-30 keV
1 Interval 2-60 keV
the next intervals every 3 keV (6-9 keV, 9-12 keV and so on).
We choose the same detectors 1F 3F 4F 8F 9F
When we`ve made spectrum in photons we have too large values for
maximum for about 10^18 photons but we see about 10^10 at the file .png
May be you can explain us why we have so different values. Using
programme OSPEX we received the too high photon flux also as you
see at the picture photon_flux_2_30.ps
and I can`t make right fitting for thermal component.
We have too much thermal EM.
How should we receive correct spectrum?
How should I send you my figures? I can`t send files as ps rar and zip.
Thank you a lot!
--
Best regards,
Natalia
****************************************************
The following are plots made to reproduce her problem of getting 10^18 photons.
We are looking at 6:46 to 6:47, this is roughly the peak of the flare.
Here is the HESSI corrected count time profile:
I select the same background as her: 6:30-6:32. Note: it is not night time selection.
This is a fit in Counts vs. Energy. This fit is over all energies during 6:46-6:47 with vth+3pow.
Binning is divided the same way Natalia had it: 2-6, 6-60 divided into 3 keV bins.
This is the photon conversion of the same fit:
This is the script for the above plots:
**************************************
; OSPEX script created Tue Nov 07 14:22:04 2006 by OSPEX writescript method.
;
; Call this script with the keyword argument, obj=obj to return the
; OSPEX object reference for use at the command line as well as in the GUI.
; For example:
; ospex_script_3kev_bins, obj=obj
;
; Note that this script simply sets parameters in the OSPEX object as they
; were when you wrote the script, and optionally restores fit results.
; To make OSPEX do anything in this script, you need to add some action commands.
; For instance, the command
; obj -> dofit, /all
; would tell OSPEX to do fits in all your fit time intervals.
; See the OSPEX methods section in the OSPEX documentation at
; http://hesperia.gsfc.nasa.gov/ssw/packages/spex/doc/ospex_explanation.htm
; for a complete list of methods and their arguments.
;
pro ospex_script_3kev_bins, obj=obj
if not is_class(obj,'SPEX',/quiet) then obj = ospex()
obj-> set, spex_specfile= 'C:\RHESSI\Natalia\hsi_spectrum_20050120_063000.fits'
obj-> set, spex_drmfile= 'C:\RHESSI\Natalia\hsi_srm_20050120_063000.fits'
obj-> set, spex_source_angle= 66.5197
obj-> set, spex_source_xy= [840.000, 264.000]
obj-> set, spex_fit_time_interval= ['20-Jan-2005 06:46:00.000', $
'20-Jan-2005 06:47:00.000']
obj-> set, spex_bk_time_interval=['20-Jan-2005 06:30:00.000', '20-Jan-2005 06:32:00.000']
obj-> set, fit_function= 'vth+3pow'
obj-> set, fit_comp_params= [1.59690, 2.91409, 1.00000, 19432.5, 1.50000, 7.51810, 5.45100, $
29.7274, 2.62267]
obj-> set, fit_comp_minima= [1.00000e-020, 0.500000, 0.0100000, 1.00000e-020, 0.500000, $
1.00000, 0.500000, 10.0000, 0.500000]
obj-> set, fit_comp_maxima= [1.00000e+020, 50.0000, 10.0000, 1.00000e+020, 12.0000, $
20000.0, 12.0000, 20000.0, 12.0000]
obj-> set, fit_comp_free_mask= [1B, 1B, 0B, 1B, 0B, 1B, 1B, 1B, 1B]
obj-> set, fit_comp_spectrum= ['full', '']
obj-> set, fit_comp_model= ['chianti', '']
obj-> set, spex_autoplot_units= 'Counts'
obj-> set, spex_fitcomp_plot_units= 'Counts'
obj-> set, spex_fitcomp_plot_photons= 1
obj-> set, spex_eband= [[3.00000, 6.00000], [6.00000, 12.0000], [12.0000, 25.0000], $
[25.0000, 50.0000], [50.0000, 100.000]]
obj-> set, spex_tband= [['20-Jan-2005 06:30:00.000', '20-Jan-2005 06:45:00.000'], $
['20-Jan-2005 06:45:00.000', '20-Jan-2005 07:00:00.000'], ['20-Jan-2005 07:00:00.000', $
'20-Jan-2005 07:15:00.000'], ['20-Jan-2005 07:15:00.000', '20-Jan-2005 07:30:00.000']]
obj -> restorefit, file='C:\RHESSI\Natalia\ospex_results_3kev_bins.fits'
end
*********************************************
Conclusion: I can not replicate Natalia's problem with 10^18 photons. My data is close to Jim's, but my fits are not. So now I have the problem of trying to find out how Jim got his fits.
Now I am trying to replicate Jim's Plots. All fits below are done with Uncertainty = 0.05
Here are Jim's plots, taken from the web:
Note: Jim's Emission measure = 11.2, Temperature = 2.24 ( ~25.9)
I have done both pile-up enabled and disabled. Here are the results for pile-up disabled.
background selected (this is the same background Natalia used):
Photon spectra:
Count vs. Energy
Script:
**********************************************
; OSPEX script created Tue Nov 07 14:52:03 2006 by OSPEX writescript method.
;
; Call this script with the keyword argument, obj=obj to return the
; OSPEX object reference for use at the command line as well as in the GUI.
; For example:
; ospex_script_jim, obj=obj
;
; Note that this script simply sets parameters in the OSPEX object as they
; were when you wrote the script, and optionally restores fit results.
; To make OSPEX do anything in this script, you need to add some action commands.
; For instance, the command
; obj -> dofit, /all
; would tell OSPEX to do fits in all your fit time intervals.
; See the OSPEX methods section in the OSPEX documentation at
; http://hesperia.gsfc.nasa.gov/ssw/packages/spex/doc/ospex_explanation.htm
; for a complete list of methods and their arguments.
;
pro ospex_script_jim, obj=obj
if not is_class(obj,'SPEX',/quiet) then obj = ospex()
obj-> set, spex_specfile= 'C:\RHESSI\Natalia\hsi_spectrum_20050120_063000_jim.fits'
obj-> set, spex_drmfile= 'C:\RHESSI\Natalia\hsi_srm_20050120_063000_jim.fits'
obj-> set, spex_source_angle= 66.5197
obj-> set, spex_source_xy= [840.000, 264.000]
obj-> set, spex_erange= [3.0000000D, 250.00000D]
obj-> set, spex_fit_time_interval= ['20-Jan-2005 06:46:04.000', $
'20-Jan-2005 06:47:04.000']
obj-> set, spex_bk_time_interval=['20-Jan-2005 07:27:44.000', '20-Jan-2005 07:29:08.000']
obj-> set, fit_function= 'vth+3pow'
obj-> set, fit_comp_params= [0.606347, 2.88523, 1.00000, 212474., 1.50000, 5.48417, $
6.10311, 29.3263, 2.86196]
obj-> set, fit_comp_minima= [1.00000e-020, 0.500000, 0.0100000, 1.00000e-020, 0.500000, $
1.00000, 0.500000, 10.0000, 0.500000]
obj-> set, fit_comp_maxima= [1.00000e+020, 50.0000, 10.0000, 1.00000e+020, 12.0000, $
20000.0, 12.0000, 20000.0, 12.0000]
obj-> set, fit_comp_free_mask= [1B, 1B, 0B, 1B, 0B, 1B, 1B, 1B, 1B]
obj-> set, fit_comp_spectrum= ['full', '']
obj-> set, fit_comp_model= ['chianti', '']
obj-> set, spex_autoplot_overlay_back= 0
obj-> set, spex_autoplot_photons= 1
obj-> set, spex_autoplot_units= 'Counts'
obj-> set, spex_fitcomp_plot_units= 'Counts'
obj-> set, spex_fitcomp_plot_photons= 1
obj-> set, spex_eband= [[3.00000, 6.00000], [6.00000, 12.0000], [12.0000, 25.0000], $
[25.0000, 50.0000], [50.0000, 100.000], [100.000, 300.000]]
obj-> set, spex_tband= [['20-Jan-2005 06:30:00.000', '20-Jan-2005 06:45:00.000'], $
['20-Jan-2005 06:45:00.000', '20-Jan-2005 07:00:00.000'], ['20-Jan-2005 07:00:00.000', $
'20-Jan-2005 07:15:00.000'], ['20-Jan-2005 07:15:00.000', '20-Jan-2005 07:30:00.000']]
obj -> restorefit, file='C:\RHESSI\Natalia\ospex_results_jim.fits'
end
**********************************************
Conclusion: I still can't get Jim's fits, so I'm going to try to change background selection to night time AND enabling pile-up.
Pile-up enabled
Night time background selected:
Plot with Jim's Parameters:
Red spectra - sum of all components
Green spectra - vth
Black spectra - data
Note: I fixed all parameters. That is why you see a chi-square of 24.26 on this plot. This is to show that Jim's best fit parameters are different from my data for the same time interval!
Red spectra - sum of all components
Green spectra - vth
Black spectra - data
Note: again DIFFERENT!!!!
Jim's Parameters and fitted:
Now, I try fitting with Jim's parameters as the starting parameters and note that the emission measure is 0.945 and temperature is 2.78.
*******************************
Recap!
Jim - EM=11.2 , T=2.24 keV
Chau - EM=0.945, T=2.78 keV
Natalia - I still can't reproduce her 10^8 photons plot.
********************************************
So we go to GOES plot for a 4th opinion!
GOES - EM ~ 30 10^49 cm ^3 , T ~ 25 K (~2 keV)
Lin says that GOES measures EM to be roughly 4 times higher than RHESSI, so we predict our emission measure to be less than 7.
So next, we see if there is a trend for RHESSI to predict the EM and Temperature at 06:46:04-06:47:04.
Parameter 1 = Temperature
Parameter 0 = Emission Measure
Temperature decreases a bit, but this follows the GOES trend.
Note: Emission measure = 3.26, Temperature = 2.97, this is finishes my story to say that temperature is about 2 which is in agreement with all trend, GOES and Jim.