[gnuastro-commits] master 39a9018 1/2: Hole filling step in NoiseChisel is now configurable

[Mohammad Akhlaghi]

branch: master
commit 39a901894d97c61e570a0045fa779d1d2a5ff92f
Author: Mohammad Akhlaghi <address@hidden>
Commit: Mohammad Akhlaghi <address@hidden>

    Hole filling step in NoiseChisel is now configurable
    
    Until now, the connectivity of the holes in NoiseChisel was hard-wired into
    the code. With this commit a new `--holengb' option is added that can be
    used to specify how a hole is defined.
---
 NEWS                                | 19 +++++++++++--------
 bin/noisechisel/args.h              | 13 +++++++++++++
 bin/noisechisel/astnoisechisel.conf |  1 +
 bin/noisechisel/detection.c         |  2 +-
 bin/noisechisel/main.h              |  1 +
 bin/noisechisel/ui.c                |  3 +++
 bin/noisechisel/ui.h                |  1 +
 doc/gnuastro.texi                   | 10 ++++++++++
 lib/binary.c                        |  6 +++++-
 9 files changed, 46 insertions(+), 10 deletions(-)

diff --git a/NEWS b/NEWS
index e26fdfb..c23d08a 100644
--- a/NEWS
+++ b/NEWS
@@ -25,6 +25,17 @@ GNU Astronomy Utilities NEWS                          -*- 
outline -*-
       configured with these options.
        --outliersclip: Sigma-clipping parameters for the process.
        --outliersigma: Multiple of sigma to define an outlier.
+    --blankasforeground: Treat blank elements as foreground (regions above
+          the threshold) in the binary processing steps: initial erosion
+          and opening as well as the filling holes and opening of
+          pseudo-detections. From this version, by default, blank elements
+          in the dataset are considered to be background, so if a
+          foreground pixel is touching it, it will be eroded. This option
+          is irrelevant if the datasets contains no blank elements, but can
+          help remove false positives that are touching blank elements.
+    --holengb: allows defining the connectivity of the holes that are
+          filled when defining pseudo-detections. Until now, this was
+          hard-wired into the code (=8) and not modifiable at run-time.
     --skyfracnoblank: Ignore blank pixels when estimating the fraction of
           undetected pixels for Sky estimation. NoiseChisel only measures
           the Sky over the tiles that have a sufficiently large fraction of
@@ -36,14 +47,6 @@ GNU Astronomy Utilities NEWS                          -*- 
outline -*-
           possible to ask for ignoring blank pixels when calculating the
           fraction. This is useful when blank/masked pixels are distributed
           across the image.
-    --blankasforeground: Treat blank elements as foreground (regions above
-          the threshold) in the binary processing steps: initial erosion
-          and opening as well as the filling holes and opening of
-          pseudo-detections. From this version, by default, blank elements
-          in the dataset are considered to be background, so if a
-          foreground pixel is touching it, it will be eroded. This option
-          is irrelevant if the datasets contains no blank elements, but can
-          help remove false positives that are touching blank elements.
 
   Statistics:
     - Sky estimation: new outlier estimation algorithm similar to NoiseChisel.
diff --git a/bin/noisechisel/args.h b/bin/noisechisel/args.h
index b17acaf..e3a46f9 100644
--- a/bin/noisechisel/args.h
+++ b/bin/noisechisel/args.h
@@ -409,6 +409,19 @@ struct argp_option program_options[] =
       GAL_OPTIONS_NOT_SET
     },
     {
+      "holengb",
+      UI_KEY_HOLENGB,
+      "INT",
+      0,
+      "4 or 8 connectivity for filling holes.",
+      UI_GROUP_DETECTION,
+      &p->holengb,
+      GAL_TYPE_SIZE_T,
+      GAL_OPTIONS_RANGE_GT_0,
+      GAL_OPTIONS_MANDATORY,
+      GAL_OPTIONS_NOT_SET
+    },
+    {
       "snminarea",
       UI_KEY_SNMINAREA,
       "INT",
diff --git a/bin/noisechisel/astnoisechisel.conf 
b/bin/noisechisel/astnoisechisel.conf
index a05915b..2362820 100644
--- a/bin/noisechisel/astnoisechisel.conf
+++ b/bin/noisechisel/astnoisechisel.conf
@@ -39,6 +39,7 @@
  minskyfrac            0.7
  sigmaclip           3,0.2
  dthresh               0.0
+ holengb                 8
  snminarea              10
  minnumfalse           100
  snquant              0.99
diff --git a/bin/noisechisel/detection.c b/bin/noisechisel/detection.c
index 4497153..5b897ca 100644
--- a/bin/noisechisel/detection.c
+++ b/bin/noisechisel/detection.c
@@ -283,7 +283,7 @@ detection_fill_holes_open(void *in_prm)
 
       /* Fill the holes in this tile: holes with maximal connectivity means
          that they are most strongly bounded. */
-      gal_binary_holes_fill(copy, copy->ndim, -1);
+      gal_binary_holes_fill(copy, p->holengb==4 ? 1 : 2, -1);
       if(fho_prm->step==1)
         {
           detection_write_in_large(tile, copy);
diff --git a/bin/noisechisel/main.h b/bin/noisechisel/main.h
index bfe485c..9ebb0ed 100644
--- a/bin/noisechisel/main.h
+++ b/bin/noisechisel/main.h
@@ -73,6 +73,7 @@ struct noisechiselparams
   double         sigmaclip[2];  /* Sigma-clipping parameters.             */
   uint8_t         checkdetsky;  /* Check pseudo-detection sky value.      */
   float               dthresh;  /* Sigma threshold for Pseudo-detections. */
+  size_t              holengb;  /* Connectivity for defining a hole.      */
   size_t            snminarea;  /* Minimum pseudo-detection area for S/N. */
   uint8_t             checksn;  /* Save pseudo-detection S/N values.      */
   size_t          minnumfalse;  /* Min No. of det/seg for true quantile.  */
diff --git a/bin/noisechisel/ui.c b/bin/noisechisel/ui.c
index 0fbe803..2ce2c70 100644
--- a/bin/noisechisel/ui.c
+++ b/bin/noisechisel/ui.c
@@ -233,6 +233,9 @@ ui_read_check_only_options(struct noisechiselparams *p)
   if(p->openingngb!=4 && p->openingngb!=8)
     error(EXIT_FAILURE, 0, "%zu not acceptable for `--openingngb'. It must "
           "be 4 or 8 (specifying the type of connectivity)", p->openingngb);
+  if(p->holengb!=4 && p->holengb!=8)
+    error(EXIT_FAILURE, 0, "%zu not acceptable for `--holengb'. It must "
+          "be 4 or 8 (specifying the type of connectivity)", p->holengb);
 
   /* Make sure that the no-erode-quantile is not smaller or equal to
      qthresh. */
diff --git a/bin/noisechisel/ui.h b/bin/noisechisel/ui.h
index 30593bc..3953935 100644
--- a/bin/noisechisel/ui.h
+++ b/bin/noisechisel/ui.h
@@ -92,6 +92,7 @@ enum option_keys_enum
   UI_KEY_OPENINGNGB,
   UI_KEY_SKYFRACNOBLANK,
   UI_KEY_CHECKDETSKY,
+  UI_KEY_HOLENGB,
   UI_KEY_CHECKSN,
   UI_KEY_DETGROWMAXHOLESIZE,
   UI_KEY_CLEANGROWNDET,
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index c38f4cc..f853761 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -16368,6 +16368,16 @@ ratio of the resulting `psudo-detections' are used to 
identify true
 vs. false detections. See Section 3.1.5 and Figure 7 in Akhlaghi and
 Ichikawa (2015) for a very complete explanation.
 
address@hidden --holengb=INT
+The connectivity (defined by the number of neighbors) to fill holes after
+applying @option{--dthresh} (above) to find pseudo-detections. For example
+in a 2D image it must be 4 (the neighbors that are most strongly connected)
+or 8 (all neighbors). The stronger the connectivity, the stronger the hole
+will be enclosed. So setting a value of 8 in a 2D image means that the
+walls of the hole are 4-connected. If standard (near Sky level) values are
+given to @option{--dthresh}, setting @option{--holengb=4}, might fill the
+complete dataset and thus not create enough pseudo-detections.
+
 @item -m INT
 @itemx --snminarea=INT
 The minimum area to calculate the Signal to noise ratio on the
diff --git a/lib/binary.c b/lib/binary.c
index 3e51fc7..556857d 100644
--- a/lib/binary.c
+++ b/lib/binary.c
@@ -725,11 +725,15 @@ gal_binary_holes_fill(gal_data_t *input, int 
connectivity, size_t maxsize)
   size_t numholes, *sizes;
   gal_data_t *inv, *tile, *holelabs=NULL;
 
-  /* A small sanity check. */
+  /* Small sanity checks. */
   if( input->type != GAL_TYPE_UINT8 )
     error(EXIT_FAILURE, 0, "%s: input must have `uint8' type, but its "
           "input dataset has `%s' type", __func__,
           gal_type_name(input->type, 1));
+  if(connectivity<1 || connectivity>input->ndim)
+    error(EXIT_FAILURE, 0, "%s: connectivity value %d is not acceptable. "
+          "It has to be between 1 and the number of input's dimensions "
+          "(%zu)", __func__, connectivity, input->ndim);
 
 
   /* Make the inverse image. */