gnuastro-commits
[Top][All Lists]
Advanced

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

[gnuastro-commits] master 37ba6dbd: Library (arithmetic.h): new index, c


From: Mohammad Akhlaghi
Subject: [gnuastro-commits] master 37ba6dbd: Library (arithmetic.h): new index, counter and swap operators
Date: Sat, 29 Oct 2022 19:57:11 -0400 (EDT)

branch: master
commit 37ba6dbd68de7bac2ca537a07f5f75bbfb4df77f
Author: Mohammad Akhlaghi <mohammad@akhlaghi.org>
Commit: Mohammad Akhlaghi <mohammad@akhlaghi.org>

    Library (arithmetic.h): new index, counter and swap operators
    
    Until now, there was no way to identify the elements within a dataset
    (pixels within an image or rows in a table column) based on their position
    using Arithmetic. But they were necessary in some scenarios within tables
    and images.
    
    With this commit, a first implementation of such operators has been made
    through the 'index' and 'counter' operators. Since the don't use the
    input's pixel values, they don't pop their input from the stack (because
    the input's pixels are usually necessary afterwards!). As a result, two
    other operators have been added (called 'indexonly' and 'counteronly').
    
    Three other related changes have also been made:
    
     - Because it is sometimes necessary to have the positions after the
       input's pixel values, a new 'swap' operator has also been added. It will
       swap the two top operands on the stack.
    
     - A new '--writeall' option has been added to Arithmetic. When it is
       called, Arithmetic will save all the operands on the stack into the
       output file and won't complain if there is more than one. This is
       usually necessary when debugging an Arithmetic command.
    
     - While working on this (in particular, the 'index' operator in table
       columns), I noticed that in columns with an integer type, rows with a
       value of zero were not printed! This was because of the newly added
       default precision operators (that had a precision of 0 for
       integers!). Therefore with this commit, the default precision of 0 is
       removed when printing integers to plain-text.
    
    This fixes bug #63285.
---
 NEWS                        |  24 ++
 bin/arithmetic/args.h       |  15 ++
 bin/arithmetic/arithmetic.c | 144 +++++++++---
 bin/arithmetic/main.h       |   1 +
 bin/arithmetic/operands.c   |   1 +
 bin/arithmetic/ui.h         |   3 +-
 bin/table/arithmetic.c      |  97 ++++----
 bin/table/table.c           |  13 +-
 doc/gnuastro.texi           | 546 ++++++++++++++++++++++++++++----------------
 lib/arithmetic.c            | 121 ++++++++--
 lib/data.c                  |  31 +++
 lib/gnuastro/arithmetic.h   |   9 +-
 lib/gnuastro/data.h         |   3 +
 lib/tableintern.c           |  23 +-
 lib/txt.c                   |   2 -
 15 files changed, 728 insertions(+), 305 deletions(-)

diff --git a/NEWS b/NEWS
index 86b6da49..c2963b45 100644
--- a/NEWS
+++ b/NEWS
@@ -8,10 +8,32 @@ See the end of the file for license conditions.
 
 ** New features
 
+   Arithmetic
+   --writeall: Write all datasets on the stack as separate HDUs in the
+     output; this is useful in debugging incomplete Arithmetic commands.
+   - New operators (also available in Table).
+     - swap: swap the top two datasets on the stack of operands.
+     - index: return dataset of same size, with pixel values that are
+       integers starting from 0 and incrementing one by one. This operator
+       doesn't pop the top dataset from the stack, it only adds one (to
+       avoid extra RAM usage in most usage scenarios).
+     - counter: similar to 'index', but integers start with 1.
+     - indexonly: similar to 'index', but pops the top stack dataset.
+     - counteronly: similar to 'counter', but pops the top stack dataset.
+
    Crop:
    --metaname: Specify the name of the cropped output HDU (value to the
      'EXTNAME' keyword in FITS).
 
+   Library:
+   - GAL_ARITHMETIC_OP_SWAP: swap the top two operands.
+   - GAL_ARITHMETIC_OP_INDEX: An index (counting from 0) for every element.
+   - GAL_ARITHMETIC_OP_INDEXONLY: Similar to 'GAL_ARITHMETIC_OP_INDEX'.
+   - GAL_ARITHMETIC_OP_COUNTER: A counter (counting from 1) for every element.
+   - GAL_ARITHMETIC_OP_COUNTERONLY: Similar to 'GAL_ARITHMETIC_OP_COUNTER'.
+   - gal_data_alloc_empty: Allocate an empty dataset with a given number of
+     dimensions.
+
 ** Removed features
 
 ** Changed features
@@ -28,6 +50,8 @@ See the end of the file for license conditions.
   bug #63270: Table aborts when called with '--descending' with complaint
               on '--txtf64format' (that was added in Gnuastro
               0.19). Reported by Giulia Golini.
+  bug #63285: Table not printing anything in rows of integer columns with
+              value of zero.
 
 
 
diff --git a/bin/arithmetic/args.h b/bin/arithmetic/args.h
index 5b3d4b26..dfaf3967 100644
--- a/bin/arithmetic/args.h
+++ b/bin/arithmetic/args.h
@@ -30,6 +30,7 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 /* Definition of program-specific options. */
 struct argp_option program_options[] =
   {
+    /* Input options. */
     {
       "globalhdu",
       UI_KEY_GLOBALHDU,
@@ -87,6 +88,7 @@ struct argp_option program_options[] =
 
 
 
+    /* Output options. */
     {
       "onedasimage",
       UI_KEY_ONEDASIMAGE,
@@ -152,6 +154,19 @@ struct argp_option program_options[] =
       GAL_OPTIONS_NOT_MANDATORY,
       GAL_OPTIONS_NOT_SET
     },
+    {
+      "writeall",
+      UI_KEY_WRITEALL,
+      0,
+      0,
+      "Write all images in stack (not just a single one).",
+      GAL_OPTIONS_GROUP_OUTPUT,
+      &p->writeall,
+      GAL_OPTIONS_NO_ARG_TYPE,
+      GAL_OPTIONS_RANGE_0_OR_1,
+      GAL_OPTIONS_NOT_MANDATORY,
+      GAL_OPTIONS_NOT_SET
+    },
 
     {0}
   };
diff --git a/bin/arithmetic/arithmetic.c b/bin/arithmetic/arithmetic.c
index 0ddafdf5..f8074caa 100644
--- a/bin/arithmetic/arithmetic.c
+++ b/bin/arithmetic/arithmetic.c
@@ -1374,6 +1374,33 @@ arithmetic_set_operator(char *string, size_t 
*num_operands, int *inlib)
 
 
 
+static gal_data_t *
+arithmetic_prepare_meta(gal_data_t *d1, gal_data_t *d2, gal_data_t *d3)
+{
+  size_t i;
+  gal_data_t *out;
+
+  /* A small sanity check (these operators only need a single
+     operator). */
+  if(d2 || d3)
+    error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at '%s' to "
+          "find and fix the problem. Meta operators should only "
+          "take a single operand, but for some reason, this "
+          "hasn't occurred here", __func__, PACKAGE_BUGREPORT);
+
+  /* Allocate the output dataset. */
+  out=gal_data_alloc_empty(d1->ndim, d1->minmapsize, d1->quietmmap);
+
+  /* Fill the output dataset and return. */
+  for(i=0;i<d1->ndim;++i) out->dsize[i]=d1->dsize[i];
+  out->size=d1->size;
+  return out;
+}
+
+
+
+
+
 static void
 arithmetic_operator_run(struct arithmeticparams *p, int operator,
                         char *operator_string, size_t num_operands,
@@ -1435,6 +1462,19 @@ arithmetic_operator_run(struct arithmeticparams *p, int 
operator,
                 PACKAGE_BUGREPORT, num_operands, operator_string);
         }
 
+      /* Meta operators (where only the information about the pixels is
+         relevant, not the pixel values themselves, like the 'index'
+         operator). For these (that only accept one operand), we should add
+         the first popped operand back on the stack, then  */
+      switch(operator)
+        {
+        case GAL_ARITHMETIC_OP_INDEX:
+        case GAL_ARITHMETIC_OP_COUNTER:
+          operands_add(p, NULL, d1);
+          d1=arithmetic_prepare_meta(d1, d2, d3);
+          break;
+        }
+
       /* Run the arithmetic operation. Note that 'gal_arithmetic'
          is a variable argument function (like printf). So the
          number of arguments it uses depend on the operator. So
@@ -1542,6 +1582,66 @@ arithmetic_set_name_used_later(void *in, char *name)
 
 
 
+static void
+arithmetic_final_read_file(struct arithmeticparams *p,
+                           struct operand *operand)
+{
+  gal_data_t *out=NULL;
+  char *hdu, *filename;
+
+  /* Read the desired image and report it if necessary. */
+  hdu=operand->hdu;
+  filename=operand->filename;
+  if( gal_fits_file_recognized(filename) )
+    {
+      /* Read the data, note that the WCS has already been set. */
+      out=gal_array_read_one_ch(filename, hdu, NULL,
+                                p->cp.minmapsize,
+                                p->cp.quietmmap);
+      out->ndim=gal_dimension_remove_extra(out->ndim, out->dsize,
+                                            NULL);
+      if(!p->cp.quiet) printf(" - %s (hdu %s) is read.\n", filename, hdu);
+    }
+  else
+    error(EXIT_FAILURE, 0, "%s: a bug! please contact us at %s to fix "
+          "the problem. While 'operands->data' is NULL, the filename "
+          "('%s') is not recognized as a FITS file", __func__,
+          PACKAGE_BUGREPORT, filename);
+
+  /* Keep the read dataset (note that we don't need to free 'filename' or
+     'hdu' since their pointers are simply copied into the operand, they
+     are not copied.*/
+  operand->data=out;
+}
+
+
+
+
+
+/* Extract all the datasets of the remaining operands. */
+static gal_data_t *
+arithmetic_final_data(struct arithmeticparams *p)
+{
+  struct operand *op;
+  gal_data_t *out=NULL;
+
+  /* Go over the operands, extract their dataset and add it to the list. */
+  for(op=p->operands; op!=NULL; op=op->next)
+    {
+      gal_list_data_add(&out, op->data);
+      out->wcs=gal_wcs_copy(p->refdata.wcs);
+    }
+
+  /* Reverse the list so it goes in the same order as the operands, then
+     return it. */
+  gal_list_data_reverse(&out);
+  return out;
+}
+
+
+
+
+
 /* This function implements the reverse polish algorithm as explained
    in the Wikipedia page.
 
@@ -1550,10 +1650,11 @@ arithmetic_set_name_used_later(void *in, char *name)
 void
 reversepolish(struct arithmeticparams *p)
 {
+  char *printnum;
+  struct operand *otmp;
   size_t num_operands=0;
   gal_list_str_t *token;
-  gal_data_t *data, *col;
-  char *hdu, *filename, *printnum;
+  gal_data_t *tmp, *data, *col;
   struct gal_options_common_params *cp=&p->cp;
   int inlib, operator=GAL_ARITHMETIC_OP_INVALID;
 
@@ -1623,7 +1724,7 @@ reversepolish(struct arithmeticparams *p)
 
   /* If there is more than one node in the operands stack then the user has
      given too many operands which is an error. */
-  if(p->operands->next!=NULL)
+  if(p->writeall==0 && p->operands->next!=NULL)
     error(EXIT_FAILURE, 0, "too many operands");
 
 
@@ -1632,37 +1733,18 @@ reversepolish(struct arithmeticparams *p)
      read the contents of the file and put the resulting dataset into the
      operands 'data' element. This can happen for example if no operators
      are called and there is only one filename as an argument (which can
-     happen in scripts). */
-  if(p->operands->data==NULL && p->operands->filename)
-    {
-      /* Read the desired image and report it if necessary. */
-      hdu=p->operands->hdu;
-      filename=p->operands->filename;
-      if( gal_fits_file_recognized(filename) )
-        {
-          /* Read the data, note that the WCS has already been set. */
-          p->operands->data=gal_array_read_one_ch(filename, hdu, NULL,
-                                                  p->cp.minmapsize,
-                                                  p->cp.quietmmap);
-          data=p->operands->data;
-          data->ndim=gal_dimension_remove_extra(data->ndim, data->dsize,
-                                                NULL);
-          if(!p->cp.quiet)
-            printf(" - %s (hdu %s) is read.\n", filename, hdu);
-        }
-      else
-        error(EXIT_FAILURE, 0, "%s: a bug! please contact us at %s to fix "
-              "the problem. While 'operands->data' is NULL, the filename "
-              "('%s') is not recognized as a FITS file", __func__,
-              PACKAGE_BUGREPORT, filename);
-    }
+     happen in scripts).*/
+  for(otmp=p->operands; otmp!=NULL; otmp=otmp->next)
+    if(otmp->data==NULL && otmp->filename)
+      arithmetic_final_read_file(p, otmp);
 
 
   /* If the final data structure has more than one element, write it as a
      FITS file. Otherwise, if the user didn't call '--output', print it in
      the standard output. */
-  data=p->operands->data;
-  if(data->size==1 && data->ndim==1 && p->outnamerequested==0)
+  data=arithmetic_final_data(p);
+  if(data->next==NULL && data->size==1 && data->ndim==1
+     && p->outnamerequested==0)
     {
       /* Make the string to print the number. */
       printnum=gal_type_to_string(data->array, data->type, 0);
@@ -1687,13 +1769,13 @@ reversepolish(struct arithmeticparams *p)
 
       /* Put a copy of the WCS structure from the reference image, it
          will be freed while freeing 'data'. */
-      data->wcs=p->refdata.wcs;
       if(data->ndim==1 && p->onedasimage==0)
         gal_table_write(data, NULL, NULL, p->cp.tableformat,
                         p->onedonstdout ? NULL : p->cp.output,
                         "ARITHMETIC", 0);
       else
-        gal_fits_img_write(data, p->cp.output, NULL, PROGRAM_NAME);
+        for(tmp=data; tmp!=NULL; tmp=tmp->next)
+          gal_fits_img_write(tmp, p->cp.output, NULL, PROGRAM_NAME);
 
       /* Let the user know that the job is done. */
       if(!p->cp.quiet)
diff --git a/bin/arithmetic/main.h b/bin/arithmetic/main.h
index c0eeddae..813b719b 100644
--- a/bin/arithmetic/main.h
+++ b/bin/arithmetic/main.h
@@ -85,6 +85,7 @@ struct arithmeticparams
   char           *metaname;  /* FITS name (EXTNAME keyword) of output.  */
   char           *metaunit;  /* FITS name (BUNIT keyword) of output.    */
   char        *metacomment;  /* FITS comment of output.                 */
+  uint8_t         writeall;  /* Write all outputs.                      */
 
   /* Operating mode: */
   int        wcs_collapsed;  /* If the internal WCS is already collapsed.*/
diff --git a/bin/arithmetic/operands.c b/bin/arithmetic/operands.c
index be563b48..ebf00b90 100644
--- a/bin/arithmetic/operands.c
+++ b/bin/arithmetic/operands.c
@@ -109,6 +109,7 @@ operands_add_list(struct arithmeticparams *p, gal_data_t 
*data)
       newnode->data=tmp;
       newnode->hdu=NULL;
       newnode->filename=NULL;
+      newnode->data->next=NULL;
 
       /* Add this dataset to the top of the stack. */
       newnode->next=p->operands;
diff --git a/bin/arithmetic/ui.h b/bin/arithmetic/ui.h
index 00e5d187..8ed295a2 100644
--- a/bin/arithmetic/ui.h
+++ b/bin/arithmetic/ui.h
@@ -32,7 +32,7 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 
 /* Available letters for short options:
 
-   a b d e f i j k l m p r t v x y z
+   b d e f i j k l m p r t v x y z
    A B C E G H J L Q R X Y
 */
 enum option_keys_enum
@@ -46,6 +46,7 @@ enum option_keys_enum
   UI_KEY_METANAME        = 'n',
   UI_KEY_METAUNIT        = 'u',
   UI_KEY_METACOMMENT     = 'c',
+  UI_KEY_WRITEALL        = 'a',
 
   /* Only with long version (start with a value 1000, the rest will be set
      automatically). */
diff --git a/bin/table/arithmetic.c b/bin/table/arithmetic.c
index 0bbec8a7..afe60ce0 100644
--- a/bin/table/arithmetic.c
+++ b/bin/table/arithmetic.c
@@ -30,6 +30,7 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 
 #include <gnuastro/wcs.h>
 #include <gnuastro/type.h>
+#include <gnuastro/pointer.h>
 
 #include <gnuastro-internal/checkset.h>
 #include <gnuastro-internal/arithmetic-set.h>
@@ -434,20 +435,64 @@ arithmetic_indexs_final(struct tableparams *p)
 /*********************************************************************/
 /********************       Low-level tools      *********************/
 /*********************************************************************/
+static void
+arithmetic_placeholder_name(gal_data_t *col)
+{
+  static size_t counter=0;
+
+  /* Increment counter next time this function is called. */
+  ++counter;
+
+  /* Free any possibly existing metadata. */
+  if(col->name)    free(col->name);
+  if(col->unit)    free(col->unit);
+  if(col->comment) free(col->comment);
+
+  /* Set the new meta-data. */
+  errno=0;
+  if( asprintf(&col->name, "ARITH_%zu", counter)==-1 )
+    error(EXIT_FAILURE, errno, "%s: asprintf error for name", __func__);
+  if( asprintf(&col->unit, "arith_unit_%zu", counter)==-1)
+    error(EXIT_FAILURE, errno, "%s: asprintf error for unit", __func__);
+  if( asprintf(&col->comment, "Column from arithmetic operation %zu",
+               counter)==-1 )
+    error(EXIT_FAILURE, errno, "%s: asprintf error for comment", __func__);
+}
+
+
+
+
+
 static gal_data_t *
 arithmetic_stack_pop(gal_data_t **stack, int operator, char *errormsg)
 {
-  gal_data_t *out=*stack;
+  size_t i;
+  gal_data_t *out, *top=*stack;
 
-  /* Update the stack. */
-  if(*stack)
-    *stack=(*stack)->next;
-  else
-    error(EXIT_FAILURE, 0, "not enough operands for '%s'%s",
-          arithmetic_operator_name(operator), errormsg?errormsg:"");
+  /* Update the stack (where necessary). */
+  switch(operator)
+    {
+    /* Operators don't need to pop the dataset from the stack (they just
+       need the metadata). */
+    case GAL_ARITHMETIC_OP_INDEX:
+    case GAL_ARITHMETIC_OP_COUNTER:
+      out=gal_data_alloc_empty(top->ndim, top->minmapsize, top->quietmmap);
+      for(i=0;i<top->ndim;++i) out->dsize[i]=top->dsize[i];
+      out->size=top->size;
+      break;
+
+    /* Operators that actually need the data. */
+    default:
+      out=*stack;
+      if(*stack)
+        *stack=(*stack)->next;
+      else
+        error(EXIT_FAILURE, 0, "not enough operands for '%s'%s",
+              arithmetic_operator_name(operator), errormsg?errormsg:"");
+    }
 
   /* Arithmetic changes the contents of a dataset, so the old name (in the
-     FITS 'EXTNAME' keyword) must not be used beyond this
+     FITS 'EXTNAME' keyword) or units must not be used beyond this
      point. Furthermore, in Arithmetic, the 'name' element is used to
      identify variables (with the 'set-' operator). */
   if(out->name)    { free(out->name);    out->name=NULL;    }
@@ -807,34 +852,6 @@ arithmetic_datetosec(struct tableparams *p, gal_data_t 
**stack,
 /*********************************************************************/
 /********************          Operations        *********************/
 /*********************************************************************/
-static void
-arithmetic_placeholder_name(gal_data_t *col)
-{
-  static size_t counter=0;
-
-  /* Increment counter next time this function is called. */
-  ++counter;
-
-  /* Free any possibly existing metadata. */
-  if(col->name)    free(col->name);
-  if(col->unit)    free(col->unit);
-  if(col->comment) free(col->comment);
-
-  /* Set the new meta-data. */
-  errno=0;
-  if( asprintf(&col->name, "ARITH_%zu", counter)==-1 )
-    error(EXIT_FAILURE, errno, "%s: asprintf error for name", __func__);
-  if( asprintf(&col->unit, "arith_unit_%zu", counter)==-1)
-    error(EXIT_FAILURE, errno, "%s: asprintf error for unit", __func__);
-  if( asprintf(&col->comment, "Column from arithmetic operation %zu",
-               counter)==-1 )
-    error(EXIT_FAILURE, errno, "%s: asprintf error for comment", __func__);
-}
-
-
-
-
-
 static void
 arithmetic_operator_run(struct tableparams *p,
                         struct arithmetic_token *token,
@@ -1066,10 +1083,10 @@ arithmetic_reverse_polish(struct tableparams *p,
       ++setprm.tokencounter;
     }
 
-  /* Put everything that remains in the stack (reversed) into the final
-     table. Just note that 'gal_list_data_add' behaves differently for
-     lists, so we'll add have to manually set the 'next' element to NULL
-     before adding the column to the final table. */
+ /* Put everything that remains in the stack (reversed) into the final
+    table. Just note that 'gal_list_data_add' behaves differently for
+    lists, so we'll add have to manually set the 'next' element to NULL
+    before adding the column to the final table. */
   gal_list_data_reverse(&stack);
   while(stack!=NULL)
     {
diff --git a/bin/table/table.c b/bin/table/table.c
index 62824aae..76be7f69 100644
--- a/bin/table/table.c
+++ b/bin/table/table.c
@@ -1264,8 +1264,13 @@ table(struct tableparams *p)
   /* When any columns with blanks should be removed. */
   if(p->noblankend) table_noblankend(p);
 
-  /* Write the output. */
-  table_txt_formats(p);
-  gal_table_write(p->table, NULL, NULL, p->cp.tableformat, p->cp.output,
-                  "TABLE", p->colinfoinstdout);
+  /* Write the output or a warning/error (it can become NULL!) */
+  if(p->table)
+    {
+      table_txt_formats(p);
+      gal_table_write(p->table, NULL, NULL, p->cp.tableformat, p->cp.output,
+                      "TABLE", p->colinfoinstdout);
+    }
+  else
+    error(EXIT_FAILURE, 0, "no output columns");
 }
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index b0f37949..c82eccc3 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -522,7 +522,8 @@ Arithmetic operators
 * Random number generators::    Random numbers can be used to add noise for 
example.
 * Elliptical shape operators::  Operations that are focused on an ellipse.
 * Loading external columns::    Read a column from a table into the stack.
-* Building new dataset::        How to construct an empty dataset from scratch.
+* Size and position operators::  Extracting image size and pixel positions.
+* Building new dataset and stack management::        How to construct an empty 
dataset from scratch.
 * Operand storage in memory or a file::  Tools for complex operations in one 
command.
 
 Convolve
@@ -4744,8 +4745,8 @@ These are very important points that you will definitely 
encouter during your ow
 #!/bin/bash
 # Script to download and keep the deep region of the XDF survey.
 #
-# Copyright (C) 2021-2022 Original Author <original@@email.server>
-# Copyright (C) 2022      Your Name <yourname@@email.server>
+# Copyright (C) 2022      Your Name <yourname@@email.company>
+# Copyright (C) 2021-2022 Initial Author <incase@@there-is.any>
 #
 # This script is free software: you can redistribute it and/or modify
 # it under the terms of the GNU General Public License as published by
@@ -16334,7 +16335,8 @@ Reading NaN as a floating point number in Gnuastro is 
not case-sensitive.
 * Random number generators::    Random numbers can be used to add noise for 
example.
 * Elliptical shape operators::  Operations that are focused on an ellipse.
 * Loading external columns::    Read a column from a table into the stack.
-* Building new dataset::        How to construct an empty dataset from scratch.
+* Size and position operators::  Extracting image size and pixel positions.
+* Building new dataset and stack management::        How to construct an empty 
dataset from scratch.
 * Operand storage in memory or a file::  Tools for complex operations in one 
command.
 @end menu
 
@@ -16825,15 +16827,6 @@ for example, with the command below, the non-blank 
pixel values of @file{cropped
 $ astarithmetic cropped.fits nonblank --onedonstdout --quiet
 @end example
 
-@item size
-Size of the dataset along a given FITS (or FORTRAN) dimension (counting from 
1).
-The desired dimension should be the first popped operand and the dataset must 
be the second popped operand.
-The output will be a single unsigned integer (dimensions cannot be negative).
-For example, the following command will produce the size of the first 
extension/HDU (the default HDU) of @file{a.fits} along the second FITS axis.
-
-@example
-$ astarithmetic a.fits 2 size
-@end example
 @end table
 
 @node Stacking operators, Filtering operators, Statistical operators, 
Arithmetic operators
@@ -17988,7 +17981,7 @@ $ asttable catalog.fits \
 @end example
 @end table
 
-@node Loading external columns, Building new dataset, Elliptical shape 
operators, Arithmetic operators
+@node Loading external columns, Size and position operators, Elliptical shape 
operators, Arithmetic operators
 @subsubsection Loading external columns
 
 In the Arithmetic program, you can always load new dataset by simply giving 
their name.
@@ -18026,8 +18019,194 @@ $ asttable tab-1.txt -c'arith $1 
load-col-1-from-tab-2.txt +'
 @end example
 @end table
 
-@node Building new dataset, Operand storage in memory or a file, Loading 
external columns, Arithmetic operators
-@subsubsection Building new dataset
+@node Size and position operators, Building new dataset and stack management, 
Loading external columns, Arithmetic operators
+@subsubsection Size and position operators
+
+With the operators below you can get metadata about the top dataset on the 
stack.
+
+@table @code
+@item index
+Add a new operand to the stack with an integer type and the same size (in all 
dimensions) as top operand on the stack (before it was called; it is not 
popped!).
+The first pixel in the returned operand is zero, and every later pixel's value 
is incremented by one.
+It is important to remember that the top operand is not popped by this 
operand, so it remains on the stack.
+After this operand is finished, it adds a new operand to the stack.
+To pop the previous operand, you can use the @code{indexonly} operator.
+
+The data type of the output is always an unsigned integer, and its width is 
determined from the number of pixels/rows in the top operand.
+For example if there are only 108 rows in a table, the returned column will 
have an unsigned 8-bit integer type (that can keep 256 separate values).
+But if the top operand is a @mymath{1000\times1000=10^6} pixel image, the 
output will be a 32-bit unsigned integer.
+For the various types of integers, see @ref{Numeric data types}.
+
+To see the index image along with the actual image, you can use the 
@option{--writeall} operator to have a multi-HDU output (without 
@option{--writeall}, Arithmetic will complain if more than one operand is left 
at the end).
+After DS9 opens with the second command, flip between the two extensions.
+
+@example
+$ astarithmetic image.fits index --writeall
+$ astscript-fits-view image_arith.fits
+@end example
+
+Below is a review some usage examples of this operator:
+
+@table @asis
+@item Image: masking margins
+With the command below, we will be masking all pixels that are 20 pixels away 
from the edges of the image (on the margin).
+Here is a description of the command below (for the basics of Arithmetic's 
notation, see @ref{Reverse polish notation}):
+@itemize
+@item
+The @code{index} operator just adds a new dataset on the stack: unlike almost 
all other operators in Arithmetic, @code{index} doesn't remove the dataset 
within @file{image.fits} from the stack (use @code{indexonly} for the 
``normal'' behavior).
+This is because @code{index} is returns the pixel metadata not data.
+As a result, after @code{index}, we have two operands on the stack: the input 
image and the index image.
+@item
+With the @code{set-i} operator, the top operand (the image containing the 
index of each pixel) is popped from the stack and associated to the name 
@code{i}.
+Therefore after this, the stack only has the input image.
+For more on the @code{set-} operator, see @ref{Operand storage in memory or a 
file}.
+@item
+We need three values from the commands before Arithmetic (for the width and 
height of the image and the size of the margin).
+To make the rest of the command easier to read/use, we'll define them in 
Arithmetic as three named operators (respectively called @code{w}, @code{h} and 
@code{m}).
+All three are integers will have a positive value lower than 
@mymath{2^{16}=65536} (for a ``normal'' image!).
+Therefore, we will store tham as 16-bit unsigned integers with the 
@code{uint16} operator (this will help optimal processing in later steps).
+For more the type changing operators, see @ref{Numerical type conversion 
operators}.
+@item
+Using the modulo @code{%} and division (@code{/}) operators on the index image 
and the width, we extract the horizontal (X) and vertical (Y) positions of each 
pixel in separately named operands called @code{X} and @code{Y}.
+The maximum value in these two will also fit within an unsigned 16-bit 
integer, so we'll also store these in that type.
+@item
+For the horizonal (X) dimension, we select pixels that are less than the 
margin (@code{X m lt}) and those that are more than the width subtracted by the 
margin (@code{X w m - gt}).
+@item
+The output of the @code{lt} and @code{gt} conditional operators above is a 
binary (0 or 1 valued) image.
+We therefore merge them into one binary image using the @code{or} operator.
+For more, see @ref{Conditional operators}.
+@item
+We repate the two steps above for the vertical (Y) dimension.
+@item
+Once the images containing the to-be-masked pixels in each dimension are made, 
we combine them into one binary image with a final @code{or} operator.
+At this point, the stack only has two operands: 1) the input image and 2) the 
binary image that has a value of 1 for all pixels whose value should be chaned.
+@item
+A single-element operand (@code{nan}) is added on the stack.
+@item
+Using the @code{where} operator, we replace all the pixels that are non-zero 
in the second operand (on the margins) to the top operand's value (NaN) in the 
third popped operand (image that was read from @code{image.fits}).
+For more on the @code{where} operator, see @ref{Conditional operators}.
+@end itemize
+
+@example
+$ margin=20
+$ width=$(astfits  image.fits --keyvalue=NAXIS1 -q)
+$ height=$(astfits image.fits --keyvalue=NAXIS2 -q)
+$ astarithmetic image.fits index       set-i \
+                $width     uint16      set-w \
+                $height    uint16      set-h \
+                $margin    uint16      set-m \
+                i w %      uint16      set-X \
+                i w /      uint16      set-Y \
+                X m lt     X w m - gt     or \
+                Y m lt     Y h m - gt     or \
+                or nan where
+@end example
+
+@item Image: Masking regions outside a circle
+As another example for usage on an image, in the command below we are using 
@code{index} to define an image where each pixel contains the distance to the 
pixel with X,Y coordinates of 345,250.
+We are then using that distance image to only keep the pixels that are within 
a 50 pixel radius of that point.
+
+The basic concept behind this process is very similar to the previous example, 
with a different mathematical definition for pixels to mask.
+The major difference is that we want the distance to a pixel within the image, 
we need to have negative values and the center coordinates can be in a 
sub-pixel positions.
+The best numeric datatype for intermediate steps is therefore floating point.
+64-bit floating point can have a precision of upto 15 digits after the decimal 
point.
+This is far too much for what we need here: in astronomical imaging, the PSF 
is usually on the scale of 1 or more pixels (see @ref{Sampling theorem}).
+So even reaching a precision of one millionth of a pixel (offered by 32-bit 
floating points) is beyond our wildest dreams (see @ref{Numeric data types}).
+We will also define the horizontal (X) and vertical (Y) operands after 
shifting to the desired central point.
+
+@example
+$ radius=50
+$ centerx=345.2
+$ centery=250.3
+$ width=$(astfits image.fits --keyvalue=NAXIS1 -q)
+$ astarithmetic image.fits index set-i \
+                $width       uint16    set-w \
+                $radius      float32   set-r \
+                $centerx     float32   set-cx \
+                $centery     float32   set-cy \
+                i w % cx -             set-X \
+                i w / cy -             set-Y \
+                X X x Y Y x + sqrt r gt \
+                nan where --output=arith-masked.fits
+@end example
+
+@cartouche
+@noindent
+@strong{Optimal data types have significant benefits:} choosing the minimum 
required datatype for your operation is very important to avoid wasting your 
CPU and RAM.
+Don't simply default to 64-bit floating points for everything!
+Integer operations are much faster than floating points, and within floating 
point types, 32-bit is faster and will use half the RAM/storage!
+For more, see @ref{Numeric data types}.
+@end cartouche
+
+The example above was just a demo for usage of the @code{index} operator and 
some important concepts.
+But it is not the easiest way to achieve the desired result above!
+An easier way for the scenario above (to keep a circle within an image and set 
everything else to NaN) is to use MakeProfiles in combination with Arithmetic, 
like below:
+
+@example
+$ radius=50
+$ centerx=345.2
+$ centery=250.3
+$ echo "1 $centerx $centery 5 $radius 0 0 1 1 1" \
+       | astmkprof --background=image.fits \
+                   --mforflatpix --clearcanvas \
+                   -omkprof-mask.fits --type=uint8
+$ astarithmetic image.fits mkprof-mask.fits not \
+                nan where -g1 -omkprof-masked.fits
+@end example
+
+@item Tables: adding new columns with row index
+Within Table, you can use this operator to add an index column like below (see 
the @code{counter} operator for starting the count from one).
+
+@example
+## The index will be the second column.
+$ asttable table.fits -c'arith $1 index'
+
+## The index will be the first column
+$ asttable table.fits -c'arith $1 index swap'
+@end example
+@end table
+
+@item indexonly
+Similar to @code{index}, except that the top operand is popped from the stack 
and is no longer available afterwards.
+
+@item counter
+Similar to @code{index}, except that counting starts from one (not zero as in 
@code{index}).
+Counting from one is usually necessary when adding row counters in tables, 
like below:
+
+@example
+$ asttable table.fits -c'arith $1 counter swap'
+@end example
+
+@item counteronly
+Similar to @code{counter}, but the top operand before it is popped (no longer 
available).
+
+@item size
+Size of the dataset along a given FITS (or FORTRAN) dimension (counting from 
1).
+The desired dimension should be the first popped operand and the dataset must 
be the second popped operand.
+The output will be a single unsigned integer (dimensions cannot be negative).
+For example, the following command will produce the size of the first 
extension/HDU (the default HDU) of @file{a.fits} along the second FITS axis.
+
+@example
+$ astarithmetic a.fits 2 size
+@end example
+
+@cartouche
+@noindent
+@strong{Not optimal:} This operator reads the top element on the stack and 
then simply reads its size along the given dimension.
+On a small dataset this won't consume much RAM, but if you want to put this in 
a pipeline or use it on large image, the extra RAM and slow operation can 
become meaningful.
+To avoid such issues, you can read the size along the given dimension using 
the @option{--keyvalue} option of @ref{Keyword inspection and manipulation}.
+For example, in the code below, the X axis position of every pixel is returned:
+
+@example
+$ width=$(astfits image.fits --keyvalue=NAXIS1 -q)
+$ astarithmetic image.fits indexonly $width % -opix-x.fits
+@end example
+@end cartouche
+@end table
+
+
+@node Building new dataset and stack management, Operand storage in memory or 
a file, Size and position operators, Arithmetic operators
+@subsubsection Building new dataset and stack management
 
 With the operator here, you can create a new dataset from scratch to start 
certain operations without any input data.
 
@@ -18048,9 +18227,25 @@ To further extend the example, you can use any of the 
noise-making operators to
 @example
 $ astarithmetic 100 200 2 makenew 5 mknoise-sigma
 @end example
+
+@item swap
+Swap the top two operands on the stack.
+For example the @code{index} operator doesn't pop with the top operand (the 
input to @code{index}), it just adds the index image to the stack.
+In case you want your next operation to be on the input to @code{index}, you 
can simply call @code{swap} and continue the operations on that image, while 
keeping the indexed pixels for later steps.
+In the example below we are using the @option{--writeall} option to write the 
full stack and if you open the outputs you will see that the stack order has 
changed.
+
+@example
+## Index image is written in HDU 1.
+$ astarithmetic image.fits index      --writeall \
+                --output=ind-first.fits
+
+## image.fits in HDU 1.
+$ astarithmetic image.fits index swap --writeall \
+                --output=img-first.fits
+@end example
 @end table
 
-@node Operand storage in memory or a file,  , Building new dataset, Arithmetic 
operators
+@node Operand storage in memory or a file,  , Building new dataset and stack 
management, Arithmetic operators
 @subsubsection Operand storage in memory or a file
 
 In your early days of using Gnuastro, to do multiple operations, it is likely 
that you will simply call Arithmetic (or Table, with column arithmetic) 
multiple times: feed the output file of the first call to the second call.
@@ -18276,6 +18471,13 @@ Note it does not matter if the final output file name 
is given with the @option{
 Arithmetic treats this option differently from its default operation in other 
Gnuastro programs (see @ref{Input output options}).
 If the output file exists, when other Gnuastro programs are called with 
@option{--dontdelete}, they simply complain and abort.
 But when Arithmetic is called with @option{--dontdelete}, it will appended the 
dataset(s) to the existing extension(s) in the file.
+
+@item -a
+@itemx --writeall
+Write all datasets on the stack as separate HDUs in the output file.
+This only affects datasets with multiple dimensions (or single-dimension 
datasets when the @option{--onedasimg} is called).
+This option is useful to debug Arithmetic calls: to check all the images on 
the stack while you are designing your operation.
+The top dataset on the stack will be on HDU number 1 of the output, the second 
dataset will be on HDU number 2 and so on.
 @end table
 
 Arithmetic accepts two kinds of input: images and numbers.
@@ -31015,20 +31217,15 @@ and some comments. To deal with any generic dataset, 
Gnuastro defines the
 @node Generic data container, Dataset allocation, Library data container, 
Library data container
 @subsubsection Generic data container (@code{gal_data_t})
 
-To be able to deal with any dataset (various dimensions, numeric data
-types, units and higher-level structures), Gnuastro defines the
-@code{gal_data_t} type which is the input/output container of choice for
-many of Gnuastro library's functions. It is defined in
-@file{gnuastro/data.h}. If you will be using (`@code{# include}'ing) those
-libraries, you do not need to include this header explicitly, it is already
-included by any library header that uses @code{gal_data_t}.
+To be able to deal with any dataset (various dimensions, numeric data types, 
units and higher-level structures), Gnuastro defines the @code{gal_data_t} type 
which is the input/output container of choice for many of Gnuastro library's 
functions.
+It is defined in @file{gnuastro/data.h}.
+If you will be using (`@code{# include}'ing) those libraries, you do not need 
to include this header explicitly, it is already included by any library header 
that uses @code{gal_data_t}.
 
 @deftp {Type (C @code{struct})} gal_data_t
-The main container for datasets in Gnuastro. It can host data of any
-dimensions, with any numeric data type. It is actually a structure, but
-@code{typedef}'d as a new type to avoid having to write the @code{struct}
-before any declaration. The actual structure is shown below which is
-followed by a description of each element.
+The main container for datasets in Gnuastro.
+It can host data of any dimensions, with any numeric data type.
+It is actually a structure, but @code{typedef}'d as a new type to avoid having 
to write the @code{struct} before any declaration.
+The actual structure is shown below which is followed by a description of each 
element.
 
 @example
 typedef struct gal_data_t
@@ -31067,68 +31264,44 @@ The list below contains a description for each 
@code{gal_data_t} element.
 @cindex @code{void *}
 @table @code
 @item void *restrict array
-This is the pointer to the main array of the dataset containing the raw
-data (values). All the other elements in this data-structure are actually
-meta-data enabling us to use/understand the series of values in this
-array. It must allow data of any type (see @ref{Numeric data types}), so it
-is defined as a @code{void *} pointer. A @code{void *} array is not
-directly usable in C, so you have to cast it to proper type before using
-it, please see @ref{Library demo - reading a image} for a demonstration.
+This is the pointer to the main array of the dataset containing the raw data 
(values).
+All the other elements in this data-structure are actually meta-data enabling 
us to use/understand the series of values in this array.
+It must allow data of any type (see @ref{Numeric data types}), so it is 
defined as a @code{void *} pointer.
+A @code{void *} array is not directly usable in C, so you have to cast it to 
proper type before using it, please see @ref{Library demo - reading a image} 
for a demonstration.
 
 @cindex @code{restrict}
 @cindex C: @code{restrict}
-The @code{restrict} keyword was formally introduced in C99 and is used to
-tell the compiler that at any moment only this pointer will modify what it
-points to (a pixel in an image for example)@footnote{Also see
-@url{https://en.wikipedia.org/wiki/Restrict}.}. This extra piece of
-information can greatly help in compiler optimizations and thus the running
-time of the program. But older compilers might not have this capability, so
-at @command{./configure} time, Gnuastro checks this feature and if the
-user's compiler does not support @code{restrict}, it will be removed from
-this definition.
+The @code{restrict} keyword was formally introduced in C99 and is used to tell 
the compiler that at any moment only this pointer will modify what it points to 
(a pixel in an image for example)@footnote{Also see 
@url{https://en.wikipedia.org/wiki/Restrict}.}.
+This extra piece of information can greatly help in compiler optimizations and 
thus the running time of the program.
+But older compilers might not have this capability, so at 
@command{./configure} time, Gnuastro checks this feature and if the user's 
compiler does not support @code{restrict}, it will be removed from this 
definition.
 
 @cindex Data type
 @item uint8_t type
-A fixed code (integer) used to identify the type of data in @code{array}
-(see @ref{Numeric data types}). For the list of acceptable values to this
-variable, please see @ref{Library data types}.
+A fixed code (integer) used to identify the type of data in @code{array} (see 
@ref{Numeric data types}).
+For the list of acceptable values to this variable, please see @ref{Library 
data types}.
 
 @item size_t ndim
 The dataset's number of dimensions.
 
-
 @cindex FORTRAN
 @cindex FITS standard
 @cindex Standard, FITS
 @item size_t *dsize
-The size of the dataset along each dimension. This is an array (with
-@code{ndim} elements), of positive integers in row-major
-order@footnote{Also see
-@url{https://en.wikipedia.org/wiki/Row-_and_column-major_order}.}  (based
-on C). When a data file is read into memory with Gnuastro's libraries, this
-array is dynamically allocated based on the number of dimensions that the
-dataset has.
-
-It is important to remember that C's row-major ordering is the opposite of
-the FITS standard which is in column-major order: in the FITS standard the
-fastest dimension's size is specified by @code{NAXIS1}, and slower
-dimensions follow. The FITS standard was defined mainly based on the
-FORTRAN language which is the opposite of C's approach to multi-dimensional
-arrays (and also starts counting from 1 not 0). Hence if a FITS image has
-@code{NAXIS1==20} and @code{NAXIS2==50}, the @code{dsize} array must be
-filled with @code{dsize[0]==50} and @code{dsize[1]==20}.
-
-The fastest dimension is the one that is contiguous in memory: to increment
-by one along that dimension, just go to the next element in the array. As
-we go to slower dimensions, the number of memory cells we have to skip for
-an increment along that dimension becomes larger.
+The size of the dataset along each dimension.
+This is an array (with @code{ndim} elements), of positive integers in 
row-major order@footnote{Also see 
@url{https://en.wikipedia.org/wiki/Row-_and_column-major_order}.}  (based on C).
+When a data file is read into memory with Gnuastro's libraries, this array is 
dynamically allocated based on the number of dimensions that the dataset has.
+
+It is important to remember that C's row-major ordering is the opposite of the 
FITS standard which is in column-major order: in the FITS standard the fastest 
dimension's size is specified by @code{NAXIS1}, and slower dimensions follow.
+The FITS standard was defined mainly based on the FORTRAN language which is 
the opposite of C's approach to multi-dimensional arrays (and also starts 
counting from 1 not 0).
+Hence if a FITS image has @code{NAXIS1==20} and @code{NAXIS2==50}, the 
@code{dsize} array must be filled with @code{dsize[0]==50} and 
@code{dsize[1]==20}.
+
+The fastest dimension is the one that is contiguous in memory: to increment by 
one along that dimension, just go to the next element in the array.
+As we go to slower dimensions, the number of memory cells we have to skip for 
an increment along that dimension becomes larger.
 
 @item size_t size
-The total number of elements in the dataset. This is actually a
-multiplication of all the values in the @code{dsize} array, so it is not an
-independent parameter. However, low-level operations with the dataset
-(irrespective of its dimensions) commonly need this number, so this element
-is designed to avoid calculating it every time.
+The total number of elements in the dataset.
+This is actually a multiplication of all the values in the @code{dsize} array, 
so it is not an independent parameter.
+However, low-level operations with the dataset (irrespective of its 
dimensions) commonly need this number, so this element is designed to avoid 
calculating it every time.
 
 @item int quietmmap
 When this value is zero, and the dataset must not be allocated in RAM (see 
@code{mmapname} and @code{minmapsize}), a warning will be printed to inform the 
user when the file is created and when it is deleted.
@@ -31173,98 +31346,70 @@ When a dataset does not have any blank values, the 
@code{GAL_DATA_FLAG_HASBLANK}
 But upon initialization, all bits also get a value of zero.
 Therefore, a checker needs this flag to see if the value in 
@code{GAL_DATA_FLAG_HASBLANK} is reliable (dataset has actually been parsed for 
a blank value) or not.
 
-Also, if it is necessary to re-check the presence of flags, you just have
-to set this flag to zero and call @code{gal_blank_present} for example, to
-parse the dataset and check for blank values. Note that for improved
-efficiency, when this flag is set, @code{gal_blank_present} will not
-actually parse the dataset, it will just use @code{GAL_DATA_FLAG_HASBLANK}.
+Also, if it is necessary to re-check the presence of flags, you just have to 
set this flag to zero and call @code{gal_blank_present} for example, to parse 
the dataset and check for blank values.
+Note that for improved efficiency, when this flag is set, 
@code{gal_blank_present} will not actually parse the dataset, it will just use 
@code{GAL_DATA_FLAG_HASBLANK}.
 
 @item GAL_DATA_FLAG_HASBLANK
-This bit has a value of @code{1} when the given dataset has blank
-values. If this bit is @code{0} and @code{GAL_DATA_FLAG_BLANK_CH} is
-@code{1}, then the dataset has been checked and it did not have any blank
-values, so there is no more need for further checks.
+This bit has a value of @code{1} when the given dataset has blank values.
+If this bit is @code{0} and @code{GAL_DATA_FLAG_BLANK_CH} is @code{1}, then 
the dataset has been checked and it did not have any blank values, so there is 
no more need for further checks.
 
 @item GAL_DATA_FLAG_SORT_CH
-Marking that the dataset is already checked for being sorted or not and
-thus that the possible @code{0} values in @code{GAL_DATA_FLAG_SORTED_I} and
-@code{GAL_DATA_FLAG_SORTED_D} are meaningful. The logic behind this is
-similar to that in @code{GAL_DATA_FLAG_BLANK_CH}.
+Marking that the dataset is already checked for being sorted or not and thus 
that the possible @code{0} values in @code{GAL_DATA_FLAG_SORTED_I} and 
@code{GAL_DATA_FLAG_SORTED_D} are meaningful.
+The logic behind this is similar to that in @code{GAL_DATA_FLAG_BLANK_CH}.
 
 @item GAL_DATA_FLAG_SORTED_I
-This bit has a value of @code{1} when the given dataset is sorted in an
-increasing manner. If this bit is @code{0} and @code{GAL_DATA_FLAG_SORT_CH}
-is @code{1}, then the dataset has been checked and was not sorted
-(increasing), so there is no more need for further checks.
+This bit has a value of @code{1} when the given dataset is sorted in an 
increasing manner.
+If this bit is @code{0} and @code{GAL_DATA_FLAG_SORT_CH} is @code{1}, then the 
dataset has been checked and was not sorted (increasing), so there is no more 
need for further checks.
 
 @item GAL_DATA_FLAG_SORTED_D
-This bit has a value of @code{1} when the given dataset is sorted in a
-decreasing manner. If this bit is @code{0} and @code{GAL_DATA_FLAG_SORT_CH}
-is @code{1}, then the dataset has been checked and was not sorted
-(decreasing), so there is no more need for further checks.
-
+This bit has a value of @code{1} when the given dataset is sorted in a 
decreasing manner.
+If this bit is @code{0} and @code{GAL_DATA_FLAG_SORT_CH} is @code{1}, then the 
dataset has been checked and was not sorted (decreasing), so there is no more 
need for further checks.
 @end table
 
-The macro @code{GAL_DATA_FLAG_MAXFLAG} contains the largest internally used
-bit-position. Higher-level flags can be defined with the bit-wise shift
-operators using this macro to define internal flags for libraries/programs
-that depend on Gnuastro without causing any possible conflict with the
-internal flags discussed above or having to check the values manually on
-every release.
+The macro @code{GAL_DATA_FLAG_MAXFLAG} contains the largest internally used 
bit-position.
+Higher-level flags can be defined with the bit-wise shift operators using this 
macro to define internal flags for libraries/programs that depend on Gnuastro 
without causing any possible conflict with the internal flags discussed above 
or having to check the values manually on every release.
 
 
 @item int status
-A context-specific status values for this data-structure. This integer will
-not be set by Gnuastro's libraries. You can use it keep some additional
-information about the dataset (with integer constants) depending on your
-applications.
+A context-specific status values for this data-structure.
+This integer will not be set by Gnuastro's libraries.
+You can use it keep some additional information about the dataset (with 
integer constants) depending on your applications.
 
 @item char *name
-The name of the dataset. If the dataset is a multi-dimensional array and
-read/written as a FITS image, this will be the value in the @code{EXTNAME}
-FITS keyword. If the dataset is a one-dimensional table column, this will
-be the column name. If it is set to @code{NULL} (by default), it will be
-ignored.
+The name of the dataset.
+If the dataset is a multi-dimensional array and read/written as a FITS image, 
this will be the value in the @code{EXTNAME} FITS keyword.
+If the dataset is a one-dimensional table column, this will be the column name.
+If it is set to @code{NULL} (by default), it will be ignored.
 
 @item char *unit
-The units of the dataset (for example, @code{BUNIT} in the standard FITS
-keywords) that will be read from or written to files/tables along with the
-dataset. If it is set to @code{NULL} (by default), it will be ignored.
+The units of the dataset (for example, @code{BUNIT} in the standard FITS 
keywords) that will be read from or written to files/tables along with 
thedataset.
+If it is set to @code{NULL} (by default), it will be ignored.
 
 @item char *comment
-Any further explanation about the dataset which will be written to any
-output file if present.
+Any further explanation about the dataset which will be written to any output 
file if present.
 
 @item disp_fmt
-Format to use for printing each element of the dataset to a plain text
-file, the acceptable values to this element are defined in @ref{Table input
-output}. Based on C's @code{printf} standards.
+Format to use for printing each element of the dataset to a plain text file, 
the acceptable values to this element are defined in @ref{Table input output}.
+Based on C's @code{printf} standards.
 
 @item disp_width
-Width of printing each element of the dataset to a plain text file, the
-acceptable values to this element are defined in @ref{Table input
-output}. Based on C's @code{printf} standards.
+Width of printing each element of the dataset to a plain text file, the 
acceptable values to this element are defined in @ref{Table input output}.
+Based on C's @code{printf} standards.
 
 @item disp_precision
-Precision of printing each element of the dataset to a plain text file, the
-acceptable values to this element are defined in @ref{Table input
-output}. Based on C's @code{printf} standards.
+Precision of printing each element of the dataset to a plain text file, the 
acceptable values to this element are defined in @ref{Table input output}.
+Based on C's @code{printf} standards.
 
 @item gal_data_t *next
-Through this pointer, you can link a @code{gal_data_t} with other datasets
-related datasets, for example, the different columns in a dataset each have
-one @code{gal_data_t} associate with them and they are linked to each other
-using this element. There are several functions described below to
-facilitate using @code{gal_data_t} as a linked list. See @ref{Linked lists}
-for more on these wonderful high-level constructs.
+Through this pointer, you can link a @code{gal_data_t} with other datasets 
related datasets, for example, the different columns in a dataset each have one 
@code{gal_data_t} associate with them and they are linked to each other using 
this element.
+There are several functions described below to facilitate using 
@code{gal_data_t} as a linked list.
+See @ref{Linked lists} for more on these wonderful high-level constructs.
 
 @item gal_data_t *block
-Pointer to the start of the complete allocated block of memory. When this
-pointer is not @code{NULL}, the dataset is not treated as a contiguous
-patch of memory. Rather, it is seen as covering only a portion of the
-larger patch of memory that @code{block} points to. See @ref{Tessellation
-library} for a more thorough explanation and functions to help work with
-tiles that are created from this pointer.
+Pointer to the start of the complete allocated block of memory.
+When this pointer is not @code{NULL}, the dataset is not treated as a 
contiguous patch of memory.
+Rather, it is seen as covering only a portion of the larger patch of memory 
that @code{block} points to.
+See @ref{Tessellation library} for a more thorough explanation and functions 
to help work with tiles that are created from this pointer.
 @end table
 
 
@@ -31276,46 +31421,36 @@ These functions are declared in 
@file{gnuastro/data.h} which is also visible fro
 
 @deftypefun {gal_data_t *} gal_data_alloc (void @code{*array}, uint8_t 
@code{type}, size_t @code{ndim}, size_t @code{*dsize}, struct wcsprm 
@code{*wcs}, int @code{clear}, size_t @code{minmapsize}, int @code{quietmmap}, 
char @code{*name}, char @code{*unit}, char @code{*comment})
 
-Dynamically allocate a @code{gal_data_t} and initialize it will all the
-given values. See the description of @code{gal_data_initialize} and
-@ref{Generic data container} for more information. This function will often
-be the most frequently used because it allocates the @code{gal_data_t}
-hosting all the values @emph{and} initializes it. Once you are done with
-the dataset, be sure to clean up all the allocated spaces with
-@code{gal_data_free}.
+Dynamically allocate a @code{gal_data_t} and initialize it will all the given 
values.
+See the description of @code{gal_data_initialize} and @ref{Generic data 
container} for more information.
+This function will often be the most frequently used because it allocates the 
@code{gal_data_t} hosting all the values @emph{and} initializes it.
+Once you are done with the dataset, be sure to clean up all the allocated 
spaces with @code{gal_data_free}.
 @end deftypefun
 
 @deftypefun void gal_data_initialize (gal_data_t @code{*data}, void 
@code{*array}, uint8_t @code{type}, size_t @code{ndim}, size_t @code{*dsize}, 
struct wcsprm @code{*wcs}, int @code{clear}, size_t @code{minmapsize}, int 
@code{quietmmap}, char @code{*name}, char @code{*unit}, char @code{*comment})
 
-Initialize the given data structure (@code{data}) with all the given
-values. Note that the raw input @code{gal_data_t} must already have been
-allocated before calling this function. For a description of each variable
-see @ref{Generic data container}. It will set the values and do the
-necessary allocations. If they are not @code{NULL}, all input arrays
-(@code{dsize}, @code{wcs}, @code{name}, @code{unit}, @code{comment}) are
-separately copied (allocated) by this function for usage in @code{data}, so
-you can safely use one value to initialize many datasets or use statically
-allocated variables in this function call. Once you are done with the
-dataset, you can free all the allocated spaces with
-@code{gal_data_free_contents}.
-
-If @code{array} is not @code{NULL}, it will be directly copied into
-@code{data->array} (based on the total number of elements calculated from
-@code{dsize}) and no new space will be allocated for the array of this
-dataset, this has many low-level advantages and can be used to work on
-regions of a dataset instead of the whole allocated array (see the
-description under @code{block} in @ref{Generic data container} for one
-example). If the given pointer is not the start of an allocated block of
-memory or it is used in multiple datasets, be sure to set it to @code{NULL}
-(with @code{data->array=NULL}) before cleaning up with
-@code{gal_data_free_contents}.
-
-@code{ndim} may be zero. In this case no allocation will occur,
-@code{data->array} and @code{data->dsize} will be set to @code{NULL} and
-@code{data->size} will be zero. However (when necessary) @code{dsize} must
-not have any zero values (a dimension of length zero is not defined).
+Initialize the given data structure (@code{data}) with all the given values.
+Note that the raw input @code{gal_data_t} must already have been allocated 
before calling this function.
+For a description of each variable see @ref{Generic data container}.
+It will set the values and do the necessary allocations.
+If they are not @code{NULL}, all input arrays (@code{dsize}, @code{wcs}, 
@code{name}, @code{unit}, @code{comment}) are separately copied (allocated) by 
this function for usage in @code{data}, so you can safely use one value to 
initialize many datasets or use statically allocated variables in this function 
call.
+Once you are done with the dataset, you can free all the allocated spaces with 
@code{gal_data_free_contents}.
+
+If @code{array} is not @code{NULL}, it will be directly copied into 
@code{data->array} (based on the total number of elements calculated from 
@code{dsize}) and no new space will be allocated for the array of this dataset, 
this has many low-level advantages and can be used to work on regions of a 
dataset instead of the whole allocated array (see the description under 
@code{block} in @ref{Generic data container} for one example).
+If the given pointer is not the start of an allocated block of memory or it is 
used in multiple datasets, be sure to set it to @code{NULL} (with 
@code{data->array=NULL}) before cleaning up with @code{gal_data_free_contents}.
+
+@code{ndim} may be zero.
+In this case no allocation will occur, @code{data->array} and 
@code{data->dsize} will be set to @code{NULL} and @code{data->size} will be 
zero.
+However (when necessary) @code{dsize} must not have any zero values (a 
dimension of length zero is not defined).
+@end deftypefun
+
+@deftypefun {gal_data_t *} gal_data_alloc_empty (size_t @code{ndim}, size_t 
@code{minmapsize}, int @code{quietmmap})
+Allocate an empty dataset with a certain number of dimensions, but no 'array' 
component.
+The @code{size} element will be set to zero and the @code{dsize} array will be 
properly allocated (based on the number of dimensions), but all elements will 
be zero.
+This is useful in scenarios where you just need a @code{gal_data_t} for meta 
data.
 @end deftypefun
 
+
 @deftypefun void gal_data_free_contents (gal_data_t @code{*data})
 Free all the non-@code{NULL} pointers in @code{gal_data_t} except for 
@code{next} and @code{block}.
 All freed arrays are set to @code{NULL}.
@@ -31324,34 +31459,25 @@ For a complete description of @code{gal_data_t} and 
its contents, see @ref{Gener
 @end deftypefun
 
 @deftypefun void gal_data_free (gal_data_t @code{*data})
-Free all the non-@code{NULL} pointers in @code{gal_data_t}, then free the
-actual data structure.
+Free all the non-@code{NULL} pointers in @code{gal_data_t}, then free the 
actual data structure.
 @end deftypefun
 
 @node Arrays of datasets, Copying datasets, Dataset allocation, Library data 
container
 @subsubsection Arrays of datasets
 
-Gnuastro's generic data container (@code{gal_data_t}) is a very versatile
-structure that can be used in many higher-level contexts. One such
-higher-level construct is an array of @code{gal_data_t} structures to
-simplify the allocation (and later cleaning) of several @code{gal_data_t}s
-that are related.
-
-For example, each column in a table is usually represented by one
-@code{gal_data_t} (so it has its own name, data type, units, etc.). A
-table (with many columns) can be seen as an array of @code{gal_data_t}s
-(when the number of columns is known a-priori). The functions below are
-defined to create a cleared array of data structures and to free them when
-none are necessary any more. These functions are declared in
-@file{gnuastro/data.h} which is also visible from the function names (see
-@ref{Gnuastro library}).
+Gnuastro's generic data container (@code{gal_data_t}) is a very versatile 
structure that can be used in many higher-level contexts.
+One such higher-level construct is an array of @code{gal_data_t} structures to 
simplify the allocation (and later cleaning) of several @code{gal_data_t}s that 
are related.
+
+For example, each column in a table is usually represented by one 
@code{gal_data_t} (so it has its own name, data type, units, etc.).
+A table (with many columns) can be seen as an array of @code{gal_data_t}s 
(when the number of columns is known a-priori).
+The functions below are defined to create a cleared array of data structures 
and to free them when none are necessary any more.
+These functions are declared in @file{gnuastro/data.h} which is also visible 
from the function names (see @ref{Gnuastro library}).
 
 @deftypefun {gal_data_t *} gal_data_array_calloc (size_t @code{size})
-Allocate an array of @code{gal_data_t} with @code{size} elements. This
-function will also initialize all the values (@code{NULL} for pointers and
-0 for other types). You can use @code{gal_data_initialize} to fill each
-element of the array afterwards. The following code snippet is one example
-of doing this.
+Allocate an array of @code{gal_data_t} with @code{size} elements.
+This function will also initialize all the values (@code{NULL} for pointers 
and 0 for other types).
+You can use @code{gal_data_initialize} to fill each element of the array 
afterwards.
+The following code snippet is one example of doing this.
 
 @example
 size_t i;
@@ -31364,10 +31490,8 @@ gal_data_array_free(dataarr, 10, 1);
 @end deftypefun
 
 @deftypefun void gal_data_array_free (gal_data_t @code{*dataarr}, size_t 
@code{num}, int @code{free_array})
-Free all the @code{num} elements within @code{dataarr} and the actual
-allocated array. If @code{free_array} is not zero, then the @code{array}
-element of all the datasets will also be freed, see @ref{Generic data
-container}.
+Free all the @code{num} elements within @code{dataarr} and the actual 
allocated array.
+If @code{free_array} is not zero, then the @code{array} element of all the 
datasets will also be freed, see @ref{Generic data container}.
 @end deftypefun
 
 @deftypefun {gal_data_t **} gal_data_array_ptr_calloc (size_t @code{size})
@@ -34538,12 +34662,6 @@ Select random values from a custom distribution 
(defined by a histogram).
 For more, see the description of the respective operators in @ref{Generating 
random numbers}.
 @end deffn
 
-@deffn Macro GAL_ARITHMETIC_OP_SIZE
-Size operator that will return a single value for datasets of any kind. When 
@code{gal_arithmetic} is called with this operator, it requires two arguments.
-The first is the dataset, and the second is a single integer value.
-The output type is a single integer.
-@end deffn
-
 @deffn Macro GAL_ARITHMETIC_OP_STITCH
 Stitch a list of input datasets along the requested dimension.
 See the description of the @code{stitch} operator in Arithmetic 
(@ref{Dimensionality changing operators}).
@@ -34632,6 +34750,28 @@ Constant components of the @command{load-col-} 
operator (see @ref{Loading extern
 These are just fixed strings (and their lengths) that are placed in between 
the various components of that operator to allow choosing a certain column of a 
certain HDU of a certain file.
 @end deffn
 
+@deffn  Macro GAL_ARITHMETIC_OP_INDEX
+@deffnx Macro GAL_ARITHMETIC_OP_COUNTER
+@deffnx Macro GAL_ARITHMETIC_OP_INDEXONLY
+@deffnx Macro GAL_ARITHMETIC_OP_COUNTERONLY
+Return a dataset with the same number of elements and dimensionality as the 
first (and only!) input dataset.
+But each output pixel's value will be replaced by its index (counting from 0) 
or counter (counting from 1).
+Note that the @code{GAL_ARITHMETIC_OP_INDEX} and 
@code{GAL_ARITHMETIC_OP_INDEXONLY} operators are identical within the library 
(same for the counter operators).
+They are given separate macros here to help the higher-level callers to manage 
their inputs separately (see @ref{Size and position operators}).
+@end deffn
+
+@deffn Macro GAL_ARITHMETIC_OP_SIZE
+Size operator that will return a single value for datasets of any kind. When 
@code{gal_arithmetic} is called with this operator, it requires two arguments.
+The first is the dataset, and the second is a single integer value.
+The output type is a single integer.
+@end deffn
+
+@deffn Macro GAL_ARITHMETIC_OP_SWAP
+Return the first dataset, but with the second dataset being placed in the 
@code{next} element of the first.
+This is useful to swap the operators on the stacks of the higher-level 
programs that call the arithmetic library.
+@end deffn
+
+
 @deftypefun {gal_data_t *} gal_arithmetic (int @code{operator}, size_t 
@code{numthreads}, int @code{flags}, ...)
 Apply the requested arithmetic operator on the operand(s).
 The @emph{operator} is identified through the macros above (that start with 
@code{GAL_ARITHMETIC_OP_}).
diff --git a/lib/arithmetic.c b/lib/arithmetic.c
index fa2b36e5..da5d33f2 100644
--- a/lib/arithmetic.c
+++ b/lib/arithmetic.c
@@ -2738,6 +2738,78 @@ arithmetic_makenew(gal_data_t *sizes)
 
 
 
+/* Build a dataset with the index of each element. */
+#define ARITHMETIC_FILL_INDEX(IT) {                                     \
+    IT i=0, *a=out->array, *af=a+out->size; do *a = i++; while(++a<af); }
+#define ARITHMETIC_FILL_COUNTER(IT) {                                   \
+    IT i=0, *a=out->array, *af=a+out->size; do *a = ++i; while(++a<af); }
+static gal_data_t *
+arithmetic_index_counter(gal_data_t *input, int operator, int flags)
+{
+  uint8_t otype;
+  gal_data_t *out;
+  size_t isize=input->size;
+
+  /* Find the best type for the output. */
+  if(      isize<UINT8_MAX  ) otype=GAL_TYPE_UINT8;
+  else if( isize<UINT16_MAX ) otype=GAL_TYPE_UINT16;
+  else if( isize<UINT32_MAX ) otype=GAL_TYPE_UINT32;
+  else                        otype=GAL_TYPE_UINT64;
+
+  /* Allocate the necessary dataset. */
+  out=gal_data_alloc(NULL, otype, input->ndim, input->dsize, NULL,
+                     0, input->minmapsize, input->quietmmap, NULL,
+                     NULL, NULL);
+
+  /* Do the respective operation. */
+  switch(operator)
+    {
+    case GAL_ARITHMETIC_OP_INDEX:
+    case GAL_ARITHMETIC_OP_INDEXONLY:
+      switch(otype)
+        {
+        case GAL_TYPE_UINT8:  ARITHMETIC_FILL_INDEX( uint8_t  ); break;
+        case GAL_TYPE_UINT16: ARITHMETIC_FILL_INDEX( uint16_t ); break;
+        case GAL_TYPE_UINT32: ARITHMETIC_FILL_INDEX( uint32_t ); break;
+        case GAL_TYPE_UINT64: ARITHMETIC_FILL_INDEX( uint64_t ); break;
+        default:
+          error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at '%s' "
+                "to find and fix the problem. The type code '%d' isn't "
+                "recognized for the 'index' operator", __func__,
+                PACKAGE_BUGREPORT, otype);
+        }
+      break;
+    case GAL_ARITHMETIC_OP_COUNTER:
+    case GAL_ARITHMETIC_OP_COUNTERONLY:
+      switch(otype)
+        {
+        case GAL_TYPE_UINT8:  ARITHMETIC_FILL_COUNTER( uint8_t  ); break;
+        case GAL_TYPE_UINT16: ARITHMETIC_FILL_COUNTER( uint16_t ); break;
+        case GAL_TYPE_UINT32: ARITHMETIC_FILL_COUNTER( uint32_t ); break;
+        case GAL_TYPE_UINT64: ARITHMETIC_FILL_COUNTER( uint64_t ); break;
+        default:
+          error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at '%s' "
+                "to find and fix the problem. The type code '%d' isn't "
+                "recognized for the 'counter' operator", __func__,
+                PACKAGE_BUGREPORT, otype);
+        }
+      break;
+    default:
+      error(EXIT_FAILURE, 0, "%s: a bug! Please contact us at '%s' to "
+            "find and fix the problem. The code '%d' isn't recognized "
+            "for the 'operator' variable", __func__, PACKAGE_BUGREPORT,
+            operator);
+    }
+
+  /* Clean up and return. */
+  if(flags & GAL_ARITHMETIC_FLAG_FREE) gal_data_free(input);
+  return out;
+}
+
+
+
+
+
 gal_data_t *
 gal_arithmetic_load_col(char *str, int searchin, int ignorecase,
                         size_t minmapsize, int quietmmap)
@@ -3021,10 +3093,6 @@ gal_arithmetic_set_operator(char *string, size_t 
*num_operands)
   else if (!strcmp(string, "random-from-hist-raw"))
     { op=GAL_ARITHMETIC_OP_RANDOM_FROM_HIST_RAW; *num_operands=3; }
 
-  /* The size operator */
-  else if (!strcmp(string, "size"))
-    { op=GAL_ARITHMETIC_OP_SIZE;              *num_operands=2;  }
-
   /* Stitching */
   else if (!strcmp(string, "stitch"))
     { op=GAL_ARITHMETIC_OP_STITCH;            *num_operands=-1; }
@@ -3113,9 +3181,21 @@ gal_arithmetic_set_operator(char *string, size_t 
*num_operands)
   else if (!strcmp(string, "box-around-ellipse"))
     { op=GAL_ARITHMETIC_OP_BOX_AROUND_ELLIPSE;*num_operands=3;  }
 
-  /* New dataset. */
+  /* Size and position operators. */
+  else if (!strcmp(string, "swap"))
+    { op=GAL_ARITHMETIC_OP_SWAP;              *num_operands=2;  }
+  else if (!strcmp(string, "index"))
+    { op=GAL_ARITHMETIC_OP_INDEX;             *num_operands=1;  }
+  else if (!strcmp(string, "indexonly"))
+    { op=GAL_ARITHMETIC_OP_INDEXONLY;         *num_operands=1;  }
+  else if (!strcmp(string, "counter"))
+    { op=GAL_ARITHMETIC_OP_COUNTER;           *num_operands=1;  }
+  else if (!strcmp(string, "counteronly"))
+    { op=GAL_ARITHMETIC_OP_COUNTERONLY;       *num_operands=1;  }
   else if (!strcmp(string, "makenew"))
     { op=GAL_ARITHMETIC_OP_MAKENEW;           *num_operands=-1;  }
+  else if (!strcmp(string, "size"))
+    { op=GAL_ARITHMETIC_OP_SIZE;              *num_operands=2;  }
 
   /* Operator not defined. */
   else
@@ -3254,6 +3334,11 @@ gal_arithmetic_operator_string(int operator)
 
     case GAL_ARITHMETIC_OP_BOX_AROUND_ELLIPSE: return "box-around-ellipse";
 
+    case GAL_ARITHMETIC_OP_SWAP:            return "swap";
+    case GAL_ARITHMETIC_OP_INDEX:           return "index";
+    case GAL_ARITHMETIC_OP_COUNTER:         return "counter";
+    case GAL_ARITHMETIC_OP_INDEXONLY:       return "indexonly";
+    case GAL_ARITHMETIC_OP_COUNTERONLY:     return "counteronly";
     case GAL_ARITHMETIC_OP_MAKENEW:         return "makenew";
 
     default:                                return NULL;
@@ -3447,13 +3532,6 @@ gal_arithmetic(int operator, size_t numthreads, int 
flags, ...)
         out=arithmetic_mknoise(operator, flags, d1, d2);
       break;
 
-    /* Size operator */
-    case GAL_ARITHMETIC_OP_SIZE:
-      d1 = va_arg(va, gal_data_t *);
-      d2 = va_arg(va, gal_data_t *);
-      out=arithmetic_size(operator, flags, d1, d2);
-      break;
-
     /* Stitch multiple datasets. */
     case GAL_ARITHMETIC_OP_STITCH:
       d1 = va_arg(va, gal_data_t *);
@@ -3498,11 +3576,28 @@ gal_arithmetic(int operator, size_t numthreads, int 
flags, ...)
       out=arithmetic_box_around_ellipse(d1, d2, d3, operator, flags);
       break;
 
-    /* Build dataset from scratch. */
+    /* Size and position operators. */
+    case GAL_ARITHMETIC_OP_SIZE:
+      d1 = va_arg(va, gal_data_t *);
+      d2 = va_arg(va, gal_data_t *);
+      out=arithmetic_size(operator, flags, d1, d2);
+      break;
     case GAL_ARITHMETIC_OP_MAKENEW:
       d1 = va_arg(va, gal_data_t *);
       out=arithmetic_makenew(d1);
       break;
+    case GAL_ARITHMETIC_OP_INDEX:
+    case GAL_ARITHMETIC_OP_COUNTER:
+    case GAL_ARITHMETIC_OP_INDEXONLY:
+    case GAL_ARITHMETIC_OP_COUNTERONLY:
+      d1 = va_arg(va, gal_data_t *);
+      out=arithmetic_index_counter(d1, operator, flags);
+      break;
+    case GAL_ARITHMETIC_OP_SWAP:
+      d1 = va_arg(va, gal_data_t *);
+      d2 = va_arg(va, gal_data_t *);
+      d1->next=d2; d2->next=NULL;  out=d1;
+      break;
 
     /* When the operator is not recognized. */
     default:
diff --git a/lib/data.c b/lib/data.c
index a720be26..267a51bb 100644
--- a/lib/data.c
+++ b/lib/data.c
@@ -213,6 +213,37 @@ gal_data_initialize(gal_data_t *data, void *array, uint8_t 
type,
 
 
 
+/* Allocate an empty (meta) dataset with a certain number of dimensions,
+   but no 'array' component, and all 'size' elements set to zero. */
+gal_data_t *
+gal_data_alloc_empty(size_t ndim, size_t minmapsize, int quietmmap)
+{
+  gal_data_t *out;
+  size_t i, *dsize=gal_pointer_allocate(GAL_TYPE_SIZE_T, ndim,
+                                        0, __func__, "dsize");
+
+  /* Fill the 'dsize' array with 1 values (so too much space isn't
+     allocated!), then allocate it. */
+  for(i=0;i<ndim;++i) dsize[i]=1;
+  out=gal_data_alloc(NULL, GAL_TYPE_UINT8, ndim, dsize, NULL, 0,
+                     minmapsize, quietmmap, NULL, NULL, NULL);
+
+  /* Update the sizes. */
+  out->size=0;
+  for(i=0;i<ndim;++i) out->dsize[i]=0;
+
+  /* Clean up the allocated space for 'out->array', and the extra 'dsize',
+     then return. */
+  free(out->array);
+  out->array=NULL;
+  free(dsize);
+  return out;
+}
+
+
+
+
+
 /* Free the allocated contents of a data structure, not the structure
    itsself. The reason that this function is separate from 'gal_data_free'
    is that the data structure might be allocated as an array (statically
diff --git a/lib/gnuastro/arithmetic.h b/lib/gnuastro/arithmetic.h
index 47fbfcdc..5304f8cb 100644
--- a/lib/gnuastro/arithmetic.h
+++ b/lib/gnuastro/arithmetic.h
@@ -192,8 +192,6 @@ enum gal_arithmetic_operators
   GAL_ARITHMETIC_OP_RANDOM_FROM_HIST,/* Randoms from a histogram (uniform).*/
   GAL_ARITHMETIC_OP_RANDOM_FROM_HIST_RAW,/* Randoms from a histogram (raw).*/
 
-  GAL_ARITHMETIC_OP_SIZE,         /* Size of the dataset along an axis     */
-
   GAL_ARITHMETIC_OP_STITCH,       /* Stitch multiple datasets together.    */
 
   GAL_ARITHMETIC_OP_TO_UINT8,     /* Convert to uint8_t.                   */
@@ -209,7 +207,14 @@ enum gal_arithmetic_operators
 
   GAL_ARITHMETIC_OP_BOX_AROUND_ELLIPSE, /* Width/Height of box over ellipse*/
 
+  /* Meta operators */
   GAL_ARITHMETIC_OP_MAKENEW,      /* Build a new dataset, containing zeros.*/
+  GAL_ARITHMETIC_OP_SIZE,         /* Size of the dataset along an axis     */
+  GAL_ARITHMETIC_OP_INDEX,        /* New with the index (counting from 0). */
+  GAL_ARITHMETIC_OP_COUNTER,      /* New with the index (counting from 0). */
+  GAL_ARITHMETIC_OP_INDEXONLY,    /* New with the index (counting from 0). */
+  GAL_ARITHMETIC_OP_COUNTERONLY,  /* New with the index (counting from 1). */
+  GAL_ARITHMETIC_OP_SWAP,         /* Swap the top two operands.            */
 
   GAL_ARITHMETIC_OP_LAST_CODE,    /* Last code of the library operands.    */
 };
diff --git a/lib/gnuastro/data.h b/lib/gnuastro/data.h
index f394a612..9dd8a00a 100644
--- a/lib/gnuastro/data.h
+++ b/lib/gnuastro/data.h
@@ -244,6 +244,9 @@ gal_data_initialize(gal_data_t *data, void *array, uint8_t 
type,
                     int clear, size_t minmapsize, int quietmmap,
                     char *name, char *unit, char *comment);
 
+gal_data_t *
+gal_data_alloc_empty(size_t ndim, size_t minmapsize, int quietmmap);
+
 void
 gal_data_free_contents(gal_data_t *data);
 
diff --git a/lib/tableintern.c b/lib/tableintern.c
index 5ff56fe8..05632ad0 100644
--- a/lib/tableintern.c
+++ b/lib/tableintern.c
@@ -318,9 +318,18 @@ gal_tableintern_col_print_info(gal_data_t *col, int 
tableformat,
         }
       else width=( col->disp_width<=0 ? GAL_TABLE_DEF_WIDTH_INT
                     : col->disp_width );
-      precision=( col->disp_precision==GAL_BLANK_INT
-                  ? GAL_TABLE_DEF_PRECISION_INT
-                  : col->disp_precision );
+
+      /* For integers, there shouldn't be any default precision. If the
+         caller didn't specify it, it should just be the full set of
+         available digits. because for 'printf' the precision of integers
+         means (according to the GNU C Library manual): "the minimum number
+         of digits to appear; leading zeros are produced if necessary. If
+         you don’t specify a precision, the number is printed with as many
+         digits as it needs.  If you convert a value of zero with an
+         explicit precision of zero, then no characters at all are
+         produced". The main problem is that integer columns can contain
+         zero and in that case '%.0u' or '%.0d' will not print anything! */
+      precision=col->disp_precision;
       break;
 
 
@@ -332,9 +341,7 @@ gal_tableintern_col_print_info(gal_data_t *col, int 
tableformat,
       fmt[0] = tableformat==GAL_TABLE_FORMAT_TXT ? 'd' : 'I';
       width = ( col->disp_width<=0 ? GAL_TABLE_DEF_WIDTH_INT
                 : col->disp_width );
-      precision=( col->disp_precision==GAL_BLANK_INT
-                  ? GAL_TABLE_DEF_PRECISION_INT
-                  : col->disp_precision );
+      precision=col->disp_precision; /* See description in unsigned types.*/
       break;
 
 
@@ -345,9 +352,7 @@ gal_tableintern_col_print_info(gal_data_t *col, int 
tableformat,
       fmt[0] = tableformat==GAL_TABLE_FORMAT_TXT ? 'd' : 'I';
       width=( col->disp_width<=0 ? GAL_TABLE_DEF_WIDTH_LINT
               : col->disp_width );
-      precision=( col->disp_precision==GAL_BLANK_INT
-                  ? GAL_TABLE_DEF_PRECISION_INT
-                  : col->disp_precision );
+      precision=col->disp_precision; /* See description in unsigned types.*/
       break;
 
 
diff --git a/lib/txt.c b/lib/txt.c
index e21e27ca..b6accb2d 100644
--- a/lib/txt.c
+++ b/lib/txt.c
@@ -1436,7 +1436,6 @@ make_fmts_for_printf(gal_data_t *datall, int leftadjust, 
size_t *len)
                              ? strlen(fmts[i*FMTS_COLS+2])
                              : data->disp_width );
 
-
       /* Print the result into the allocated string and add its length to
          the final length of the overall format statement. The space in the
          end of 'fmts[i*2]' is to ensure that the columns don't merge, even
@@ -1461,7 +1460,6 @@ make_fmts_for_printf(gal_data_t *datall, int leftadjust, 
size_t *len)
                                 data->disp_precision, lng, fmt);
         }
 
-
       /* Set the string for the Gnuastro type. For strings, we also need to
          write the maximum number of characters.*/
       if(data->type==GAL_TYPE_STRING)



reply via email to

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