Last modified: 30 Jan 2024

URL: https://cxc.cfa.harvard.edu/ciao/threads/expmap_hrci/

HRC-I Exposure Map and Exposure-corrected Image

CIAO 4.16 Science Threads


Overview

Synopsis:

mkexpmap generates an exposure map which may be used to convert a counts image of a source to an image in flux units. The computed exposure map is essentially an image of the effective area at each sky position, accounting for the effects of dither motion which are especially important near the edges of the detector. The fluximage script automates the creation of an exposure-corrected image for a Chandra observation.

Purpose:

To build an exposure map for an HRC-I observation, create an exposure-corrected image, and find an approximation for the source flux.

Related Links:

Last Update: 30 Jan 2024 - noted csh syntax for setting variable and added bash equivalent.


Contents


Get Started

Download the sample data: 6202 (HRC-I, M31*)

unix% cd 6202
unix% punlearn chandra_repro
unix% chandra_repro mode=h

Using the fluximage Script

How fluximage works

When running fluximage, you are only required to provide an input event file. The script will read the related data product filenames - bad pixels, aspect solution, mask (ACIS), and dead time correction (HRC) - and look for them in the working directory. If the files are in a different location or you wish to be explicit in what files are used, all of the input filenames may be set in the parameter file.

The fluximage script runs the following tools:

  • dmcopy: to create an event image at with the specified binning factor
  • hrc_bkgrnd_lookup, reproject_events, and dmimgcalc: to subtract the particle background (HRC-I only)
  • asphist: to build the aspect histogram(s)
  • mkinstmap: to calculate the instrument map(s) for the center of each energy band
  • mkexpmap: to calculate the exposure map(s) in each energy band
  • dmimgcalc: to combine the exposure maps (multi-chip/plate case only)
  • dmimgthresh: to make a "threshold cut" before dividing the image by the exposure map, removing the hot pixels at the edges (optional)
  • dmimgcalc: to normalize the image by the exposure map

For the multi-chip ACIS or multi-plate HRC-S cases, the tools are run once per chip or plate and are then combined into an image of the full detector.

By default, the intermediate per-chip data products are removed after the script has completed running. To save these (potentially numerous) files, set the cleanup parameter to no.


Run the script

In this example, the input event file is provided, and the supporting data filenames are read from the header.

Since energy is not explicitly resolved in HRC observations, the center of the energy band is chosen at the user's discretion. Here, the full energy range is set in the bands parameter and a center energy of 1.1 keV is chosen. We chose to use a binning scale of 8, which creates pixels just over 1 arcsecond in size, using the xygrid parameter to restrict the analysis to a region around the target (in this case M31).

unix% punlearn fluximage
unix% fluximage repro/ m31 bands=::1.1 bin=8 xygrid=13176.5:20248.5:8,13240.5:18584.5:8
Running fluximage
Version: 04 November 2021

Found repro/hrcf06202_repro_evt2.fits
Using event file repro/hrcf06202_repro_evt2.fits
Using PI=: with a monochromatic energy of 1.1 keV.
Aspect solution repro/pcadf06202_001N001_asol1.fits found.
Bad-pixel file repro/hrcf06202_repro_bpix1.fits found.
Mask file repro/hrcf06202_001N006_msk1.fits found.
Dtf file repro/hrcf06202_001N006_dtf1.fits found.
Background file for repro/hrcf06202_repro_evt2.fits found: /data/chandra_caldb/ciao4.12/data/chandra/hrc/bkgrnd/hrciD2005-01-01bkgrndN0002.fits

The output images will have 884 by 668 pixels,
    and cover x=13176.5:20248.5:8,y=13240.5:18584.5:8.

Running tasks in parallel with 4 processors.
Creating aspect histogram for obsid 6202
# asphist (CIAO): WARNING: skipping 489 livetime correction records (from time: 223263168.067467 to time: 223264168.467512)


Creating instrument map for obsid 6202
Setting up the HRC-I background for obsid 6202
# dmhistory (CIAO): WARNING: Found and corrected "pixlib" library parameters

# dmhistory (CIAO): WARNING: Found and corrected "pixlib" library parameters

Background filter (obsid 6202): [status=xxxxxx00xxxx0xxx00000000x0000000]
Subtracting HRC-I background for obsid 6202
Creating exposure map for obsid 6202
Thresholding data for obsid 6202
Exposure-correcting image for obsid 6202

The following files were created:

 The clipped counts image is:
     m31_all_thresh.img

 The observation FOV is:
     m31_6202.fov

 The clipped exposure map is:
     m31_all_thresh.expmap

 The exposure-corrected image is:
     m31_all_flux.img

You can check the parameter file that was used with dmhistory:

unix% dmhistory m31_all_flux.img fluximage
fluximage infile="repro/" outroot="m31" bands="::1.1" xygrid="13176.5:20248.5:8,13240.5:18584.5:8" 
  binsize="INDEF" asolfile="" badpixfile="" maskfile="" dtffile="" units="default" 
  expmapthresh="1.5%" background="default" parallel="yes" nproc="INDEF" tmpdir="/tmp" cleanup="yes"
  clobber="no" verbose="1" 

Now proceed to the Calculate the Source Flux section.


WARNING about creating a large image

If the chosen binning factor is small enough (e.g. binsize=1 for an HRC-I observation) then the following warning will be seen:

# DMCOPY (CIAO): WARNING: Creating large image: 840 MB. Current max set at 500 MB.
Increase maximum using [opt mem=n] or increase blocking to reduce size.

The output files are completely valid; the warning is just there to let you know that the files that are being created are large. At present there is no way to hide this message.


Step-by-Step Guide

Please ensure that you have set up ardlib to use the bad pixel file for your observation before following this thread; see the Setting the Observation-specific Bad Pixel Files thread for more information.

1. Select a background events file

An HRC-I background events file that has been reprojected to match the data is required for the "background-subtract the image" step. Follow the HRC-I Background Event Files thread to create the background file.

The background event file used in this thread is called background.evt2.


2. Create an Image

First, we need to create the image which will ultimately be normalized by the exposure map. Here we decided to block the image by a factor of 8, restricting attention to a region around M31:

unix% dmcopy "hrcf06202_repro_evt2.fits[sky=box(16712.5,15912.5,7072,5344)][bin sky=::8]" img.bin8
unix% get_sky_limits img.bin8
Running: get_sky_limits
  version: 07 October 2016
Checking binning of image: img.bin8
  Image has 884 x 668 pixels
  Pixel size is 8.0 by 8.0
  Lower left (0.5,0.5) corner is x,y= 13176.5, 13240.5
  Upper right (884.5,668.5) corner is x,y= 20248.5, 18584.5
  DM filter is:
    x=13176.5:20248.5:#884,y=13240.5:18584.5:#668
  mkexpmap xygrid value is:
    13176.5:20248.5:#884,13240.5:18584.5:#668

3. Background-subtract the Image

The final flux map may contain ring-shaped artifacts due to applying a HRMA vignetting correction to the particle background. Subtracting an estimate for the particle background from the observation before exposure correcting eliminates these artifacts.

Bin the HRC-I background events file to create an image that matches the observation:

unix% dmcopy "background.evt2[bin x=13176.5:20248.5:#884,y=13240.5:18584.5:#668]" bgimg.bin8

The dmimgcalc tool is used to subtract the background from the data. The operation (op) parameter is set to use the exposure times to weight the background file:

unix% dmimgcalc img.bin8,bgimg.bin8 none \
      bgsub.bin8 \
      op="imgout=img1-((img1_exposure/img2_exposure)*img2)"
arning: ASPTYPE has different value...Merged...
BTIMDRFT not present in all input files...FAIL...
BTIMNULL not present in all input files...FAIL...
BTIMRATE not present in all input files...FAIL...
warning: DATAMODE has different value...Merged...
warning: DS_IDENT has different value...Merged...
omit - FOC_LEN values different more than 1.000000
OBI_NUM values are different...FAIL...
warning: OBJECT has different value...Merged...
warning: OBSERVER has different value...Merged...
omit - ROLL_NOM values different more than 1.000000
warning: SEQ_NUM has different value...Merged...
warning: TITLE has different value...Merged...

The resulting subtracted image is named bgsub.bin8.


4. Compute Exposure Map

Compute the Aspect Histogram

With the aspect offsets file, we can create a binned histogram, detailing the aspect history of the observation.

In some cases there will be more than one aspect solution file (pcadXXX_asol1.fits) for an observation. All the files must be input, either as a list or as a stack. If you used chandra_repro to re-process the data then it has created a stack file for you, called hrcf<obsid>_asol1.lis, which we use in this case (although as we only have a single aspect solution we could also have just used it directly):

unix% cat hrcf06202_asol1.lis 
/data/ciao/6202/primary/pcadf06202_001N001_asol1.fits

unix% punlearn asphist
unix% pset asphist evtfile=hrcf06202_repro_evt2.fits
unix% pset asphist dtffile=hrcf06202_001N006_dtf1.fits
unix% asphist @hrcf06202_asol1.lis asphist.fits 
Event List Files (hrcf06202_repro_evt2.fits): 
Live Time Correction List Files for HRC (hrcf06202_001N006_dtf1.fits): 
# asphist (CIAO): WARNING: skipping 489 livetime correction records (from time: 223263168.067467 to time: 223264168.467512)

Calculate the Instrument Map

Since the mirror effective area is used to create the instrument map, and that area is energy dependent, it is necessary to decide at what energy to perform the calculation (or whether to use a spectrum as weights). Since energy is not explicitly resolved in HRC observations, the monoenergy parameter is determined at the discretion of the observer (the default value is 1 keV).

Note that it is not necessary for the instrument map to be congruent with the exposure map. We set the pixelgrid parameter to create a 2048 by 2048 pixel image that covers the entire detector area.

Users are reminded to be sure that they have setup the badpixel file in their ardlib.par parameter file. By default, chandra_repro will take care of this but users should double check to make sure that the correct file is being used for the obsevation being analyzed.

unix% plist ardlib | grep HRC-I_BAD
AXAF_HRC-I_BADPIX_FILE = NONE             Enter HRC-I Badpix file

unix% pset ardlib AXAF_HRC-I_BADPIX_FILE=hrcf06202_repro_bpix1.fits
unix% punlearn mkinstmap
unix% pset mkinstmap obsfile=hrcf06202_repro_evt2.fits
unix% pset mkinstmap detsubsys=HRC-I
unix% pset mkinstmap pixelgrid="1:16384:8,1:16384:8"
unix% pset mkinstmap maskfile=hrcf06202_001N004_msk1.fits
unix% mkinstmap instmap.img NONE 1.1
Pixel grid specification x0:x1:#nx,y0:y1:#ny (1:16384:8,1:16384:8): 
Name of fits file + extension with obs info (hrcf06202_repro_evt2.fits): 
Detector Name (HRC-I): 
Grating for zeroth order ARF (NONE|LETG|HETG) (NONE): 
NONE, or name of ACIS window mask file (NONE): 

Be sure to include the maskfile especially for observations where HRC edge-blanking was used.

Calculate the Exposure Map

Now we use mkexpmap and the aspect information stored in the histogram to project the instrument map onto the sky. We set the xygrid parameter to produce a 2048 x 2048 pixel output map; this corresponds to a bin size of 16 for the sky axes, which was used when creating an image from the event list. The get_sky_limits script can be used to easily calculate this information from the existing image:

unix% get_sky_limits bgsub.bin8 
Running: get_sky_limits
  version: 07 October 2016
Checking binning of image: bgsub.bin8
  Image has 884 x 668 pixels
  Pixel size is 8.0 by 8.0
  Lower left (0.5,0.5) corner is x,y= 13176.5, 13240.5
  Upper right (884.5,668.5) corner is x,y= 20248.5, 18584.5
  DM filter is:
    x=13176.5:20248.5:#884,y=13240.5:18584.5:#668
  mkexpmap xygrid value is:
    13176.5:20248.5:#884,13240.5:18584.5:#668

You can then set the xygrid parameter using the information provided by the script, either manually or via:

unix% pset mkexpmap xygrid=")get_sky_limits.xygrid"

(if the latter, do not run get_sky_limits again until after running mkexmap).

If you are computing a low-resolution exposure map and speed is more important than accuracy, set useavgaspect=yes. In doing so, only the average aspect pointing will be used to derive the exposure map; otherwise all points in the aspect histogram will be used. The time required to compute the exposure map is proportional to the number of bins in the aspect histogram; if the aspect histogram contains 100 bins, then the use of this option reduces the run time by a factor of 100, approximately (you may also want to set verbose to 2, since this causes mkexpmap to output percentage-completed information). Using the full aspect solution will help accurately account for chip edges, bad pixels, etc.

(t)csh
unix% set xygrid=`pget get_sky_limits xygrid`

bash/zsh
unix% xygrid=`pget get_sky_limits xygrid`
unix% echo $xygrid
13176.5:20248.5:#884,13240.5:18584.5:#668
unix% punlearn mkexpmap
unix% pset mkexpmap xygrid=$xygrid normalize=no
unix% mkexpmap asphist.fits expmap.img instmap.img mode=h

The exposure map can be displayed in ds9 (Figure 1).

Figure 1: Exposure map

[Thumbnail image: The exposure map of the chip shows variations as one moves outward on the detector.]

[Version: full-size]

[Print media version: The exposure map of the chip shows variations as one moves outward on the detector.]

Figure 1: Exposure map

unix% ds9 expmap.img -scale log -cmap a -zoom 0.5

Since we set the normalize parameter to no, the exposure map has units of [cm2*s*counts/photon]. This allows us to simply divide the image by the exposure map to derive an image in units of flux ([photons/cm2/s/pixel]). If the setting had been left as yes (the default), the units of the exposure map would be [cm2*counts/photon]; see the help file for mkexpmap for more details on this.


5. Normalize the Image by the Exposure Map

The strongly variable exposure near the edge of a dithered field may produce "hot" pixels when divided into an image. While technically proper, these hot pixels can be an eyesore, drawing attention to a noisy, uninteresting portion of the image. The dmimgthresh tool may be used to make a "threshold cut" before dividing the image by the exposure map, thus removing the hot pixels:

unix% punlearn dmimgthresh
unix% dmimgthresh bgsub.bin8 bgsub.thresh.bin8 expfile=expmap.img cut=1.5%
unix% dmimgthresh expmap.img expmap.thresh.img  cut=1.5%

Here we set our threshold at 1.5% of the maximum value of the exposure map. All image pixels with values of exposure less than this value will be set to 0.0 in the output file. You may want to adjust these values for your own observation.

Note: in this example the thresholding is not going to remove any pixels since we have chosen the analysis region such that it does not include the chip edges, and so avoids areas where this might be an issue.

The exposure map is in units of [cm2*s*counts/photon] since it was created by projecting the instrument map (in [cm2*counts/photon]) onto the tangent plane of the observation. To create an image in units of [photon/cm2/s/pixel], we simply need to divide by the exposure map. This can be performed in one step with dmimgcalc:

unix% punlearn dmimgcalc
unix% dmimgcalc bgsub.thresh.bin8 expmap.thresh.img norm.img div
warning: ASPTYPE has different value...Merged...
BTIMDRFT not present in all input files...FAIL...
BTIMNULL not present in all input files...FAIL...
BTIMRATE not present in all input files...FAIL...
warning: CONTENT has 1 different values.
warning: DATAMODE has different value...Merged...
warning: DS_IDENT has different value...Merged...
omit - keyword FOC_LEN not present in all input files
OBI_NUM not present in all input files...FAIL...
warning: OBJECT has different value...Merged...
warning: OBSERVER has different value...Merged...
omit - keyword ROLL_NOM not present in all input files
warning: SEQ_NUM has different value...Merged...
warning: TITLE has different value...Merged...

The messages are related to how the tool merges the header information in the input files. The merging_rules ahelp file explains the rules and how they affect the output file header.

The units of norm.img (Figure 2) are [photon/cm2/s/pixel]; the image is shown in Figure 2.

Figure 2: Exposure-corrected Image

[Thumbnail image: The exposure-corrected image of the detector is displayed in ds9.]

[Version: full-size]

[Print media version: The exposure-corrected image of the detector is displayed in ds9.]

Figure 2: Exposure-corrected Image

unix% ds9 norm.img -scale mode 99.5 -smooth -cmap b

Calculate the Source Flux

Since the units of the exposure-corrected image are [photon/cm2/s/pixel], adding up the pixel values around a source results in the source flux in [photon/cm2/s]. Note that this flux is an approximation - as discussed in An Introduction to Exposure Maps (PS, 12pp) - since a spectral shape was assumed when using mkinstmap (in this example, a monochromatic source).

Using the source region "flux.reg" (which was chosen to contain several sources):

unix% cat core.reg
# Region file format: CIAO version 1.0
circle(16822.2,16447.4,120)

The flux can be calculated in several ways:

  1. From the CIAO analysis menu in ds9. Load the data and region file, then run "Analysis → CIAO → Statistics → All (no centroid)".

  2. with dmstat:

    unix% dmstat "m31_all_flux.img[sky=region(core.reg)]" centroid=no
    PFLUX_IMAGE[/cm**2 /s]
        min:        -3.4538932734e-07             @:        ( 16316.5 16564.5 )
        max:        2.4539862411e-05              @:        ( 16404.5 16484.5 )
       mean:        3.7787575081e-07 
      sigma:        1.5950733808e-06 
        sum:        0.00143972889
       good:        697 
       null:        264 
    
  3. with dmextract:

    unix% punlearn dmextract
    unix% dmextract "m31_all_flux.img[bin sky=@core.reg]" source.flux opt=generic
    unix% dmlist "source.flux[cols counts]" data,clean
    #  COUNTS
          0.00143972889004
    

    Since the input to dmextract was an image, not an event list, the COUNTS column actually reports the total flux (in [photon/cm2/s]) for the source region. While slightly more involved, the dmextract method can be used on multiple sources in a single command, and the results are conveniently stored in a table.

To compute robust source intensity quantities (net counts, source rate, photon flux, energy flux) and the related confidence intervals, use the aprates tool. The Compute Net Counts, Rate, or Flux for Point Sources thread shows how to run aprates.

If you are working with event lists, the eff2evt tool can be used to compute the approximate flux, and calculate the QE and Effective Area for sources. The Calculate the Flux for a Position thread describes how to use this tool.

The srcflux script can be used to compute net counts rates, photon, and energy fluxes automatically.


Data Caveat: "rings" in the exposure-corrected image

If the background subtraction step is skipped, the normalized image contains a pattern of rings. Viewing the image with ds9, using color map i8 makes these rings (Figure 3) most evident. The rings are the result of applying the HRMA vignetting correction in the exposure map to the un-vignetted particle background.

The effect can be mitigated by background-subtracting the data image.

Figure 3: "Rings" in the exposure-corrected image

[Thumbnail image: The image is displayed in ds9 with color map i8.]

[Version: full-size]

[Print media version: The image is displayed in ds9 with color map i8.]

Figure 3: "Rings" in the exposure-corrected image

The rings are the result of applying the HRMA vignetting correction in the exposure map to the un-vignetted particle background.

This image was created using:

unix% fluximage repro/ rings background=none bands=::1.1
...
unix% ds9 rings/all_flux.img -scale log -cmap i8

History

09 Jan 2012 reviewed for CIAO 4.4: added the option of using the CIAO analysis menu in ds9 to calculate the source flux
06 Feb 2012 fluximage updates were released in the 06 Feb 2012 scripts package: HRC-I background is now correctly filtered before subtraction (the application of the STATUS filter was lost in the update of the script for CIAO 4.4).

Additional changes are: setting badpixfile=CALDB uses the bad pixel file from the CALDB rather than the per-observation version; the HRC-I instrument map is now created with the same binning as the exposure map rather than the native HRC pixel scale; the HRC-I background subtraction now only occurs if the HRC background caldb package has been installed.

16 Feb 2012 fluximage updates were released in the 16 Feb 2012 scripts package: setting badpixfile=NONE uses no bad pixel file when creating the instrument, and hence exposure, maps.
15 Oct 2012 The fluximage script has been updated in the 15 Oct 2012 scripts package: changes include an updated parameter file; output file names are different; and support for spectrally-weighted exposure maps. It is now also possible to turn off the background subtraction with the new background parameter. The observation used as an example has been changed to ObsId 6202, one of the (many) observations of M31.
03 Dec 2012 Review for CIAO 4.5; remove mkexpmap chatter
04 Dec 2013 Review for CIAO 4.6; an event file, rather than aspect histogram, should be used for the obsfile parameter of mkinstmap.
17 Dec 2014 Reviewed for CIAO 4.7; minor edits.
30 Jan 2017 Reviewed for CIAO 4.9. Added information about setting badpixel file and mask file for edge-blanking.
18 Jan 2022 Reviewed for CIAO 4.14. Updated for Repro5/CALDB 4.9.6.
30 Jan 2024 noted csh syntax for setting variable and added bash equivalent.