[gnuastro-commits] master be3b565: Library (arithmetic.h): mag-to-counts is a new operator

[Mohammad Akhlaghi]

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

    Library (arithmetic.h): mag-to-counts is a new operator
    
    Until now we had the 'counts-to-mag' operator in the arithmetic library
    (and thus the Arithmetic program or column arithmetic in Table). But in
    some scenarios (like adding Poisson noise to a mock image with a background
    value in magnitudes), we need the inverse operator.
    
    With this commit, the inverse operator has been added for this job and
    besides its own documentation and example, an example of using it has been
    added to the 'mknoise-poisson' operator.
---
 NEWS                      |  6 +++---
 doc/gnuastro.texi         | 31 +++++++++++++++++++++++++++----
 lib/arithmetic.c          |  6 ++++++
 lib/gnuastro/arithmetic.h |  1 +
 lib/gnuastro/units.h      |  3 +++
 lib/units.c               | 10 ++++++++++
 6 files changed, 50 insertions(+), 7 deletions(-)

diff --git a/NEWS b/NEWS
index 648cb3f..7e1580f 100644
--- a/NEWS
+++ b/NEWS
@@ -67,6 +67,7 @@ See the end of the file for license conditions.
      - mknoise-poisson: Add Poisson noise with the given background.
      - mknoise-uniform: Add uniform noise around existing value.
      - counts-to-mag: Convert counts to magnitudes with given zero point.
+     - mag-to-counts: Convert magnitudes to counts with given zero point.
      - counts-to-jy: Convert counts to Janskys through a zero point based
           on AB magnitudes.
      Trigonometric operators proposed by Zahra Sharbaf and Samane Raji.
@@ -178,9 +179,6 @@ See the end of the file for license conditions.
      the help of Joseph Mazzarella.
 
   Library:
-   - gal_ds9_reg_read_polygon: Parse the polygon from an SAO DS9 region file.
-   - gal_units_counts_to_mag: Convert counts to magnitudes.
-   - gal_units_counts_to_jy: Convert counts to Janskys.
    - New arithmetic operator macros (for the 'gal_arithmetic' function):
      - GAL_ARITHMETIC_OP_SIN: sine (input in deg).
      - GAL_ARITHMETIC_OP_COS: cosine (input in deg).
@@ -209,8 +207,10 @@ See the end of the file for license conditions.
      - GAL_WCS_COORDSYS_GALACTIC: Galactic coordinates.
      - GAL_WCS_COORDSYS_SUPERGALACTIC: Supergalactic coordinates.
    - gal_array_file_recognized: For FITS images, check contents also.
+   - gal_ds9_reg_read_polygon: Parse the polygon from an SAO DS9 region file.
    - gal_fits_file_recognized: Check the file contents when suffix isn't FITS.
    - gal_units_counts_to_mag: Convert counts to magnitudes.
+   - gal_units_mag_to_counts: Convert magnitudes to counts.
    - gal_units_counts_to_jy: Convert counts to Janskys.
    - gal_wcs_coordsys_from_string: WCS coordinate system from string.
    - gal_wcs_coordsys_identify: Parse WCS struct to find  coordinate system.
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 13fee77..c41e903 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -12225,15 +12225,25 @@ These operators take a single operand.
 Convert counts (usually CCD outputs) to magnitudes using the given zeropoint.
 The zero point is the first popped operand and the count value is the second.
 For example assume you have measured the standard deviation of the noise in an 
image to be @code{0.1}, and the image's zero point is @code{22.5}.
-You can therefore measure the @emph{per-pixel} surface brightness limit of the 
dataset (which is the magnitude of the noise standrard deviation) with this 
simple command below.
-Note that because the output is a simple number, we are using @option{--quiet} 
to avoid printing extra information.
+You can therefore measure the @emph{per-pixel} surface brightness limit of the 
dataset (which is the magnitude of the noise standrard deviation) with this 
simple command below (for more on surface brightness see @ref{Brightness flux 
magnitude}).
+Note that because the output is a single number, we are using @option{--quiet} 
to avoid printing extra information.
 
 @example
-astarithmetic 0.1 22.5 counts-to-mag --quiet
+$ astarithmetic 0.1 22.5 counts-to-mag --quiet
 @end example
 
 Of course, you can also convert every pixel in an image (or table column in 
Table's @ref{Column arithmetic}) with this operator if you replace the second 
popped operand with an image/column.
 
+@item mag-to-counts
+Convert magnitudes to counts (usually CCD outputs) using the given zeropoint.
+The zero point is the first popped operand and the magnitude value is the 
second.
+For example if an object has a magnitude of 20, you can estimate the counts 
corresponding to it (when the image has a zeropoint of 24.8) with this command:
+Note that because the output is a single number, we are using @option{--quiet} 
to avoid printing extra information.
+
+@example
+$ astarithmetic 20 24.8 mag-to-counts --quiet
+@end example
+
 @item counts-to-jy
 @cindex AB magnitude
 @cindex Magnitude, AB
@@ -12964,12 +12974,25 @@ But this feature has a side-effect: that if the order 
of calling the @code{mknoi
 If you need this feature please send us an email at 
@code{bug-gnuastro@@gnu.org} (to motivate us in its implementation).}.
 
 @item mknoise-poisson
+@cindex Poisson noise
 Add Poisson noise to each element of the input dataset (see @ref{Photon 
counting noise}).
-This operator takes two arguments: the top/first popped operand is the 
background, the next popped operand is the dataset that the noise should be 
added to.
+This operator takes two arguments: the top/first popped operand is the 
background value (in units of electron counts), the next popped operand is the 
dataset that the noise should be added to.
 
 Except for the noise-model, this operator is very similar to 
@code{mknoise-sigma} and the examples there apply here too.
 The main difference with @code{mknoise-sigma} is that in a Poisson 
distribution the scatter/sigma will depend on each element's value.
 
+For example, let's assume you have made a mock image called @file{mock.fits} 
with @ref{MakeProfiles} and its assumed zeropoint is 22.5 (for more on the zero 
point, see @ref{Brightness flux magnitude}).
+Let's assume the background level for the Poisson noise has a value of 19 
magnitudes.
+You can first use the @code{mag-to-counts} operator to convert this background 
magnitude into counts, then feed the background value in counts to 
@code{mknoise-poisson} operator:
+
+@example
+$ astarithmetic mock.fits 19 22.5 mag-to-counts \
+                mknoise-poisson
+@end example
+
+Try changing the background value from 19 to 10 to see the effect!
+Recall that the tutorial @ref{Sufi simulates a detection} shows how you can 
use MakeProfiles to build mock images.
+
 @item mknoise-uniform
 Add uniform noise to each element of the input dataset.
 This operator takes two arguments: the top/first popped operand is the width 
of the interval, the second popped operand is the dataset that the noise should 
be added to (each element will be the center of the interval).
diff --git a/lib/arithmetic.c b/lib/arithmetic.c
index f68bb83..4be99df 100644
--- a/lib/arithmetic.c
+++ b/lib/arithmetic.c
@@ -1908,6 +1908,8 @@ arithmetic_function_binary_flt(int operator, int flags, 
gal_data_t *il,
       BINFUNC_F_OPERATOR_SET( atan2, *180.0f/pi ); break;
     case GAL_ARITHMETIC_OP_COUNTS_TO_MAG:
       BINFUNC_F_OPERATOR_SET( gal_units_counts_to_mag, +0 ); break;
+    case GAL_ARITHMETIC_OP_MAG_TO_COUNTS:
+      BINFUNC_F_OPERATOR_SET( gal_units_mag_to_counts, +0 ); break;
     case GAL_ARITHMETIC_OP_COUNTS_TO_JY:
       BINFUNC_F_OPERATOR_SET( gal_units_counts_to_jy, +0 ); break;
     default:
@@ -2095,6 +2097,8 @@ gal_arithmetic_set_operator(char *string, size_t 
*num_operands)
     { op=GAL_ARITHMETIC_OP_DEGREE_TO_DEC;     *num_operands=1;  }
   else if (!strcmp(string, "counts-to-mag"))
     { op=GAL_ARITHMETIC_OP_COUNTS_TO_MAG;     *num_operands=2;  }
+  else if (!strcmp(string, "mag-to-counts"))
+    { op=GAL_ARITHMETIC_OP_MAG_TO_COUNTS;     *num_operands=2;  }
   else if (!strcmp(string, "counts-to-jy"))
     { op=GAL_ARITHMETIC_OP_COUNTS_TO_JY;      *num_operands=2;  }
 
@@ -2279,6 +2283,7 @@ gal_arithmetic_operator_string(int operator)
     case GAL_ARITHMETIC_OP_DEGREE_TO_RA:    return "degree-to-ra";
     case GAL_ARITHMETIC_OP_DEGREE_TO_DEC:   return "degree-to-dec";
     case GAL_ARITHMETIC_OP_COUNTS_TO_MAG:   return "counts-to-mag";
+    case GAL_ARITHMETIC_OP_MAG_TO_COUNTS:   return "mag-to-counts";
     case GAL_ARITHMETIC_OP_COUNTS_TO_JY:    return "counts-to-jy";
 
     case GAL_ARITHMETIC_OP_MINVAL:          return "minvalue";
@@ -2407,6 +2412,7 @@ gal_arithmetic(int operator, size_t numthreads, int 
flags, ...)
     case GAL_ARITHMETIC_OP_POW:
     case GAL_ARITHMETIC_OP_ATAN2:
     case GAL_ARITHMETIC_OP_COUNTS_TO_MAG:
+    case GAL_ARITHMETIC_OP_MAG_TO_COUNTS:
     case GAL_ARITHMETIC_OP_COUNTS_TO_JY:
       d1 = va_arg(va, gal_data_t *);
       d2 = va_arg(va, gal_data_t *);
diff --git a/lib/gnuastro/arithmetic.h b/lib/gnuastro/arithmetic.h
index f6c2712..560601c 100644
--- a/lib/gnuastro/arithmetic.h
+++ b/lib/gnuastro/arithmetic.h
@@ -127,6 +127,7 @@ enum gal_arithmetic_operators
   GAL_ARITHMETIC_OP_DEGREE_TO_RA, /* right ascension to decimal. */
   GAL_ARITHMETIC_OP_DEGREE_TO_DEC,/* declination to decimal. */
   GAL_ARITHMETIC_OP_COUNTS_TO_MAG,/* Counts to magnitude. */
+  GAL_ARITHMETIC_OP_MAG_TO_COUNTS,/* Magnitude to counts. */
   GAL_ARITHMETIC_OP_COUNTS_TO_JY, /* Counts to Janskys with AB-mag zeropoint. 
*/
 
   GAL_ARITHMETIC_OP_MINVAL,       /* Minimum value of array.               */
diff --git a/lib/gnuastro/units.h b/lib/gnuastro/units.h
index 3243f1f..229c610 100644
--- a/lib/gnuastro/units.h
+++ b/lib/gnuastro/units.h
@@ -77,6 +77,9 @@ double
 gal_units_counts_to_mag(double counts, double zeropoint);
 
 double
+gal_units_mag_to_counts(double mag, double zeropoint);
+
+double
 gal_units_counts_to_jy(double counts, double zeropoint_ab);
 
 __END_C_DECLS    /* From C++ preparations */
diff --git a/lib/units.c b/lib/units.c
index ce9e873..d93dcbf 100644
--- a/lib/units.c
+++ b/lib/units.c
@@ -368,6 +368,16 @@ gal_units_counts_to_mag(double counts, double zeropoint)
 
 
 
+/* Convert magnitude to counts using the given zeropoint. */
+double
+gal_units_mag_to_counts(double mag, double zeropoint)
+{
+  return pow(10, (mag - zeropoint)/(-2.5f));
+}
+
+
+
+
 /* Convert Pixel values to Janskys with an AB-magnitude based
    zero-point. See the "Brightness, Flux, Magnitude and Surface
    brightness". */