Repatch (from linux) of OSX IRAF

1 files changed, 772 insertions, 0 deletions

diff --git a/vendor/cfitsio/cookbook.f b/vendor/cfitsio/cookbook.f
new file mode 100644
index 00000000..8becfdd2
--- /dev/null
+++ b/vendor/cfitsio/cookbook.f
@@ -0,0 +1,772 @@
+      program main
+
+C  This is the FITSIO cookbook program that contains an annotated listing of
+C  various computer programs that read and write files in FITS format
+C  using the FITSIO subroutine interface.  These examples are
+C  working programs which users may adapt and modify for their own
+C  purposes.  This Cookbook serves as a companion to the FITSIO User's
+C  Guide that provides more complete documentation on all the
+C  available FITSIO subroutines.
+
+C  Call each subroutine in turn:
+
+      call writeimage
+      call writeascii
+      call writebintable
+      call copyhdu
+      call selectrows
+      call readheader
+      call readimage
+      call readtable
+      print *
+      print *,"All the fitsio cookbook routines ran successfully."
+
+      end
+C *************************************************************************
+      subroutine writeimage
+
+C  Create a FITS primary array containing a 2-D image
+
+      integer status,unit,blocksize,bitpix,naxis,naxes(2)
+      integer i,j,group,fpixel,nelements,array(300,200)
+      character filename*80
+      logical simple,extend
+
+C  The STATUS parameter must be initialized before using FITSIO.  A
+C  positive value of STATUS is returned whenever a serious error occurs.
+C  FITSIO uses an `inherited status' convention, which means that if a
+C  subroutine is called with a positive input value of STATUS, then the
+C  subroutine will exit immediately, preserving the status value. For 
+C  simplicity, this program only checks the status value at the end of 
+C  the program, but it is usually better practice to check the status 
+C  value more frequently.
+
+      status=0
+
+C  Name of the FITS file to be created:
+      filename='ATESTFILEZ.FITS'
+
+C  Delete the file if it already exists, so we can then recreate it.
+C  The deletefile subroutine is listed at the end of this file.
+      call deletefile(filename,status)
+
+C  Get an unused Logical Unit Number to use to open the FITS file.
+C  This routine is not required;  programmers can choose any unused
+C  unit number to open the file.
+      call ftgiou(unit,status)
+
+C  Create the new empty FITS file.  The blocksize parameter is a
+C  historical artifact and the value is ignored by FITSIO.
+      blocksize=1
+      call ftinit(unit,filename,blocksize,status)
+
+C  Initialize parameters about the FITS image.
+C  BITPIX = 16 means that the image pixels will consist of 16-bit
+C  integers.  The size of the image is given by the NAXES values. 
+C  The EXTEND = TRUE parameter indicates that the FITS file
+C  may contain extensions following the primary array.
+      simple=.true.
+      bitpix=16
+      naxis=2
+      naxes(1)=300
+      naxes(2)=200
+      extend=.true.
+
+C  Write the required header keywords to the file
+      call ftphpr(unit,simple,bitpix,naxis,naxes,0,1,extend,status)
+
+C  Initialize the values in the image with a linear ramp function
+      do j=1,naxes(2)
+          do i=1,naxes(1)
+              array(i,j)=i - 1 +j - 1
+          end do
+      end do
+
+C  Write the array to the FITS file.
+C  The last letter of the subroutine name defines the datatype of the
+C  array argument; in this case the 'J' indicates that the array has an
+C  integer*4 datatype. ('I' = I*2, 'E' = Real*4, 'D' = Real*8).
+C  The 2D array is treated as a single 1-D array with NAXIS1 * NAXIS2
+C  total number of pixels.  GROUP is seldom used parameter that should
+C  almost always be set = 1.
+      group=1
+      fpixel=1
+      nelements=naxes(1)*naxes(2)
+      call ftpprj(unit,group,fpixel,nelements,array,status)
+
+C  Write another optional keyword to the header
+C  The keyword record will look like this in the FITS file:
+C
+C  EXPOSURE=                 1500 / Total Exposure Time
+C
+      call ftpkyj(unit,'EXPOSURE',1500,'Total Exposure Time',status)
+
+C  The FITS file must always be closed before exiting the program. 
+C  Any unit numbers allocated with FTGIOU must be freed with FTFIOU.
+      call ftclos(unit, status)
+      call ftfiou(unit, status)
+
+C  Check for any errors, and if so print out error messages.
+C  The PRINTERROR subroutine is listed near the end of this file.
+      if (status .gt. 0)call printerror(status)
+      end
+C *************************************************************************
+      subroutine writeascii
+
+C  Create an ASCII table containing 3 columns and 6 rows.  For convenience,
+C  the ASCII table extension is appended to the FITS image file created 
+C  previously by the WRITEIMAGE subroutine.
+
+      integer status,unit,readwrite,blocksize,tfields,nrows,rowlen
+      integer nspace,tbcol(3),diameter(6), colnum,frow,felem
+      real density(6)
+      character filename*40,extname*16
+      character*16 ttype(3),tform(3),tunit(3),name(6)
+      data ttype/'Planet','Diameter','Density'/
+      data tform/'A8','I6','F4.2'/
+      data tunit/' ','km','g/cm'/
+      data name/'Mercury','Venus','Earth','Mars','Jupiter','Saturn'/
+      data diameter/4880,12112,12742,6800,143000,121000/
+      data density/5.1,5.3,5.52,3.94,1.33,0.69/
+
+C  The STATUS parameter must always be initialized.
+      status=0
+
+C  Name of the FITS file to append the ASCII table to:
+      filename='ATESTFILEZ.FITS'
+
+C  Get an unused Logical Unit Number to use to open the FITS file.
+      call ftgiou(unit,status)
+
+C  Open the FITS file with write access.
+C  (readwrite = 0 would open the file with readonly access).
+      readwrite=1
+      call ftopen(unit,filename,readwrite,blocksize,status)
+
+C  FTCRHD creates a new empty FITS extension following the current
+C  extension and moves to it.  In this case, FITSIO was initially
+C  positioned on the primary array when the FITS file was first opened, so
+C  FTCRHD appends an empty extension and moves to it.  All future FITSIO
+C  calls then operate on the new extension (which will be an ASCII
+C  table).
+      call ftcrhd(unit,status)
+
+C  define parameters for the ASCII table (see the above data statements)
+      tfields=3
+      nrows=6
+      extname='PLANETS_ASCII'
+      
+C  FTGABC is a convenient subroutine for calculating the total width of
+C  the table and the starting position of each column in an ASCII table.
+C  Any number of blank spaces (including zero)  may be inserted between
+C  each column of the table, as specified by the NSPACE parameter.
+      nspace=1
+      call ftgabc(tfields,tform,nspace,rowlen,tbcol,status)
+
+C  FTPHTB writes all the required header keywords which define the
+C  structure of the ASCII table. NROWS and TFIELDS give the number of
+C  rows and columns in the table, and the TTYPE, TBCOL, TFORM, and TUNIT
+C  arrays give the column name, starting position, format, and units,
+C  respectively of each column. The values of the ROWLEN and TBCOL parameters
+C  were previously calculated by the FTGABC routine.
+      call ftphtb(unit,rowlen,nrows,tfields,ttype,tbcol,tform,tunit,
+     &            extname,status)
+
+C  Write names to the first column, diameters to 2nd col., and density to 3rd
+C  FTPCLS writes the string values to the NAME column (column 1) of the
+C  table.  The FTPCLJ and FTPCLE routines write the diameter (integer) and
+C  density (real) value to the 2nd and 3rd columns.  The FITSIO routines
+C  are column oriented, so it is usually easier to read or write data in a
+C  table in a column by column order rather than row by row.
+      frow=1
+      felem=1
+      colnum=1
+      call ftpcls(unit,colnum,frow,felem,nrows,name,status)
+      colnum=2
+      call ftpclj(unit,colnum,frow,felem,nrows,diameter,status)  
+      colnum=3
+      call ftpcle(unit,colnum,frow,felem,nrows,density,status)  
+
+C  The FITS file must always be closed before exiting the program. 
+C  Any unit numbers allocated with FTGIOU must be freed with FTFIOU.
+      call ftclos(unit, status)
+      call ftfiou(unit, status)
+
+C  Check for any error, and if so print out error messages.
+C  The PRINTERROR subroutine is listed near the end of this file.
+      if (status .gt. 0)call printerror(status)
+      end
+C *************************************************************************
+      subroutine writebintable
+
+C  This routine creates a FITS binary table, or BINTABLE, containing
+C  3 columns and 6 rows.  This routine is nearly identical to the
+C  previous WRITEASCII routine, except that the call to FTGABC is not
+C  needed, and FTPHBN is called rather than FTPHTB to write the
+C  required header keywords.
+
+      integer status,unit,readwrite,blocksize,hdutype,tfields,nrows
+      integer varidat,diameter(6), colnum,frow,felem
+      real density(6)
+      character filename*40,extname*16
+      character*16 ttype(3),tform(3),tunit(3),name(6)
+      data ttype/'Planet','Diameter','Density'/
+      data tform/'8A','1J','1E'/
+      data tunit/' ','km','g/cm'/
+      data name/'Mercury','Venus','Earth','Mars','Jupiter','Saturn'/
+      data diameter/4880,12112,12742,6800,143000,121000/
+      data density/5.1,5.3,5.52,3.94,1.33,0.69/
+
+C  The STATUS parameter must always be initialized.
+      status=0
+
+C  Name of the FITS file to append the ASCII table to:
+      filename='ATESTFILEZ.FITS'
+
+C  Get an unused Logical Unit Number to use to open the FITS file.
+      call ftgiou(unit,status)
+
+C  Open the FITS file, with write access.
+      readwrite=1
+      call ftopen(unit,filename,readwrite,blocksize,status)
+
+C  Move to the last (2nd) HDU in the file (the ASCII table).
+      call ftmahd(unit,2,hdutype,status)
+
+C  Append/create a new empty HDU onto the end of the file and move to it.
+      call ftcrhd(unit,status)
+
+C  Define parameters for the binary table (see the above data statements)
+      tfields=3
+      nrows=6
+      extname='PLANETS_BINARY'
+      varidat=0
+      
+C  FTPHBN writes all the required header keywords which define the
+C  structure of the binary table. NROWS and TFIELDS gives the number of
+C  rows and columns in the table, and the TTYPE, TFORM, and TUNIT arrays
+C  give the column name, format, and units, respectively of each column.
+      call ftphbn(unit,nrows,tfields,ttype,tform,tunit,
+     &            extname,varidat,status)
+
+C  Write names to the first column, diameters to 2nd col., and density to 3rd
+C  FTPCLS writes the string values to the NAME column (column 1) of the
+C  table.  The FTPCLJ and FTPCLE routines write the diameter (integer) and
+C  density (real) value to the 2nd and 3rd columns.  The FITSIO routines
+C  are column oriented, so it is usually easier to read or write data in a
+C  table in a column by column order rather than row by row.  Note that
+C  the identical subroutine calls are used to write to either ASCII or
+C  binary FITS tables.
+      frow=1
+      felem=1
+      colnum=1
+      call ftpcls(unit,colnum,frow,felem,nrows,name,status)
+      colnum=2
+      call ftpclj(unit,colnum,frow,felem,nrows,diameter,status)  
+      colnum=3
+      call ftpcle(unit,colnum,frow,felem,nrows,density,status)  
+
+C  The FITS file must always be closed before exiting the program. 
+C  Any unit numbers allocated with FTGIOU must be freed with FTFIOU.
+      call ftclos(unit, status)
+      call ftfiou(unit, status)
+
+C  Check for any error, and if so print out error messages.
+C  The PRINTERROR subroutine is listed near the end of this file.
+      if (status .gt. 0)call printerror(status)
+      end
+C *************************************************************************
+      subroutine copyhdu
+
+C  Copy the 1st and 3rd HDUs from the input file to a new FITS file
+
+      integer status,inunit,outunit,readwrite,blocksize,morekeys,hdutype
+      character infilename*40,outfilename*40
+
+C  The STATUS parameter must always be initialized.
+      status=0
+
+C     Name of the FITS files:
+      infilename='ATESTFILEZ.FITS'
+      outfilename='BTESTFILEZ.FITS'
+
+C  Delete the file if it already exists, so we can then recreate it
+C  The deletefile subroutine is listed at the end of this file.
+      call deletefile(outfilename,status)
+
+C  Get  unused Logical Unit Numbers to use to open the FITS files.
+      call ftgiou(inunit,status)
+      call ftgiou(outunit,status)
+
+C  Open the input FITS file, with readonly access
+      readwrite=0
+      call ftopen(inunit,infilename,readwrite,blocksize,status)
+
+C  Create the new empty FITS file (value of blocksize is ignored)
+      blocksize=1
+      call ftinit(outunit,outfilename,blocksize,status)
+
+C  FTCOPY copies the current HDU from the input FITS file to the output
+C  file.  The MOREKEY parameter allows one to reserve space for additional
+C  header keywords when the HDU is created.   FITSIO will automatically
+C  insert more header space if required, so programmers do not have to
+C  reserve space ahead of time, although it is more efficient to do so if
+C  it is known that more keywords will be appended to the header.
+      morekeys=0
+      call ftcopy(inunit,outunit,morekeys,status)
+
+C  Append/create a new empty extension on the end of the output file
+      call ftcrhd(outunit,status)
+
+C  Skip to the 3rd extension in the input file which in this case
+C  is the binary table created by the previous WRITEBINARY routine.
+      call ftmahd(inunit,3,hdutype,status)
+
+C  FTCOPY now copies the binary table from the input FITS file
+C  to the output file.
+      call ftcopy(inunit,outunit,morekeys,status)  
+
+C  The FITS files must always be closed before exiting the program. 
+C  Any unit numbers allocated with FTGIOU must be freed with FTFIOU.
+C  Giving -1 for the value of the first argument causes all previously
+C  allocated unit numbers to be released.
+
+      call ftclos(inunit, status)
+      call ftclos(outunit, status)
+      call ftfiou(-1, status)
+
+C  Check for any error, and if so print out error messages.
+C  The PRINTERROR subroutine is listed near the end of this file.
+      if (status .gt. 0)call printerror(status)
+      end
+C *************************************************************************
+      subroutine selectrows
+
+C  This routine copies selected rows from an input table into a new output
+C  FITS table.  In this example all the rows in the input table that have
+C  a value of the DENSITY column less that 3.0 are copied to the output
+C  table.  This program illustrates several generally useful techniques,
+C  including:
+C      how to locate the end of a FITS file
+C      how to create a table when the total number of rows in the table
+C      is not known until the table is completed
+C      how to efficiently copy entire rows from one table to another.
+
+      integer status,inunit,outunit,readwrite,blocksize,hdutype
+      integer nkeys,nspace,naxes(2),nfound,colnum,frow,felem
+      integer noutrows,irow,temp(100),i
+      real nullval,density(6)
+      character infilename*40,outfilename*40,record*80
+      logical exact,anynulls
+
+C  The STATUS parameter must always be initialized.
+      status=0
+
+C     Names of the FITS files:
+      infilename='ATESTFILEZ.FITS'
+      outfilename='BTESTFILEZ.FITS'
+
+C  Get  unused Logical Unit Numbers to use to open the FITS files.
+      call ftgiou(inunit,status)
+      call ftgiou(outunit,status)
+
+C  The input FITS file is opened with READONLY access, and the output
+C  FITS file is opened with WRITE access.
+      readwrite=0
+      call ftopen(inunit,infilename,readwrite,blocksize,status)
+      readwrite=1
+      call ftopen(outunit,outfilename,readwrite,blocksize,status)
+
+C  move to the 3rd HDU in the input file (a binary table in this case)
+      call ftmahd(inunit,3,hdutype,status)
+
+C  This do-loop illustrates how to move to the last extension in any FITS
+C  file.  The call to FTMRHD moves one extension at a time through the
+C  FITS file until an `End-of-file' status value (= 107) is returned.
+      do while (status .eq. 0)
+          call ftmrhd(outunit,1,hdutype,status)
+      end do
+
+C  After locating the end of the FITS file, it is necessary to reset the
+C  status value to zero and also clear the internal error message stack
+C  in FITSIO.  The previous `End-of-file' error will have produced
+C  an unimportant message on the error stack which can be cleared with
+C  the call to the FTCMSG routine (which has no arguments).
+
+      if (status .eq. 107)then
+          status=0
+          call ftcmsg
+      end if
+
+C  Create a new empty extension in the output file.
+      call ftcrhd(outunit,status)
+
+C  Find the number of keywords in the input table header.
+      call ftghsp(inunit,nkeys,nspace,status)
+
+C  This do-loop of calls to FTGREC and FTPREC copies all the keywords from
+C  the input to the output FITS file.  Notice that the specified number
+C  of rows in the output table, as given by the NAXIS2 keyword, will be
+C  incorrect.  This value will be modified later after it is known how many
+C  rows will be in the table, so it does not matter how many rows are specified
+C  initially.
+      do i=1,nkeys
+          call ftgrec(inunit,i,record,status)
+          call ftprec(outunit,record,status)
+      end do
+
+C  FTGKNJ is used to get the value of the NAXIS1 and NAXIS2 keywords,
+C  which define the width of the table in bytes, and the number of
+C  rows in the table.
+      call ftgknj(inunit,'NAXIS',1,2,naxes,nfound,status)
+
+C  FTGCNO gets the column number of the `DENSITY' column; the column
+C  number is needed when reading the data in the column.  The EXACT
+C  parameter determines whether or not the match to the column names
+C  will be case sensitive.
+      exact=.false.
+      call ftgcno(inunit,exact,'DENSITY',colnum,status)
+
+C  FTGCVE reads all 6 rows of data in the `DENSITY' column.  The number
+C  of rows in the table is given by NAXES(2). Any null values in the
+C  table will be returned with the corresponding value set to -99
+C  (= the value of NULLVAL).  The ANYNULLS parameter will be set to TRUE
+C  if any null values were found while reading the data values in the table.
+      frow=1
+      felem=1
+      nullval=-99.
+      call ftgcve(inunit,colnum,frow,felem,naxes(2),nullval,
+     &            density,anynulls,status)
+
+C  If the density is less than 3.0, copy the row to the output table.
+C  FTGTBB and FTPTBB are low-level routines to read and write, respectively,
+C  a specified number of bytes in the table, starting at the specified
+C  row number and beginning byte within the row.  These routines do
+C  not do any interpretation of the bytes, and simply pass them to or
+C  from the FITS file without any modification.  This is a faster
+C  way of transferring large chunks of data from one FITS file to another,
+C  than reading and then writing each column of data individually.
+C  In this case an entire row of bytes (the row length is specified
+C  by the naxes(1) parameter) is transferred.  The datatype of the 
+C  buffer array (TEMP in this case) is immaterial so long as it is
+C  declared large enough to hold the required number of bytes.
+      noutrows=0
+      do irow=1,naxes(2)
+          if (density(irow) .lt. 3.0)then
+              noutrows=noutrows+1
+              call ftgtbb(inunit,irow,1,naxes(1),temp,status)
+              call ftptbb(outunit,noutrows,1,naxes(1),temp,status)
+          end if
+      end do
+
+C  Update the NAXIS2 keyword with the correct no. of rows in the output file.
+C  After all the rows have been written to the output table, the
+C  FTMKYJ routine is used to overwrite the NAXIS2 keyword value with
+C  the correct number of rows.  Specifying `\&' for the comment string
+C  tells FITSIO to keep the current comment string in the keyword and
+C  only modify the value.  Because the total number of rows in the table
+C  was unknown when the table was first created, any value (including 0)
+C  could have been used for the initial NAXIS2 keyword value.
+      call ftmkyj(outunit,'NAXIS2',noutrows,'&',status)
+
+C  The FITS files must always be closed before exiting the program. 
+C  Any unit numbers allocated with FTGIOU must be freed with FTFIOU.
+      call ftclos(inunit, status)
+      call ftclos(outunit, status)
+      call ftfiou(-1, status)
+
+C  Check for any error, and if so print out error messages.
+C  The PRINTERROR subroutine is listed near the end of this file.
+      if (status .gt. 0)call printerror(status)
+      end
+C *************************************************************************
+      subroutine readheader
+
+C  Print out all the header keywords in all extensions of a FITS file
+
+      integer status,unit,readwrite,blocksize,nkeys,nspace,hdutype,i,j
+      character filename*80,record*80
+
+C  The STATUS parameter must always be initialized.
+      status=0
+
+C  Get an unused Logical Unit Number to use to open the FITS file.
+      call ftgiou(unit,status)
+
+C     name of FITS file 
+      filename='ATESTFILEZ.FITS'
+
+C     open the FITS file, with read-only access.  The returned BLOCKSIZE
+C     parameter is obsolete and should be ignored. 
+      readwrite=0
+      call ftopen(unit,filename,readwrite,blocksize,status)
+
+      j = 0
+100   continue
+      j = j + 1
+
+      print *,'Header listing for HDU', j
+
+C  The FTGHSP subroutine returns the number of existing keywords in the
+C  current header data unit (CHDU), not counting the required END keyword,
+      call ftghsp(unit,nkeys,nspace,status)
+
+C  Read each 80-character keyword record, and print it out.
+      do i = 1, nkeys
+          call ftgrec(unit,i,record,status)
+          print *,record
+      end do
+
+C  Print out an END record, and a blank line to mark the end of the header.
+      if (status .eq. 0)then
+          print *,'END'
+          print *,' '
+      end if
+
+C  Try moving to the next extension in the FITS file, if it exists.
+C  The FTMRHD subroutine attempts to move to the next HDU, as specified by
+C  the second parameter.   This subroutine moves by a relative number of
+C  HDUs from the current HDU.  The related FTMAHD routine may be used to
+C  move to an absolute HDU number in the FITS file.  If the end-of-file is
+C  encountered when trying to move to the specified extension, then a
+C  status = 107 is returned.
+      call ftmrhd(unit,1,hdutype,status)
+
+      if (status .eq. 0)then
+C         success, so jump back and print out keywords in this extension
+          go to 100
+
+      else if (status .eq. 107)then
+C         hit end of file, so quit
+          status=0
+      end if
+
+C  The FITS file must always be closed before exiting the program. 
+C  Any unit numbers allocated with FTGIOU must be freed with FTFIOU.
+      call ftclos(unit, status)
+      call ftfiou(unit, status)
+
+C  Check for any error, and if so print out error messages.
+C  The PRINTERROR subroutine is listed near the end of this file.
+      if (status .gt. 0)call printerror(status)
+      end
+C *************************************************************************
+      subroutine readimage
+
+C  Read a FITS image and determine the minimum and maximum pixel value.
+C  Rather than reading the entire image in
+C  at once (which could require a very large array), the image is read
+C  in pieces, 100 pixels at a time.  
+
+      integer status,unit,readwrite,blocksize,naxes(2),nfound
+      integer group,firstpix,nbuffer,npixels,i
+      real datamin,datamax,nullval,buffer(100)
+      logical anynull
+      character filename*80
+
+C  The STATUS parameter must always be initialized.
+      status=0
+
+C  Get an unused Logical Unit Number to use to open the FITS file.
+      call ftgiou(unit,status)
+
+C  Open the FITS file previously created by WRITEIMAGE
+      filename='ATESTFILEZ.FITS'
+      readwrite=0
+      call ftopen(unit,filename,readwrite,blocksize,status)
+
+C  Determine the size of the image.
+      call ftgknj(unit,'NAXIS',1,2,naxes,nfound,status)
+
+C  Check that it found both NAXIS1 and NAXIS2 keywords.
+      if (nfound .ne. 2)then
+          print *,'READIMAGE failed to read the NAXISn keywords.'
+          return
+       end if
+
+C  Initialize variables
+      npixels=naxes(1)*naxes(2)
+      group=1
+      firstpix=1
+      nullval=-999
+      datamin=1.0E30
+      datamax=-1.0E30
+
+      do while (npixels .gt. 0)
+C         read up to 100 pixels at a time 
+          nbuffer=min(100,npixels)
+      
+          call ftgpve(unit,group,firstpix,nbuffer,nullval,
+     &            buffer,anynull,status)
+
+C         find the min and max values
+          do i=1,nbuffer
+              datamin=min(datamin,buffer(i))
+              datamax=max(datamax,buffer(i))
+          end do
+
+C         increment pointers and loop back to read the next group of pixels
+          npixels=npixels-nbuffer
+          firstpix=firstpix+nbuffer
+      end do
+
+      print *
+      print *,'Min and max image pixels = ',datamin,datamax
+
+C  The FITS file must always be closed before exiting the program. 
+C  Any unit numbers allocated with FTGIOU must be freed with FTFIOU.
+      call ftclos(unit, status)
+      call ftfiou(unit, status)
+
+C  Check for any error, and if so print out error messages.
+C  The PRINTERROR subroutine is listed near the end of this file.
+      if (status .gt. 0)call printerror(status)
+      end
+C *************************************************************************
+      subroutine readtable
+
+C  Read and print data values from an ASCII or binary table
+C  This example reads and prints out all the data in the ASCII and
+C  the binary tables that were previously created by WRITEASCII and
+C  WRITEBINTABLE.  Note that the exact same FITSIO routines are
+C  used to read both types of tables.
+
+      integer status,unit,readwrite,blocksize,hdutype,ntable
+      integer felem,nelems,nullj,diameter,nfound,irow,colnum
+      real nulle,density
+      character filename*40,nullstr*1,name*8,ttype(3)*10
+      logical anynull
+
+C  The STATUS parameter must always be initialized.
+      status=0
+
+C  Get an unused Logical Unit Number to use to open the FITS file.
+      call ftgiou(unit,status)
+
+C  Open the FITS file previously created by WRITEIMAGE
+      filename='ATESTFILEZ.FITS'
+      readwrite=0
+      call ftopen(unit,filename,readwrite,blocksize,status)
+
+C  Loop twice, first reading the ASCII table, then the binary table
+      do ntable=2,3
+
+C  Move to the next extension
+          call ftmahd(unit,ntable,hdutype,status)
+
+          print *,' '
+          if (hdutype .eq. 1)then
+              print *,'Reading ASCII table in HDU ',ntable
+          else if (hdutype .eq. 2)then
+              print *,'Reading binary table in HDU ',ntable
+          end if
+
+C  Read the TTYPEn keywords, which give the names of the columns
+          call ftgkns(unit,'TTYPE',1,3,ttype,nfound,status)
+          write(*,2000)ttype
+2000      format(2x,"Row   ",3a10)
+
+C  Read the data, one row at a time, and print them out
+          felem=1
+          nelems=1
+          nullstr=' '
+          nullj=0
+          nulle=0.
+          do irow=1,6
+C             FTGCVS reads the NAMES from the first column of the table.
+              colnum=1
+              call ftgcvs(unit,colnum,irow,felem,nelems,nullstr,name,
+     &                    anynull,status)
+
+C             FTGCVJ reads the DIAMETER values from the second column.
+              colnum=2
+              call ftgcvj(unit,colnum,irow,felem,nelems,nullj,diameter,
+     &                    anynull,status)
+
+C             FTGCVE reads the DENSITY values from the third column.
+              colnum=3
+              call ftgcve(unit,colnum,irow,felem,nelems,nulle,density,
+     &                    anynull,status)
+              write(*,2001)irow,name,diameter,density
+2001          format(i5,a10,i10,f10.2)
+          end do
+      end do
+
+C  The FITS file must always be closed before exiting the program. 
+C  Any unit numbers allocated with FTGIOU must be freed with FTFIOU.
+      call ftclos(unit, status)
+      call ftfiou(unit, status)
+
+C  Check for any error, and if so print out error messages.
+C  The PRINTERROR subroutine is listed near the end of this file.
+      if (status .gt. 0)call printerror(status)
+      end
+C *************************************************************************
+      subroutine printerror(status)
+
+C  This subroutine prints out the descriptive text corresponding to the
+C  error status value and prints out the contents of the internal
+C  error message stack generated by FITSIO whenever an error occurs.
+
+      integer status
+      character errtext*30,errmessage*80
+
+C  Check if status is OK (no error); if so, simply return
+      if (status .le. 0)return
+
+C  The FTGERR subroutine returns a descriptive 30-character text string that
+C  corresponds to the integer error status number.  A complete list of all
+C  the error numbers can be found in the back of the FITSIO User's Guide.
+      call ftgerr(status,errtext)
+      print *,'FITSIO Error Status =',status,': ',errtext
+
+C  FITSIO usually generates an internal stack of error messages whenever
+C  an error occurs.  These messages provide much more information on the
+C  cause of the problem than can be provided by the single integer error
+C  status value.  The FTGMSG subroutine retrieves the oldest message from
+C  the stack and shifts any remaining messages on the stack down one
+C  position.  FTGMSG is called repeatedly until a blank message is
+C  returned, which indicates that the stack is empty.  Each error message
+C  may be up to 80 characters in length.  Another subroutine, called
+C  FTCMSG, is available to simply clear the whole error message stack in
+C  cases where one is not interested in the contents.
+      call ftgmsg(errmessage)
+      do while (errmessage .ne. ' ')
+          print *,errmessage
+          call ftgmsg(errmessage)
+      end do
+      end
+C *************************************************************************
+      subroutine deletefile(filename,status)
+
+C  A simple little routine to delete a FITS file
+
+      integer status,unit,blocksize
+      character*(*) filename
+
+C  Simply return if status is greater than zero
+      if (status .gt. 0)return
+
+C  Get an unused Logical Unit Number to use to open the FITS file
+      call ftgiou(unit,status)
+
+C  Try to open the file, to see if it exists
+      call ftopen(unit,filename,1,blocksize,status)
+
+      if (status .eq. 0)then
+C         file was opened;  so now delete it 
+          call ftdelt(unit,status)
+      else if (status .eq. 103)then
+C         file doesn't exist, so just reset status to zero and clear errors
+          status=0
+          call ftcmsg
+      else
+C         there was some other error opening the file; delete the file anyway
+          status=0
+          call ftcmsg
+          call ftdelt(unit,status)
+      end if
+
+C  Free the unit number for later reuse
+      call ftfiou(unit, status)
+      end