Simulating an LETG/ACIS-S observation

In this example, we simulate an LETG/ACIS-S observation with marx. The primary goal of this exercise is to see how well marx can simulate this instrument combination by performing a spectral analysis of the result. For simplicity, we assume an on-axis point source and an absorbed powerlaw spectrum with a column density of 10^19 atoms/cm^2, a spectral index of 1.8, and a normalization of 0.1 photons/keV/cm^2/sec at 1 keV, Rather than specifying an exposure time for the observation as in the other examples, here we specify the desired number of detected events. In this case, the simulation will be run until 10,000,000 (1e7) events have been detected.

The first step is to create a 2-column text file that tabulates the absorbed powerlaw flux [photons/sec/keV/cm^2] (second column) as a function of energy [keV] (first column). The easiest way to create such a file is to make use of a spectral modeling program. Here as in the other examples, we use isis,

In isis the absorbed powerlaw model is specified using:

isis> fit_fun("phabs(1)*powerlaw(1)");
isis> set_par(1, 0.001);
isis> set_par(2, 0.1);
isis> set_par(3, 1.8);
isis> list_par;
phabs(1)*powerlaw (1)
 idx  param             tie-to  freeze         value         min         max
  1  phabs(1).nH            0     0            0.001           0      100000  10^22
  2  powerlaw(1).norm       0     0              0.1           0        0.01
  3  powerlaw(1).PhoIndex   0     0              1.8           1           3
isis> save_par ("letgplaw.p");

The next step is to convert the parameter file letgplaw.p to the spectrum file that marx expects. The marxflux script may be used to create a file called letgplawflux.tbl in the appropriate format via

 unix% ./marxflux -e '[0.03:12.0:0.001]' letgplaw.p letgplawflux.tbl

The letgplawflux.tbl file is input to marx using the following marx.par parameters:

   SpectrumType=FILE
   SpectrumFile=letgplawflux.tbl
   SourceFlux=-1

    marx SourceFlux=-1 SpectrumType="FILE" SpectrumFile="letgplawflux.tbl" \
       NumRays=-10000000 ExposureTime=0 TStart=2012.5 \
       OutputDir=letgplaw GratingType="LETG" DetectorType="ACIS-S" \
       DitherModel="INTERNAL" RA_Nom=30 Dec_Nom=40 Roll_Nom=50 \
       SourceRA=30 SourceDEC=40 \
       Verbose=yes mode=h

    marx2fits letgplaw letgplaw_evt1.fits

The fits file (letgplaw_evt1.fits) can be further processed with standard CIAO tools. As some of these tools require the aspect history, the marxasp program is used to create an aspect solution file that matches the simulation:

    marxasp MarxDir="letgplaw" OutputFile="letgplaw_asol1.fits"

It is interesting to look at the event file with a viewer such as ds9. Here we use evt2img, which displays an RGB image of the events projected to the sky plane.

.

The energy bands were color-coded as follows: red: 0.0 to 0.5 keV, green: 0.5 to 2.0 keV, blue: 2.0 to 12 keV.

This image shows what the pattern of events looks like near zeroth order. There are a number of features in this image that are worth pointing out. The thick blue line in the middle of the fan-like structure that extends from the top left of the figure to the bottom right corresponds to diffracted photons from the primary grating bars of the LEG. The events in this line are the subject of the spatial extraction step when creating the event histograms for spectral analysis, as described below. The other blue lines in the fan-like structure correspond to photons that have also been diffracted by the fine support structure of the LEG. This fine support structure forms a grating that diffracts in a direction orthogonal to the primary diffraction direction. The events in the multicolored line that runs from the bottom left to the upper right consists of two contributions. The first is from photons that arrived at the detector while the image was being read out between frames. The second is from the undiffracted photons from the primary grating that have been diffracted from the LEG fine support structure. The red star-like pattern in the center of the figure was produced by photons that have diffracted from the LEG's coarse support structure, whose 2 mm period diffracts only the very lowest energy photons through an appreciable angle. Finally, if one looks closely enough, the shadowing of the mirror support struts may be seen near zeroth order.

For this example, we wish to verify that the marx simulation is consistent with the input spectrum. To this end, we use isis to fit an absorbed powerlaw to the pha spectra. The spectral fits to the plus and minus first order histograms can be carried out using the following isis commands:

    m1 = load_data ("letgplaw_pha2.fits", 1);
    p1 = load_data ("letgplaw_pha2.fits", 2);
    assign_arf (load_arf ("letgplaw_m1.garf"), m1);
    assign_arf (load_arf ("letgplaw_p1.garf"), p1);
    assign_rmf (load_rmf ("letgplaw_m1.grmf"), m1);
    assign_rmf (load_rmf ("letgplaw_p1.grmf"), p1);
    fit_fun ("phabs(1)*powerlaw(1)");
    set_par (1, 0.001); freeze (1);
    group_data ([m1,p1], 2);
    fit_counts();
    list_par();

 Parameters[Variable] = 3[2]
            Data bins = 8192
           Chi-square = 6987.855
   Reduced chi-square = 0.8532179
phabs(1)*powerlaw(1)
 idx  param             tie-to  freeze         value         min         max
  1  phabs(1).nH            0     1            0.001           0      100000  10^22
  2  powerlaw(1).norm       0     0        0.1011472           0       1e+10
  3  powerlaw(1).PhoIndex   0     0         1.795557          -2           9

Plots of the spectral fit to count histograms are constructed using the rplot_counts command applied to the minus first order and plus first order datasets.

isis> rplot_counts(m1);

isis> rplot_counts(p1);

(Since we are working with grating data, it is most natural to work with wavelength units. The corresponding plots in energy units may be found below.)

The main systematic features in the residuals of the plots are due to the clipping of events by too narrow of a window in order sorting space. This may be seen in the next figure, the so-called order-sorting plot, which gives the number of events as a function of dispersion coordinate (abscissa) and effective order (ordinate). Events with effective orders below 0.7 and above 1.2 were excluded from the count histograms. In this simulation, such events appear mainly at wavelengths less than about 7 angstroms in both orders, and just above 20 angstroms in the positive first order.

For grating data, isis defaults to using wavelength units. We use the plot_unit command to tell isis to use energy units for plotting:

    plot_unit("keV");
    plot_bin_density;
    xrange(0.2,5.0);    % starts ~ 60A
    xlog; ylog;
    group ([m1,p1]; min_chan=2, min_sn=sqrt(10.0));
    eval_counts();