gnuastro-commits
[Top][All Lists]
Advanced

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

[gnuastro-commits] master 8a8ae1e: MakeCatalog: new columns related to m


From: Mohammad Akhlaghi
Subject: [gnuastro-commits] master 8a8ae1e: MakeCatalog: new columns related to maximum value and half of it
Date: Sun, 18 Oct 2020 18:04:12 -0400 (EDT)

branch: master
commit 8a8ae1e6aecc00fab7332c915bae6c0437c9cbbc
Author: Mohammad Akhlaghi <mohammad@akhlaghi.org>
Commit: Mohammad Akhlaghi <mohammad@akhlaghi.org>

    MakeCatalog: new columns related to maximum value and half of it
    
    Recently the '--fwhmobs' option was added, which can be a very good measure
    of the flux distribution within a labeled region. However, we didn't have
    any method to access its reliability.
    
    Therefore with this commit there is a new column called '--maximum', which
    returns the exact maximum used by the FWHM measurement. By comparing the
    value in '--maximum' with the local standard deviation we can see how
    reliable the result is. This is described fully in the book under
    '--fwhmobs'.
    
    Also, since the maximum value is much less affected by noise or deblending
    issues (compared to the total flux), the area of the pixels above this
    value can also shed some light on the profile shape and surface
    brightness. So three more columns were added for such measurements:
    '--halfmaxarea', '--halfmaxsum' and '--halfmaxsb'.
---
 NEWS                      |  6 ++-
 bin/mkcatalog/args.h      | 56 +++++++++++++++++++++++++++
 bin/mkcatalog/columns.c   | 98 +++++++++++++++++++++++++++++++++++++++++++++--
 bin/mkcatalog/main.h      |  6 ++-
 bin/mkcatalog/mkcatalog.c |  3 ++
 bin/mkcatalog/parse.c     | 66 ++++++++++++++++++++++---------
 bin/mkcatalog/ui.c        |  4 ++
 bin/mkcatalog/ui.h        |  6 ++-
 doc/gnuastro.texi         | 33 ++++++++++++++--
 9 files changed, 249 insertions(+), 29 deletions(-)

diff --git a/NEWS b/NEWS
index 6dd8ad6..efbfd17 100644
--- a/NEWS
+++ b/NEWS
@@ -42,8 +42,12 @@ See the end of the file for license conditions.
      within the image.
 
   MakeCatalog:
+   --maximum: maximum value of labeled regions pixels (clump/object).
    --areaarcsec2: area of labeled region (clump/object) in arcsec^2.
-   --fwhmobs: The observed FWHM in pixels, along the major axis.
+   --fwhmobs: observed FWHM in pixels, along the major axis.
+   --halfmaxarea: number of pixels with a value larger than half the maximum.
+   --halfmaxsum: sum of pixels with a value larger than half the maximum.
+   --halfmaxsb: surf. brightness within half max (--halfmaxsum/--halfmaxarea).
    --fracsum: fractions to use in '--fracsumarea1' or '--fracsumarea2'.
    --fracsumarea1: area of given fraction of the summed object or clump values.
    --fracsumarea2: area of given fraction of the summed object or clump values.
diff --git a/bin/mkcatalog/args.h b/bin/mkcatalog/args.h
index 35ee8f0..f69c419 100644
--- a/bin/mkcatalog/args.h
+++ b/bin/mkcatalog/args.h
@@ -1145,6 +1145,20 @@ struct argp_option program_options[] =
       ui_column_codes_ll
     },
     {
+      "maximum",
+      UI_KEY_MAXIMUM,
+      0,
+      0,
+      "Maximum value (mean of top three pixels)",
+      UI_GROUP_COLUMNS_BRIGHTNESS,
+      0,
+      GAL_TYPE_INVALID,
+      GAL_OPTIONS_RANGE_ANY,
+      GAL_OPTIONS_NOT_MANDATORY,
+      GAL_OPTIONS_NOT_SET,
+      ui_column_codes_ll
+    },
+    {
       "magnitude",
       UI_KEY_MAGNITUDE,
       0,
@@ -1645,6 +1659,48 @@ struct argp_option program_options[] =
       ui_column_codes_ll
     },
     {
+      "halfmaxarea",
+      UI_KEY_HALFMAXAREA,
+      0,
+      0,
+      "No. pixels valued above half the max.",
+      UI_GROUP_COLUMNS_MORPHOLOGY,
+      0,
+      GAL_TYPE_INVALID,
+      GAL_OPTIONS_RANGE_ANY,
+      GAL_OPTIONS_NOT_MANDATORY,
+      GAL_OPTIONS_NOT_SET,
+      ui_column_codes_ll
+    },
+    {
+      "halfmaxsum",
+      UI_KEY_HALFMAXSUM,
+      0,
+      0,
+      "Sum of pixels above half the maxium.",
+      UI_GROUP_COLUMNS_MORPHOLOGY,
+      0,
+      GAL_TYPE_INVALID,
+      GAL_OPTIONS_RANGE_ANY,
+      GAL_OPTIONS_NOT_MANDATORY,
+      GAL_OPTIONS_NOT_SET,
+      ui_column_codes_ll
+    },
+    {
+      "halfmaxsb",
+      UI_KEY_HALFMAXSB,
+      0,
+      0,
+      "Surface brightness within half the maximum.",
+      UI_GROUP_COLUMNS_MORPHOLOGY,
+      0,
+      GAL_TYPE_INVALID,
+      GAL_OPTIONS_RANGE_ANY,
+      GAL_OPTIONS_NOT_MANDATORY,
+      GAL_OPTIONS_NOT_SET,
+      ui_column_codes_ll
+    },
+    {
       "halfsumarea",
       UI_KEY_HALFSUMAREA,
       0,
diff --git a/bin/mkcatalog/columns.c b/bin/mkcatalog/columns.c
index 932bdb5..218e592 100644
--- a/bin/mkcatalog/columns.c
+++ b/bin/mkcatalog/columns.c
@@ -1202,6 +1202,20 @@ columns_define_alloc(struct mkcatalogparams *p)
                                    ciflag[ CCOL_RIV_SUM ] = 1;
           break;
 
+        case UI_KEY_MAXIMUM:
+          name           = "MAXIMUM";
+          unit           = MKCATALOG_NO_UNIT;
+          ocomment       = "Maximum of sky subtracted values.";
+          ccomment       = "Maximum of pixels subtracted by rivers.";
+          otype          = GAL_TYPE_FLOAT32;
+          ctype          = GAL_TYPE_FLOAT32;
+          disp_fmt       = GAL_TABLE_DISPLAY_FMT_GENERAL;
+          disp_width     = 10;
+          disp_precision = 5;
+          oiflag[ OCOL_NUM     ] = ciflag[ CCOL_NUM     ] = 1;
+          oiflag[ OCOL_MAXIMUM ] = ciflag[ CCOL_MAXIMUM ] = 1;
+          break;
+
         case UI_KEY_SIGCLIPNUMBER:
           name           = "SIGCLIP-NUMBER";
           unit           = "counter";
@@ -1665,10 +1679,53 @@ columns_define_alloc(struct mkcatalogparams *p)
           oiflag[ OCOL_NUMHALFSUM ] = ciflag[ CCOL_NUMHALFSUM ] = 1;
           break;
 
+        case UI_KEY_HALFMAXAREA:
+          name           = "HALF_MAX_AREA";
+          unit           = "counter";
+          ocomment       = "Number of pixels with a value larger than half the 
maximum.";
+          ccomment       = ocomment;
+          otype          = GAL_TYPE_INT32;
+          ctype          = GAL_TYPE_INT32;
+          disp_fmt       = 0;
+          disp_width     = 6;
+          disp_precision = 0;
+          oiflag[ OCOL_NUM         ] = ciflag[ CCOL_NUM         ] = 1;
+          oiflag[ OCOL_NUMHALFMAX  ] = ciflag[ CCOL_NUMHALFMAX  ] = 1;
+          break;
+
+        case UI_KEY_HALFMAXSUM:
+          name           = "HALF_MAX_SUM";
+          unit           = MKCATALOG_NO_UNIT;
+          ocomment       = "Sum of pixels with a value larger than half the 
maximum.";
+          ccomment       = ocomment;
+          otype          = GAL_TYPE_FLOAT32;
+          ctype          = GAL_TYPE_FLOAT32;
+          disp_fmt       = 0;
+          disp_width     = 6;
+          disp_precision = 0;
+          oiflag[ OCOL_NUM         ] = ciflag[ CCOL_NUM         ] = 1;
+          oiflag[ OCOL_SUMHALFMAX  ] = ciflag[ CCOL_SUMHALFMAX  ] = 1;
+          break;
+
+        case UI_KEY_HALFMAXSB:
+          name           = "HALF_MAX_SB";
+          unit           = MKCATALOG_NO_UNIT;
+          ocomment       = "Brightness within half the maximum, divided by its 
area.";
+          ccomment       = ocomment;
+          otype          = GAL_TYPE_FLOAT32;
+          ctype          = GAL_TYPE_FLOAT32;
+          disp_fmt       = 0;
+          disp_width     = 6;
+          disp_precision = 0;
+          oiflag[ OCOL_NUM         ] = ciflag[ CCOL_NUM         ] = 1;
+          oiflag[ OCOL_NUMHALFMAX  ] = ciflag[ CCOL_NUMHALFMAX  ] = 1;
+          oiflag[ OCOL_SUMHALFMAX  ] = ciflag[ CCOL_SUMHALFMAX  ] = 1;
+          break;
+
         case UI_KEY_HALFSUMSB:
           name           = "HALF_SUM_SB";
           unit           = MKCATALOG_NO_UNIT;
-          ocomment       = "Surface brightness within area containing half 
total sum.";
+          ocomment       = "Half the brightness, divided by its area.";
           otype          = GAL_TYPE_FLOAT32;
           ctype          = GAL_TYPE_FLOAT32;
           disp_fmt       = GAL_TABLE_DISPLAY_FMT_GENERAL;
@@ -2329,6 +2386,10 @@ columns_fill(struct mkcatalog_passparams *pp)
                                       : NAN );
           break;
 
+        case UI_KEY_MAXIMUM:
+          ((float *)colarr)[oind] = oi[ OCOL_MAXIMUM ];
+          break;
+
         case UI_KEY_SIGCLIPNUMBER:
           ((int32_t *)colarr)[oind] = oi[ OCOL_SIGCLIPNUM ];
           break;
@@ -2438,6 +2499,18 @@ columns_fill(struct mkcatalog_passparams *pp)
           ((int32_t *)colarr)[oind] = oi[OCOL_NUMHALFSUM];
           break;
 
+        case UI_KEY_HALFMAXAREA:
+          ((int32_t *)colarr)[oind] = oi[OCOL_NUMHALFMAX];
+          break;
+
+        case UI_KEY_HALFMAXSUM:
+          ((float *)colarr)[oind] = oi[OCOL_SUMHALFMAX];
+          break;
+
+        case UI_KEY_HALFMAXSB:
+          ((float *)colarr)[oind] = oi[OCOL_SUMHALFMAX]/oi[OCOL_NUMHALFMAX];
+          break;
+
         case UI_KEY_FRACSUMAREA1:
           ((int32_t *)colarr)[oind] = oi[OCOL_NUMFRACSUM1];
           break;
@@ -2472,7 +2545,7 @@ columns_fill(struct mkcatalog_passparams *pp)
             }
           tmp = sqrt( oi[tmpind]/(tmp*M_PI) );
           if(key==UI_KEY_FWHMOBS)
-            ((float *)colarr)[oind] = (tmp<1e-6 ? NAN : tmp)*2;
+            ((float *)colarr)[oind] = tmp<1e-6 ? NAN : (tmp*2);
           else
             ((float *)colarr)[oind] = tmp<1e-6 ? NAN : tmp;
           break;
@@ -2623,6 +2696,10 @@ columns_fill(struct mkcatalog_passparams *pp)
                                         ? ci[ CCOL_MEDIAN ] : NAN );
             break;
 
+          case UI_KEY_MAXIMUM:
+            ((float *)colarr)[cind] = ci[ CCOL_MAXIMUM ];
+            break;
+
           case UI_KEY_SIGCLIPNUMBER:
             ((int32_t *)colarr)[cind] = ci[ CCOL_SIGCLIPNUM ];
             break;
@@ -2739,6 +2816,18 @@ columns_fill(struct mkcatalog_passparams *pp)
             ((int32_t *)colarr)[cind] = ci[CCOL_NUMHALFSUM];
             break;
 
+          case UI_KEY_HALFMAXAREA:
+            ((int32_t *)colarr)[cind] = ci[CCOL_NUMHALFMAX];
+            break;
+
+          case UI_KEY_HALFMAXSUM:
+            ((float *)colarr)[cind] = ci[CCOL_SUMHALFMAX];
+            break;
+
+          case UI_KEY_HALFMAXSB:
+            ((float *)colarr)[cind] = ci[CCOL_SUMHALFMAX]/ci[CCOL_NUMHALFMAX];
+            break;
+
           case UI_KEY_FRACSUMAREA1:
             ((int32_t *)colarr)[cind] = ci[CCOL_NUMFRACSUM1];
             break;
@@ -2765,7 +2854,10 @@ columns_fill(struct mkcatalog_passparams *pp)
               case UI_KEY_FRACRADIUSOBS2: tmpind=CCOL_NUMFRACSUM2; break;
               }
             tmp = sqrt( ci[tmpind]/(tmp*M_PI) );
-            ((float *)colarr)[cind] = tmp<1e-6 ? NAN : tmp;
+            if(key==UI_KEY_FWHMOBS)
+              ((float *)colarr)[cind] = tmp<1e-6 ? NAN : (tmp*2);
+            else
+              ((float *)colarr)[cind] = tmp<1e-6 ? NAN : tmp;
             break;
 
           default:
diff --git a/bin/mkcatalog/main.h b/bin/mkcatalog/main.h
index 7801679..1e5b68e 100644
--- a/bin/mkcatalog/main.h
+++ b/bin/mkcatalog/main.h
@@ -76,7 +76,8 @@ enum objectcols
     OCOL_NUMXY,          /* Number of values in the first two dims.   */
     OCOL_SUM,            /* Sum of (value-sky) in object.             */
     OCOL_SUM_VAR,        /* Variance including values (not just sky). */
-    OCOL_MEDIAN,         /* Median of value in object.                */
+    OCOL_MEDIAN,         /* Median of values in object.               */
+    OCOL_MAXIMUM,        /* Maximum value in object.                  */
     OCOL_SIGCLIPNUM,     /* Sigma-clipped mean of this object.        */
     OCOL_SIGCLIPSTD,     /* Sigma-clipped mean of this object.        */
     OCOL_SIGCLIPMEAN,    /* Sigma-clipped mean of this object.        */
@@ -110,6 +111,7 @@ enum objectcols
     OCOL_UPPERLIMIT_Q,   /* Quantile of object in random distribution.*/
     OCOL_UPPERLIMIT_SKEW,/* (Mean-Median)/STD of random distribution. */
     OCOL_NUMHALFMAX,     /* Area/Number of pixels above half of max.  */
+    OCOL_SUMHALFMAX,     /* Sum of pixels above half of max.          */
     OCOL_NUMHALFSUM,     /* Area/Number containing half of total sum. */
     OCOL_NUMFRACSUM1,    /* Area/Number containing frac of total sum. */
     OCOL_NUMFRACSUM2,    /* Area/Number containing frac of total sum. */
@@ -137,6 +139,7 @@ enum clumpcols
     CCOL_SUM,            /* River subtracted brightness.              */
     CCOL_SUM_VAR,        /* Variance including values (not just sky). */
     CCOL_MEDIAN,         /* Median of values in clump.                */
+    CCOL_MAXIMUM,        /* Maximum value in clump.                   */
     CCOL_SIGCLIPNUM,     /* Sigma-clipped mean of this clump.         */
     CCOL_SIGCLIPSTD,     /* Sigma-clipped mean of this clump.         */
     CCOL_SIGCLIPMEAN,    /* Sigma-clipped mean of this clump.         */
@@ -179,6 +182,7 @@ enum clumpcols
     CCOL_UPPERLIMIT_Q,   /* Quantile of object in random distribution.*/
     CCOL_UPPERLIMIT_SKEW,/* (Mean-Median)/STD of random distribution. */
     CCOL_NUMHALFMAX,     /* Area/Number of pixels above half of max.  */
+    CCOL_SUMHALFMAX,     /* Sum of pixels above half of max.          */
     CCOL_NUMHALFSUM,     /* Area/Number containing half of total sum. */
     CCOL_NUMFRACSUM1,    /* Area/Number containing frac of total sum. */
     CCOL_NUMFRACSUM2,    /* Area/Number containing frac of total sum. */
diff --git a/bin/mkcatalog/mkcatalog.c b/bin/mkcatalog/mkcatalog.c
index d05f80d..c8186ec 100644
--- a/bin/mkcatalog/mkcatalog.c
+++ b/bin/mkcatalog/mkcatalog.c
@@ -182,6 +182,9 @@ mkcatalog_single_object(void *in_prm)
       /* If an order-based calculation is requested, another pass is
          necessary. */
       if( p->oiflag[ OCOL_MEDIAN ]
+          || p->oiflag[ OCOL_MAXIMUM ]
+          || p->oiflag[ OCOL_SUMHALFMAX ]
+          || p->oiflag[ OCOL_NUMHALFMAX ]
           || p->oiflag[ OCOL_NUMHALFSUM ]
           || p->oiflag[ OCOL_SIGCLIPNUM ]
           || p->oiflag[ OCOL_SIGCLIPSTD ]
diff --git a/bin/mkcatalog/parse.c b/bin/mkcatalog/parse.c
index 5bd36e8..ab7744a 100644
--- a/bin/mkcatalog/parse.c
+++ b/bin/mkcatalog/parse.c
@@ -1182,7 +1182,7 @@ parse_clumps(struct mkcatalog_passparams *pp)
 
 
 
-static double
+static size_t
 parse_frac_find(gal_data_t *sorted_d, double value, double frac, int dosum)
 {
   size_t i;
@@ -1211,11 +1211,12 @@ parse_area_of_frac_sum(struct mkcatalog_passparams *pp, 
gal_data_t *values,
 {
   struct mkcatalogparams *p=pp->p;
 
-  double max, *sorted;
+  size_t i, ind;
   gal_data_t *sorted_d;
+  double sum=0.0f, max, *sorted;
   uint8_t *flag = o1c0 ? p->oiflag : p->ciflag;
-  double sum = o1c0 ? outarr[OCOL_SUM] : outarr[CCOL_SUM];
   double *fracsum = p->fracsum ? p->fracsum->array : NULL;
+  double sumlab = o1c0 ? outarr[OCOL_SUM] : outarr[CCOL_SUM];
 
   /* Allocate the array to use. */
   sorted_d = ( values->type==GAL_TYPE_FLOAT64
@@ -1229,26 +1230,45 @@ parse_area_of_frac_sum(struct mkcatalog_passparams *pp, 
gal_data_t *values,
   /* Set the required fractions. */
   if(flag[ o1c0 ? OCOL_NUMHALFSUM : CCOL_NUMHALFSUM ])
     outarr[ o1c0 ? OCOL_NUMHALFSUM : CCOL_NUMHALFSUM ]
-      = parse_frac_find(sorted_d, sum, 0.5f, 1);
+      = parse_frac_find(sorted_d, sumlab, 0.5f, 1);
 
   if(flag[ o1c0 ? OCOL_NUMFRACSUM1 : CCOL_NUMFRACSUM1 ])
     outarr[ o1c0 ? OCOL_NUMFRACSUM1 : CCOL_NUMFRACSUM1 ]
-      = parse_frac_find(sorted_d, sum, fracsum[0], 1);
+      = parse_frac_find(sorted_d, sumlab, fracsum[0], 1);
 
   if(flag[ o1c0 ? OCOL_NUMFRACSUM2 : CCOL_NUMFRACSUM2 ])
     outarr[ o1c0 ? OCOL_NUMFRACSUM2 : CCOL_NUMFRACSUM2 ]
-      = parse_frac_find(sorted_d, sum, fracsum[1], 1);
+      = parse_frac_find(sorted_d, sumlab, fracsum[1], 1);
 
-  /* For the FWHM, we'll use the median of the top three pixels for the
-     maximum (to avoid noise). */
-  if(flag[ o1c0 ? OCOL_NUMHALFMAX : CCOL_NUMHALFMAX ])
+  /* Values related to the maximum. */
+  if( flag[    o1c0 ? OCOL_MAXIMUM    : CCOL_MAXIMUM ]
+      || flag[ o1c0 ? OCOL_NUMHALFMAX : CCOL_NUMHALFMAX ]
+      || flag[ o1c0 ? OCOL_SUMHALFMAX : CCOL_SUMHALFMAX ])
     {
+      /* Set the array and maximum value. We'll use the median of the top
+         three pixels for the maximum (to avoid noise) */
       sorted=sorted_d->array;
       max = ( sorted_d->size>3
               ? (sorted[0]+sorted[1]+sorted[2])/3
               : sorted[0] );
-      outarr[ o1c0 ? OCOL_NUMHALFMAX : CCOL_NUMHALFMAX ]
-        = parse_frac_find(sorted_d, max, 0.5f, 0);
+
+      /* If we want the maximum value, then write it in. */
+      if(flag[ o1c0 ? OCOL_MAXIMUM : CCOL_MAXIMUM ])
+        outarr[ o1c0 ? OCOL_MAXIMUM : CCOL_MAXIMUM ] = max;
+
+      /* The number of pixels within half the maximum. */
+      if(flag[ o1c0 ? OCOL_NUMHALFMAX : CCOL_NUMHALFMAX ])
+        outarr[ o1c0 ? OCOL_NUMHALFMAX : CCOL_NUMHALFMAX ]
+          = parse_frac_find(sorted_d, max, 0.5f, 0);
+
+      /* The sum of the pixels within half the maximum. */
+      if(flag[ o1c0 ? OCOL_SUMHALFMAX : CCOL_SUMHALFMAX ])
+        {
+          sum=0.0f;
+          ind=parse_frac_find(sorted_d, max, 0.5f, 0);
+          for(i=0;i<ind;++i) sum+=sorted[i];
+          outarr[ o1c0 ? OCOL_SUMHALFMAX : CCOL_SUMHALFMAX ]=sum;
+        }
     }
 
   /* For a check.
@@ -1286,6 +1306,8 @@ parse_order_based(struct mkcatalog_passparams *pp)
   if(tmpsize==0)
     {
       if(p->oiflag[ OCOL_MEDIAN        ]) pp->oi[ OCOL_MEDIAN       ] = NAN;
+      if(p->oiflag[ OCOL_MAXIMUM       ]) pp->oi[ OCOL_MAXIMUM      ] = NAN;
+      if(p->oiflag[ OCOL_SUMHALFMAX    ]) pp->oi[ OCOL_SUMHALFMAX   ] = NAN;
       if(p->oiflag[ OCOL_NUMHALFMAX    ]) pp->oi[ OCOL_NUMHALFMAX   ] = 0;
       if(p->oiflag[ OCOL_NUMHALFSUM    ]) pp->oi[ OCOL_NUMHALFSUM   ] = 0;
       if(p->oiflag[ OCOL_NUMFRACSUM1   ]) pp->oi[ OCOL_NUMFRACSUM1  ] = 0;
@@ -1299,14 +1321,16 @@ parse_order_based(struct mkcatalog_passparams *pp)
           {
             ci=&pp->ci[ i * CCOL_NUMCOLS ];
             if(p->ciflag[ CCOL_MEDIAN        ]) ci[ CCOL_MEDIAN      ] = NAN;
-            if(p->oiflag[ OCOL_NUMHALFMAX    ]) ci[ OCOL_NUMHALFMAX  ] = 0;
-            if(p->ciflag[ OCOL_NUMHALFSUM    ]) ci[ CCOL_NUMHALFSUM  ] = 0;
-            if(p->oiflag[ OCOL_NUMFRACSUM1   ]) ci[ OCOL_NUMFRACSUM1 ] = 0;
-            if(p->oiflag[ OCOL_NUMFRACSUM2   ]) ci[ OCOL_NUMFRACSUM2 ] = 0;
+            if(p->ciflag[ CCOL_MAXIMUM       ]) ci[ CCOL_MAXIMUM     ] = NAN;
+            if(p->ciflag[ CCOL_SUMHALFMAX    ]) ci[ CCOL_SUMHALFMAX  ] = NAN;
+            if(p->ciflag[ CCOL_NUMHALFMAX    ]) ci[ CCOL_NUMHALFMAX  ] = 0;
+            if(p->ciflag[ CCOL_NUMHALFSUM    ]) ci[ CCOL_NUMHALFSUM  ] = 0;
+            if(p->ciflag[ CCOL_NUMFRACSUM1   ]) ci[ CCOL_NUMFRACSUM1 ] = 0;
+            if(p->ciflag[ CCOL_NUMFRACSUM2   ]) ci[ CCOL_NUMFRACSUM2 ] = 0;
             if(p->ciflag[ CCOL_SIGCLIPNUM    ]) ci[ CCOL_SIGCLIPNUM  ] = 0;
             if(p->ciflag[ CCOL_SIGCLIPSTD    ]) ci[ CCOL_SIGCLIPSTD  ] = 0;
             if(p->ciflag[ CCOL_SIGCLIPMEAN   ]) ci[ CCOL_SIGCLIPMEAN ] = NAN;
-            if(p->ciflag[ CCOL_SIGCLIPMEDIAN ]) ci[CCOL_SIGCLIPMEDIAN] = NAN;
+            if(p->ciflag[ CCOL_SIGCLIPMEDIAN ]) ci[ CCOL_SIGCLIPMEDIAN] = NAN;
           }
       return;
     }
@@ -1414,7 +1438,9 @@ parse_order_based(struct mkcatalog_passparams *pp)
     }
 
   /* Fractional values. */
-  if( p->oiflag[    OCOL_NUMHALFMAX  ]
+  if( p->oiflag[    OCOL_MAXIMUM     ]
+      || p->oiflag[ OCOL_NUMHALFMAX  ]
+      || p->oiflag[ OCOL_SUMHALFMAX  ]
       || p->oiflag[ OCOL_NUMHALFSUM  ]
       || p->oiflag[ OCOL_NUMFRACSUM1 ]
       || p->oiflag[ OCOL_NUMFRACSUM2 ] )
@@ -1468,8 +1494,10 @@ parse_order_based(struct mkcatalog_passparams *pp)
             }
 
           /* Estimate half of the total sum. */
-          if( p->ciflag[    CCOL_NUMHALFMAX ]
-              || p->ciflag[ CCOL_NUMHALFSUM ]
+          if( p->ciflag[    CCOL_MAXIMUM     ]
+              || p->ciflag[ CCOL_NUMHALFMAX  ]
+              || p->ciflag[ CCOL_SUMHALFMAX  ]
+              || p->ciflag[ CCOL_NUMHALFSUM  ]
               || p->ciflag[ CCOL_NUMFRACSUM1 ]
               || p->ciflag[ CCOL_NUMFRACSUM2 ] )
             parse_area_of_frac_sum(pp, clumpsvals[i], ci, 0);
diff --git a/bin/mkcatalog/ui.c b/bin/mkcatalog/ui.c
index 0fdb289..f87770f 100644
--- a/bin/mkcatalog/ui.c
+++ b/bin/mkcatalog/ui.c
@@ -967,6 +967,7 @@ ui_necessary_inputs(struct mkcatalogparams *p, int *values, 
int *sky,
         case OCOL_SUM:                *values        = 1;          break;
         case OCOL_SUM_VAR:            *values = *std = 1;          break;
         case OCOL_MEDIAN:             *values        = 1;          break;
+        case OCOL_MAXIMUM:            *values        = 1;          break;
         case OCOL_SIGCLIPNUM:         *values        = 1;          break;
         case OCOL_SIGCLIPMEDIAN:      *values        = 1;          break;
         case OCOL_SIGCLIPMEAN:        *values        = 1;          break;
@@ -998,6 +999,7 @@ ui_necessary_inputs(struct mkcatalogparams *p, int *values, 
int *sky,
         case OCOL_UPPERLIMIT_Q:       *values        = 1;          break;
         case OCOL_UPPERLIMIT_SKEW:    *values        = 1;          break;
         case OCOL_NUMHALFMAX:         *values        = 1;          break;
+        case OCOL_SUMHALFMAX:         *values        = 1;          break;
         case OCOL_NUMHALFSUM:         *values        = 1;          break;
         case OCOL_NUMFRACSUM1:        *values        = 1;          break;
         case OCOL_NUMFRACSUM2:        *values        = 1;          break;
@@ -1032,6 +1034,7 @@ ui_necessary_inputs(struct mkcatalogparams *p, int 
*values, int *sky,
           case CCOL_SUM:              *values        = 1;          break;
           case CCOL_SUM_VAR:          *values = *std = 1;          break;
           case CCOL_MEDIAN:           *values        = 1;          break;
+          case CCOL_MAXIMUM:          *values        = 1;          break;
           case CCOL_SIGCLIPNUM:       *values        = 1;          break;
           case CCOL_SIGCLIPMEDIAN:    *values        = 1;          break;
           case CCOL_SIGCLIPMEAN:      *values        = 1;          break;
@@ -1072,6 +1075,7 @@ ui_necessary_inputs(struct mkcatalogparams *p, int 
*values, int *sky,
           case CCOL_UPPERLIMIT_Q:     *values        = 1;          break;
           case CCOL_UPPERLIMIT_SKEW:  *values        = 1;          break;
           case CCOL_NUMHALFMAX:       *values        = 1;          break;
+          case CCOL_SUMHALFMAX:       *values        = 1;          break;
           case CCOL_NUMHALFSUM:       *values        = 1;          break;
           case CCOL_NUMFRACSUM1:      *values        = 1;          break;
           case CCOL_NUMFRACSUM2:      *values        = 1;          break;
diff --git a/bin/mkcatalog/ui.h b/bin/mkcatalog/ui.h
index d139016..8c8af27 100644
--- a/bin/mkcatalog/ui.h
+++ b/bin/mkcatalog/ui.h
@@ -105,10 +105,10 @@ enum option_keys_enum
   UI_KEY_CHECKUPLIM,
   UI_KEY_NOCLUMPSORT,
   UI_KEY_FRACSUM,
-  UI_KEY_AREAARCSEC2,
 
   UI_KEY_OBJID,                         /* Catalog columns. */
   UI_KEY_IDINHOSTOBJ,
+  UI_KEY_AREAARCSEC2,
   UI_KEY_AREAXY,
   UI_KEY_CLUMPSAREA,
   UI_KEY_WEIGHTAREA,
@@ -152,6 +152,7 @@ enum option_keys_enum
   UI_KEY_BRIGHTNESSNORIVER,
   UI_KEY_MEAN,
   UI_KEY_MEDIAN,
+  UI_KEY_MAXIMUM,
   UI_KEY_CLUMPSMAGNITUDE,
   UI_KEY_UPPERLIMIT,
   UI_KEY_UPPERLIMITONESIGMA,
@@ -171,6 +172,9 @@ enum option_keys_enum
   UI_KEY_GEOAXISRATIO,
   UI_KEY_GEOPOSITIONANGLE,
   UI_KEY_FWHMOBS,
+  UI_KEY_HALFMAXAREA,
+  UI_KEY_HALFMAXSUM,
+  UI_KEY_HALFMAXSB,
   UI_KEY_HALFSUMAREA,
   UI_KEY_HALFSUMSB,
   UI_KEY_HALFRADIUSOBS,
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index bb14c01..728a328 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -16163,6 +16163,12 @@ For clumps, the average river flux is subtracted from 
the sky subtracted mean.
 The median sky subtracted value of pixels within the object or clump.
 For clumps, the average river flux is subtracted from the sky subtracted 
median.
 
+@item --maximum
+The maximum value of pixels within the object or clump.
+When the label (object or clump) is larger than three pixels, the maximum is 
actually derived by the mean of the brightest three pixels, not the largest 
pixel value of the same label.
+This is because noise fluctuations can be very strong in the extreme values of 
the objects/clumps due to Poisson noise (which gets stronger as the mean gets 
higher).
+Simply using the maximum pixel value will create a strong scatter in results 
that depend on the maximum (for example the @option{--fwhmobs} option also uses 
this value internally).
+
 @item --sigclip-number
 The number of elements/pixels in the dataset after sigma-clipping the object 
or clump.
 The sigma-clipping parameters can be set with the @option{--sigmaclip} option 
described in @ref{MakeCatalog inputs and basic settings}.
@@ -16278,12 +16284,31 @@ projection onto the first two dimensions would be a 
narrow-band image.
 @item --fwhmobs
 @cindex FWHM
 The full width at half maximum (in units of pixels, along the semi-major axis) 
of the labeled region (object or clump).
-The maximum value is estimated from the mean of the top-three pixels with the 
highest values (to avoid noise fluctuations which can be large due to Poisson 
noise).
-The number of pixels that have half the value of that maximum are then found 
and a radius is estimated from the area.
+The maximum value is estimated from the mean of the top-three pixels with the 
highest values, see the description under @option{--maximum}.
+The number of pixels that have half the value of that maximum are then found 
(value in the @option{--halfmaxarea} column) and a radius is estimated from the 
area.
 See the description under @option{--halfradiusobs} for more on converting area 
to radius along major axis.
 
-It is therefore important to highlight that this FWHM is not estimated from a 
parametric fit to the labeled region: it does not account for the PSF, and it 
will be strongly affected by noise.
-So when half the maximum value is too close to the noise level, the value 
returned by this function is meaningless (its just noise peaks which are 
randomly distributed over the area).
+Because of its non-parametric nature, this column is most reliable on clumps 
and should only be used in objects with great causion.
+This is because objects can have more than one clump (peak with true signal) 
and multiple peaks are not treated separately in objects, so the result of this 
column will be biased.
+
+Also, because of its non-parametric nature, this FWHM it does not account for 
the PSF, and it will be strongly affected by noise if the object is faint.
+So when half the maximum value (which can be requested using the 
@option{--maximum} column) is too close to the local noise level (which can be 
requested using the @option{--std} column), the value returned in this column 
is meaningless (its just noise peaks which are randomly distributed over the 
area).
+You can therefore use the @option{--maximum} and @option{--std} columns to 
remove unreliable FWHMs.
+For example if a labeled region's maximum is less than 2 times the sky 
standard deviation, the value will certainly be unreliable (half of that is 
@mymath{1\sigma}!).
+For a more reliable value, this fraction should be around 4 (so half the 
maximum is 2@mymath{\sigma}).
+
+@item --halfmaxarea
+The number of pixels containing half the maximum flux within the labeled 
region.
+This option is used to estimate @option{--fwhmobs}, so please read the notes 
there for the cavaeats and necessary precautions.
+
+@item --halfmaxsum
+The sum of the pixel values containing half the maximum flux within the 
labeled region (or those that are counted in @option{--halfmaxarea}).
+This option uses the pixels within @option{--fwhmobs}, so please read the 
notes there for the cavaeats and necessary precautions.
+
+@item --halfmaxsb
+The surface brightness within the region that contains half the maximum value 
of the labeled region.
+This column is just the division of the @option{--halfmaxsum} column by 
@option{--halfmaxarea} column.
+So please see their documentation for the caveats and necessary precautions.
 
 @item --halfsumarea
 The number of pixels that contain half the object or clump's total sum of 
pixels (half the value in the @option{--brightness} column).



reply via email to

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