.help extrpsf Feb.96 ftools.caltools
.ih
NAME

extrpsf -->  ExtractS radial PSF (RPSF) dataset from an event file.
.ih
USAGE:

evtfil rpsfil coordsys x_in y_in equinox chanmin chanmax rad_length \
       bkgrnd_pres bck_rad nbins error areawgt area_wgt_option wgt_rad
 
.ih
DESCRIPTION:

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.

.ih
CALCULATION & PROGATION OF ERRORS

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:
.ls ERROR = 'Gauss'
whereby the errors are calculated using SQRT(N).
.le
.ls ERROR = 'POISS-1'(the default) 
whereby the algorithm of Gehrels (1986 ApJ, 303, 336)
eqn 7 with S=1 is used:
.br
.nj
   error = 1.0 + SQRT(N + 0.75)
.br
.ju
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).
.le
.ls ERROR = 'POISS-2'
whereby the algorithm of Gehrels (1986 ApJ, 303, 336) eqn 11,
with S=1 is used.
.br
.nj
   error = SQRT(N - 0.25)
.br
.ju
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.
 
.le
.ls ERROR = 'POISS-3'
whereby the mean of the errors given by POISS-1 and POISS-2 above is
used.
.le
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:
.br
.nj
 N   true 1-sigma        errors calcd using   errors calcd
.br
.nj
     Poisson errors      Gehrels approx       using SQRT(N)
.br
.nj
  0  +1.84   -0.00       +1.87   -0.00        +0.00    -0.00
.br
.nj
  1  +2.30   -0.83       +2.32   -0.67        +1.00    -1.00
.br
.nj
  2  +2.63   -1.92       +2.66   -1.33        +1.41    -1.41
.br
.nj
  3  +2.92   -1.63       +2.94   -1.66        +1.73    -1.73
.br
.nj
  4  +3.16   -1.91       +3.18   -1.94        +2.00    -2.00
.br
.nj
  5  +3.38   -2.16       +3.40   -2.18        +2.24    -2.24
.br
.nj
 10  +4.27   -3.11       +4.28   -3.12        +3.16    -3.16
.br
.nj
 50  +8.12   -7.05       +8.12   -7.05        +7.07    -7.07
.br
.nj
100 +11.00   -9.98      +11.00   -9.99       +10.00   -10.00
.ju


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:
.br
.nj
For example, adding two PHA datasets each with N=5 counts in a given
channel will give a value of 10, and errors of:
.br
.nj
   +4.24   -4.24   (using error='POISS-1')
.br
.nj
   +3.08   -3.08   (using error='POISS-2')
.br
.nj
   +3.95   -3.95   (using error='POISS-3')
.br
compared to:
.br
.nj
   +4.27   -3.11   (statistically correct values)
.br
.nj
   +3.16   -3.16   (propagating 'Gaussian' errors)
.br
.ju

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.

.ih
CALCULATION OF AREA WEIGHTING FACTOR 

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.
.ih
WARNINGS ON USAGE:
.ls None 
.ih
PARAMETERS:
.ls evtfil [character string]
The name of the input event file.
.le 
.ls rpsfil [character string]
The name of the output FITS file containing the Point Spread Function 
for the observational data.
.le
.ls 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'.  
.le
.ls 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.
.le
.ls 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.
.le
.ls 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.
.le
.ls 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.
.le
.ls rad_length [real]
Total Radial length through which the RPSF is requested. This length 
must be entered in unit of arc-min.   
.le
.ls 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. 
.le
.ls 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.
.le
.ls nbins [integer]
Total number of bins user wants to be divide the whole region into. 
This number should be less than 1000.  
.le
.ls error [character string]
Error to be chosen by the user. Three different options of error
calculation is provided. 
.br
.le
.ls 
Gaussian error:
.br
1. error = SQRT(N) 
.br
.le
.ls 
Poisson error:  
.br
2. error = 1.0 + SQRT(N + 0.75)
.br
.le
.ls
3. error = SQRT(N - 0.25)
.br
.le
.ls
4. mean of the errors given by 2. and 3. above is used, i.e.,
          (1.0+SQRT(N+0.75))+(SQRT(N-0.25))      
         -----------------------------------                         
                        2                                            

.le
.ls (properr=true) [boolean]
This flag is for calculation of propagation error which is set to
true by default.
.le
.ls 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.   
.le
.ls 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.
.le
.ls 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 
.le
.ls region_pres [boolean]
Whether calculation to be done with region file included
.le
.ls 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.
.le
.ls (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.
.le
.ls (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.
.le
.ih
BUGS:
None known 
.ih
SEE ALSO

Help on MATHPHA for error calculation.

.ih
LOG OF SIGNIFICANT CHANGES:
.br
v1.0.0 (Feb, 1996) created
.br
.le
v1.4.0 (Feb, 1997)Modified. The maximum of number of bins (nbins) 
is increased from 100 to 1000.  
.ls

.ih
PRIMARY AUTHOR:
.nj
Banashree Mitra Seifert
.br
.nj
HEASARC, NASA/GSFC
.br
.nj
http://heasarc.gsfc.nasa.gov/cgi-bin/ftoolshelp
.ju
.endhelp
