PAPFIP: Python and Pgplot Fits Image Plotter

What is PAPFIP?

PAPFIP is a FITS image plotter based on the python bindings to the pgplot graphics library. PAPFIP is an acronym for "Python and PGPLOT-based FITS image plotter". It is intended as a non-interactive (command-line or script-driven) astronomical image plotter intended for the automated or repeated production of finalized publication-quality postscipt files, rather than a method of inspecting or visualizing the data.

Obtaining and download PAPFIP

The current version of PAPFIP (version 1.1, released 2008-10-24) can be downloaded here. Create a directory for it, untar it, and read the INSTALL document for instructions on how to install and run PAPFIP. PAPFIP's software requirments are detailed below.

Features and limitations

PAPFIP is still undergoing development and its features continue to expand and evolve. The following is a reasonably up-to-date and accurate description of its capabilities, but looking at the examples will also give you a better idea of what PAPFIP can do. Viewing the command line options is also informative.

• Plot a FITS primary image in
  • greyscale using a linear, square-root, asinh or log10 intensity scaling, with user-defined image minima and maxima. (-g or --plotgrey, --scale, --min, --max, --asinh options).
  • and/or in contour form, using a comma-separated list of contour levels (-c or --plotcontour, -l or --levels options).
  • Contour levels may be labelled with their numerical values with a user-defined format (--contourlabel, --contourlab* options).
• Over-plotting of sources or regions within a ds9-format region file in phys or RA/DEC (sexagesimal or decimal degrees) coords (-r or --regions options). A non-official WCS format may also be used.
  • Multiple region files (of the same format) may be used simply by supplying a comma-separated list of files to the --regions option.
  • Supported region types are circle(), box(), ellipse(), polygon(), cross, x and line style markers.
  • The line style, line width and line color of the shapes/symbols from the region files can be controlled (--regionsformat option).
  • Included/excluded regions are plotted with different line styles (-e or --excludestyle).
• Image overlays are fully supported.
  • This overlays a second fits image in contour form over the first image, where the first image is presented using style defined with the grey-scale image options. The over-lay is handled using WCS coordinate keywords that must be present in the headers of each fits file. For astronomical images with axes of right ascension verse declination the two images must use the same RADECSYS and be north-aligned.
• The field of view, i.e. axis minima and maxima, may be user defined. By default the region fully covering the first input fits image is plotted (-y or --axisregion).
  • Axes are plotted either in terms of raw pixel number (-x or --axes=pixels), WCS coordinates (--axes=wcs), or units of angular distance (arcseconds, arminutes or degrees) offset from either the exact pixel center of the image or a user-defined "central" RA,DEC or WCS coordinate. (--axes=offset, -C or --center ).
  • User-defined x and y axis labels, and optionally a data label, may be specified should the user desire it. By default PAPFIP makes use of the axis information in the fits file to construct axis labels (--datalabel, --xlabel, --ylabel).
  • The plot region typically has the same numerical aspect ratio as the data, but we can force the axis box to be exactly square (--windowscale=noscale).
• An optional scale bar can be plotted (-b or --bar).

The following is a list of what PAPFIP can not currently do, but which are features I intend to add or fix in future versions.

• The latest version of ds9 (version 4.0 and above) produces region files that have a different format to versions 3.x. PAPFIP currently only supports ds9 format region files from version 3.x and lower.
• Color plots with user-specified color tables, in addition to the current grey-scale plotting option.
• The intensity-scale bar currently labels the value in the scaled units, not in the original units (for example if the maximum value is 100 and you choose log scaling, the scale bar shows the maximum as 2.0). It would be nice to fix this.
• Coordinate axes with tick marks and labeled RA and DEC coordinates. The tricky part is dealing with cases with declinations near the poles.
• The ability to over-lay and plot multiple files (where N>2), extending the current ability to over-plot a second image in contour form over the grey-scale first image.
• 2-D plots of arbitrary slices out of 3-D fits files.
• 2-D plots of arbitrary slices out of 3-D HDF4 files.
• It uses Numeric rather than numpy. The python pgplot bindings PAPFIP uses have not been updated to use numpy.

Command line options

Run papfip with "--help" appended to the command line will show you the full set of command line options available.

ryo-ohki ~/Fun/Python/Papfip-1.1 13:46>python papfip1.py --help
Usage: papfip1.py [options] -i file -d device -g or -c

Options:
  -h, --help            show this help message and exit
  -i INPUT, --input=INPUT
                        (required) INPUT fits file to plot. For plotting
                        overlays input a comma-separated list of two file
                        names.
  -d DEVICE, --device=DEVICE
                        (required) PGPLOT plotting device.
  -g, --plotgrey        (required) Plot greyscale image. Grey-scale and/or
                        contour plotting must be defined. For overlays the
                        first file is plotted in grey-scale and the second
                        with contours.
  -l L1,L2..,LN, --levels=L1,L2..,LN
                        Comma separated list of contour levels to plot.
                        Required if contour option selected.
  -c, --plotcontour     Plot image in contour form with user-specified levels.
  -r REGIONFILES, --regions=REGIONFILES
                        Comma separated list of ds9 region files
  -e ESTYLE, --excludestyle=ESTYLE
                        Excluded region line style.
  -s SCALE, --scale=SCALE
                        Grey-scale image scaling (lin|sqrt|log|asinh). A
                        scaling value is required for asinh scaling.
  -m SCALEPAR1, --min=SCALEPAR1
                        Scale parameter, min data value to plot
  -n SCALEPAR2, --max=SCALEPAR2
                        Scale parameter, max data value to plot.
  -a SCALEPAR3, --asinh=SCALEPAR3
                        Scale parameter for asinh plotting. The data will be
                        divided by this value before plotting.
  -b, --bar             Plot a scale bar on the RHS of the image.
  -x AXES, --axes=AXES  Axis plotting options (offset|pixels|wcs).
  -C CENTER, --center=CENTER
                        Offset axes from this RA,DEC position (decimal
                        degrees).
  -u UNITS, --axisunits=UNITS
                        Axis units.
  -y XMIN,XMAX,YMIN,YMAX, --axisregion=XMIN,XMAX,YMIN,YMAX
                        Comma separated list of xmin,xmax,ymin,ymax for
                        plotting a specific axis region rather than the entire
                        image. Units assumed are the same as the axisunits
                        specified.
  --datalabel=DATALABEL
                        Data label, e.g. 'T (K)'
  --contourlabel        Labels contour levels.
  --contlabspacing=INT  Scale parameter, min data value to plot
  --contlabmincell=INT  Scale parameter, min data value to plot
  --contlabformat=FORMAT
                        C-style format statement.
  --xlabel=LABEL        For example '\gb' gives Greek beta
  --ylabel=LABEL        For example 'log\d10\u v\d7\u (km/s)'
  --regionsformat=LS,LW,LC
                        Comma-separted integer line style, width and color
  --windowscale=scale|noscale
                        Scale axis window to follow data or force square axes?

Examples of use and screenshots

The files necessary to reproduce the examples shown here are provided as part of the PAPFIP distribution, and are in the Examples/ directory along with simple csh scripts that run the examples.

The images below are from examples 3 and 4 respectively.

Example 1: Running papfip without a file to get the command line usage and options list. See the command line options shown above.

The next two examples of PAPFIP use two of astronomical data sets provided as fits files in the Examples/ directory of the PAPFIP distrbution (nucl_m82_dks_v4_im_e0320_b1.fits and nucl_m82_dks_v4_im_e0320_b1_sm3as.fits). These are X-ray images (in the E=0.3-2.0 keV energy band) of the nucleus of M82 as detected by the Chandra X-ray Observatory from Sequence 600270 / ObsID 2933. The examples shown below uses the command line, but if you're in the Examples/ directory you can just execute "csh example2.csh" to achieve the same result.

Example 2: A simple grey-scale plot of the fits file to an X-window, adopting the default values for all other parameters. This is really only useful for inspecting the image. The axes are displayed in the default form, which is offset in arcseconds from the center of the image. This example assumes you are in the Examples/ directory of the PAPFIP distribution.

python ../papfip1.py \
     -i nucl_m82_dks_v4_im_e0320_b1_sm3as.fits -d /xserve -g
  readImg: Stripping NaN values from data.
  Axis plotting will show offset from image center.
  Offset will be in arcsec.
  Good. Image has a radecsys keyword> 'ICRS    ' <
  True data minimum and maximum values:   0.00025014980929 15.6706542969
  Scaled data minimum and maximum values: 0.00025014980929 15.6706542969
  Scaled data min and max values: 0.00025014980929 15.6706542969
Finished.

Example 3: Using the same image as before, we now apply an asinh intensity scaling between the limits of 0.1 and 10 counts/pixel, with the pivot point at 3 couts/pixel. A intensity scale bar is shown at the top of the plot. In addition to the grey-scale display, contours of equal intensity are plotted at levels of 0.5, 1, 2 and 3 counts/pixel. The axes are shown in units of arcminutes offset from the center of the image. The contents of a ds9-format region file (technically a ds9 version 3 region file) containing the location of the detected soft X-ray point sources is shown in green. Finally, instead of plotting to an X-window the plot is drawn directly to a color postsript file with the file name example3.ps (PGPLOT device name "example3.ps/vcps").

python papfip1.py -i nucl_m82_dks_v4_im_e0320_b1_sm3as.fits \
     -d example3.ps/vcps \
     -g -l 0.5,1,2,3 -s asinh -m 0.1 -n 10 -a 3 \
     -r accept_srclist_src_m82_dks_s3_im_e0320_b2_radec.reg \
     -b -u arcmin
  readImg: Stripping NaN values from data.
  Axis plotting will show offset from image center.
  Offset will be in arcmin.
  Good. Image has a radecsys keyword> 'ICRS    ' <
  scaleData: asinh new data min and max  0.165904550269 1.9188964721
  True data minimum and maximum values:   0.00025014980929 15.6706542969
  Scaled data minimum and maximum values: 0.165904550269 1.9188964721
  Scaled data min and max values: 0.165904550269 1.9188964721
[ 0.5  1.   2.   3. ]
  Preparing to plot regions on image.
    Regions plotted will note offset from image center.
  Good. Image has a radecsys keyword> 'ICRS    ' <
Finished.

Example 4: Finally an example overlaying two different images, both of the starburst galaxy NGC 6810, with the first (optical H-alpha emission) shown in grey-scale and the second (soft X-ray emission) overlaid in contour form. A region file is used to draw an ellipse showing the optical extent of the galaxy (the D25 ellipse). In addition to also using many of the same options used in example 3, here the coordinate system is centered using the -C flag on the user-specified RA,DEC pair corresponding to the center of the galaxy (instead of the center of the fits image), and the axis region shown is a subset of the full image (using the -y flag). A postscript version of the resulting image is plot_n6810_halphadeep.eps.

python ../papfip1.py -g -i n6810_halpha_goodwcs_sm2as.fits,n6810_v2_m12_bscim_e0311_b50_asm6.fits -d /xserve -s asinh -m 2 -n 100 -a 10 -l 2.5e-6,5.0e-6,10.0e-6,20.0e-6,40.0e-6,80.0e-6,160.0e-6 -r reg_n6810_d25i.reg -y -79,79,-79,79 -C 19:43:34.4,-58:39:20.6
  readImg: Stripping NaN values from data.
  Axis plotting will show offset from image center.
  Offset will be in arcsec.
  Good. Image has a radecsys keyword> 'FK5     ' <
  stringRaDecToRaDec: RA: 295.893333333   Dec: -58.6557222222
  scaleData: asinh new data min and max  0.198690110349 2.9982229503
  True data minimum and maximum values:   -154.741409302 7209.13378906
  Scaled plotting minimum and maximum values: 0.198690110349 2.9982229503
  Scaled data min and max values: 0.198690110349 2.9982229503
  User-specified axis region selected.
  Loading second data file for overlay: n6810_v2_m12_bscim_e0311_b50_asm6.fits
  readImg: Stripping NaN values from data.
  Axis plotting will show offset from image center.
  Offset will be in arcsec.
  Good. Image has a radecsys keyword> 'FK5     ' <
  stringRaDecToRaDec: RA: 295.893333333   Dec: -58.6557222222
[  2.49999994e-06   4.99999987e-06   9.99999975e-06   1.99999995e-05
        3.99999990e-05   7.99999980e-05   1.59999996e-04]
  Preparing to plot regions on image.
    Regions plotted will note offset from image center.
  Good. Image has a radecsys keyword> 'FK5     ' <
  stringRaDecToRaDec: RA: 295.893333333   Dec: -58.6557222222
Finished.

Example 5: This example, while not visually stunning in itself, highlights some of the new functionality in version 1.1 of PAPFIP, specifically dealing with (a) non-astromical fits images that nonetheless have axis coordinates specified using WCS coordinate keywords, and (b) over-plotting of polygonal regions using a modified version of the ds9 region file syntax. The --regionsformat option is used to set the region file plotting line style (2=dashed), line width (3) and line color (1=black). The --windowscale=none option is used to force the plot to appear square even though numerically its axis ratio is 0.9 by 3.6. The default axis labeling is replaced with user-defined greek symbols using the --xlabel and --ylabel options.

python ../papfip1.py -i 2d_modb_combined_matches.fits -g -d /xserve --axes wcs --min=0.0 --max=3 --xlabel=\ge --ylabel=\gb --regions=reg_2d_scaledcc_fine_combined_matches.reg --regionsformat=2,3,1 --windowscale=none
  readImg: Stripping NaN values from data.
  Axis plotting will use fits header WCS coordinates.
  headerInfo: Image header keyword size and coordinates.
    Image is 19 x 26 pixels in size.
    RADECSYS
    CTYPE1 'epsilon '   CUNIT1 'epsilon '
    CRVAL1 1.0   CRPIX1 19.0   CDELT1 0.05
    CTYPE2 'beta    '   CUNIT2 'beta    '
    CRVAL2 1.0   CRPIX2 1.0   CDELT2 0.1
    Image size is 3420.0 x 9360.0 arcsec.
  Note: No RADECSYS keyword found in FITS header.
  No RADECSYS keyword found. Proceeding with caution.
  True data minimum and maximum values:   0.0 1.0
  Scaled plotting minimum and maximum values: 0.0 3.0
  Scaled data min and max values: 0.0 1.0
  Preparing to plot regions on image.
    Regions plotted will use image WCS coordinate system.
  convertRegionstoWCS Warning: This routine current assumes that:
    DS9-format region files are in version 3.0, not version 4.
    Regions are in the exact units needed, no conversion needed.
  No RADECSYS keyword found. Proceeding with caution.
Finished.

PAPFIP software requirements

• Python versions 2.4 or higher.
• PGPLOT version 5.2.0 or higher.
• The Numeric package. Numeric, rather than NumPy or Numarray, is used purely because I need to use pyhdf for another project and want to be able reuse ppgplot-based functions in both projects.
• Papfip uses (and comes with) a modified version of Andrew Williams's fits.py code for python FITS I/O.
• Papfip makes use of (and is supplied with) Janko Hauser's NumAdd.py module for Numeric.
• The python PGPLOT bindings ppgplot, version 1.3.

PAPFIP development history and major changes by release

PAPFIP version 1.1: Released 2008-10-24. Fits WCS coordinates are now accepted in input images and in ds9 (version 3) format region files. Contour level plotting has been improved, including numerical labelling of individual contour levels. Command line axis and data labelling options have been added. Some other minor options have been added and bugs fixed.

PAPFIP version 1.0: Initial freshmeat announcement on 2008-05-01, although the code had not been altered since July 2007.

Pre-release PAPFIP: Development began in February 2005 on an informal basis, and as an excuse to mess around with Python. The primary coding is complete by July 2006, associated with the plotting requirements for Strickland 2007 and Strickland & Heckman 2007.

Why bother writing this when IDL, pdl, ds9, fv, etc. exist?

Many small reasons.

• ds9 is a great data visualization tool - interactively its one of the best tools available for exploring an astronomical images. The xpa interface allows some command-line or scripted control of ds9, but you still require some human interaction to produce a postscipt file. For publication-quality plots the setting of contours, RA/DEC axis labeling, image intensity scaling, source region marking - all of these require ineractive fine tuning. Its OK to do once, but its a real pain to remake a plot, e.g. after remaking the raw fits file with a new calibration.
• IDL is expensive, syntactically and structurally ugly, and burdened with licensing issues managers - I don't use it. A purely free-software project allows use on any (unix) machine any where any time. GDL is an truly-open gnu clean-room reimplementation of IDL which might be great for people with a large pre-existing investment in using IDL.
• PDL and pgperl can do many of the things I require, but I haven't done much perl programming and didn't like the arcana of its syntax. I had trouble extending and modifying the online examples, problems I didn't have using python and ppgplot.
• Python has a nice simple syntax, and its becoming quite necessary to know python as it becomes more important as the scripting language in PyRAF and probably XSPEC too (and now as python bindings to the s-lang functions of CIAO have appeared). I've been meaning to learn the language for a few years now. A small project like this is a great way of learning a language in a practical way.
• I've long wanted a image plotting method that could reproduce the fully command line compatable plotting capabilities of Asterix, itself based on PGPLOT. Image-on-image overlays, RA/DEC-aware plotting, plots with mixed images, line-plot and points, multi-panel plots, Asterix did it all AND remembered what it had done - the data file stored the graphics operations you had set up so all you had to do even years later to remake the postscipt plot was say "gdraw m82.sdf epsf_p" and "voila!" it was done. Didn't like the way you'd contoured it, you could also alter or reset the graphics meta data stored in the SDF. If you'd remade the data file with a new reduction, and wanted to apply all the labeling, contouring, source-region marking, image stretch and colour table operations to the new file you could even "hcopy" the the graphics meta data from within the old data file to the new one and it'd work. Or you could just rerun some csh-script you'd written to set up the graphics meta data in the first place. You simply can not do all of that with ds9, or fv, interactively or not interactively. Some small aspects of that yes, but not all, and none of them have the meta-data approach. Of course, Asterix is sadly not longer under development or maintainence, and worked best on starlink SDFs rather than FITS files. Eventually the old libc5 binary distribution of Asterix is not going to work on newer linux distributions...