|
This is a part of the UK Mirror of the
HEASARC HEAsoft pages,
maintained by LEDAS.
All HEAsoft-related links on this page are to HTML files or data stored at Leicester, unless indicated by a  .
Please refer any problems with these pages to Duncan Law-Green (ledas-help@star.le.ac.uk). |
These are working programs written in ANSI C which illustrate how easy it is to 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.
fitscopy infile[ext][filters] outfileThis 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 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' extensionThis 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 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 arrayThis 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 infile[ext] keyword newvalueThis 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 keywordThis 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 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 - arithmetically combine 2 images
imarith image1[ext] image2[ext] operation outimageAdd, subtract, multiply, or divide the 2 input images 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 a out.fit - add the 2 images
imarith 'in1.fit[1:20,1:20]' 'in2.fit[1:20,1:20]' a \!out.fit
- add 2 image sections; also overwrite out.fit if it exists.
imarith data.fit[1] data.fit[2] m out.fit - multiply the images
in the 1st and 2nd extensions of the file data.fit
This program opens the 2 input images and 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 from each of
the input images. The program computes the output image values by
reading one row at a time from the input images, 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 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 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 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 intable expression colname outtableCompute 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 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).