extrpsf -- ExtractS radial PSF (RPSF) dataset from an event file.


        evtfil rpsfil coordsys x_in y_in equinox chanmin chanmax rad_length 
            bkgrnd_pres bck_rad nbins error areawgt area_wgt_option wgt_rad


This task extracts a Radial Point-Spread Function (RPSF) dataset from an event file. An RPSF dataset consists of the number of counts per unit area as a function of radius (using series of concentric annulus). This is calculated using the X & Y (ie RA & dec) columns in the event file centered on a point defined by the user (in either RA & dec or pixel coordinates) for a user-defined number of annuli and outer-radius. If desired, the user can also specify the radius (within the outer radius) beyond which will be used to estimate the number of counts per unit area corresponding to the background. The user can also specify any of several prescriptions to be used to calculate the statistical errors associated with the RPSF data. Another option is possible that user can define region file to exclude region(s). Six shapes are supported at present: CIRCLE, BOX, POLYGON, POINT, ELLIPSE, and ANNULUS.

NOTE: region file is used ONLY to exclude regions from input file.


There are several options open to users regarding the calculation or propagation of errors, controlled by the parameters error and properr.

The error parameter controls which prescription is used to calculate errors, should the task need to do so. The value of this keyword is therefore important if the errors are to be propagated.

The hidden parameter (properr='yes' by default) is used to know if Poissonian statistics should be used to calculate the statistical errors, or if errors are not to be propagated (properr='no').

Currently, the following prescriptions are available. If N is the number of counts observed in a given radial bin, then:

ERROR = 'Gauss'
whereby the errors are calculated using SQRT(N).

ERROR = 'POISS-1'(the default)
whereby the algorithm of Gehrels (1986 ApJ, 303, 336) eqn 7 with S=1 is used:

                error = 1.0 + SQRT(N + 0.75)

The value is statistically that of the (larger) +ve error of a Poissonian distribution, but within this task this value is assigned to both the +ve and -ve error on the counts in that channel. For small N, the errors created using this prescription is significantly GREATER (and hence more conservative) than both the true -ve error, and that obtained by simply using SQRT(N). Thus, this prescription is recommended unless a user fully understands the implications of using a different prescription. Both differences quickly reduce as one moves to larger N (and the Poissonian distribution becomes more symmetric/Gaussian).

whereby the algorithm of Gehrels (1986 ApJ, 303, 336) eqn 11, with S=1 is used.

                error = SQRT(N - 0.25)

The value is statistically that of the (smaller) -ve error of a Poissonian distribution, but within this task this value is assigned to both the +ve and -ve error on the counts in that channel. This error prescription UNDERESTIMATES the error for small N.

whereby the mean of the errors given by POISS-1 and POISS-2 above is used.

Caution is urged, particularly when using ERRMETH = 'Gauss', as unexpected and/or misleading results can be produced. See OGIP/95-008 for further details.

The following tables enables direct comparisons to be made between these approximations:

      N   true 1-sigma        errors calcd using   errors calcd
          Poisson errors      Gehrels approx       using SQRT(N)
       0  +1.84   -0.00       +1.87   -0.00        +0.00    -0.00
       1  +2.30   -0.83       +2.32   -0.67        +1.00    -1.00
       2  +2.63   -1.92       +2.66   -1.33        +1.41    -1.41
       3  +2.92   -1.63       +2.94   -1.66        +1.73    -1.73
       4  +3.16   -1.91       +3.18   -1.94        +2.00    -2.00
       5  +3.38   -2.16       +3.40   -2.18        +2.24    -2.24
      10  +4.27   -3.11       +4.28   -3.12        +3.16    -3.16
      50  +8.12   -7.05       +8.12   -7.05        +7.07    -7.07
     100 +11.00   -9.98      +11.00   -9.99       +10.00   -10.00

The parameter properr (hidden) controls whether the errors are to be propagated during the algebra or (if properr='no') whether the errors are simply calculated from the resultant PHA dataset. In the former case, the errors are propagated in the normal manner (err3 = SQRT[err1**2 + err2**2]). This is an approximation when in the Poissonian regime (ie for low N). Whilst it is true that the variances of the Poissonian distribution are combined in this normal way, confidence limit error bars are not simply related to the variance like they are for Gaussian statistics.

HOWEVER, these approximations work well for all but the smallest N, and is certainly superior to either assuming SQRT(N) errors, or neglecting errors altogther:

For example, adding two PHA datasets each with N=5 counts in a given channel will give a value of 10, and errors of:

        +4.24   -4.24   (using error='POISS-1')
        +3.08   -3.08   (using error='POISS-2')
        +3.95   -3.95   (using error='POISS-3')

compared to:

        +4.27   -3.11   (statistically correct values)
        +3.16   -3.16   (propagating 'Gaussian' errors)

Severely misleading and/or incorrect results are possible if the errors are not propagated (ie if properr=no). Turning off error propagation is only reccomended when one is adding non background-subtracted datasets, and one fully understands the risks.


The radial bin to which each event is assigned is determined by the annulus in which the center of that pixel lies (ie events are NOT "split" between annuli according to the fraction of the pixel area in each). This results in an "active" area in a given radial bin running from R_lo to R_hi slightly different from the pi*(R_hi^2 - R_lo^2) expected from purely geometrical considerations. This difference is clearly greatest in the innermost radial bins. The so-called "Area weighting factor", defined as the ratio of the active area within each annulus to the geometrical area, is therefore calculated within the task, used in the calculation of the RPSF (counts per unit area) and written to the output file. Whilst strictly the area weighting factor should be calculated for all radial bins, its determination can be rather CPU-intensive under many circumstances. However, assuming that the radial bins are larger than the pixel size, the area weighting factor tends to unity towards the outer radial bins. Thus the run-time of the task can be greatly reduced by setting area_wgt_option=2, calculating the area weighting factor within a user-defined radius wgt_rad, and assuming it is unity at larger radii. This option is only suggested for "quick-looks" at the RPSF profile, but it is recommended an RPSF is re-extracted with area_wgt_option=1 should the RPSF warrant further detailed analysis.


See above


evtfil [character string]
The name of the input event file.

rpsfil [character string]
The name of the output FITS file containing the Point Spread Function for the observational data.

coordsys [character string]
The name of the coordinate system the input centroid is supposed to be. If this is in pixel coordinate then it is 'IMAGE' system and if the input is in RA and DEC (fractional degrees) coordinate system is 'CELESTIAL'.

x_in [real]
The centroid of the event. If 'IMAGE', then it is the x-coordinate of the pixel, and if 'CELESTIAL' it is the Right Ascension (in fraction of degrees) of the event.

y_in [real]
The centroid of the event. If 'IMAGE', then it is the y-coordinate of the pixel, and if 'CELESTIAL' it is the Declination (in fraction of degrees) of the event.

chanmin [integer]
Minimum PI channel number. This is required for the ROSAT PSPC and gives the lowest PI channel to be used to construct the radial profile.

chanmax [integer]
Maximum PI channel number. This is required for the ROSAT PSPC and gives the highest PI channel to be used to construct the radial profile.

rad_length [real]
Total Radial length through which the RPSF is requested. This length must be entered in unit of arc-min.

bkgrnd_pres [boolean]
This parameter determines whether the user would like to calculate back-ground counts. If YES/yes then it asks for inner radius for background calculation.

bck_rad [real]
Background radius in arcmin. The program will use the total number of counts in a circular annulus from bck_rad to rad_length as an estimate of the number of background counts.

nbins [integer]
Total number of bins user wants to be divide the whole region into. This number should be less than 1000.

error [character string]
Error to be chosen by the user. Three different options of error calculation is provided.

Gaussian error:
1. error = SQRT(N)

Poisson error:
2. error = 1.0 + SQRT(N + 0.75)

3. error = SQRT(N - 0.25)

4. mean of the errors given by 2. and 3. above is used, i.e.,


(properr=true) [boolean]
This flag is for calculation of propagation error which is set to true by default.

areawgt [boolean]
This flag is for calculation of the area weighting factor. If this is set to YES/yes, the area weighting factor will be calculated and otherwise skip this calculation.

area_wgt_option [string]
Two options are provided for calculation of area weighting factor. area weighting for the entire region(option=1) or upto a certain distance from the center of the event region only (option=2). In case of option=2, it then asks for the outer radius for this calculation.

wgt_rad [integer]
If area weighting is wanted upto a certain distance from the center of a event, then it asks for the radial distance (in arc-min) through which it is required to be calculated

region_pres [boolean]
Whether calculation to be done with region file included

regfil [string]
If region_pres is YES/yes, then is asked for the name of the ascii region file. This is region file extracted from saoimage. The region file is used to define EXCLUDED region, if any. It is immaterial of the sign(-,+ or !) preceeding the region. Region file is used for only excluding the specified region.

(chatter=9) [integer]
Flag to set the chattyness at execution. Default value is set to 9, which gives the user the task version and few warnings. Lower/higher values produces quieter/verbose output on the screen.

(clobber = false) [boolean]
Flag specifying whether or not a pre-existing file with the same name as that requested as the output file from this task will be overwritten.


None known


Help on MATHPHA for error calculation.


v1.0.0 (Feb, 1996) created

v1.4.0 (Feb, 1997)Modified. The maximum of number of bins (nbins) is increased from 100 to 1000.


          Banashree Mitra Seifert


Feb.96 ftools.caltools