HRC-I Exposure Map and Exposure-corrected Image
CIAO 4.12 Science Threads
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.
To build an exposure map for an HRC-I observation, create an exposure-corrected image, and find an approximation for the source flux.
- Analysis Guide: HRC Imaging
Last Update: 30 Jan 2017 - Reviewed for CIAO 4.9. Added information about setting badpixel file and mask file for edge-blanking.
- Get Started
- Using the fluximage Script
- Step-by-Step Guide
- Calculate the Source Flux
- Data Caveat: "rings" in the exposure-corrected image
Download the sample data: 6202 (HRC-I, M31*)
unix% download_chandra_obsid 6202
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: 13 December 2013 Found repro/hrcf06202_repro_evt2.fits Using event file repro/hrcf06202_repro_evt2.fits Using a monochromatic energy of 1.1 keV. Aspect solution repro/pcadf223264168N003_asol1.fits found. Bad-pixel file repro/hrcf06202_repro_bpix1.fits found. Mask file repro/hrcf06202_001N004_msk1.fits found. Dtf file repro/hrcf06202_001N004_dtf1.fits found. Background file for repro/hrcf06202_repro_evt2.fits found: $CALDB/data/chandra/hrc/bkgrnd/hrciD2005-01-01bkgrndN0001.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) Setting up the HRC-I background for obsid 6202 Creating instrument map for obsid 6202 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 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.
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: 12 September 2012 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
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/pcadf223264168N002_asol1.fits unix% punlearn asphist unix% pset asphist evtfile=hrcf06202_repro_evt2.fits unix% pset asphist dtffile=hrcf06202_001N004_dtf1.fits unix% asphist @hrcf06202_asol1.lis asphist.fits Event List Files (hrcf06202_repro_evt2.fits): Live Time Correction List Files for HRC (hrcf06202_001N004_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: 12 September 2012 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.
unix% set 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
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.
Figure 2: Exposure-corrected Image
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(16380.5,16516.5,119.71632)
The flux can be calculated in several ways:
From the CIAO analysis menu in ds9. Load the data and region file, then run "Analysis → CIAO → Statistics → All (no centroid)".
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.00026337939831 good: 697 null: 264
unix% punlearn dmextract unix% dmextract "m31_all_flux.img[bin firstname.lastname@example.org]" source.flux opt=generic unix% dmlist "source.flux[cols counts]" data,clean # COUNTS 0.0002633793983
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
|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.|