Simulation for Suzaku: Evaluate HXD systematic errors
OverviewLast Update: 1 Dec 2006 - updated for CIAO 3.4: removed "AO1" from thread title Synopsis: This thread describes how to simulate Suzaku/HXD (Hard X-ray Detector) observations and evaluate the HXD background using Sherpa. It is based on the Cookbook for HXD Simulations described by the HXD Team, but uses a somewhat different approach that leads to similar but not identical results. It is presented here to show how HXD simulations can be done with Sherpa, and to show how systematic errors can be treated within Sherpa. Related Links:
|
Contents
- Getting Started
- Simulate the Model Spectra
- Examine the HXD Background Level
- Effect of Background Systematic Error
- Save the Results
- Commands Used in this thread
- History
- Images
Getting Started
The sensitivity of the HXD is dominated by the systematic error of the instrumental background estimates. Only objects brighter than the systematic error of the background estimate at the energy band in question are detectable. The HXD team strongly suggests that users check the effect of background systematic error (BGD-sys-err) in the simulation.
The PHA and response files used in this thread may be obtained from the Cookbook for HXD Simulations. The following files are needed:
ae_hxd_gso_20051019.rsp ae_hxd_gsobkg_20051105.pha ae_hxd_pinbkg_20051105.pha ae_hxd_pinhxnom_20051104.rsp
Note that we use the HXD nominal postion from the file ae_hxd_pinhxnom_20051104.rsp in these simulations. If you choose to use the default XIS nominal position instead, you will need to replace that file with ae_hxd_pinxinom_20051104.rsp in your analysis.
Simulate the Model Spectra
In order to create a simulated dataset, it is necessary to define a source model and a grid over which to evaluate the model stack. We will use the datasets ae_hxd_pinbkg_20051105.pha and ae_hxd_gsobkg_20051105.pha to create the grid; the input datasets will be overwritten by the simulated data created by FAKEIT. Note that we read each file twice, first as a data set (with the "data" command) and then as a background ("back" command). After reading the data and the background files, we input the response files and set the instrument models for the two data sets.
sherpa> paramprompt off Model parameter prompting is off sherpa> sherpa.dataplot.x_log = 1 sherpa> sherpa.dataplot.y_log = 1 sherpa> data 1 ae_hxd_pinbkg_20051105.pha The inferred file type is PHA. If this is not what you want, please specify the type explicitly in the data command. sherpa> back 1 ae_hxd_pinbkg_20051105.pha The inferred file type is PHA. If this is not what you want, please specify the type explicitly in the data command. sherpa> data 2 ae_hxd_gsobkg_20051105.pha The inferred file type is PHA. If this is not what you want, please specify the type explicitly in the data command. sherpa> back 2 ae_hxd_gsobkg_20051105.pha The inferred file type is PHA. If this is not what you want, please specify the type explicitly in the data command. sherpa> rsp[pin] sherpa> rsp[gso] sherpa> pin.rmf = ae_hxd_pinhxnom_20051104.rsp sherpa> gso.rmf = ae_hxd_gso_20051019.rsp sherpa> instrument 1 = pin sherpa> instrument 2 = gso
Next, define the source model and set the inital value of various model parameters. Here we use a power law model multiplied by a constant:
sherpa> source 1,2 = const1d[c1]*powlaw1d[p1] sherpa> c1 integrate off sherpa> c1.c0 = 0.1 sherpa> p1.gamma = 2.1 sherpa> p1.ref = 1.0 sherpa> p1.ampl.min = 0.0 sherpa> p1.ampl.max = 100.0 sherpa> p1.ampl = 9.7 sherpa> freeze p1.ampl
The simulation time is set to 1e5 seconds for each of the datasets before running fakeit to create the spectra:
sherpa> fakeit 1 time = 1e5 sherpa> fakeit 1 FAKEIT: The current background data have been added to the faked spectrum. sherpa> fakeit 2 time = 1e5 sherpa> fakeit 2 FAKEIT: The current background data have been added to the faked spectrum.
Finally, restrict the energy range of each file to the default energy band of the Hard X-ray Detector (HXD):
sherpa> ignore 1 energy 60.0: sherpa> ignore 2 energy :30.0,600.0:
Examine the HXD Background Level
Now we can examine the HXD background levels and plot the data and the background. First we subtract the background from the simulated PHA files, then overplot the data and background for both instruments. The color of the curves is changed to distinguish between them easily.
sherpa> subtract 1,2 sherpa> oplot data 1 data 2 back 1 back 2 Warning: negative and zero values ignored in log scale sherpa> limits x 10 600 sherpa> limits y 0.001 1 sherpa> c 2 blue sherpa> c 3 green sherpa> c 4 red sherpa> redraw
Figure 1 ![[Link to Image 1: Background-subtracted data and backgrounds]](../../imgs/imageicon.gif){kind=small}
shows the resulting
plot. Another filter is used to ignore the undetectable band:
sherpa> ignore 2 energy 250.0:
Effect of Background Systematic Error
We can evaluate the background systematic error using the backscale in Sherpa. First we unsubtract the data and set the backscale parameter to 1 for both simulated data sets and the background files. Then we fit the data to get the best fit parameters for the photon index and the constant, and obtain the confidence range for these parameters (projection). Note that the power law amplitude remains frozen at the simulated value.
sherpa> unsubtract 1,2 sherpa> setdata 1 backscale = 1.0 sherpa> setdata 2 backscale = 1.0 sherpa> setback 1 backscale = 1.0 sherpa> setback 2 backscale = 1.0 sherpa> subtract 1,2 sherpa> fit LVMQT: V2.0 LVMQT: initial statistic value = 76.9445 LVMQT: final statistic value = 76.5725 at iteration 9 c1.c0 0.0997023 p1.gamma 2.09921 sherpa> projection Projection complete for parameter: c1.c0 Projection complete for parameter: p1.gamma Computed for sherpa.proj.sigma = 1 -------------------------------------------------------- Parameter Name Best-Fit Lower Bound Upper Bound -------------------------------------------------------- c1.c0 0.0997023 -0.000915563 +0.000925169 p1.gamma 2.09921 -0.00321998 +0.00322589 sherpa> oplot data 1 data 2 back 1 back 2 Warning: negative and zero values ignored in log scale sherpa> limits x 10 600 sherpa> limits y 0.001 1 sherpa> c 2 blue sherpa> c 3 green sherpa> c 4 red sherpa> redraw
Figure 2 shows the results with
all backscales set to 1.0. The simulated values of photon
index and constant parameters are recovered with a 1 sigma
lower and upper bound as indicated by projection. This range
needs to be compared to the systematic background error,
which is 10%.
Next the background is overestimated by setting the background backscales to 0.9:
sherpa> unsubtract 1,2
sherpa> setback 1 backscale = 0.9
sherpa> setback 2 backscale = 0.9
sherpa> subtract 1,2
sherpa> fit
LVMQT: V2.0
LVMQT: initial statistic value = 6215.34
LVMQT: final statistic value = 4071.34 at iteration 9
c1.c0 0.144673
p1.gamma 2.24062
sherpa> projection
Projection complete for parameter: c1.c0
Projection complete for parameter: p1.gamma
Computed for sherpa.proj.sigma = 1
--------------------------------------------------------
Parameter Name Best-Fit Lower Bound Upper Bound
--------------------------------------------------------
c1.c0 0.144673 -0.00140331 +0.00142631
p1.gamma 2.24062 -0.00346758 +0.00350342
sherpa> oplot data 1 data 2 back 1 back 2
Warning: negative and zero values ignored in log scale
sherpa> limits x 10 600
sherpa> limits y 0.001 1
sherpa> c 2 blue
sherpa> c 3 green
sherpa> c 4 red
sherpa> redraw
Figure 3 shows the results with
background backscales set to 0.9.
Finally we underestimate the background by setting the background backscales to 1.1:
sherpa> unsubtract 1,2
sherpa> setback 1 backscale = 1.1
sherpa> setback 2 backscale = 1.1
sherpa> subtract 1,2
sherpa> fit
LVMQT: V2.0
LVMQT: initial statistic value = 11579.9
LVMQT: final statistic value = 2763.36 at iteration 12
c1.c0 0.0679399
p1.gamma 1.95657
sherpa> projection
Projection complete for parameter: c1.c0
Projection complete for parameter: p1.gamma
Computed for sherpa.proj.sigma = 1
--------------------------------------------------------
Parameter Name Best-Fit Lower Bound Upper Bound
--------------------------------------------------------
c1.c0 0.0679399 -0.000576908 +0.000579968
p1.gamma 1.95657 -0.00290537 +0.00289905
sherpa> oplot data 1 data 2 back 1 back 2
Warning: negative and zero values ignored in log scale
sherpa> limits x 10 600
sherpa> limits y 0.001 1
sherpa> c 2 blue
sherpa> c 3 green
sherpa> c 4 red
sherpa> redraw
Figure 4 shows the results with
background backscales set to 1.1.
The three p1.gamma values obtained in this thread for different background backscale indicate a range of the systematic error in the simulation. The systematic error is larger than the statistical error (found by projection in each case). We assumed p1.gamma=2.1 and c1.c0=0.1 in the simulation. The systematic error in this simulation indicates p1.gamma values between 1.95-2.24 and c1.c0 values between 0.06-0.14.
Save the Results
First we write out the simulated datasets. The save and write mdl commands are then used to record the rest of the session information.
sherpa> write data fake1.pha PHA sherpa> write data 2 fake2.pha PHA sherpa> save all suzaku.shp sherpa> write mdl "suzaku_mdl.fits"
suzaku.shp is text file with all the information needed to restore this session. suzaku_mdl.fits (the MDL file) stores the data, model, and other information in FITS format.
Commands Used in this thread
The text file suzaku_commands.txt contains all the commands issued in this thread:
paramprompt off sherpa.dataplot.x_log = 1 sherpa.dataplot.y_log = 1 data 1 ae_hxd_pinbkg_20051105.pha back 1 ae_hxd_pinbkg_20051105.pha data 2 ae_hxd_gsobkg_20051105.pha back 2 ae_hxd_gsobkg_20051105.pha rsp[pin] rsp[gso] pin.rmf = ae_hxd_pinhxnom_20051104.rsp gso.rmf = ae_hxd_gso_20051019.rsp instrument 1 = pin instrument 2 = gso source 1,2 = const1d[c1]*powlaw1d[p1] c1 integrate off c1.c0 = 0.1 p1.gamma = 2.1 p1.ref = 1.0 p1.ampl.min = 0.0 p1.ampl.max = 100.0 p1.ampl = 9.7 freeze p1.ampl fakeit 1 time = 1e5 fakeit 1 fakeit 2 time = 1e5 fakeit 2 ignore 1 energy 60.0: ignore 2 energy :30.0,600.0: subtract 1,2 oplot data 1 data 2 back 1 back 2 limits x 10 600 limits y 0.001 1 c 2 blue c 3 green c 4 red redraw ignore 2 energy 250.0: unsubtract 1,2 setdata 1 backscale = 1.0 setdata 2 backscale = 1.0 setback 1 backscale = 1.0 setback 2 backscale = 1.0 subtract 1,2 fit projection oplot data 1 data 2 back 1 back 2 limits x 10 600 limits y 0.001 1 c 2 blue c 3 green c 4 red redraw unsubtract 1,2 setback 1 backscale = 0.9 setback 2 backscale = 0.9 subtract 1,2 fit projection oplot data 1 data 2 back 1 back 2 limits x 10 600 limits y 0.001 1 c 2 blue c 3 green c 4 red redraw unsubtract 1,2 setback 1 backscale = 1.1 setback 2 backscale = 1.1 subtract 1,2 fit projection oplot data 1 data 2 back 1 back 2 limits x 10 600 limits y 0.001 1 c 2 blue c 3 green c 4 red redraw write data fake1.pha PHA write data 2 fake2.pha PHA save all suzaku.shp write mdl "suzaku_mdl.fits"
History
| 14 Jan 2005 | reviewed for CIAO 3.2: no changes |
| 14 Dec 2005 | updated for CIAO 3.3: Astro-E2 is now named Suzaku; new HXD data files, fits and plots updated accordingly |
| 01 Dec 2006 | updated for CIAO 3.4: removed "AO1" from thread title |
