gnuastro-commits
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[gnuastro-commits] master 6d69bf8: Convolve program can now convolve 1D,


From: Mohammad Akhlaghi
Subject: [gnuastro-commits] master 6d69bf8: Convolve program can now convolve 1D, for example spectra
Date: Thu, 20 Dec 2018 13:25:58 -0500 (EST)

branch: master
commit 6d69bf8807f1a1b7e9aec29dd280d18ed630d6ae
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>

    Convolve program can now convolve 1D, for example spectra
    
    Gnuastro's spatial domain convolution function is
    dimension-agnostic. However, until now, the convolution program's user
    interface was highly intertwined with a multi-dimensional (2D to be
    particular) dataset. With this commit, that dependency has been broken,
    allowing convolve to smooth 1D datasets also.
---
 NEWS                            |  19 +-
 bin/convertt/ui.c               |   2 +-
 bin/convolve/args.h             |  28 ++-
 bin/convolve/astconvolve.conf   |   1 -
 bin/convolve/convolve.c         |  31 ++-
 bin/convolve/main.h             |   7 +-
 bin/convolve/ui.c               | 424 +++++++++++++++++++++++++++++-----------
 bin/convolve/ui.h               |   8 +-
 bin/match/ui.c                  |   3 +-
 bin/mkprof/ui.c                 |   2 +-
 bin/statistics/ui.c             |   2 +-
 bin/table/ui.c                  |   4 +-
 doc/gnuastro.texi               |  30 ++-
 lib/fits.c                      |   2 +
 lib/gnuastro-internal/options.h |   4 +-
 lib/options.c                   |  26 +--
 16 files changed, 430 insertions(+), 163 deletions(-)

diff --git a/NEWS b/NEWS
index 1043f92..a46813a 100644
--- a/NEWS
+++ b/NEWS
@@ -9,13 +9,20 @@ GNU Astronomy Utilities NEWS                          -*- 
outline -*-
      of output when it is FITS. This is only relevant for some programs,
      (for example not the Fits or Table programs).
    - Standard input (for example from pipes) is now available to feed input
-     to all programs that accept plain text input (ConvertType, Match,
-     MakeProfiles, Statistics, Table).
+     to all programs that accept plain text input (ConvertType, Convolve,
+     Match, MakeProfiles, Statistics, Table).
    - Updated acknowledgement statement (output of `--cite' option).
 
   Arithmetic:
     --onedasimage: write output as an image if it has one dimension, not table.
 
+  Convolve:
+    - Convolves 1D arrays (table columns, for example spectra)
+      also. Therefore two new options have been added to it: `--column'
+      (`-c', similar to other programs that read table columns), and
+      `--kernelcolumn' (to specify the column of the kernel in its own
+      file/extension).
+
   Fits:
     --numhdus: prints the number of HDUs in the given FITS file.
 
@@ -90,6 +97,11 @@ GNU Astronomy Utilities NEWS                          -*- 
outline -*-
       FITS image/array. This can be changed with the new `--onedasimage'
       option.
 
+  Convolve:
+    - The short option for `--numchannels' is now `-n'. Until now, it was
+      `-c', but that would conflict with the short option used for
+      `--column' in all the other programs that also read from a table.
+
   MakeProfiles:
     --mergedsize: new name for the old `--naxis' option. Since the option
           names and values are now written into the FITS header of the
@@ -134,6 +146,9 @@ GNU Astronomy Utilities NEWS                          -*- 
outline -*-
   bug #55157: No sanity check on values given to Crop's --section.
 
 
+
+
+
 * Noteworthy changes in release 0.7 (library 5.0.0) (2018-08-08) [stable]
 
 ** New features
diff --git a/bin/convertt/ui.c b/bin/convertt/ui.c
index ee45765..e9a3afd 100644
--- a/bin/convertt/ui.c
+++ b/bin/convertt/ui.c
@@ -514,7 +514,7 @@ ui_make_channels_ll(struct converttparams *p)
   /* If there weren't any channels, abort with an error. */
   if(p->numch==0)
     error(EXIT_FAILURE, 0, "%s",
-          gal_options_stdin_error(p->cp.stdintimeout, 0));
+          gal_options_stdin_error(p->cp.stdintimeout, 0, "input"));
 
   /* Reverse the list of channels into the input order. */
   gal_list_data_reverse(&p->chll);
diff --git a/bin/convolve/args.h b/bin/convolve/args.h
index f8a17dd..8e03ea6 100644
--- a/bin/convolve/args.h
+++ b/bin/convolve/args.h
@@ -42,7 +42,20 @@ struct argp_option program_options[] =
       &p->kernelname,
       GAL_TYPE_STRING,
       GAL_OPTIONS_RANGE_ANY,
-      GAL_OPTIONS_MANDATORY,
+      GAL_OPTIONS_NOT_MANDATORY,
+      GAL_OPTIONS_NOT_SET
+    },
+    {
+      "column",
+      UI_KEY_COLUMN,
+      "STR",
+      0,
+      "Column name or number if input is a table.",
+      GAL_OPTIONS_GROUP_INPUT,
+      &p->column,
+      GAL_TYPE_STRING,
+      GAL_OPTIONS_RANGE_ANY,
+      GAL_OPTIONS_NOT_MANDATORY,
       GAL_OPTIONS_NOT_SET
     },
     {
@@ -59,6 +72,19 @@ struct argp_option program_options[] =
       GAL_OPTIONS_NOT_SET
     },
     {
+      "kernelcolumn",
+      UI_KEY_KERNELCOLUMN,
+      "STR",
+      0,
+      "Column name or number if kernel is a table.",
+      GAL_OPTIONS_GROUP_INPUT,
+      &p->kernelcolumn,
+      GAL_TYPE_STRING,
+      GAL_OPTIONS_RANGE_ANY,
+      GAL_OPTIONS_NOT_MANDATORY,
+      GAL_OPTIONS_NOT_SET
+    },
+    {
       "nokernelflip",
       UI_KEY_NOKERNELFLIP,
       0,
diff --git a/bin/convolve/astconvolve.conf b/bin/convolve/astconvolve.conf
index 46f5965..1ff2e7d 100644
--- a/bin/convolve/astconvolve.conf
+++ b/bin/convolve/astconvolve.conf
@@ -18,7 +18,6 @@
 # warranty.
 
 # Input:
- kernel             kernel.fits
  khdu               1
 
 # Output:
diff --git a/bin/convolve/convolve.c b/bin/convolve/convolve.c
index fae4a18..4e05635 100644
--- a/bin/convolve/convolve.c
+++ b/bin/convolve/convolve.c
@@ -761,6 +761,7 @@ void
 convolve(struct convolveparams *p)
 {
   gal_data_t *out, *check;
+  int multidim=p->input->ndim>1;
   struct gal_options_common_params *cp=&p->cp;
 
 
@@ -768,10 +769,10 @@ convolve(struct convolveparams *p)
   if(p->domain==CONVOLVE_DOMAIN_SPATIAL)
     {
       /* Prepare the mesh structure. */
-      gal_tile_full_two_layers(p->input, &cp->tl);
+      if(multidim) gal_tile_full_two_layers(p->input, &cp->tl);
 
       /* Save the tile IDs if they are requested. */
-      if(cp->tl.tilecheckname)
+      if(multidim && cp->tl.tilecheckname)
         {
           check=gal_tile_block_check_tiles(cp->tl.tiles);
           gal_fits_img_write(check, cp->tl.tilecheckname, NULL, PROGRAM_NAME);
@@ -782,12 +783,14 @@ convolve(struct convolveparams *p)
          want to do spatial domain convolution with this Convolve program
          is edge correction. So by default we assume it and will only
          ignore it if the user asks.*/
-      out=gal_convolve_spatial(cp->tl.tiles, p->kernel, cp->numthreads,
-                               !p->noedgecorrection, cp->tl.workoverch);
+      out=gal_convolve_spatial(multidim ? cp->tl.tiles : p->input, p->kernel,
+                               cp->numthreads,
+                               multidim ? !p->noedgecorrection : 1,
+                               multidim ? cp->tl.workoverch : 1 );
 
       /* Clean up: free the actual input and replace it's pointer with the
          convolved dataset to save as output. */
-      gal_tile_full_free_contents(&cp->tl);
+      if(&cp->tl) gal_tile_full_free_contents(&cp->tl);
       gal_data_free(p->input);
       p->input=out;
     }
@@ -795,12 +798,20 @@ convolve(struct convolveparams *p)
     convolve_frequency(p);
 
   /* Save the output (which is in p->input) array. */
-  gal_fits_img_write_to_type(p->input, cp->output, NULL, PROGRAM_NAME,
-                             cp->type);
+  if(p->input->ndim==1)
+    gal_table_write(p->input, NULL, p->cp.tableformat, p->cp.output,
+                    "CONVOLVED", 0);
+  else
+    gal_fits_img_write_to_type(p->input, cp->output, NULL, PROGRAM_NAME,
+                               cp->type);
 
   /* Write Convolve's parameters as keywords into the first extension of
      the output. */
-  gal_fits_key_write_filename("input", p->filename, &cp->okeys, 1);
-  gal_fits_key_write_config(&cp->okeys, "Convolve configuration",
-                            "CONVOLVE-CONFIG", cp->output, "0");
+  if( gal_fits_name_is_fits(p->cp.output) )
+    {
+      gal_fits_key_write_filename("input", p->filename, &cp->okeys, 1);
+      gal_fits_key_write_config(&cp->okeys, "Convolve configuration",
+                                "CONVOLVE-CONFIG", cp->output, "0");
+    }
+  printf("  - Output: %s\n", p->cp.output);
 }
diff --git a/bin/convolve/main.h b/bin/convolve/main.h
index 11c0f53..f3b8812 100644
--- a/bin/convolve/main.h
+++ b/bin/convolve/main.h
@@ -39,7 +39,7 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 
 /* Macros */
 #define CONVFLOATINGPOINTERR 1e-10
-
+#define INPUT_USE_TYPE       GAL_TYPE_FLOAT32
 
 
 
@@ -54,7 +54,6 @@ enum complex_to_real
   COMPLEX_TO_REAL_REAL,
 };
 
-
 enum domain_codes
 {
   CONVOLVE_DOMAIN_INVALID,           /* ==0 by C standard. */
@@ -73,8 +72,10 @@ struct convolveparams
   /* From command-line */
   struct gal_options_common_params cp; /* Common parameters.              */
   char             *filename;  /* Name of input file.                     */
+  char               *column;  /* Name of column if input is a table.     */
   char           *kernelname;  /* File name of kernel.                    */
   char                 *khdu;  /* HDU of kernel.                          */
+  char         *kernelcolumn;  /* Column to read the input kernel.        */
   uint8_t       nokernelflip;  /* Do not flip the kernel.                 */
   uint8_t       nokernelnorm;  /* Do not normalize the kernel.            */
   double        minsharpspec;  /* Deconvolution: min spect. of sharp img. */
@@ -84,6 +85,8 @@ struct convolveparams
   uint8_t   noedgecorrection;  /* Do not correct spatial edge effects.    */
 
   /* Internal */
+  int                 isfits;  /* Input is a FITS file.                   */
+  int               hdu_type;  /* Type of HDU (image or table).           */
   int                 domain;  /* Frequency or spatial domain conv.       */
   gal_data_t          *input;  /* Input image array.                      */
   gal_data_t         *kernel;  /* Input Kernel array.                     */
diff --git a/bin/convolve/ui.c b/bin/convolve/ui.c
index f8ed0e6..5495b36 100644
--- a/bin/convolve/ui.c
+++ b/bin/convolve/ui.c
@@ -129,10 +129,7 @@ ui_initialize_options(struct convolveparams *p,
         break;
 
       case GAL_OPTIONS_KEY_LOG:
-      case GAL_OPTIONS_KEY_SEARCHIN:
       case GAL_OPTIONS_KEY_IGNORECASE:
-      case GAL_OPTIONS_KEY_TABLEFORMAT:
-      case GAL_OPTIONS_KEY_STDINTIMEOUT:
       case GAL_OPTIONS_KEY_INTERPNUMNGB:
       case GAL_OPTIONS_KEY_INTERPONLYBLANK:
         cp->coptions[i].flags=OPTION_HIDDEN;
@@ -213,11 +210,6 @@ ui_read_check_only_options(struct convolveparams *p)
 {
   struct gal_options_common_params *cp=&p->cp;
 
-  /* Make sure the kernel name is a FITS file and a HDU is given. */
-  if( gal_fits_name_is_fits(p->kernelname)==0 )
-    error(EXIT_FAILURE, 0, "`%s' is not a recognized FITS file name",
-          p->kernelname);
-
 
   /* Read the domain from a string into an integer. */
   if( !strcmp("spatial", p->domainstr) )
@@ -255,15 +247,53 @@ ui_read_check_only_options(struct convolveparams *p)
 static void
 ui_check_options_and_arguments(struct convolveparams *p)
 {
-  /* Make sure an input file name was given and if it was a FITS file, that
-     a HDU is also given. */
-  if(p->filename==NULL)
-    error(EXIT_FAILURE, 0, "no input file is specified");
-
-  /* Make sure the input name is a FITS file name. */
-  if( gal_fits_name_is_fits(p->filename)==0 )
-    error(EXIT_FAILURE, 0, "`%s' is not a recognized FITS file name",
-          p->filename);
+  int kernel_type;
+
+  if(p->filename)
+    {
+      /* If input is FITS. */
+      if( (p->isfits=gal_fits_name_is_fits(p->filename)) )
+        {
+          /* Check if a HDU is given. */
+          if( p->cp.hdu==NULL )
+            error(EXIT_FAILURE, 0, "no HDU specified. When the input is a "
+                  "FITS file, a HDU must also be specified, you can use "
+                  "the `--hdu' (`-h') option and give it the HDU number "
+                  "(starting from zero), extension name, or anything "
+                  "acceptable by CFITSIO");
+
+          /* If its an image, make sure column isn't given (in case the
+             user confuses an image with a table). */
+          p->hdu_type=gal_fits_hdu_format(p->filename, p->cp.hdu);
+          if(p->hdu_type==IMAGE_HDU && p->column)
+            error(EXIT_FAILURE, 0, "%s (hdu: %s): is a FITS image "
+                  "extension. The `--column' option is only applicable "
+                  "to tables.", p->filename, p->cp.hdu);
+        }
+    }
+
+  if(p->kernelname)
+    {
+      /* If input is FITS. */
+      if( gal_fits_name_is_fits(p->kernelname) )
+        {
+          /* Check if a HDU is given. */
+          if( p->khdu==NULL )
+            error(EXIT_FAILURE, 0, "no HDU specified. When the kernel is a "
+                  "FITS file, a HDU must also be specified, you can use "
+                  "the `--khdu' (`-u') option and give it the HDU number "
+                  "(starting from zero), extension name, or anything "
+                  "acceptable by CFITSIO");
+
+          /* If its an image, make sure column isn't given (in case the
+             user confuses an image with a table). */
+          kernel_type=gal_fits_hdu_format(p->kernelname, p->khdu);
+          if(kernel_type==IMAGE_HDU && p->kernelcolumn)
+            error(EXIT_FAILURE, 0, "%s (hdu: %s): is a FITS image "
+                  "extension. The `--kernelcolumn' option is only "
+                  "applicable to tables.", p->kernelname, p->khdu);
+        }
+    }
 }
 
 
@@ -288,27 +318,158 @@ ui_check_options_and_arguments(struct convolveparams *p)
 /**************************************************************/
 /***************       Preparations         *******************/
 /**************************************************************/
+static gal_data_t *
+ui_read_column(struct convolveparams *p, int i0k1)
+{
+  int tformat;
+  char *source;
+  gal_data_t *out, *cinfo;
+  gal_list_str_t *column=NULL;
+  size_t ncols, nrows, counter=0;
+  char *hdu        = i0k1==0 ? p->cp.hdu   : p->khdu;
+  char *name       = i0k1==0 ? "input"     : "kernel";
+  char *filename   = i0k1==0 ? p->filename : p->kernelname;
+  char *columnname = i0k1==0 ? p->column   : p->kernelcolumn;
+  gal_list_str_t *lines = gal_options_check_stdin(filename,
+                                                  p->cp.stdintimeout, name);
+
+  /* If no column is specified, Convolve will abort and an error will be
+     printed when the table has more than one column. If there is only one
+     column, there is no need to specify any, so Convolve will use it. */
+  if(columnname==NULL)
+    {
+      /* Get the basic table information. */
+      cinfo=gal_table_info(filename, hdu, lines, &ncols, &nrows, &tformat);
+      gal_data_array_free(cinfo, ncols, 1);
+
+      /* See how many columns it has and take the proper action. */
+      switch(ncols)
+        {
+        case 0:
+          error(EXIT_FAILURE, 0, "%s contains no usable information",
+                ( filename
+                  ? gal_checkset_dataset_name(filename, hdu)
+                  : "Standard input" ));
+        case 1:
+          gal_checkset_allocate_copy("1", &columnname);
+          break;
+        default:
+          error(EXIT_FAILURE, 0, "%s is a table containing more than one "
+                "column. However, the specific column to work on isn't "
+                "specified.\n\n"
+                "Please use the `--column' (`-c') or `--kernelcolumn' "
+                "options (depending on which dataset it is) to specify a "
+                "column. You can either give it the column number "
+                "(couting from 1), or a match/search in its meta-data (e.g., "
+                "column names).\n\n"
+                "For more information, please run the following command "
+                "(press the `SPACE' key to go down and `q' to return to the "
+                "command-line):\n\n"
+                "    $ info gnuastro \"Selecting table columns\"\n",
+                ( filename
+                  ? gal_checkset_dataset_name(filename, hdu)
+                  : "Standard input" ));
+        }
+
+    }
+  gal_list_str_add(&column, columnname, 0);
+
+  /* Read the desired column(s). */
+  out=gal_table_read(filename, hdu, lines, column, p->cp.searchin,
+                     p->cp.ignorecase, p->cp.minmapsize, NULL);
+  gal_list_str_free(lines, 1);
+
+  /* Confirm if only one column was read (it is possible to match more than
+     one column). */
+  if(out->next!=NULL)
+    {
+      if(filename)
+        gal_checkset_dataset_name(filename, hdu);
+      else
+        source="standard-input";
+      error(EXIT_FAILURE, 0, "%s: more than one column in input table mached "
+            "the search criteria. Please limit the match by specifying the "
+            "exact name (if its unique) or column number", source);
+    }
+
+  /* Make sure it is a usable datatype. */
+  switch(out->type)
+    {
+    case GAL_TYPE_BIT:
+    case GAL_TYPE_STRLL:
+    case GAL_TYPE_STRING:
+    case GAL_TYPE_COMPLEX32:
+    case GAL_TYPE_COMPLEX64:
+      error(EXIT_FAILURE, 0, " read column number %zu has a %s type, "
+            "which is not currently supported by %s", counter,
+            gal_type_name(out->type, 1), PROGRAM_NAME);
+    }
+  out=gal_data_copy_to_new_type_free(out, INPUT_USE_TYPE);
+
+  /* If the input was from standard input, we can actually write this into
+     it (for future reporting). */
+  if(filename==NULL)
+    gal_checkset_allocate_copy("standard-input",
+                               ( i0k1==0 ? &p->filename : &p->kernelname ));
+
+  /* Clean up and return. */
+  gal_list_str_free(column, 0);
+  return out;
+}
+
+
+
+
+
+/* Read the input dataset. */
+static void
+ui_read_input(struct convolveparams *p)
+{
+  /* To see if we should read it as a table. */
+  p->input=NULL;
+
+  /* If the input is a FITS image or any recognized array file format, then
+     read it as an array, otherwise, as a table. */
+  if( p->filename && gal_array_name_recognized(p->filename) )
+    if (p->isfits && p->hdu_type==IMAGE_HDU)
+      {
+        p->input=gal_array_read_one_ch_to_type(p->filename, p->cp.hdu, NULL,
+                                               INPUT_USE_TYPE,
+                                               p->cp.minmapsize);
+        p->input->wcs=gal_wcs_read(p->filename, p->cp.hdu, 0, 0,
+                                   &p->input->nwcs);
+      }
+
+  /* The input isn't an image (wasn't read yet), so we'll read it as a
+     column. */
+  if(p->input==NULL)
+    p->input=ui_read_column(p, 0);
+}
+
+
+
+
+
 /* Read the kernel. VERY IMPORTANT: We can't use the `fits_img_read_kernel'
    because the Convolve program also does de-convolution. */
 static void
 ui_read_kernel(struct convolveparams *p)
 {
-  float *f, *ff;
-
   /* Read the image into file. */
-  p->kernel = gal_array_read_one_ch_to_type(p->kernelname, p->khdu,
-                                            NULL, GAL_TYPE_FLOAT32,
-                                            p->cp.minmapsize);
+  if( p->kernelname
+      && p->input->ndim>1
+      && gal_array_name_recognized(p->kernelname)  )
+    p->kernel = gal_array_read_one_ch_to_type(p->kernelname, p->khdu,
+                                              NULL, INPUT_USE_TYPE,
+                                              p->cp.minmapsize);
+  else
+    p->kernel=ui_read_column(p, 1);
 
-  /* Convert all the NaN pixels to zero if the kernel contains blank
-     pixels, also update the flags so it is not checked any more. */
-  if(gal_blank_present(p->kernel, 1))
-    {
-      ff = (f=p->kernel->array) + p->kernel->size;
-      do *f = isnan(*f) ? 0.0f : *f; while(++f<ff);
-      p->kernel->flag |= GAL_DATA_FLAG_BLANK_CH;
-      p->kernel->flag &= ~GAL_DATA_FLAG_HASBLANK;
-    }
+  /* Make sure that the kernel and input have the same number of
+     dimensions. */
+  if(p->kernel->ndim!=p->input->ndim)
+    error(EXIT_FAILURE, 0, "input datasets must have the same number of "
+          "dimensions");
 }
 
 
@@ -318,45 +479,27 @@ ui_read_kernel(struct convolveparams *p)
 static void
 ui_preparations(struct convolveparams *p)
 {
+  int check=0;
+  double sumv=0;
   size_t i, size;
   gal_data_t *sum;
-  float *kernel, tmp;
+  float *f, *fp, tmp, *kernel;
   struct gal_options_common_params *cp=&p->cp;
   char *outsuffix = p->makekernel ? "_kernel.fits" : "_convolved.fits";
 
-  /* Set the output name if the user hasn't set it. */
-  if(cp->output==NULL)
-    cp->output=gal_checkset_automatic_output(cp, p->filename, outsuffix);
-  gal_checkset_writable_remove(cp->output, 0, cp->dontdelete);
-  if(p->checkfreqsteps)
-    {
-      p->freqstepsname=gal_checkset_automatic_output(cp, p->filename,
-                                                     "_freqsteps.fits");
-      gal_checkset_writable_remove(p->freqstepsname, 0, cp->dontdelete);
-    }
-  if(cp->tl.checktiles)
-    {
-      cp->tl.tilecheckname=gal_checkset_automatic_output(cp, p->filename,
-                                                         "_tiled.fits");
-      gal_checkset_writable_remove(cp->tl.tilecheckname, 0,
-                                   cp->dontdelete);
-    }
-
 
-  /* Read the input image as a float64 array. */
-  p->input=gal_array_read_one_ch_to_type(p->filename, cp->hdu, NULL,
-                                         GAL_TYPE_FLOAT32, cp->minmapsize);
-  p->input->wcs=gal_wcs_read(p->filename, cp->hdu, 0, 0, &p->input->nwcs);
+  /* Read the input dataset. */
+  ui_read_input(p);
 
 
-  /* Currently Convolve only works on 2D images. */
-  if(p->input->ndim!=2)
+  /* Currently Convolve only works on 1D and 2D datasets. */
+  if(p->input->ndim>2)
     error(EXIT_FAILURE, 0, "%s (hdu %s) has %zu dimensions. Currently "
-          "Convolve only operates on 2D images", p->filename, cp->hdu,
-          p->input->ndim);
+          "Convolve only operates on 1D (table column) and 2D (image) "
+          "datasets", p->filename, cp->hdu, p->input->ndim);
 
 
-  /* See if there are any blank values. */
+  /* Domain-specific checks. */
   if(p->domain==CONVOLVE_DOMAIN_FREQUENCY)
     {
       if( gal_blank_present(p->input, 1) )
@@ -371,10 +514,18 @@ ui_preparations(struct convolveparams *p)
                 "----------------------------------------\n\n",
                 PROGRAM_NAME, p->filename, cp->hdu, cp->output,
                 PROGRAM_NAME);
+
+      /* Frequency domain is only implemented in 2D. */
+      if( p->input->ndim==1 )
+        error(EXIT_FAILURE, 0, "Frequency domain convolution is currently "
+              "not implemented on 1D datasets. Please use `--domain=spatial' "
+              "to convolve this dataset");
     }
   else
-    gal_tile_full_sanity_check(p->filename, cp->hdu, p->input, &cp->tl);
-
+    {
+      if(p->input->ndim>1)
+        gal_tile_full_sanity_check(p->filename, cp->hdu, p->input, &cp->tl);
+    }
 
 
   /* Read the file specified by --kernel. If makekernel is specified, then
@@ -382,31 +533,39 @@ ui_preparations(struct convolveparams *p)
      argument) is the blurry image. */
   if(p->makekernel)
     {
-      /* Read in the kernel array. */
-      ui_read_kernel(p);
+      /* Currently this is not implemented in 1D. */
+      if(p->kernel->ndim==1)
+        error(EXIT_FAILURE, 0, "`--makekernel' is currently not available "
+              "on 1D datasets");
+      else
+        {
+          /* Read in the kernel array. */
+          ui_read_kernel(p);
 
-      /* Make sure the size of the kernel is the same as the input */
-      if( p->input->dsize[0]!=p->kernel->dsize[0]
-          || p->input->dsize[1]!=p->kernel->dsize[1] )
-        error(EXIT_FAILURE, 0, "with the `--makekernel' (`-m') option, "
-              "the input image and the image specified with the `--kernel' "
-              "(`-k') option should have the same size. The lower resolution "
-              "input image (%s) has %zux%zu pixels while the sharper image "
-              "(%s) specified with the kernel option has %zux%zu pixels",
-              p->filename, p->input->dsize[1], p->input->dsize[0],
-              p->kernelname, p->kernel->dsize[1], p->kernel->dsize[0]);
-
-      /* Divide both images by their sum so their lowest frequency becomes
-         1 and their division (in the frequency domain) would be
-         meaningful. */
-      sum=gal_statistics_sum(p->input);
-      sum=gal_data_copy_to_new_type_free(sum, GAL_TYPE_FLOAT32);
-      p->input = gal_arithmetic(GAL_ARITHMETIC_OP_DIVIDE,
-                                GAL_ARITHMETIC_FLAGS_ALL, p->input, sum);
-      sum=gal_statistics_sum(p->kernel);
-      sum=gal_data_copy_to_new_type_free(sum, GAL_TYPE_FLOAT32);
-      p->kernel = gal_arithmetic(GAL_ARITHMETIC_OP_DIVIDE,
-                                GAL_ARITHMETIC_FLAGS_ALL, p->kernel, sum);
+          /* Make sure the size of the kernel is the same as the input */
+          if( p->input->dsize[0]!=p->kernel->dsize[0]
+              || p->input->dsize[1]!=p->kernel->dsize[1] )
+            error(EXIT_FAILURE, 0, "with the `--makekernel' (`-m') option, "
+                  "the input image and the image specified with the "
+                  "`--kernel' (`-k') option should have the same size. The "
+                  "lower resolution input image (%s) has %zux%zu pixels "
+                  "while the sharper image (%s) specified with the kernel "
+                  "option has %zux%zu pixels", p->filename,
+                  p->input->dsize[1], p->input->dsize[0], p->kernelname,
+                  p->kernel->dsize[1], p->kernel->dsize[0]);
+
+          /* Divide both images by their sum so their lowest frequency becomes
+             1 and their division (in the frequency domain) would be
+             meaningful. */
+          sum=gal_statistics_sum(p->input);
+          sum=gal_data_copy_to_new_type_free(sum, GAL_TYPE_FLOAT32);
+          p->input = gal_arithmetic(GAL_ARITHMETIC_OP_DIVIDE,
+                                    GAL_ARITHMETIC_FLAGS_ALL, p->input, sum);
+          sum=gal_statistics_sum(p->kernel);
+          sum=gal_data_copy_to_new_type_free(sum, GAL_TYPE_FLOAT32);
+          p->kernel = gal_arithmetic(GAL_ARITHMETIC_OP_DIVIDE,
+                                     GAL_ARITHMETIC_FLAGS_ALL, p->kernel, sum);
+        }
     }
 
   /* Read the kernel. If there is anything particular to Convolve, then
@@ -416,44 +575,71 @@ ui_preparations(struct convolveparams *p)
      into one loop. */
   else
     {
-      if(p->nokernelnorm || p->nokernelflip)
-        {
-          /* Read in the kernel array: */
-          ui_read_kernel(p);
-
-          /* Check its size (must be odd). */
-          if(p->kernel->dsize[0]%2==0 || p->kernel->dsize[1]%2==0)
-            error(EXIT_FAILURE, 0, "the kernel image has to have an odd "
-                  "number of pixels on both sides (there has to be on pixel "
-                  "in the center). %s (hdu: %s) is %zu by %zu",
-                  p->kernelname, p->khdu, p->kernel->dsize[1],
-                  p->kernel->dsize[0]);
+      /* Read in the kernel array: */
+      ui_read_kernel(p);
 
-          /* Normalize the kernel: */
-          if( !p->nokernelnorm )
+      /* Check if the size along each dimension of the kernel is an odd
+         number. If they are all an odd number, then the for each dimension,
+         check will be incremented once. */
+      for(i=0;i<p->kernel->ndim;++i)
+        check += p->kernel->dsize[i]%2;
+      if(check!=p->kernel->ndim)
+        error(EXIT_FAILURE, 0, "%s: the kernel has to have an odd number of "
+              "elements in all dimensions (there has to be one element/pixel "
+              "in the center). At least one of the dimensions of doesn't "
+              "have an odd number of pixels",
+              gal_checkset_dataset_name(p->kernelname, p->khdu));
+
+
+      /* If there are any NaN pixels, set them to zero and normalize it. A
+         blank pixel in a kernel is going to make a completely blank
+         output.*/
+      if( !p->nokernelnorm )
+        {
+          sumv=0;
+          fp=(f=p->kernel->array)+p->kernel->size;
+          do
             {
-              sum=gal_statistics_sum(p->kernel);
-              p->kernel = gal_arithmetic(GAL_ARITHMETIC_OP_DIVIDE,
-                                         GAL_ARITHMETIC_FLAGS_ALL,
-                                         p->kernel, sum);
+              if(isnan(*f)) *f=0.0f;
+              else          sumv+=*f;
             }
+          while(++f<fp);
+          p->kernel->flag |= GAL_DATA_FLAG_BLANK_CH;
+          p->kernel->flag &= ~GAL_DATA_FLAG_HASBLANK;
+          f=p->kernel->array; do *f++ *= 1/sumv; while(f<fp);
+        }
 
-          /* Flip the kernel: */
-          if( !p->nokernelflip )
+      /* Flip the kernel: */
+      if( !p->nokernelflip )
+        {
+          size=p->kernel->size;
+          kernel=p->kernel->array;
+          for(i=0;i<p->kernel->size/2;++i)
             {
-              size=p->kernel->size;
-              kernel=p->kernel->array;
-              for(i=0;i<p->kernel->size/2;++i)
-                {
-                  tmp                     = kernel[ i            ];
-                  kernel[ i             ] = kernel[ size - i - 1 ];
-                  kernel[ size -  i - 1 ] = tmp;
-                }
+              tmp                     = kernel[ i            ];
+              kernel[ i             ] = kernel[ size - i - 1 ];
+              kernel[ size -  i - 1 ] = tmp;
             }
         }
-      else
-        p->kernel = gal_fits_img_read_kernel(p->kernelname, p->khdu,
-                                             cp->minmapsize);
+    }
+
+
+  /* Set the output name if the user hasn't set it. */
+  if(cp->output==NULL)
+    cp->output=gal_checkset_automatic_output(cp, p->filename, outsuffix);
+  gal_checkset_writable_remove(cp->output, 0, cp->dontdelete);
+  if(p->checkfreqsteps)
+    {
+      p->freqstepsname=gal_checkset_automatic_output(cp, p->filename,
+                                                     "_freqsteps.fits");
+      gal_checkset_writable_remove(p->freqstepsname, 0, cp->dontdelete);
+    }
+  if(cp->tl.checktiles)
+    {
+      cp->tl.tilecheckname=gal_checkset_automatic_output(cp, p->filename,
+                                                         "_tiled.fits");
+      gal_checkset_writable_remove(cp->tl.tilecheckname, 0,
+                                   cp->dontdelete);
     }
 }
 
@@ -483,8 +669,10 @@ ui_print_intro(struct convolveparams *p)
 {
   printf("%s started on %s", PROGRAM_NAME, ctime(&p->rawtime));
   printf("  - Using %zu CPU threads.\n", p->cp.numthreads);
-  printf("  - Input: %s (hdu: %s)\n", p->filename, p->cp.hdu);
-  printf("  - Kernel: %s (hdu: %s)\n", p->kernelname, p->khdu);
+  printf("  - Input: %s\n",
+         gal_checkset_dataset_name(p->filename, p->cp.hdu));
+  printf("  - Kernel: %s\n",
+         gal_checkset_dataset_name(p->kernelname, p->khdu));
 }
 
 
diff --git a/bin/convolve/ui.h b/bin/convolve/ui.h
index 947e1a3..81e67f9 100644
--- a/bin/convolve/ui.h
+++ b/bin/convolve/ui.h
@@ -42,7 +42,7 @@ enum program_args_groups
 
 /* Available letters for short options:
 
-   a b e f g i j l n p s v w x y z
+   a b e f g i j l p s v w x y z
    A B E G J L O Q R W X Y
 */
 enum option_keys_enum
@@ -53,14 +53,16 @@ enum option_keys_enum
   UI_KEY_MINSHARPSPEC   = 'H',
   UI_KEY_CHECKFREQSTEPS = 'C',
   UI_KEY_TILESIZE       = 't',
-  UI_KEY_NUMCHANNELS    = 'c',
+  UI_KEY_COLUMN         = 'c',
+  UI_KEY_NUMCHANNELS    = 'n',
   UI_KEY_REMAINDERFRAC  = 'r',
   UI_KEY_DOMAIN         = 'd',
   UI_KEY_MAKEKERNEL     = 'm',
 
   /* Only with long version (start with a value 1000, the rest will be set
      automatically). */
-  UI_KEY_NOKERNELFLIP = 1000,
+  UI_KEY_KERNELCOLUMN = 1000,
+  UI_KEY_NOKERNELFLIP,
   UI_KEY_NOKERNELNORM,
   UI_KEY_NOEDGECORRECTION,
 };
diff --git a/bin/match/ui.c b/bin/match/ui.c
index c7a81e6..dae5ca8 100644
--- a/bin/match/ui.c
+++ b/bin/match/ui.c
@@ -423,7 +423,8 @@ ui_read_columns_to_double(struct matchparams *p, char 
*filename, char *hdu,
      once, we don't want to write a blank list over it (the Standard input
      will be empty after being read). */
   if(p->stdinlines==NULL)
-    p->stdinlines=gal_options_check_stdin(filename, p->cp.stdintimeout);
+    p->stdinlines=gal_options_check_stdin(filename, p->cp.stdintimeout,
+                                          "input");
   tout=gal_table_read(filename, hdu, filename ? NULL : p->stdinlines,
                       cols, cp->searchin, cp->ignorecase, cp->minmapsize,
                       NULL);
diff --git a/bin/mkprof/ui.c b/bin/mkprof/ui.c
index 590eda2..16c3571 100644
--- a/bin/mkprof/ui.c
+++ b/bin/mkprof/ui.c
@@ -620,7 +620,7 @@ ui_read_cols(struct mkprofparams *p)
   gal_list_str_reverse(&colstrs);
 
   /* Read the desired columns from the file. */
-  lines=gal_options_check_stdin(p->catname, p->cp.stdintimeout);
+  lines=gal_options_check_stdin(p->catname, p->cp.stdintimeout, "input");
   cols=gal_table_read(p->catname, p->cp.hdu, lines, colstrs, p->cp.searchin,
                       p->cp.ignorecase, p->cp.minmapsize, NULL);
   gal_list_str_free(lines, 1);
diff --git a/bin/statistics/ui.c b/bin/statistics/ui.c
index 295753e..2774022 100644
--- a/bin/statistics/ui.c
+++ b/bin/statistics/ui.c
@@ -701,7 +701,7 @@ ui_read_columns(struct statisticsparams *p)
   gal_data_t *cols, *tmp, *cinfo;
   size_t size, ncols, nrows, counter=0;
   gal_list_str_t *lines=gal_options_check_stdin(p->inputname,
-                                                p->cp.stdintimeout);
+                                                p->cp.stdintimeout, "input");
 
   /* If a reference column is also given, add it to the list of columns to
      read. */
diff --git a/bin/table/ui.c b/bin/table/ui.c
index b94cc96..1612aa2 100644
--- a/bin/table/ui.c
+++ b/bin/table/ui.c
@@ -276,7 +276,7 @@ ui_print_info_exit(struct tableparams *p)
   size_t i, numcols, numrows;
 
   /* Read the table information for the number of columns and rows. */
-  lines=gal_options_check_stdin(p->filename, p->cp.stdintimeout);
+  lines=gal_options_check_stdin(p->filename, p->cp.stdintimeout, "input");
   allcols=gal_table_info(p->filename, p->cp.hdu, lines, &numcols,
                          &numrows, &tableformat);
   if(p->filename==NULL) p->filename="Standard-input";
@@ -374,7 +374,7 @@ ui_preparations(struct tableparams *p)
   ui_columns_prepare(p);
 
   /* Read in the table columns. */
-  lines=gal_options_check_stdin(p->filename, p->cp.stdintimeout);
+  lines=gal_options_check_stdin(p->filename, p->cp.stdintimeout, "input");
   p->table=gal_table_read(p->filename, cp->hdu, lines, p->columns,
                           cp->searchin, cp->ignorecase, cp->minmapsize,
                           NULL);
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 8491f25..2feed22 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -13933,8 +13933,9 @@ value.
 @node Invoking astconvolve,  , Convolution kernel, Convolve
 @subsection Invoking Convolve
 
-Convolve an input image with a known kernel or make the kernel necessary to
-match two PSFs. The general template for Convolve is:
+Convolve an input dataset (2D image or 1D spectrum for example) with a
+known kernel, or make the kernel necessary to match two PSFs. The general
+template for Convolve is:
 
 @example
 $ astconvolve [OPTION...] ASTRdata
@@ -13953,15 +13954,27 @@ $ astconvolve observedimg.fits --kernel=psf.fits 
--domain=spatial
 ## Find the kernel to match sharper and blurry PSF images:
 $ astconvolve --kernel=sharperimage.fits --makekernel=10           \
               blurryimage.fits
+
+## Convolve a Spectrum (column 14 in the FITS table below) with a
+## custom kernel (the kernel will be normalized internally, so only
+## the ratios are important). Sed is used to replace the spaces with
+## new line characters so Convolve sees them as values in one column.
+$ echo "1 3 10 3 1" | sed 's/ /\n/g' | astconvolve spectra.fits -c14
 @end example
 
 The only argument accepted by Convolve is an input image file. Some of the
 options are the same between Convolve and some other Gnuastro
 programs. Therefore, to avoid repetition, they will not be repeated
 here. For the full list of options shared by all Gnuastro programs, please
-see @ref{Common options}. In particular, in the spatial domain convolve
-uses Gnuastro's tessellation, see @ref{Tessellation} and the common options
-related to that in @ref{Processing options}.
+see @ref{Common options}. In particular, in the spatial domain, on a
+multi-dimensional datasets, convolve uses Gnuastro's tessellation to speed
+up the run, see @ref{Tessellation}. Common options related to tessellation
+are described in in @ref{Processing options}.
+
+1-dimensional datasets (for example spectra) are only read as columns
+within a table (see @ref{Tables} for more on how Gnuastro programs read
+tables). Note that currently 1D convolution is only implemented in the
+spatial domain and thus kernel-matching is also not supported.
 
 Here we will only explain the options particular to Convolve. Run Convolve
 with @option{--help} in order to see the full list of options Convolve
@@ -13969,6 +13982,11 @@ accepts, irrespective of where they are explained in 
this book.
 
 @table @option
 
address@hidden --kernelcolumn
+Column containing the 1D kernel. When the input dataset is a 1-dimensional
+column, and the host table has more than one column, use this option to
+specify which column should be used.
+
 @item --nokernelflip
 Do not flip the kernel after reading it the spatial domain
 convolution. This can be useful if the flipping has already been
@@ -14085,6 +14103,8 @@ direction, so crop the central pixel of the warped 
image to have a
 clean image for the de-convolution.
 @end itemize
 
+Note that this feature is not yet supported in 1-dimensional datasets.
+
 @item -c
 @itemx --minsharpspec
 (@option{=FLT}) The minimum frequency spectrum (or coefficient, or pixel
diff --git a/lib/fits.c b/lib/fits.c
index 9e6f77e..5f724c0 100644
--- a/lib/fits.c
+++ b/lib/fits.c
@@ -1781,6 +1781,8 @@ gal_fits_img_read_kernel(char *filename, char *hdu, 
size_t minmapsize)
       else          sum+=*f;
     }
   while(++f<fp);
+  kernel->flag |= GAL_DATA_FLAG_BLANK_CH;
+  kernel->flag &= ~GAL_DATA_FLAG_HASBLANK;
   f=kernel->array; do *f++ *= 1/sum; while(f<fp);
 
   /* Flip the kernel about the center (necessary for non-symmetric
diff --git a/lib/gnuastro-internal/options.h b/lib/gnuastro-internal/options.h
index d00f8ed..1476642 100644
--- a/lib/gnuastro-internal/options.h
+++ b/lib/gnuastro-internal/options.h
@@ -304,10 +304,10 @@ error_t
 gal_options_common_argp_parse(int key, char *arg, struct argp_state *state);
 
 char *
-gal_options_stdin_error(long stdintimeout, int precedence);
+gal_options_stdin_error(long stdintimeout, int precedence, char *name);
 
 gal_list_str_t *
-gal_options_check_stdin(char *inputname, long stdintimeout);
+gal_options_check_stdin(char *inputname, long stdintimeout, char *name);
 
 
 
diff --git a/lib/options.c b/lib/options.c
index 5d66667..b0b580d 100644
--- a/lib/options.c
+++ b/lib/options.c
@@ -1364,20 +1364,19 @@ gal_options_common_argp_parse(int key, char *arg, 
struct argp_state *state)
 
 
 char *
-gal_options_stdin_error(long stdintimeout, int precedence)
+gal_options_stdin_error(long stdintimeout, int precedence, char *name)
 {
   char *out;
 
-  if( asprintf(&out, "no input!\n\n"
-               "The (first) input dataset can be read from a file "
-               "(specified as an argument), or the standard input.%s "
-               "Standard input can come from a pipe (output of another "
-               "program) or typed on the command-line before %ld "
-               "micro-seconds (configurable with the `--stdintimeout' "
-               "option).", ( precedence
-                             ? "If both are provided, a file takes precedence"
-                             : " " ),
-               stdintimeout )<0 )
+  if( asprintf(&out, "no %s!\n\n"
+               "The %s can be read from a file (specified as an argument), "
+               "or the standard input.%s Standard input can come from a "
+               "pipe (output of another program) or typed on the "
+               "command-line before %ld micro-seconds (configurable with "
+               "the `--stdintimeout' option).", name, name,
+               ( precedence
+                 ? " If both are provided, a file takes precedence."
+                 : "" ), stdintimeout )<0 )
     error(EXIT_FAILURE, 0, "%s: `asprintf' allocation error", __func__);
 
   return out;
@@ -1390,13 +1389,14 @@ gal_options_stdin_error(long stdintimeout, int 
precedence)
 /* Make the notice that is printed before program terminates, when no input
    is given and Standard input is also available. */
 gal_list_str_t *
-gal_options_check_stdin(char *inputname, long stdintimeout)
+gal_options_check_stdin(char *inputname, long stdintimeout, char *name)
 {
   gal_list_str_t *lines=inputname ? NULL : gal_txt_stdin_read(stdintimeout);
 
   /* See if atleast one of the two inputs is given. */
   if(inputname==NULL && lines==NULL)
-    error( EXIT_FAILURE, 0, "%s", gal_options_stdin_error(stdintimeout,1));
+    error( EXIT_FAILURE, 0, "%s", gal_options_stdin_error(stdintimeout,
+                                                          1, name));
 
   /* Return the output. */
   return lines;



reply via email to

[Prev in Thread] Current Thread [Next in Thread]