National Aeronautics and Space Administration

These are working programs written in ANSI C that illustrate how one can easily read, write, and modify FITS files using the CFITSIO library. Most of these programs are very short, containing only a few 10s of lines of executable code or less, yet they perform quite useful operations on FITS files. Copy the programs to your local machine, then compile, and link them with the CFITSIO library. A short description of how to use each program can be displayed by executing the program without any command line arguments.

You may freely modify, reuse, and redistribute these programs as you wish. It is often easier to use one of these programs as a template when writing a new program, rather than coding the new program completely from scratch.

Available FITS Utility Programs

• fitscopy - copy a file with optional filtering
• listhead - list header keywords
• liststruc - list the structure of a FITS file.
• modhead - write or modify a header keyword
• imcopy - copy an image with optional compression/decompression
• imarith - add, subtract, multiply, or divide an image by a constant or another image
• imlist - list pixel values in an image
• imstat - compute statistics of image pixels
• tablist - display the contents of a FITS table
• tabcalc - general table calculator
• tabmerge - merge rows from one table with another table.
• cookbook - additional example routines.

fitscopy - copy a file with optional filtering

fitscopy infile[ext][filters] outfile

This seemingly simple program can apply complex filtering methods which transform the input file as it is copied. For example, it can select a subimage out of a larger image, select rows in a table based on a logical expression, create new table columns or modify the values in existing columns, and filter tables based on time (using GTI extensions) or on spatial position (using spatial region files). See the Extended File Name syntax page for more examples of the filtering syntax, and refer to the CFITSIO User's Guide for a complete description.

Note that in these examples it is often necessary to enclose the input filename in single quote characters on the Unix command line if the name contains special characters such as '[' or '*'.

EXAMPLES
fitscopy in.fit out.fit                - simple file copy
fitscopy in.fit \!out.fit              - overwrite out.fit file
  (In Unix,  the '!' must be preceded by a '\' on the command line)
fitscopy 'in.fit[11:50,21:50]' out.fit - copy 40x30 pixel subimage
fitscopy 'in.fit[-*,*]' out.fit        - mirror reverse the image
fitscopy 'stdin[11:50,21:60]' stdout   - piped subimage 
fitscopy iniraf.imh out.fit            - convert IRAF image to FITS
fitscopy 'in.dat[i512,512]' out.fit    - binary file to FITS
fitscopy 'in.fit[evt][pi>35]' out.fit  - copy rows which have pi>35
fitscopy 'in.fit[2][bin X,Y]' out.fit  - bin X,Y cols to make image
fitscopy 'in.fit[evt][col Xc=3.*X]' out.fit - create new Xc column
fitscopy 'in.fit[evt][gtifilter()]' out.fit - apply GTI time filter
fitscopy 'in.fit[2][regfilter("pow.reg")]' out.fit - region filter

The program itself is almost trivial, with just 6 lines of code required to open the input file, create the empty output file, copy each HDU from the input file to the output file, then close the 2 files. (See a listing of fitscopy.c).

listhead - list header keywords

listhead infile[ext]

This program will list the header keywords in the specified HDU (Header Data Unit) of a file. If a HDU name or number is not appended to the input root file name, then the program will list the keywords in every HDU in the file.

EXAMPLES;
listhead file.fit        - list all the headers in the file
listhead 'file.fit[0]'   - list the primary array header
listhead 'file.fit[2]'   - list the header of 2nd extension
listhead file.fit+2      - same as above
listhead 'file.fit[GTI]' - list header of the 'GTI' extension

This program is very short. It simply opens the input file and then reads and prints out every keyword in the current HDU in sequence. If a specific HDU was not specified as part of the input file name, it then trys to move to the next HDU in the file and list its keywords. This continues until it reaches the end of the file. (See a listing of listhead).

liststruc - list file structure

liststruc infile[ext]

This program will list the organizational structure of the specified HDU of a file. If a HDU name or number is not appended to the input root file name, then the program will list the structure of every HDU in the file.

If the HDU is the primary array or IMAGE extension then it lists the dimensions and datatype of the image. Otherwise, if the HDU is a table then it lists the name and data format of all the columns in the table.

EXAMPLES:
liststruc file.fit      - list structure of every HDU in the file
liststruc 'file.fit[0]' - list structure of the primary array

This program opens the file and determines the type (i.e., image or table) of the current HDU. The relevant structural parameters for that HDU are then printed out. If a specific HDU is not given as part of the input file name, then the program moves to successive HDUs in the file, if any, and repeats the process until it reaches the end of the file. (See a listing of liststruc).

modhead - write or modify a header keyword

modhead infile[ext] keyword newvalue

This program can write or modify the value of a header keyword. If the 4th 'newvalue' parameter is not specified then the program just prints out the current value of the keyword. Otherwise, it replaces the value of the keyword if it exists, or writes a new keyword if it doesn't exist.

Examples: 
  modhead file.fits dec      - list the value of the DEC keyword 
  modhead file.fits dec 30.0 - write or modify the DEC keyword

This program opens the input file and attempts to read the specified keyword. If the keyword exists, the current value is printed out. If the 4th 'newvalue' argument is specified then the program first checks to make sure that we are not attempting to modify a structural keyword (like NAXIS) which could corrupt the file format. It then constructs a new keyword record string with the new value and overwrites or appends that record to the header. (See a listing of modhead).

imcopy - copy a FITS image with optional data compression

imcopy inimage[ext] outimage[compress]

Copy the input FITS image to the new output file with optional data compression. The output image maybe externally compressed with the gzip or Unix compress algorithm, or it may be internally compressed with the new tile compression format. The output image may be written in the tile-compressed FITS image format by adding the compression qualifiers in square brackets following the output file name.

Examples: 

  imcopy infile.fit 'outfile.fit[compress]' 

       This will use the default compression algorithm (Rice) and the
       default tile size (row by row)

  imcopy infile.fit 'outfile.fit[compress GZIP]' 

       This will use the GZIP compression algorithm and the default
       tile size (row by row).  The allowed compression algorithms are
       Rice, GZIP, and PLIO.  Only the first letter of the algorithm
       name needs to be specified.

  imcopy infile.fit 'outfile.fit[compress G 100,100]' 

       This will use the GZIP compression algorithm and 100 X 100 pixel
       tiles.

  imcopy infile.fit 'outfile.fit[compress R 100,100; 4]' 

       This will use the Rice compression algorithm, 100 X 100 pixel
       tiles, and noise_bits = 4 (assuming the input image has a
       floating point data type).  Decreasing the value of noisebits
       will improve the overall compression efficiency at the expense
       of losing more information.

  imcopy infile.fit outfile.fit

       If the input file is in tile-compressed format, then it will be
       uncompressed to the output file.  Otherwise, it simply copies
       the input image to the output image.

  imcopy 'infile.fit[1001:1500,2001:2500]'  outfile.fit

       This extracts a 500 X 500 pixel section of the much larger
       input image (which may be in tile-compressed format).  The
       output is a normal uncompressed FITS image.

  imcopy 'infile.fit[1001:1500,2001:2500]'  outfile.fit.gz

       Same as above, except the output file is externally compressed
       using the gzip algorithm.

(See a listing of imcopy).

imarith - perform arithmetic on 1 or 2 images

imarith image1[ext] image2[ext] operation outimage  (2 images)
or
imarith image[ext] value operation outimage         (1 image)

Add, subtract, multiply, or divide the pixels in one image by another image or a constant and write the result to a new output image. Supported arithmetic operations are 'add', 'sub', 'mul', or 'div' (only the first character need be supplied)

Examples: 
  imarith in1.fit in2.fit add out.fit     - add the 2 images
  imarith in.fit 1.2 mul out.fit         - multiply image by 1.1
  imarith 'in1.fit[1:20,1:20]' 'in2.fit[1:20,1:20]' add \!out.fit
     - add 2 image sections; also overwrite out.fit if it exists.
  imarith data.fit[1] data.fit[2] mul out.fit - multiply the images
               in the 1st and 2nd extensions of the file data.fit

This program first opens the input images. If 2 images, it checks that they have the same dimensions. It then creates the empty output file and copies the header of the first image into it (thus duplicating the size and data type of that image). Memory is allocated to hold one row of the images. The program computes the output image values by reading one row at a time from the input image(s), performing the desired arithmetic operation, and then writing the resulting row to the output file. The program reads and writes the row of data in double precision format, regardless of the intrinsic datatype of the images; CFITSIO transparently converts the data format if necessary as the rows are read and written. This program also supports 3D data cubes by looping through each plane of the cube. (See a listing of imarith).

imlist - list pixel values

imlist infile[ext][section filter]

List the pixel values in all or part of an image. In order to reduce the amount of output, it is usually more useful to print only a small section of the image by specifying the minimum and maximum pixel in each dimension to be printed, enclosed in square brackets following the root file name and optional extension specifier.

Examples: 
  imlist infile.fit                - print the whole image
  imlist 'in.fit[*:10,*:10]'       - print every 10th pixel
  imlist 'in.fit[20:29,10:19]'     - print 10x10 section
  imlist 'in.fit[20:29:2,10:19:2]' - print every other pixel
  imlist 'in.fit[3][20:29,10:19]'  - print 3rd IMAGE extension
  imlist 'intab.fit[1][bin (X,Y)=32] - print image generated by
     binning the X and Y table columns with a binning size = 32

This program opens the file (which may have been filtered) and checks that the HDU is a 1D or 2D image. It then allocates memory to hold 1 row of pixels, and loops through the image from top to bottom reading and printing out each row. The data format of the image (integer or floating point) determines the output display format that is used for the pixel values. Note that the program reads the image into an array of doubles, regardless of the intrinsic datatype of the image. (See a listing of imlist).

imstat - compute statistics of an image

imstat infile[ext]

Compute simple statistics for the pixels in an input image. It computes the sum of all the pixels, the mean pixel value, and the minimum and maximum pixel values.

Examples: 
  imstat infile.fit                - stats of the whole image
  imstat 'in.fit[20:29,10:19]'     - stats of 10x10 pixel section
  imstat 'in.fit[3][20:29,10:19]'  - stats of 3rd IMAGE extension
  imstat 'intab.fit[1][bin (X,Y)=32] - stats of an image generated
     by binning the X and Y table columns with a binning size = 32

This program opens the file (which may have been filtered) and checks that the HDU is a 2D image. (It could easily be extended to support 3D data cubes as well, as is done in the imarith program). It then allocates memory to hold 1 row of pixels, and loops through the image, accumulating the statistics row by row. Note that the program reads the image into an array of doubles, regardless of the intrinsic datatype of the image. (See a listing of imstat).

tablist - list table contents

tablist infile[ext][col filter][rowfilter]

List the values in all the rows and columns of the input table. In order to reduce the amount of output, it is sometimes useful to select only certain rows or certain columns in the table to be printed by filtering the input file as shown in the examples.

Examples: 
  tablist 'intab.fit[EVENTS]'       - print 'EVENTS' table
  tablist 'intab.fit[1]'            - print 1st extension
  tablist 'intab.fit[1][col X;Y]'   - print X and Y columns
  tablist 'intab.fit[1][PHA > 24]'  - print rows with PHA > 24
  tablist 'intab.fit[1][col X;Y][PHA > 24]' - as above
  tablist 'intab.fit[1][col X;Y;Rad=sqrt(X**2+Y**2)] - print  
     X, Y, and a new 'Rad' column that is computed on the fly.

This program opens the file (which may have been filtered to select only certain rows or columns in the table) and checks that the HDU is a table. It then calculates the number of columns that can be printed on an 80-column output line. The program reads and prints out the values in these columns row by row. Once all the rows have been printed out, it repeats the process for the next set of columns, if any, until the entire table has been printed. All the column values are read as formatted character strings, regardless of the intrinsic datatype of the column. CFITSIO uses the TDISPn keyword if it exists to format the string (e.g., TDISP1 = 'F8.2', or TDISP2 = 'E15.6'), otherwise it uses a default display format for each column datatype. It should be cautioned that reading a table this way, one element at a time column by column and row by row, is not very efficient. For large tables, it is more efficient to read each column over least 100 rows at a time into temporary arrays. After processing all these rows, then read the next 100 rows from the table, and so on, until all the rows have been processed. (See a listing of tablist).

tabcalc - general table calculator

tabcalc intable expression colname outtable

Compute new values for the specified table column using the input arithmetic expression, which may be a function of the values in other table columns or keyword values. The input file is first copied to the output file, then the output file is updated with the new column values. If the column doesn't already exist, then a new column will be appended to the table. Instead of supplying the arithmetic expression directly on the command line, you can also put the expression in an ASCII text file, and supply the name of the text file prepended with the '@' character.

The following arithmetic operators and functions may be used in the expression:

    "addition"         +          "subtraction"        -
    "multiplication"   *          "division"           /
    "negation"         -          "exponentiation"     **   ^
    "absolute value"   abs(x)     "cosine"             cos(x)
    "sine"             sin(x)     "tangent"            tan(x)
    "arc cosine"       arccos(x)  "arc sine"           arcsin(x)
    "arc tangent"      arctan(x)  "arc tangent"        arctan2(x,y)
    "exponential"      exp(x)     "square root"        sqrt(x)
    "natural log"      log(x)     "common log"         log10(x)
    "modulus"          i % j      "random # [0.0,1.0)" random()
    "minimum"          min(x,y)   "maximum"            max(x,y)
    "if-then-else"     b?x:y

Note that the fitscopy program can perform the same operations as tabcalc, using a slightly different command line syntax.

Examples: 
  tabcalc in.fit+1 'X+0.5' Xc out.fit - calculate new Xc column
  fitscopy 'in.fit+1[col Xc=X+0.5]' out.fit - same as above
  tabcalc in.fit+1 @expr.txt Xc out.fit - read expression from file

Other example expressions:
  'counts/#exposure' - divide counts column by EXPOSURE keyword

  'counts > 0 ? counts/#exposure : -99' - if counts column is
  greater than 0 then divide counts by the exposure keyword,
  otherwise set the value equal to -99.

  '(counts{-1} * 0.25 + counts * 0.5 + counts{+1} * 0.25) / 3.' 
   - compute weighted running mean of counts by averaging the
     values in the previous row, the current row, and the next
     row.  The row offset is enclosed in curly brackets.

This program, like the fitscopy program, is very short and simple because almost all the work is done by CFITSIO. Only 6 lines of code are needed to open the input file, create the empty output file, copy all HDUs from the input to output file, and finally calculate the new column with a call to the 'fits_calculator' routine. (See a listing of tabcalc).

tabmerge - merge rows from 2 tables

tabmerge intable[ext][filters] outtable[ext]

Merge 2 tables by appending all the rows from the 1st table onto the 2nd table. The 2 tables must have identical structure, with the same number of columns with the same datatypes.

Examples: 
1.  tabmerge in.fit+1 out.fit+2
    - append the rows from the 1st extension of the input file
      into the table in the 2nd extension of the output file.

2.  tabmerge 'in.fit+1[[PI > 45]' out.fit+2
    - Same as the 1st example, except only rows that have a PI
      column value > 45 will be merged into the output table.

In this program, a series of 'if' statements are used to verify that the input files exist and that they point to valid tables with identical sets of columns. If all these checks are satisfied correctly, then additional empty rows are appended to the end of the output table, and the rows are copied one by one from the input table to the output table. Since the tables have identical structures, it is possible here to use the low-level CFITSIO routines that read and write each row of the table as a raw string of bytes. This is much faster than having to read and write the numerical values in each column individually.

Note that in this example the output file is modified in place, rather than creating a whole new output file. It is often easier to modify the existing file in this way rather than create a new one, but there is a slight risk that the output file could be corrupted if the program crashes before the modifications have been completed (e.g., due to a power failure, or the user exits the program with Cntl-C, or there is insufficient disk space to write the modified file).

(See a listing of tabmerge).

cookbook - more example routines

This is a somewhat longer program that gives additional examples of reading and writing images and tables. (See a listing of cookbook).