[gnuastro-commits] master 0a32a82 027/125: Any number of searched columns from FITS are read

[Mohammad Akhlaghi]

branch: master
commit 0a32a829f5572d8bf21eb728f0a086dbe35a660e
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>

    Any number of searched columns from FITS are read
    
    The high level `gal_table_read_cols' function can now read any number of
    columns from a given FITS table. Also the display format FITS header
    keywords are now read into the data structure when they are present. The
    great advantage of this is that since it is now implemented as a library,
    it can be used by any of the programs (and even users of the libraries) to
    greatly simplify the process of reading the input data into the program for
    further analysis.
---
 bin/table/main.h     |   7 +-
 bin/table/ui.c       |  17 ++---
 lib/data.c           |  23 ++++---
 lib/fits.c           | 189 +++++++++++++++++++++++++++++++++++++++++++++++----
 lib/gnuastro/data.h  |  11 ++-
 lib/gnuastro/fits.h  |   4 ++
 lib/gnuastro/table.h |  27 +++++++-
 lib/table.c          |  61 ++++++++++++-----
 8 files changed, 281 insertions(+), 58 deletions(-)

diff --git a/bin/table/main.h b/bin/table/main.h
index b3943ff..9225956 100644
--- a/bin/table/main.h
+++ b/bin/table/main.h
@@ -63,12 +63,11 @@ struct tableparams
   struct gal_commonparams  cp;  /* Common parameters.                   */
 
   /* Input */
-  struct gal_linkedlist_stll *columns; /* List of given columns. */
+  struct gal_linkedlist_stll *columns; /* List of given columns.        */
 
   /* Output: */
-  size_t                nrows;  /* Number of rows in table.             */
-  size_t               nocols;  /* Number of output columns.            */
-  int                 outtype;  /* Type of output table.                */
+  int            outtabletype;  /* Type of output table (FITS, txt).    */
+  gal_data_t           *table;  /* Linked list of output table columns. */
 
   /* Operating modes */
   int             information;  /* ==1, only print FITS information.    */
diff --git a/bin/table/ui.c b/bin/table/ui.c
index 825020f..11d5e4d 100644
--- a/bin/table/ui.c
+++ b/bin/table/ui.c
@@ -266,9 +266,9 @@ sanitycheck(struct tableparams *p)
   if(gal_fits_name_is_fits(up->filename))
     {
       if( !strcmp(up->fitstabletype, "ascii") )
-        p->outtype=GAL_TABLE_TYPE_AFITS;
+        p->outtabletype=GAL_TABLE_TYPE_AFITS;
       else if( !strcmp(up->fitstabletype, "binary") )
-        p->outtype=GAL_TABLE_TYPE_BFITS;
+        p->outtabletype=GAL_TABLE_TYPE_BFITS;
       else
         error(EXIT_FAILURE, 0, "the value to the `--fitstabletype' "
               "option on the command line or `fitstabletype' variable in "
@@ -276,7 +276,7 @@ sanitycheck(struct tableparams *p)
               "`binary'. You have given `%s'", up->fitstabletype);
     }
   else
-    p->outtype=GAL_TABLE_TYPE_TXT;
+    p->outtabletype=GAL_TABLE_TYPE_TXT;
 
   /* Set the searchin integer value. */
   p->searchin=gal_table_searchin_from_str(p->up.searchin);
@@ -370,12 +370,13 @@ print_information_exit(struct tableparams *p)
 void
 preparearrays(struct tableparams *p)
 {
-  /* Pull out the columns that we are looking for. */
+  /* Reverse the list of column search criteria that we are looking for
+     (since this is a last-in-first-out linked list, the order that
+     elements were added to the list is the reverse of the order that they
+     will be popped). */
   gal_linkedlist_reverse_stll(&p->columns);
-  gal_table_read_cols(p->up.filename, p->cp.hdu, p->columns,
-                      p->searchin, p->ignorecase);
-
-
+  p->table=gal_table_read_cols(p->up.filename, p->cp.hdu, p->columns,
+                               p->searchin, p->ignorecase, p->cp.minmapsize);
 }
 
 
diff --git a/lib/data.c b/lib/data.c
index e3c0c2f..548f70b 100644
--- a/lib/data.c
+++ b/lib/data.c
@@ -458,7 +458,7 @@ gal_data_initialize(gal_data_t *data, void *array, int type,
 gal_data_t *
 gal_data_alloc(void *array, int type, size_t ndim, long *dsize,
                struct wcsprm *wcs, int clear, size_t minmapsize,
-               char *title, char *unit, char *comment)
+               char *name, char *unit, char *comment)
 {
   gal_data_t *out;
 
@@ -473,7 +473,7 @@ gal_data_alloc(void *array, int type, size_t ndim, long 
*dsize,
 
   /* Initialize the allocated array. */
   gal_data_initialize(out, array, type, ndim, dsize, wcs, clear, minmapsize,
-                      title, unit, comment);
+                      name, unit, comment);
 
 
   /* Return the final structure. */
@@ -499,21 +499,22 @@ gal_data_alloc(void *array, int type, size_t ndim, long 
*dsize,
 void
 gal_data_free(gal_data_t *data, int only_contents)
 {
-  char **strarray=data->array;
-  char *str=strarray[0], *strf=str+data->size;
-
   /* Free all the possible allocations. */
-  free(data->dsize);
-  if(data->name) free(data->name);
-  if(data->unit) free(data->unit);
-  if(data->wcs) wcsfree(data->wcs);
+  if(data->dsize)   free(data->dsize);
+  if(data->name)    free(data->name);
+  if(data->unit)    free(data->unit);
+  if(data->wcs)     wcsfree(data->wcs);
   if(data->comment) free(data->comment);
 
   /* If the data type is string, then each element in the array is actually
      a pointer to the array of characters, so free them before freeing the
      actual array. */
-  if(data->type==GAL_DATA_TYPE_STRING)
-    do if(str) free(str); while(++str<strf);
+  if(data->type==GAL_DATA_TYPE_STRING && data->array)
+    {
+      char **strarray=data->array;
+      char *str=strarray[0], *strf=str+data->size;
+      do if(str) free(str); while(++str<strf);
+    }
 
   /* Free the array. */
   if(data->mmapname)
diff --git a/lib/fits.c b/lib/fits.c
index 255f273..f7166ad 100644
--- a/lib/fits.c
+++ b/lib/fits.c
@@ -1305,7 +1305,7 @@ gal_fits_read_img_hdu(char *filename, char *hdu, char 
*maskname,
   free(blank);
 
 
-  /* Close the FITS file, and return the data pointer */
+  /* Close the input FITS file. */
   fits_close_file(fptr, &status);
   gal_fits_io_error(status, NULL);
 
@@ -1640,17 +1640,105 @@ remove_trailing_space(char *str)
 
 
 
+/* The general format of the TDISPn keywords in FITS is like this: `Tw.p',
+   where `T' specifies the general format, `w' is the width to be given to
+   this column and `p' is the precision. For integer types, percision is
+   actually the minimum number of integers, for floats, it is the number of
+   decimal digits beyond the decimal point.
+
+ */
+static void
+set_display_format(char *tdisp, gal_data_t *data, char *filename, char *hdu,
+                   char *keyname)
+{
+  int isanint;
+  char *tailptr;
+
+  /* First, set the general display format */
+  switch(tdisp[0])
+    {
+    case 'I':
+      isanint=1;
+      data->disp_fmt=GAL_TABLE_DISPLAY_FMT_DECIMAL;
+      break;
+
+    case 'O':
+      isanint=1;
+      data->disp_fmt=GAL_TABLE_DISPLAY_FMT_OCTAL;
+      break;
+
+    case 'Z':
+      isanint=1;
+      data->disp_fmt=GAL_TABLE_DISPLAY_FMT_HEX;
+      break;
+
+    case 'F':
+      isanint=0;
+      data->disp_fmt=GAL_TABLE_DISPLAY_FMT_FLOAT;
+      break;
+
+    case 'E':
+    case 'D':
+      isanint=0;
+      data->disp_fmt=GAL_TABLE_DISPLAY_FMT_EXP;
+      break;
+
+    case 'G':
+      isanint=0;
+      data->disp_fmt=GAL_TABLE_DISPLAY_FMT_GENERAL;
+      break;
+
+    default:
+      error(EXIT_FAILURE, 0, "%s (hdu: %s): Format character `%c' in the "
+            "value (%s) of the keywork %s not recognized", filename, hdu,
+            tdisp[0], tdisp, keyname);
+    }
+
+  /* Parse the rest of the string to see if a width and precision are given
+     or not. */
+  data->disp_width=strtol(&tdisp[1], &tailptr, 0);
+  switch(*tailptr)
+    {
+    case '.':      /* Width is set, go onto finding the precision. */
+      data->disp_precision = strtol(&tailptr[1], &tailptr, 0);
+      if(*tailptr!='\0')
+        error(EXIT_FAILURE, 0, "%s (hdu: %s): The value `%s' of the "
+              "`%s' keyword could not recognized (it doesn't finish after "
+              "the precision)", filename, hdu, tdisp, keyname);
+      break;
+
+    case '\0':     /* No precision given, use a default value.     */
+      data->disp_precision = ( isanint
+                               ? GAL_TABLE_DEF_FLT_PRECISION
+                               : GAL_TABLE_DEF_FLT_PRECISION );
+      break;
+
+    default:
+      error(EXIT_FAILURE, 0, "%s (hdu: %s): The value `%s' of the "
+            "`%s' keyword could not recognized (it doesn't have a `.', or "
+            "finish, after the width)", filename, hdu, tdisp,
+            keyname);
+    }
+
+
+}
+
+
+
+
+
 /* See the descriptions of `gal_table_info'. */
 gal_data_t *
 gal_fits_table_info(char *filename, char *hdu, size_t *numcols,
                     int *tabletype)
 {
+  long repeat;
   int tfields;        /* The maximum number of fields in FITS is 999 */
-  int status=0;
   size_t index;
   fitsfile *fptr;
   size_t i, numrows;
   gal_data_t *cols=NULL;
+  int status=0, datatype;
   char *tailptr, keyname[FLEN_KEYWORD]="XXXXXXXXXXXXX", value[FLEN_VALUE];
 
 
@@ -1670,6 +1758,11 @@ gal_fits_table_info(char *filename, char *hdu, size_t 
*numcols,
           tfields*sizeof *cols);
 
 
+  /* Save the number of rows as the data structure size and also length
+     along the first (and only) dimension. */
+  for(i=0;i<*numcols;++i) cols[i].size=numrows;
+
+
   /* Read all the keywords one by one and if they match, then put them in
      the correct value. Note that we are starting from keyword 9 because
      according to the FITS standard, the first 8 keys in a FITS table are
@@ -1689,20 +1782,25 @@ gal_fits_table_info(char *filename, char *hdu, size_t 
*numcols,
          value[0] is a single quote.*/
       if(strncmp(keyname, "TFORM", 5)==0)
         {
+          /* Remove the ending trailing space and quotation sign. */
+          remove_trailing_space(value);
+          if(*tabletype==GAL_TABLE_TYPE_AFITS)
+            fits_ascii_tform(&value[1], &datatype, NULL, NULL, &status);
+          else
+            fits_binary_tform(&value[1], &datatype, &repeat, NULL, &status);
+
           /* Small sanity check. */
-          if(isdigit(value[1] && value[1]!='1'))
-             error(EXIT_FAILURE, 0, "The repeat value of %s is "
-                   "%c, currently Table can only use columns with a repeat "
-                   "of 1.", keyname, value[1]);
-
-          /* If they don't start with a `1', the values to TFORM are only a
-             single character. so start the pointer to copy at 1 and put
-             the string terminator at 3. */
-          if(value[1]=='1') value[1]=value[2];
-          value[2]='\0';
+          if(repeat>1)
+             error(EXIT_FAILURE, 0, "The repeat value of %s is %ld, "
+                   "currently we can only use columns with a repeat "
+                   "of 1. Please get in touch with us at %s to add this "
+                   "feature", keyname, repeat, PACKAGE_BUGREPORT);
+
+          /* See which column this information was for and add it. In the
+             meantime, also do a sanity check. */
           index = strtoul(&keyname[5], &tailptr, 10) - 1;
           if(index<tfields)     /* Counting from zero was corrected above. */
-            cols[index].type=gal_fits_tform_to_type(value[1]);
+            cols[index].type=gal_fits_datatype_to_type(datatype);
         }
 
       /* COLUMN NAME. All strings in CFITSIO start and finish with single
@@ -1737,11 +1835,74 @@ gal_fits_table_info(char *filename, char *hdu, size_t 
*numcols,
           if(index<tfields)
             gal_checkset_allocate_copy(&value[1], &cols[index].comment);
         }
-    }
 
+      /* COLUMN DISPLAY FORMAT */
+      else if(strncmp(keyname, "TDISP", 5)==0)
+        {
+          /* similar to tname, see above.*/
+          remove_trailing_space(value);
+          index = strtoul(&keyname[5], &tailptr, 10) - 1;
+          if(index<tfields)
+            set_display_format(&value[1], &cols[index], filename, hdu,
+                               keyname);
+        }
+    }
 
   /* Close the FITS file and report an error if we had any. */
   fits_close_file(fptr, &status);
   gal_fits_io_error(status, NULL);
   return cols;
 }
+
+
+
+
+
+/* Read the column indexs given in the `indexll' linked list from a FITS
+   table into a linked list of data structures, note that this is a
+   low-level function, so the output data linked list is the inverse of the
+   input indexs linked list. You can use */
+gal_data_t *
+gal_fits_read_cols(char *filename, char *hdu, gal_data_t *colinfo,
+                   struct gal_linkedlist_sll *indexll, int minmapsize)
+{
+  size_t ind;
+  void *blank;
+  long dsize[1];
+  fitsfile *fptr;
+  int status=0, anynul;
+  gal_data_t *out=NULL, *col;
+
+  /* Open the FITS file */
+  fptr=gal_fits_read_hdu(filename, hdu, 1);
+
+  /* Pop each index and read/store the array. */
+  while(indexll!=NULL)
+    {
+      /* Pop the index. */
+      gal_linkedlist_pop_from_sll(&indexll, &ind);
+
+      /* Allocate the necessary data structure (including the array) for
+         this column. */
+      dsize[0]=colinfo[ind].size;
+      col=gal_data_alloc(NULL, colinfo[ind].type, 1, dsize, NULL, 0,
+                         minmapsize, colinfo[ind].name, colinfo[ind].unit,
+                         colinfo[ind].comment);
+
+      /* Allocate a blank value for the give type and read/store the
+         column using CFITSIO. Afterwards, free the blank value. */
+      blank=gal_data_alloc_blank(col->type);
+      fits_read_col(fptr, gal_fits_type_to_datatype(col->type), ind+1, 1, 1,
+                    col->size, blank, col->array, &anynul, &status);
+      free(blank);
+
+      /* Add the column to the final list of data structures. */
+      col->next=out;
+      out=col;
+    }
+
+  /* Close the FITS file */
+  fits_close_file(fptr, &status);
+  gal_fits_io_error(status, NULL);
+  return out;
+}
diff --git a/lib/gnuastro/data.h b/lib/gnuastro/data.h
index fd30a96..deeea19 100644
--- a/lib/gnuastro/data.h
+++ b/lib/gnuastro/data.h
@@ -200,9 +200,9 @@ enum gal_data_operators
       data structures to follow the same number and decide if they should
       be mmap'd or allocated.
 
-      - `minmapsize' ==0:  array is definitely mmap'd.
+      - `minmapsize' ==  0:  array is definitely mmap'd.
 
-      - `minmapsize' ==-1: array is definitely in RAM.
+      - `minmapsize' == -1: array is definitely in RAM.
 
     - The `dsize' array is in the `long' type because CFITSIO uses the long
       type and this will make it easier to call CFITSIO functions.*/
@@ -227,6 +227,11 @@ typedef struct gal_data_t
   char              *unit;  /* Units of the data.                          */
   char           *comment;  /* A more detailed description of the data.    */
 
+  /* For printing */
+  int            disp_fmt;  /* See `gal_table_diplay_formats'.             */
+  int          disp_width;  /* Width of space to print in ASCII.           */
+  int      disp_precision;  /* Precision to print in ASCII.                */
+
   /* As linked list. */
   struct gal_data_t *next;  /* To use it as a linked list if necessary.    */
 } gal_data_t;
@@ -265,7 +270,7 @@ gal_data_initialize(gal_data_t *data, void *array, int type,
 gal_data_t *
 gal_data_alloc(void *array, int type, size_t ndim, long *dsize,
                struct wcsprm *wcs, int clear, size_t minmapsize,
-               char *title, char *unit, char *comment);
+               char *name, char *unit, char *comment);
 
 void
 gal_data_free(gal_data_t *data, int only_contents);
diff --git a/lib/gnuastro/fits.h b/lib/gnuastro/fits.h
index e9faf32..8ccd595 100644
--- a/lib/gnuastro/fits.h
+++ b/lib/gnuastro/fits.h
@@ -257,6 +257,10 @@ gal_data_t *
 gal_fits_table_info(char *filename, char *hdu, size_t *numcols,
                     int *tabletype);
 
+gal_data_t *
+gal_fits_read_cols(char *filename, char *hdu, gal_data_t *colinfo,
+                   struct gal_linkedlist_sll *indexll, int minmapsize);
+
 
 
 
diff --git a/lib/gnuastro/table.h b/lib/gnuastro/table.h
index 007f8ff..ab620b6 100644
--- a/lib/gnuastro/table.h
+++ b/lib/gnuastro/table.h
@@ -54,6 +54,13 @@ __BEGIN_C_DECLS  /* From C++ preparations */
 
 
 
+/* Macros: */
+
+#define GAL_TABLE_DEF_INT_PRECISION 0
+#define GAL_TABLE_DEF_FLT_PRECISION 8
+
+
+
 /* Types of table storage for input or output. */
 enum gal_table_types
 {
@@ -80,6 +87,24 @@ enum gal_table_where_to_search
 
 
 
+/* Particular display formats for different types: Integers or floating
+   point numbers can be printed in various display formats. The values here
+   are stored in `gal_data_t' when necessary to help in printing of the
+   data.*/
+enum gal_table_diplay_formats
+{
+  GAL_TABLE_DISPLAY_FMT_DECIMAL,        /* Integers: decimal.          */
+  GAL_TABLE_DISPLAY_FMT_OCTAL,          /* Integers: octal.            */
+  GAL_TABLE_DISPLAY_FMT_HEX,            /* Integers: hexadecimal.      */
+  GAL_TABLE_DISPLAY_FMT_FLOAT,          /* Floats: with decimal point. */
+  GAL_TABLE_DISPLAY_FMT_EXP,            /* Floats: as exponential.     */
+  GAL_TABLE_DISPLAY_FMT_GENERAL,        /* Floats: general (%g in C).  */
+};
+
+
+
+
+
 /* Functions */
 gal_data_t *
 gal_table_info(char *filename, char *hdu, size_t *numcols, int *tabletype);
@@ -90,7 +115,7 @@ gal_table_searchin_from_str(char *searchin_str);
 gal_data_t *
 gal_table_read_cols(char *filename, char *hdu,
                     struct gal_linkedlist_stll *cols, int searchin,
-                    int ignorecase);
+                    int ignorecase, int minmapsize);
 
 
 __END_C_DECLS    /* From C++ preparations */
diff --git a/lib/table.c b/lib/table.c
index 3ed5781..b81402c 100644
--- a/lib/table.c
+++ b/lib/table.c
@@ -38,13 +38,14 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 /************************************************************************/
 /* Store the information of each column in a table (either as a text file
    or as a FITS table) into an array of data structures with `numcols'
-   structures (one data structure for each column). Note that the arrays in
-   the data structures will be empty.
-
-   The array of data structures will be returned, the number of columns
-   will be put in the `numcols' argument and the type of the table (e.g.,
-   ascii text file, or FITS binary or ASCII table) will be put in
-   `tabletype' (macros defined in `gnuastro/table.h'. */
+   structures (one data structure for each column). The number of rows is
+   stored as the `size' element of each data structure. The type of the
+   table (e.g., ascii text file, or FITS binary or ASCII table) will be put
+   in `tabletype' (macros defined in `gnuastro/table.h'.
+
+   Note that other than the character strings (column name, units and
+   comments), nothing in the data structure(s) will be allocated by this
+   function for the actual data (e.g., the `array' or `dsize' elements). */
 gal_data_t *
 gal_table_info(char *filename, char *hdu, size_t *numcols, int *tabletype)
 {
@@ -298,10 +299,10 @@ make_list_of_indexs(struct gal_linkedlist_stll *cols, 
gal_data_t *allcols,
 gal_data_t *
 gal_table_read_cols(char *filename, char *hdu,
                     struct gal_linkedlist_stll *cols, int searchin,
-                    int ignorecase)
+                    int ignorecase, int minmapsize)
 {
   int tabletype;
-  size_t numcols;
+  size_t i, numcols;
   gal_data_t *allcols, *out=NULL;
   struct gal_linkedlist_sll *indexll;
 
@@ -312,16 +313,42 @@ gal_table_read_cols(char *filename, char *hdu,
   indexll=make_list_of_indexs(cols, allcols, numcols, searchin,
                               ignorecase, filename, hdu);
 
-  gal_linkedlist_print_sll(indexll);
-
-  /* Reverse the list of indexs so the data structure linked list will be
-     in the same order as the input list of columns.
-  while(indexll!=NULL)
+  /* Depending on the table type, read the columns into the output
+     structure. Note that the functions here pop each index, read/store the
+     desired column and pop the next, so after these functions, the output
+     linked list will have the opposite order of its input `indexll'
+     list. So before calling any of them, we will first reverse the
+     `indexll' list, so the output data structure list will have the same
+     order as the input list of desired columns. Also note that after these
+     functions, the `indexll' will be all freed (each popped element is
+     actually freed).*/
+  gal_linkedlist_reverse_sll(&indexll);
+  switch(tabletype)
     {
+    case GAL_TABLE_TYPE_TXT:
+      error(EXIT_FAILURE, 0, "reading columns not yet implemented for "
+            "plain text files");
+      break;
+
+    case GAL_TABLE_TYPE_AFITS:
+    case GAL_TABLE_TYPE_BFITS:
+      out=gal_fits_read_cols(filename, hdu, allcols, indexll, minmapsize);
+      break;
+
+    default:
+      error(EXIT_FAILURE, 0, "table type code %d not recognized for "
+            "`tabletype' in `gal_table_read_cols'", tabletype);
+    }
 
+  /* Clean up. */
+  for(i=0;i<numcols;++i)
+    {
+      allcols[i].wcs=NULL;
+      allcols[i].dsize=NULL;
+      gal_data_free(&allcols[i], 1);
     }
-  */
-  printf("\n\n-----\ntmp stop.\n\n");
-  exit(0);
+  free(allcols);
+
+  /* Return the final linked list. */
   return out;
 }