[gnuastro-commits] master 6bc0984 1/3: Library (units.h): function to convert counts to Janskys

[Mohammad Akhlaghi]

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

    Library (units.h): function to convert counts to Janskys
    
    Until now there was no easy way to convert counts to Janskys. An AWK
    command was given in the section on Brightness and magnitudes, but it is
    annoying to go and copy-paste it into a command!
    
    With this command, a new function has been added to the units library:
    'gal_units_counts_to_jy'. With this function in place, an operator in the
    Arithmetic library has also been defined for this task and through the
    library, it is possible to do this conversion on images or tables directly
    on the command-line.
---
 NEWS                      |  9 +++++----
 doc/gnuastro.texi         | 37 ++++++++++++++++++++++++++++++++++++-
 lib/arithmetic.c          |  6 ++++++
 lib/gnuastro/arithmetic.h |  9 +++++----
 lib/gnuastro/units.h      |  3 +++
 lib/units.c               | 32 ++++++++++++++++++++++++++++++++
 6 files changed, 87 insertions(+), 9 deletions(-)

diff --git a/NEWS b/NEWS
index a1f0ebc..b373a9d 100644
--- a/NEWS
+++ b/NEWS
@@ -23,10 +23,7 @@ See the end of the file for license conditions.
      the other columns (any MakeCatalog measurement is possible).
 
   Arithmetic:
-   - New operators (the trigonometric/hyperbolic functions were previously
-     only avaialble in Table's column arithmetic, but they have been moved
-     into the Gnuastro library and are thus now available on images within
-     Arithmetic also):
+   - New operators (all also available in Table's column arithmetic):
      - sin: Trigonometric sine (input in degrees).
      - cos: Trigonometric cosine (input in degrees).
      - tan: Trigonometric tangent (input in degrees).
@@ -39,6 +36,8 @@ See the end of the file for license conditions.
      - asinh: Inverse of hyperbolic sine.
      - acosh: Inverse of hyperbolic cosine.
      - atabh: Inverse of hyperbolic tangent.
+     - counts-to-jy: Convert counts to Janskys through a zero point based
+          on AB magnitudes.
 
   Table:
    - When given a value of '_all', the '--noblank' option (that will remove
@@ -73,6 +72,7 @@ See the end of the file for license conditions.
      them for example.
 
   Library:
+   - 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).
@@ -87,6 +87,7 @@ See the end of the file for license conditions.
      - GAL_ARITHMETIC_OP_ASINH: Inverse hyperbolic sine.
      - GAL_ARITHMETIC_OP_ACOSH: Inverse hyperbolic cosine.
      - GAL_ARITHMETIC_OP_ATANH: Inverse hyperbolic tangent.
+     - GAL_ARITHMETIC_OP_COUNTS_TO_JY: Convert counts to Janskys.
 
 ** Removed features
 
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 8a8b412..de3e345 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -11658,6 +11658,21 @@ These operators take a single operand.
 Inverse Hyperbolic sine, cosine, and tangent.
 These operators take a single operand.
 
+@item counts-to-jy
+@cindex AB magnitude
+@cindex Magnitude, AB
+Convert counts (CCD outputs) to Janskys through an AB-magnitude based 
zeropoint.
+The top-popped operand is assumed to be the AB-magnitude zero point and the 
second-popped operand is assumed to be a dataset in units of counts (an image 
in Arithmetic, and a column in Table's @ref{Column arithmetic}).
+For the full equation and basic definitions, see @ref{Brightness flux 
magnitude}.
+
+@cindex SDSS
+For example SDSS images are calibrated in units of nano-maggies, with a fixed 
zero point magnitude of 22.5.
+Therefore you can convert the units of SDSS image pixels to Janskys with the 
command below:
+
+@example
+$ astarithmetic sdss-image.fits 22.5 counts-to-jy
+@end example
+
 @item minvalue
 Minimum value in the first popped operand, so ``@command{a.fits minvalue}'' 
will push the minimum pixel value in this image onto the stack.
 When this operator acts on a single image, the output (operand that is put 
back on the stack) will no longer be an image, but a number.
@@ -18441,6 +18456,9 @@ Using the zero point magnitude (@mymath{Z}), we can 
write the magnitude relation
 @cindex AB magnitude
 @cindex Magnitude, AB
 Having the zero point of an image, you can convert its pixel values to 
physical units of microJanskys (or @mymath{\mu{}Jy}) to enable direct 
pixel-based comparisons with images from other instruments (just note that this 
assumes instrument and observation signatures are corrected, things like the 
flat-field or the Sky).
+Jansky is used for measuring spectral flux density.
+One Jansky is equivalent to @mymath{10^{-26} W/m^2/Hz} (watts per square meter 
per hertz).
+
 This conversion can be done with the fact that in the AB magnitude 
standard@footnote{@url{https://en.wikipedia.org/wiki/AB_magnitude}}, 
@mymath{3631Jy} corresponds to the zero-th magnitude, therefore 
@mymath{B\equiv3631\times10^{6}\mu{Jy}} and @mymath{m\equiv0}.
 We can therefore estimate the brightness (@mymath{B_z}, in @mymath{\mu{Jy}}) 
corresponding to the image zero point (@mymath{Z}) using this equation:
 
@@ -18452,7 +18470,7 @@ We can therefore estimate the brightness (@mymath{B_z}, 
in @mymath{\mu{Jy}}) cor
 Because the image zero point corresponds to a pixel value of @mymath{1}, the 
@mymath{B_z} value calculated above also corresponds to a pixel value of 
@mymath{1}.
 Therefore you simply have to multiply your image by @mymath{B_z} to convert it 
to @mymath{\mu{Jy}}.
 Don't forget that this only applies when your zero point was also estimated in 
the AB magnitude system.
-On the command-line, you can easily estimate this value for a certain zero 
point with AWK, then multiply it to all the pixels in the image with 
@ref{Arithmetic}.
+On the command-line, you can estimate this value for a certain zero point with 
AWK, then multiply it to all the pixels in the image with @ref{Arithmetic}.
 For example let's assume you are using an SDSS image with a zero point of 22.5:
 
 @example
@@ -18460,6 +18478,15 @@ bz=$(echo 22.5 | awk '@{print 3631 * 10^(6-$1/2.5)@}')
 astarithmetic sdss.fits $bz x --output=sdss-in-muJy.fits
 @end example
 
+@noindent
+But in Gnuastro, it gets even easier: Arithmetic has an operator called 
@code{counts-to-jy}.
+This will directly convert your image pixels (in units of counts) to Janskys 
though a provided AB Magnitude-based zero point like below.
+See @ref{Arithmetic operators} for more.
+
+@example
+$ astarithmetic sdss.fits 22.5 counts-to-jy
+@end example
+
 @cindex Steradian
 @cindex Angular coverage
 @cindex Celestial sphere
@@ -28326,6 +28353,14 @@ Convert the input Declination (Dec) degree (a single 
floating point number) to o
 If @code{usecolon!=0}, then the delimiters between the components will be 
colons: @code{_:_:_}.
 @end deftypefun
 
+@deftypefun double gal_units_counts_to_jy (double @code{counts}, double 
@code{zeropoint_ab})
+@cindex Jansky (Jy)
+@cindex AB Magnitude
+@cindex Magnitude, AB
+Convert counts to Janskys through an AB magnitude-based zero point.
+For more on the equation, see @ref{Brightness flux magnitude}.
+@end deftypefun
+
 @node Spectral lines library, Cosmology library, Unit conversion library 
(@file{units.h}), Gnuastro library
 @subsection Spectral lines library (@file{speclines.h})
 
diff --git a/lib/arithmetic.c b/lib/arithmetic.c
index e7f68d9..aed186c 100644
--- a/lib/arithmetic.c
+++ b/lib/arithmetic.c
@@ -1775,6 +1775,8 @@ arithmetic_function_binary_flt(int operator, int flags, 
gal_data_t *il,
       BINFUNC_F_OPERATOR_SET( pow,   +0 );         break;
     case GAL_ARITHMETIC_OP_ATAN2:
       BINFUNC_F_OPERATOR_SET( atan2, *180.0f/pi ); break;
+    case GAL_ARITHMETIC_OP_COUNTS_TO_JY:
+      BINFUNC_F_OPERATOR_SET( gal_units_counts_to_jy, +0 ); break;
     default:
       error(EXIT_FAILURE, 0, "%s: operator code %d not recognized",
             __func__, operator);
@@ -1958,6 +1960,8 @@ gal_arithmetic_set_operator(char *string, size_t 
*num_operands)
     { op=GAL_ARITHMETIC_OP_DEGREE_TO_RA;      *num_operands=1;  }
   else if (!strcmp(string, "degree-to-dec"))
     { op=GAL_ARITHMETIC_OP_DEGREE_TO_DEC;     *num_operands=1;  }
+  else if (!strcmp(string, "counts-to-jy"))
+    { op=GAL_ARITHMETIC_OP_COUNTS_TO_JY;      *num_operands=2;  }
 
   /* Statistical/higher-level operators. */
   else if (!strcmp(string, "minvalue"))
@@ -2131,6 +2135,7 @@ gal_arithmetic_operator_string(int operator)
     case GAL_ARITHMETIC_OP_DEC_TO_DEGREE:   return "dec-to-degree";
     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_JY:    return "counts-to-jy";
 
     case GAL_ARITHMETIC_OP_MINVAL:          return "minvalue";
     case GAL_ARITHMETIC_OP_MAXVAL:          return "maxvalue";
@@ -2253,6 +2258,7 @@ gal_arithmetic(int operator, size_t numthreads, int 
flags, ...)
     /* Binary function operators. */
     case GAL_ARITHMETIC_OP_POW:
     case GAL_ARITHMETIC_OP_ATAN2:
+    case GAL_ARITHMETIC_OP_COUNTS_TO_JY:
       d1 = va_arg(va, gal_data_t *);
       d2 = va_arg(va, gal_data_t *);
       out=arithmetic_function_binary_flt(operator, flags, d1, d2);
diff --git a/lib/gnuastro/arithmetic.h b/lib/gnuastro/arithmetic.h
index b1719df..6ad36b0 100644
--- a/lib/gnuastro/arithmetic.h
+++ b/lib/gnuastro/arithmetic.h
@@ -120,10 +120,11 @@ enum gal_arithmetic_operators
   GAL_ARITHMETIC_OP_ACOSH,        /* Inverse hyperbolic cosine. */
   GAL_ARITHMETIC_OP_ATANH,        /* Inverse hyperbolic tangent. */
 
-  GAL_ARITHMETIC_OP_RA_TO_DEGREE, /* right ascension to decimal      */
-  GAL_ARITHMETIC_OP_DEC_TO_DEGREE,/* declination to decimal          */
-  GAL_ARITHMETIC_OP_DEGREE_TO_RA, /* right ascension to decimal      */
-  GAL_ARITHMETIC_OP_DEGREE_TO_DEC,/* declination to decimal          */
+  GAL_ARITHMETIC_OP_RA_TO_DEGREE, /* right ascension to decimal. */
+  GAL_ARITHMETIC_OP_DEC_TO_DEGREE,/* declination to decimal. */
+  GAL_ARITHMETIC_OP_DEGREE_TO_RA, /* right ascension to decimal. */
+  GAL_ARITHMETIC_OP_DEGREE_TO_DEC,/* declination to decimal. */
+  GAL_ARITHMETIC_OP_COUNTS_TO_JY, /* Counts to Janskys with AB-mag zeropoint. 
*/
 
   GAL_ARITHMETIC_OP_MINVAL,       /* Minimum value of array.               */
   GAL_ARITHMETIC_OP_MAXVAL,       /* Maximum value of array.               */
diff --git a/lib/gnuastro/units.h b/lib/gnuastro/units.h
index 0543e67..52ee3b8 100644
--- a/lib/gnuastro/units.h
+++ b/lib/gnuastro/units.h
@@ -73,6 +73,9 @@ gal_units_degree_to_ra (double decimal, int usecolon);
 char *
 gal_units_degree_to_dec (double decimal, int usecolon);
 
+double
+gal_units_counts_to_jy(double counts, double zeropoint_ab);
+
 __END_C_DECLS    /* From C++ preparations */
 
 #endif           /* __GAL_UNITS_H__ */
diff --git a/lib/units.c b/lib/units.c
index 20146e2..b20cf65 100644
--- a/lib/units.c
+++ b/lib/units.c
@@ -331,3 +331,35 @@ gal_units_degree_to_dec(double decimal, int usecolon)
   /* Return the final string. */
   return dec;
 }
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+/**********************************************************************/
+/****************          Flux conversions           *****************/
+/**********************************************************************/
+
+/* Convert Pixel values to Janskys with an AB-magnitude based
+   zero-point. See the "Brightness, Flux, Magnitude and Surface
+   brightness". */
+double
+gal_units_counts_to_jy(double counts, double zeropoint_ab)
+{
+  return counts * 3631 * pow(10, -1 * zeropoint_ab / 2.5);
+}