[gnuastro-commits] master 857c66f0: psf-stmap: nan stamp if star is completely outside of field of view

[Mohammad Akhlaghi]

branch: master
commit 857c66f0c7882fd987e877a79903ae070c434298
Author: Sepideh Eskandarlou <sepideh.eskandarlou@gmail.com>
Commit: Mohammad Akhlaghi <mohammad@akhlaghi.org>

    psf-stmap: nan stamp if star is completely outside of field of view
    
    Until now, if a star was completely outside of input image and did not have
    any overlap with the input image, the script would abort and print
    error. However, later, if a normalization value was found to be NaN, it
    would produce a fully NaN-valued image. This inconsistant exit behavior
    makes it hard to use the output in pipelines.
    
    With this commit, a unified output is always created upon aborting: an
    output image is going to be produced, but an image that is fully
    NaN-valued. An informative warning is also printed just before the program
    ends in both cases.
    
    While inspecting Sepideh's commit, Mohammad noticed the following two
    issues that are also implemented in this commit:
    
     - The Fits program can't write NaN-valued FITS keywords (which was why we
       were not including the 'NORMVAL' keyword when the normalization value
       was NaN until now). Therefore, the proper check was added within the
       FITS program so when a given value is NaN, it prints an empty-valued
       keyword (permitted by the FITS standard, and CFITSIO).
    
     - The Library's 'gal_blank_is' function wasn't prepared for checking
       strings that aren't "n/a" (the blank string value)! If the string was
       something else, it would go to the next 'case' element which was for
       complex numbers and produce confusing error messages.
    
    This fixes bug #62794.
---
 NEWS                    |  2 ++
 bin/fits/keywords.c     |  6 ++++--
 bin/script/psf-stamp.in | 46 ++++++++++++++++++++++++++++++----------------
 doc/gnuastro.texi       | 20 ++++++++++++++------
 lib/blank.c             |  8 +++++++-
 5 files changed, 57 insertions(+), 25 deletions(-)

diff --git a/NEWS b/NEWS
index 142a13bf..acd348b9 100644
--- a/NEWS
+++ b/NEWS
@@ -246,6 +246,8 @@ See the end of the file for license conditions.
   bug #62784: Conversion of sexagesimal RA or Dec ignores sign when first
               digit is zero (for example '-00d12m34'). Found by Manuel
               Sánchez-Benavente.
+  bug #62794: psf-stamp scripts produces inconsistant outputs on
+              failure. Found by Sepideh Eskandarlou and Nafise Sedighi.
 
 
 
diff --git a/bin/fits/keywords.c b/bin/fits/keywords.c
index 6f3bca2d..21d00da9 100644
--- a/bin/fits/keywords.c
+++ b/bin/fits/keywords.c
@@ -203,7 +203,8 @@ keywords_write_update(struct fitsparams *p, fitsfile **fptr,
         {
           if(u1w2==1)
             {
-              if(keyll->value)
+              if( keyll->value
+                  && gal_blank_is(keyll->value, keyll->type)==0 )
                 {
                   if( fits_update_key(*fptr,
                                       gal_fits_type_to_datatype(keyll->type),
@@ -220,7 +221,8 @@ keywords_write_update(struct fitsparams *p, fitsfile **fptr,
             }
           else if (u1w2==2)
             {
-              if(keyll->value)
+              if( keyll->value
+                  && gal_blank_is(keyll->value, keyll->type)==0 )
                 {
                   if( fits_write_key(*fptr,
                                      gal_fits_type_to_datatype(keyll->type),
diff --git a/bin/script/psf-stamp.in b/bin/script/psf-stamp.in
index ec038126..0f4ee685 100644
--- a/bin/script/psf-stamp.in
+++ b/bin/script/psf-stamp.in
@@ -246,6 +246,12 @@ EOF
     fi
 }
 
+all_nan_warning() {
+    cat <<EOF
+$scriptname: WARNING: all the pixels in the requested normalization radius are 
NaN. Therefore the output stamp will be all NaN. If you later use this stamp 
(among others) with any stack operator of 'astarithmetic' this image will 
effectively be ignored, so you can safely continue with stacking all the files
+EOF
+}
+
 
 
 
@@ -547,17 +553,33 @@ fi
 # Crop the object around its center with the given stamp size width. It may
 # happen that the given coordinate is fully outside of the image (within
 # the requested stamp-width). In this case Crop won't generate any output,
-# so we are checking the output name here.
+# so we are checking the existance of the '$cropped' file. If not created
+# we will create a fully NaN-valued image that can be ignored in stacks.
+#
+# We are setting '--checkcenter=0' in Crop so it doesn't check if the
+# central pixel is covered or not (we may be interested in the outer parts
+# of the PSF of a star that is centered outside of the image).
 cropped=$tmpdir/cropped-$objectid.fits
-astcrop $inputs --hdu=$hdu --mode=img \
+astcrop $inputs --hdu=$hdu --mode=img --checkcenter=0 \
         --width=$xstampwidth,$ystampwidth \
         --center=$xcenter,$ycenter \
         --output=$cropped $quiet
 if ! [ -f $cropped ]; then
-    cat <<EOF
-$scriptname: ERROR: the coordinate '$center' (in '$mode' mode) is not within 
the image (within $xstampwidth x $ystampwidth pixels from the edges of the 
input image).
-EOF
-    exit 1
+
+    # 'makenew' will generate a fully 0-valued image. Adding (with '+') any
+    # number by NaN will be NaN. Therefore the output of the Arithmetic
+    # command below will be a fully NaN-valued image.
+    astarithmetic $xstampwidth $ystampwidth 2 makenew float32 \
+                  nan + $quiet --output=$output
+
+    # Write an empty 'NORMVAL' keyword (as is expected from this program).
+    astfits $output --write=NORMVAL,nan,"Normalization value"
+
+    # Let the user know what happened.
+    if [ x"$quiet" = x ]; then all_nan_warning; fi
+
+    # Return successfully and don't continue
+    exit 0
 fi
 
 
@@ -708,9 +730,7 @@ astarithmetic $cropped_masked --hdu=1 $normvalue float32 / \
 # and clump that were not masked and the normalization value. If not masks
 # were applied, then set the clump and object labels to 'none'.
 if [ x"$clab" = x ]; then
-    if ! [ x"$normvalue" = xnan ]; then
-        astfits $output --write=NORMVAL,$normvalue,"Normalization value"
-    fi
+    astfits $output --write=NORMVAL,$normvalue,"Normalization value"
 else
     astfits $output --write=NORMVAL,$normvalue,"Normalization value" \
             --write=CLABEL,$clab,"CLUMP label that was not masked" \
@@ -724,13 +744,7 @@ fi
 # Print warning. We are printing this warning here (at the end of the
 # script) because we dont want it to be mixed with the outputs of the
 # previus commands.
-if [ $normvalue = nan ]; then
-    if [ x"$quiet" = x ]; then
-        cat <<EOF
-$scriptname: WARNING: all the pixels in the requested normalization radius are 
NaN. Therefore the output stamp will be fully NaN. When you use this stamp with 
any stack operator of 'astsrithmetic' (for making the PSF) this image will 
effectively be ignored
-EOF
-    fi
-fi
+if [ x"$quiet" = x ] && [ $normvalue = nan ]; then all_nan_warning; fi
 
 
 
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 5776c6b9..d02d6e31 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -26205,13 +26205,21 @@ The output will be an image with a size specified by 
@option{--stampwidth}, cent
 If no center is specified, then it is assumed that the object of interest is 
already in the center of the image.
 In the same way, if no normalization ring is considered, the output will not 
be normalized.
 
-If the given coordinate is not within the image at all (such that the full 
requested stamp is outside the input image), then this script will abort with a 
message and a non-zero return value.
-If all the pixels within the normalization radii are NaN, the output stamp 
will be fully NaN.
-This can happen when the stamp is only partially out of the original image 
field of view, or there is a bad pixel region over the star.
-A fully NaN image can safely be used with the Stacking operators of Arithmetic 
(see @ref{Stacking operators}): in effect it will be ignored.
+In the following cases, this script will produce a fully NaN-valued stamp (of 
the size given to @option{--stampwidth}).
+A fully NaN image can safely be used with the stacking operators of Arithmetic 
(see @ref{Stacking operators}) because they will be ignored.
+In case you don't want to waste storage with fully NaN images, you can 
compress them with @code{gzip --best output.fits}, and give the resulting 
@file{.fits.gz} file to Arithmetic.
+@itemize
+@item
+The requested box (center coordinate with desired width) is not within the 
input image at all.
+@item
+If a normalization radius is requested, and all the pixels within the 
normalization radii are NaN.
+Here are some scenarios that this can happen:
+1) You have a saturated star (where the saturated pixels are NaN), and your 
normalization radius falls within the saturated region.
+2) The star is outside the image by more than your larger normalization radius 
(so there are no pixels for doing normalization), but the full stamp width 
still overlaps part of the image.
+@end itemize
 
-More options are available with the goal of obtaining better stamps.
-A full description of each option is given below.
+@noindent
+The full set of options are listed below for optimal customization in 
different scenarios:
 
 @table @option
 
diff --git a/lib/blank.c b/lib/blank.c
index 187d4af3..8444d10e 100644
--- a/lib/blank.c
+++ b/lib/blank.c
@@ -332,6 +332,11 @@ gal_blank_as_string(uint8_t type, int width)
 int
 gal_blank_is(void *pointer, uint8_t type)
 {
+  /* A small sanity check. */
+  if(pointer==NULL)
+    error(EXIT_FAILURE, 0, "%s: input pointer is NULL", __func__);
+
+  /* Do the checks based on the type. */
   switch(type)
     {
     /* Numeric types */
@@ -347,7 +352,8 @@ gal_blank_is(void *pointer, uint8_t type)
     case GAL_TYPE_FLOAT64:   return isnan( *(double *)(pointer) );
 
     /* String. */
-    case GAL_TYPE_STRING:    if(!strcmp(pointer,GAL_BLANK_STRING)) return 1;
+    case GAL_TYPE_STRING:
+      if(!strcmp(pointer,GAL_BLANK_STRING)) return 1; else return 0;
 
     /* Complex types */
     case GAL_TYPE_COMPLEX32: