Search https://cxc.harvard.edu/ciao/
Contact the CXC HelpDesk
Last modified: 22 Feb 2022

URL: https://cxc.cfa.harvard.edu/ciao/threads/psf_image/
/ ciao / threads / psf_image /
Return to Threads Page: Top | All | psf

Creating an Image of the PSF

CIAO 4.16 Science Threads

Overview

Synopsis:

The output of MARX is a pseudo event file containing a list of rays that have been projected onto the detector (recall that MARX includes the instrumental response when doing this). This thread shows how to bin the output file into an image that (optionally) matches that of your data. This image can then be used in your scientific analysis.

It is important to read the caveats before using the ChaRT PSF simulations in analysis.

Related Links:

Last Update: 22 Feb 2022 - Review for CIAO 4.14. Updated for Repro5/CALDB 4.9.6

Contents

Get Started

This example uses marx_output_i0000.fits, which was created in the Using MARX to Create an Event File thread. Also:

Download the sample data: 942 (ACIS-S, NGC 4244)

unix% download_chandra_obsid 942

Examine the Data

The following sections describe how to create an image of the PSF that is either centered on the PSF or matches an existing image of the dataset.

The projected rayfile (marx_output_i0000.fits) can be examined with dmlist: 

unix% dmlist marx_output_i0000.fits blocks
 
--------------------------------------------------------------------------------
Dataset: marx_output_i0000.fits
--------------------------------------------------------------------------------
 
     Block Name                          Type         Dimensions
--------------------------------------------------------------------------------
Block    1: NULL                           Null        
Block    2: EVENTS                         Table        22 cols x 1563     rows
Block    3: GTI4                           Table         2 cols x 1        rows
Block    4: GTI5                           Table         2 cols x 1        rows
Block    5: GTI6                           Table         2 cols x 1        rows
Block    6: GTI7                           Table         2 cols x 1        rows
Block    7: GTI8                           Table         2 cols x 1        rows
Block    8: GTI9                           Table         2 cols x 1        rows
Block    9: MARX_PAR                       Table         1 cols x 522      rows


unix% dmlist marx_output_i0000.fits cols
 
--------------------------------------------------------------------------------
Columns for Table Block EVENTS
--------------------------------------------------------------------------------
 
ColNo  Name                 Unit        Type             Range
   1   TIME                 s            Real8          362912400.0:362962714.1465607882 time since observation start
   2   CCD_ID                            Int2           0:9                  CCD id number
   3   NODE_ID                           Int2           -                    0-4
   4   EXPNO                             Int4           0:2147483647         Exposure number
   5   chip(CHIPX,CHIPY)    pixel        Int2           2:1023               CHIP X
   6   tdet(TDETX,TDETY)    pixel        Int4           2:8191               Detector X
   7   det(DETX,DETY)       pixel        Real4          0.50:     8192.50    Focal Plane X
   8   sky(X,Y)             pixel        Real4          0.50:     8192.50    sky X pixel
   9   PHA                  adu          Int4           0:36855              Total PHA for event
  10   ENERGY               eV           Real4          0:  1000000.0        Nominal energy of event
  11   PI                   Chan         Int4           -                    pulse invariant energy of event
  12   FLTGRADE                          Int2           0:255                Event Grade Code
  13   GRADE                             Int2           -                    ACIS grade code
  14   STATUS[4]                         Bit(4)                              status flags
  15   SHELL                             Int2           -                    Mirror Shell (0=1,1=3,2=4,3=6)
  16   ZCOS                              Real4          -Inf:+Inf            Z direction cosine
  17   YCOS                              Real4          -Inf:+Inf            Y direction cosine
  18   XCOS                              Real4          -Inf:+Inf            X direction cosine
  19   ZPOS                 mm           Real4          -Inf:+Inf            Z position in MARX coords
  20   YPOS                 mm           Real4          -Inf:+Inf            Y position in MARX coord system
  21   XPOS                 mm           Real4          -Inf:+Inf            X position in MARX coord system
  22   MARX_ENERGY          keV          Real4          -Inf:+Inf            Energy of ray
 
--------------------------------------------------------------------------------
World Coord Transforms for Columns in Table Block EVENTS
--------------------------------------------------------------------------------
 
ColNo    Name
8:    EQPOS(RA ) = (+184.3419) +TAN[(-0.000136667)* (sky(X)-(+4096.50))]
           (DEC)   (+37.7819 )      (+0.000136667)  (   (Y) (+4096.50))

Note that the table contains coordinates for the events in a number of coordinate systems. We shall assume that we want an image of the PSF in SKY coordinates.

Centering the PSF

Get the coordinates

The center of the PSF can be found by using the dmcoords tool to convert the input RA and Dec values used to create the PSF into SKY coordinates. Note that "centering" may be a bit of a misnomer; far off-axis, it is likely that the "center" of the PSF is not coincident with either the centroid and/or the sky pixel matching the off-axis angle.

The necessary values are recorded in the header of the ChaRT output (HRMA_ra184.23746_dec37.72658_source_flux_chart.dat_dithered_i0000_rays.fits): 

unix% dmlist HRMA_ra184.23746_dec37.72658_source_flux_chart.dat_dithered_i0000_rays.fits header | egrep 'SRC_(RA|DEC)'
0010 SRC_RA               184.2374583333333 [deg]        String       Input source Right Ascension
0011 SRC_DEC              37.72658055555556 [deg]        String       Input source Declination

dmcoords can now be used to convert between this celestial position and physical sky pixels 

unix% punlearn dmcoords
unix% dmcoords acisf00942_repro_evt2.fits op=cel ra=184.237458333 dec=37.7265805556 celfmt=deg
unix% pget dmcoords x y
4700.906879836509
3692.224597519996

Calculate the binning factor

We choose to create a 512 x 512 pixel image around the central position (with a binning factor of 1) in order to encompass most of the PSF signal. To calculate a 512 x 512 region centered at (4700.9,3692.2) 

512/2 = 256

X direction: 4700.9-256 = 4444.9
             4700.0+256 = 4956.9

Y direction: 3692.2-256 = 3436.2
             3692.2+256 = 3948.2

The binning expression then is [bin x=4444.9:4956.9:1,y=3436.2:3948.2:1]

Create the image

Now we can use dmcopy with the binning syntax to create an image: 

unix% dmcopy "marx_output_i0000.fits[bin x=4444.9:4956.9:1,y=3436.2:3948.2:1]" psf_512_i0000.fits

Figure 1: Image of the PSF

[Thumbnail image: 512 by 512 pixel image of PSF]

[Version: full-size]

[Print media version: 512 by 512 pixel image of PSF]

Figure 1: Image of the PSF

The PSF image is 512 x 512 pixels.

Binning the Rayfile to Match an Image

An alternative approach is to bin the PSF image to match an existing image of your data. The image of the data in this example was made by creating an event file of chip S2 (ccd_id=6) and simultaneously binning it into an image: 

unix% apply_fov_limits "acisf00942_repro_evt2.fits[ccd_id=6,energy=400:6000]" ccd6_bin1.fits bin=1
Running: apply_fov_limits
  version: 01 April 2020
Observation: ObsId 942 - ACIS-235678
Using ccd_id=6 from ccd6_evt.fits[energy=400:6000]
Calculating FOV file using:
  Aspect solution pcadf00942_001N001_asol1.fits
  Mask file       acisf00942_001N005_msk1.fits

The output image will have 1482 by 1481 pixels, pixel size of 0.492 arcsec,
    and cover x=3814.5:5296.5:1,y=2760.5:4241.5:1.

Created: ccd6_bin1.fits

Figure 2: Image of the event data

[Thumbnail image: Image of the event data (ObsID 942)]

[Version: full-size]

[Print media version: Image of the event data (ObsID 942)]

Figure 2: Image of the event data

The event file for ObsID 942 is filtered on ccd_id=6, energy from 0.4 to 6.0 keV, and sky filtered using the field-of-view (FOV) file. The data are then binned by a factor of 1 to create an image.

The get_sky_limits script (part of the CIAO Scripts distribution) can be used to find out the necessary DM binning specification to match the PSF to this image: 

unix% get_sky_limits ccd6_bin1.fits
Running: get_sky_limits
  version: 07 October 2016
Checking binning of image: ccd6_bin1.fits
  Image has 1482 x 1481 pixels
  Pixel size is 1.0 by 1.0
  Lower left (0.5,0.5) corner is x,y= 3814.5, 2760.5
  Upper right (1482.5,1481.5) corner is x,y= 5296.5, 4241.5
  DM filter is:
    x=3814.5:5296.5:#1482,y=2760.5:4241.5:#1481
  mkexpmap xygrid value is:
    3814.5:5296.5:#1482,2760.5:4241.5:#1481

Finally, use dmcopy to create the image: 

unix% dmcopy "marx_output_i0000.fits[bin x=3814.5:5296.5:#1482,y=2760.5:4241.5:#1481]" psf_match_i0000.fits

Displaying the two images side-by-side: 

unix% ds9 psf_match_i0001.fits ccd6_bin1.fits &

Figure 3: Image of the PSF matched to the data

[Thumbnail image: Image of the PSF matched to the data in ds9]

[Version: full-size]

[Print media version: Image of the PSF matched to the data in ds9]

Figure 3: Image of the PSF matched to the data

The PSF image is on the left and the event data is on the right. It is easy to see that the PSF lines up with the chosen source (see the Preparing to Run ChaRT thread).

Create a Sub-pixelated PSF Image

To create the sub-pixelated image with dmcopy, use a similar filter to the one given in the Centering the PSF section. 

unix% dmcopy "marx_output_i0000.fits[bin x=3820.5:5302.5:0.1,y=2768.5:4248.5:0.1]" psf_subpix.fits

Figure 4: Sub-pixelated PSF image created with dmcopy

[Thumbnail image: Sub-pixelated PSF image created with dmcopy]

[Version: full-size]

[Print media version: Sub-pixelated PSF image created with dmcopy]

Figure 4: Sub-pixelated PSF image created with dmcopy

The image of the PSF created at a binning of 0.1. The green square in the lower right hand size is 1 full resolution pixel (0.492" x 0.492").

The MARX output can also be binned at subpixel resolution directly in ds9. Since the projected event file contains X and Y values without any pixel quantization, it can be binned into an image at any desired sub-pixel scale. The commands to bin and recenter the event file can be specified on the command line: 

unix% ds9 marx_output_i0000.fits -pan to 4704.01 3699.93 -bin factor 0.1

or by going through the Binning menu. The resulting image looks similar to Figure 4.

Combine multiple realizations

To combine the multiple ChaRT realizations, they much each be individually projected using MARX, made into events files. Those event files can then be merged and binned into an image. The tcsh and bash shell commands to loop over the remaining iterations (i0001,i0002,i0003,i0004) are shown below. 

tcsh shell:

unix% foreach d ( i0001 i0002 i0003 i0004 )
  marx @@./marx OutputDir=marx_output_${d}.dir SAOSACFile=HRMA_ra184.23746_dec37.72658_source_flux_chart.dat_dithered_${d}_rays.fits
  marx2fits --pixadj=EDSER marx_output_${d}.dir marx_output_${d}.fits
end

bash shell:

unix$ for d in i0001 i0002 i0003 i0004 
do
  marx @@./marx OutputDir=marx_output_${d}.dir SAOSACFile=HRMA_ra184.23746_dec37.72658_source_flux_chart.dat_dithered_${d}_rays.fits
  marx2fits --pixadj=EDSER marx_output_${d}.dir marx_output_${d}.fits
done

The marx.par parameters were already pset in the Using MARX to Create an Event File from ChaRT Rays thread, therefore only the OutputDir and SAOSACFile names need to be specified here.

Finally the iterations are combined together into a single pseudo event file using dmmerge 

unix% dmmerge "marx_output_i00*.fits" marx_output.fits
unix% set cts = `dmlist marx_output.fits counts`
unix% dmimgcalc "marx_output.fits[bin x=3821.5:5302.5:#1481,y=2768.5:4248.5:#1480]" none norm_psf.fits op="imgout=img1/((float)${cts})" clob+

Or dmcopy can be used to create an unnormalized image 

unix% dmcopy "marx_output.fits[bin x=3821.5:5302.5:#1481,y=2768.5:4248.5:#1480]" cts_psf.fits clob+

The final combined PSF from the 5 realizations is shown in Figure 5.

Figure 5: PSF combined from multiple simulations

[Thumbnail image: PSF combined from multiple simulations]

[Version: full-size]

[Print media version: PSF combined from multiple simulations]

Figure 5: PSF combined from multiple simulations

Users can check the number of events in the final combined PSF image using dmlist 

unix% dmlist "marx_output.fits" counts
15658

showing a total of 10803 events in this example.

Summary

The image of the PSF can now be used with CIAO.

History

27 Jun 2003 original version, updated for CIAO 3.0: layout
16 Feb 2005 reviewed for CIAO 3.3: no changes
18 Aug 2008 updated for CIAO 4.0: version N003 event file, minor changes to screen output; converted images to inline
17 Feb 2010 reviewed for CIAO 4.2: minor updates to screen output
15 Dec 2010 reviewed for CIAO 4.3: no changes
29 Mar 2013 updated file names from N003 to N004
11 May 2015 Updated with ChaRT v2 example. Added a section on combining mulitple iterations into an averaged PSF.
22 Feb 2022 Review for CIAO 4.14. Updated for Repro5/CALDB 4.9.6
Return to Threads Page: Top | All | psf
/ ciao / threads / psf_image /