[gnuastro-commits] master 94acc6b5 2/5: PSF subtract: adding --modelonly to not subtract the PSF

[Mohammad Akhlaghi]

branch: master
commit 94acc6b569467d570d79eec7f1cd3e4331498cec
Author: Raul Infante-Sainz <infantesainz@gmail.com>
Commit: Mohammad Akhlaghi <mohammad@akhlaghi.org>

    PSF subtract: adding --modelonly to not subtract the PSF
    
    Until now, the script in charge of making the subtraction of the PSF was
    not able to provide the modeled PSF (for each star). However, in some
    situations it is good to be able to do it. For example if the modeling of
    the stars within one image can be done independently (because the flux of
    the stars don't affect each other). In this case, it would be worth to be
    able to model all stars in paralel to speed up the process.
    
    With this commit, a new option has been added in order to obtain the
    scattered light field instead the default subtraction of the PSF. It adds
    more flexibility and power when modeling the stars with the PSF.
---
 bin/script/psf-subtract.in | 53 +++++++++++++++++++++++++++++++---------------
 doc/gnuastro.texi          | 46 ++++++++++++++++++++++++++++------------
 2 files changed, 69 insertions(+), 30 deletions(-)

diff --git a/bin/script/psf-subtract.in b/bin/script/psf-subtract.in
index fb6189f8..cb1949a8 100644
--- a/bin/script/psf-subtract.in
+++ b/bin/script/psf-subtract.in
@@ -44,6 +44,7 @@ center=""
 keeptmp=0
 output=""
 tmpdir=""
+modelonly=0
 version=@VERSION@
 scriptname=@SCRIPT_NAME@
 
@@ -93,6 +94,7 @@ $scriptname options:
   -o, --output            Output table with the radial profile.
   -t, --tmpdir            Directory to keep temporary files.
   -k, --keeptmp           Keep temporal/auxiliar files.
+  -m, --modelonly         Give the model as output (don't subtract it).
 
  Operating mode:
   -?, --help              Print this help list.
@@ -202,14 +204,16 @@ do
         -s*)                 scale=$(echo "$1"  | sed -e's/-f//');      
check_v "$1" "$scale";  shift;;
 
         # Output parameters
-        -k|--keeptmp)     keeptmp=1; shift;;
-        -k*|--keeptmp=*)  on_off_option_error --keeptmp -k;;
-        -t|--tmpdir)      tmpdir="$2";                          check_v "$1" 
"$tmpdir";  shift;shift;;
-        -t=*|--tmpdir=*)  tmpdir="${1#*=}";                     check_v "$1" 
"$tmpdir";  shift;;
-        -t*)              tmpdir=$(echo "$1" | sed -e's/-t//'); check_v "$1" 
"$tmpdir";  shift;;
-        -o|--output)      output="$2";                          check_v "$1" 
"$output"; shift;shift;;
-        -o=*|--output=*)  output="${1#*=}";                     check_v "$1" 
"$output"; shift;;
-        -o*)              output=$(echo "$1" | sed -e's/-o//'); check_v "$1" 
"$output"; shift;;
+        -k|--keeptmp)        keeptmp=1; shift;;
+        -k*|--keeptmp=*)     on_off_option_error --keeptmp -k;;
+        -m|--modelonly)      modelonly=1; shift;;
+        -m*|--modelonly=*)   on_off_option_error --modelonly -m;;
+        -t|--tmpdir)         tmpdir="$2";                          check_v 
"$1" "$tmpdir";  shift;shift;;
+        -t=*|--tmpdir=*)     tmpdir="${1#*=}";                     check_v 
"$1" "$tmpdir";  shift;;
+        -t*)                 tmpdir=$(echo "$1" | sed -e's/-t//'); check_v 
"$1" "$tmpdir";  shift;;
+        -o|--output)         output="$2";                          check_v 
"$1" "$output"; shift;shift;;
+        -o=*|--output=*)     output="${1#*=}";                     check_v 
"$1" "$output"; shift;;
+        -o*)                 output=$(echo "$1" | sed -e's/-o//'); check_v 
"$1" "$output"; shift;;
 
         # Non-operating options.
         -q|--quiet)       quiet="--quiet"; shift;;
@@ -329,10 +333,12 @@ objectid="$xcoord"_"$ycoord"
 # an explicit name. If the user has defined a specific path/name for the
 # output, it will be used for saving the output file. If the user does not
 # specify a output name, then a default value containing the center and
-# mode will be generated.
+# mode will be generated. Moreover, if the user specify the option
+# --modelonly, then the output will be the modeled PSF (not subtracted from
+# the original image).
 bname_prefix=$(basename $inputs | sed 's/\.fits/ /' | awk '{print $1}')
 if [ x"$tmpdir" = x ]; then \
-  tmpdir=$(pwd)/"$bname_prefix"_psfmodelscatteredlight
+  tmpdir=$(pwd)/"$bname_prefix"_subtracted
 fi
 
 if [ -d "$tmpdir" ]; then
@@ -341,9 +347,13 @@ else
   mkdir -p "$tmpdir"
 fi
 
-# Output
+# Output, subtract or model only
 if [ x"$output" = x ]; then
-  output="$bname_prefix"_psfmodelscatteredlight_$objectid.fits
+  if [ x"$modelonly" = x1 ]; then
+    output="$bname_prefix"_model_$objectid.fits
+  else
+    output="$bname_prefix"_subtracted_$objectid.fits
+  fi
 fi
 
 
@@ -503,9 +513,12 @@ astcrop $psftranslated --mode=img 
--section=$xrange,$yrange \
 
 
 
-# Put nan values as zero values
-# -----------------------------
+# Final output: PSF-subtracted or PSF-model star
+# ----------------------------------------------
 #
+# Here, the final output is generated. By default, the final result is the
+# PSF-subtracted image. But, if the user uses the option --modelonly, then
+# the output will be the modeled PSF (not subtraction of the PSF is done).
 # In the vast majority of the cases, having zeros in the outer parts of the
 # astronomical images is not a good idea. It is better to have them as nan
 # values. However, for the scattered light field modeling, it is better to
@@ -513,9 +526,15 @@ astcrop $psftranslated --mode=img 
--section=$xrange,$yrange \
 # zero values, it is possible to make the addition of different modeled
 # stars (using the PSF) or the subtraction without having nan values in the
 # final subtracted image.
-astarithmetic $psfcropped --hdu=1 set-i \
-              i i isblank 0 where set-psf \
-              $inputs --hdu=$hdu psf - --output=$output $quiet
+if [ x"$modelonly" = x0 ]; then
+  astarithmetic $psfcropped --hdu=1 set-i \
+                i i isblank 0 where set-psf \
+                $inputs --hdu=$hdu psf - --output=$output $quiet
+else
+  astarithmetic $psfcropped --hdu=1 set-i \
+                i i isblank 0 where --output=$output $quiet
+
+fi
 
 
 
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 98a40a09..ea9a7844 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -6487,6 +6487,16 @@ $ astscript-fits-view label/67510-seg.fits 
subtracted/67510.fits \
 $ rm scattered-light.fits
 @end example
 
+You will probably have noticed that in the scattered light field there are 
some patches that correspond to the saturation of the stars.
+Since we obtained the scattered light field by subtracting PSF-subtracted 
image from the original image, it is natural that we have such saturated 
regions.
+To solve such inconvinience, this script also has an option to not make the 
subtraction of the PSF but to give as output the modeled star.
+For doing that, it is necessary to run the script with the option 
@option{--modelonly}.
+We encourage the reader to obtain such scattered light field model.
+In some scenarios it could be interesting having such way of correcting the 
PSF.
+For example, if there are many faint stars that can be modeled at the same 
time because their flux don't affect each other.
+In such situation, the task could be easily paralelized without having to wait 
to model the brighter stars before the fainter ones.
+At the end, once all stars have been modeled, a simple Arithmetic command 
could be used to sum the different modeled-PSF stamps to obtain the entire 
scattered light field.
+
 In general you see that the subtraction has been done nicely and almost all 
the extended wings of the PSF have been subtracted.
 The central regions of the stars aren't perfectly subtracted:
 
@@ -25508,8 +25518,9 @@ For example, to check that the intermediate images have 
the desired center, they
 
 @node Invoking astscript-psf-subtract,  , Invoking astscript-psf-scale-factor, 
PSF construction and subtraction
 @subsection Invoking astscript-psf-subtract
-This installed script will put the PSF provided into a given position of the 
sky within the input image.
-It is aimed to model the scattered light field of an input image.
+This installed script will put the provided PSF into a given position of the 
sky within the input image, and then it will subtract it.
+It is aimed to model and subtract the scattered light field of an input image.
+It is also possible to obtain the modeled star with the PSF (and not make the 
subtraction of it from the original image).
 
 This script can be used with the following general template:
 
@@ -25521,33 +25532,37 @@ $ astscript-psf-subtract [OPTION...] FITS-file
 Examples:
 
 @example
-## Multiply the the PSF (psf.fits) by 3 and put it into the the
-## pixel position (x,y)=(53,69) of the input image (image.fits).
+## Multiply the the PSF (psf.fits) by 3 and subtract it from the
+## input image (image.fits) at the pixel position (x,y)=(53,69).
 $ astscript-psf-subtract image.fits \
     --psf=psf.fits \
     --mode=img \
-    --fluxfactor=3 \
+    --scale=3 \
     --center=53,69 \
     --output=star-53_69.fits
 
 ## Iterate over a catalog with positions of stars that are
 ## in the the input image. Use WCS coordinates.
 $ asttable catalog.fits | while read -r ra dec mag; do
-         fluxfactor=$(cat fluxfactor-"$ra"_"$dec".txt)
+         scale=$(cat fluxfactor-"$ra"_"$dec".txt)
          astscript-psf-subtract image.fits \
              --mode=wcs \
              --psf=psf.fits \
-             --fluxfactor=$fluxfactor \
+             --scale=$scale \
              --center=$ra,$dec; done
 
 @end example
 
 The input is an image from which the star is considered.
-The result is an image with the same dimensions than the input, but with the 
PSF (@option{--psf}) situated in the position specified by @option{--center}, 
and with the flux level specified by @option{--fluxfactor}.
-The region not covered by the PSF are set to zero values.
+The result is the same image but with the star subtracted (modeled by the PSF).
+The modeling of the star is done with the PSF image specified with the option 
@option{--psf}, and flux-scaled with the option @option{--scale} at the 
position defined by @option{--center}.
+Instead of obtaining the PSF-subtracted image, it is also possible to obtain 
the modeled star by the PSF.
+To to that, use the option @option{--modelonly}.
+With this option, the output will an image with the same size as the original 
one with the PSF situated in the center and flux-scaled.
+In this case, the region not covered by the PSF are set to zero values.
+
 Note that this script works over individual objects.
 As a consequence, to generate a scattered light field of many stars it is 
necessary to make multiple calls and then combine all the individual modeled 
stars.
-
 A full description of each option is given below.
 
 @table @option
@@ -25578,11 +25593,11 @@ The positions along each dimension must be separated 
by a comma (@key{,}).
 The number of values given to this option must be the same as the dimensions 
of the input dataset.
 The units of the coordinates are read based on the value to the 
@option{--mode} option, see above.
 
-@item -f FLT
-@itemx --fluxfactor=FLT
+@item -s FLT
+@itemx --scale=FLT
 Factor by which the PSF (@option{--psf}) is multiplied.
 This factor is necessary to put the PSF with the desired flux level.
-A convenient way of obtaining this value is by using the script 
@file{astscript-model-scale-factor}, see @ref{Invoking 
astscript-psf-scale-factor}.
+A convenient way of obtaining this value is by using the script 
@file{astscript-scale-factor}, see @ref{Invoking astscript-psf-scale-factor}.
 
 @item -t STR
 @itemx --tmpdir=STR
@@ -25599,6 +25614,11 @@ You can disable the deletion of the temporary 
directory with the @option{--keept
 Don't delete the temporary directory (see description of @option{--tmpdir} 
above).
 This option is useful for debugging.
 For example, to check that the intermediate images have the desired center, 
flux level, etc.
+
+@item -m
+@itemx --modelonly
+Don't make the subtraction of the modeled star by the PSF.
+This option is useful when the user wants to obtain the scattered light field 
given by the PSF modeled star.
 @end table