[gnuastro-commits] master 03cd54d: Query: added NED's Galactic extinction calculator

[Mohammad Akhlaghi]

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

    Query: added NED's Galactic extinction calculator
    
    Until now, only one of NED's very rich dataset was available within
    Query. This was because only that one table followed the general Table
    Access Protocol (TAP) of the Virtual Observatory standard.
    
    With this commit, NED's galactic extinction calculator for a certain
    coordinate is now also available. For this feature, we use NED's custom CGI
    interface, not the standard TAP. The format of the CGI URL was kindly
    provided by Joseph Mazzarella from NED.
    
    The only issue is that NED provides a VOTable which Gnuastro doesn't yet
    support. So the necessary shell command to convert the VOTable data into a
    FITS table is also given in the manual.
    
    This commit was done with the help of Joseph Mazzarella from NED.
    
    This completes task #15911.
---
 NEWS               |   7 ++++
 THANKS             |   1 +
 bin/query/astron.c |   5 ++-
 bin/query/gaia.c   |   3 ++
 bin/query/main.h   |   1 +
 bin/query/ned.c    | 119 ++++++++++++++++++++++++++++++++++++++++++++++-------
 bin/query/query.c  |  48 ++++++++++++---------
 bin/query/ui.c     |  67 +++++++++++++++++++++++-------
 bin/query/vizier.c |   3 ++
 doc/gnuastro.texi  |  36 ++++++++++++++--
 10 files changed, 238 insertions(+), 52 deletions(-)

diff --git a/NEWS b/NEWS
index b311d55..6d05a81 100644
--- a/NEWS
+++ b/NEWS
@@ -71,6 +71,13 @@ See the end of the file for license conditions.
      program for selecting a certain sub-set of your FITS files, or sorting
      them for example.
 
+  Query:
+   - NED's Galactic extinction calculator is now available for any
+     coordinate with a command like below. See the description of the
+     'extinction' dataset of NED in the "Available datasets" section of the
+     manual.
+         astquery ned --dataset=extinction --center=49.9507,41.5116
+
   Library:
    - gal_units_counts_to_jy: Convert counts to Janskys.
    - New arithmetic operator macros (for the 'gal_arithmetic' function):
diff --git a/THANKS b/THANKS
index 6d36246..876200e 100644
--- a/THANKS
+++ b/THANKS
@@ -64,6 +64,7 @@ support in Gnuastro. The list is ordered alphabetically (by 
family name).
     Sebastián Luna Valero                sluna@iaa.es
     Alberto Madrigal                     brt.madrigal@gmail.com
     Guillaume Mahler                     guillaume.mahler@univ-lyon1.fr
+    Joseph Mazzarella                    mazz@ipac.caltech.edu
     Alireza Molaeinezhad                 amolaei@gmail.com
     Javier Moldon                        jmoldon@iaa.es
     Juan Molina Tobar                    juan.a.molina.t@gmail.com
diff --git a/bin/query/astron.c b/bin/query/astron.c
index c3f41ad..01dd25d 100644
--- a/bin/query/astron.c
+++ b/bin/query/astron.c
@@ -50,6 +50,9 @@ astron_sanity_checks(struct queryparams *p)
           gal_checkset_allocate_copy("tgssadr.main", &p->datasetstr);
         }
     }
+
+  /* Currently we assume ASTRON only uses TAP. */
+  p->usetap=1;
 }
 
 
@@ -59,7 +62,7 @@ astron_sanity_checks(struct queryparams *p)
 void
 astron_prepare(struct queryparams *p)
 {
-  /* NED-specific. */
+  /* ASTRON-specific. */
   astron_sanity_checks(p);
 
   /* Set the URLs, note that this is a simply-linked list, so we need to
diff --git a/bin/query/gaia.c b/bin/query/gaia.c
index 135e3ca..a7dd61f 100644
--- a/bin/query/gaia.c
+++ b/bin/query/gaia.c
@@ -95,6 +95,9 @@ gaia_sanity_checks(struct queryparams *p)
           gal_checkset_allocate_copy("public.tycho2", &p->datasetstr);
         }
     }
+
+  /* Currently we assume GAIA only uses TAP. */
+  p->usetap=1;
 }
 
 
diff --git a/bin/query/main.h b/bin/query/main.h
index 377d41c..62c70f3 100644
--- a/bin/query/main.h
+++ b/bin/query/main.h
@@ -70,6 +70,7 @@ struct queryparams
   char                *ra_name;  /* Name of RA column.                 */
   char               *dec_name;  /* Name of Dec columns.               */
   char           *finalcommand;  /* The final command used.            */
+  int                   usetap;  /* If a TAP-download should be used.  */
 
   /* Output: */
   time_t               rawtime;  /* Starting time of the program.      */
diff --git a/bin/query/ned.c b/bin/query/ned.c
index f2fd24f..80a9bae 100644
--- a/bin/query/ned.c
+++ b/bin/query/ned.c
@@ -37,7 +37,7 @@ along with Gnuastro. If not, see 
<http://www.gnu.org/licenses/>.
 
 
 
-
+/* Basic sanity checks. */
 static void
 ned_sanity_checks(struct queryparams *p)
 {
@@ -51,9 +51,35 @@ ned_sanity_checks(struct queryparams *p)
         }
     }
 
+  /* Database-specific checks. For example, if we should use TAP or
+     not. Note that the user may give 'NEDTAP.objdir', so we can't use the
+     'if' above (for expanding summarized names). */
+  if( !strcmp(p->datasetstr, "NEDTAP.objdir") )
+    p->usetap=1;
+  else if( !strcmp(p->datasetstr, "extinction") )
+    {
+      /* Crash for options that are not compatible with extinction. */
+      if( p->radius || p->width || p->range || p->noblank || p->columns
+          || p->head!=GAL_BLANK_SIZE_T || p->sort )
+        error(EXIT_FAILURE, 0, "NED's extinction calculator returns "
+              "the galactic extinction for a single point (in multiple "
+              "filters), therefore the following options are not "
+              "acceptable with it: '--radius', '--width', '--range', "
+              "'--noblank', '--column', '--head' and '--sort'");
+
+      /* Make sure that '--center' is given. */
+      if(p->center==NULL)
+        error(EXIT_FAILURE, 0, "no coordinate specified! Please use "
+              "'--center' to specify the RA and Dec (in J2000) of "
+              "your desired coordinate, for example "
+              "--center=10.68458,41.269166");
+    }
+
   /* Currently NED only has a single table for TAP access, so warn the
      users about this if they ask for any other table. */
-  if( p->datasetstr==NULL || strcmp(p->datasetstr, "NEDTAP.objdir") )
+  if( p->usetap
+      && ( p->datasetstr==NULL
+           || strcmp(p->datasetstr, "NEDTAP.objdir") ) )
     error(EXIT_FAILURE, 0, "NED currently only supports a single "
           "dataset with the TAP protocol called 'NEDTAP.objdir' "
           "(which you can also call in Query with '--dataset=objdir'). "
@@ -67,22 +93,87 @@ ned_sanity_checks(struct queryparams *p)
 
 
 
+/* Extinction with NED */
+void
+ned_extinction(struct queryparams *p)
+{
+  double *darr;
+  char *command;
+
+  /* If the user wants information, we'll specify a (0,0) center coordinate
+     and continue. In the end, instead of saving the file, we'll just
+     report the metadata. */
+  if(p->information)
+    error(EXIT_FAILURE, 0, "'--information' is not yet supported for "
+          "NED's extinction calculator");
+
+  /* Build the calling command. Note that the query quotes are
+     included by the function building it. */
+  darr=p->center->array;
+  if( asprintf(&command, "curl%s -o%s 
'https://ned.ipac.caltech.edu/cgi-bin/calc?in_csys=Equatorial&out_csys=Equatorial&in_equinox=J2000.0&out_equinox=J2000.0&obs_epoch=2000.0&lon=%fd&lat=%fd&of=xml_main&ext=1'",
 p->cp.quiet ? " -s" : "",
+               p->downloadname, darr[0], darr[1])<0 )
+    error(EXIT_FAILURE, 0, "%s: asprintf allocation ('command')",
+          __func__);
+
+  /* Print the calling command for the user to know. */
+  if(p->dryrun==1 || p->cp.quiet==0)
+    {
+      if(p->dryrun==0) printf("\n");
+      error(EXIT_SUCCESS, 0, "%s: %s",
+            p->dryrun ? "would run" : "running", command);
+      if(p->dryrun==0) printf("\nDownload status:\n");
+    }
+
+  /* Run the command if '--dryrun' isn't called: if the command succeeds
+     'system' returns 'EXIT_SUCCESS'. */
+  if(p->dryrun==0)
+    {
+      if(system(command)!=EXIT_SUCCESS)
+        error(EXIT_FAILURE, 0, "the query download command %sfailed%s\n",
+              p->cp.quiet==0 ? "printed above " : "",
+              p->cp.quiet==0 ? "" : " (the command can be printed "
+              "if you don't use the option '--quiet', or '-q')");
+    }
+}
+
+
+
+
+
+/* For NED's non-TAP queries. */
+void
+ned_non_tap(struct queryparams *p)
+{
+  if( !strcmp(p->datasetstr, "extinction") )
+    ned_extinction(p);
+}
+
+
+
+
+
 void
 ned_prepare(struct queryparams *p)
 {
   /* NED-specific. */
   ned_sanity_checks(p);
 
-  /* Set the URLs, note that this is a simply-linked list, so we need to
-     reverse it in the end (with 'gal_list_str_reverse') to have the same
-     order here. */
-  gal_list_str_add(&p->urls,
-                   "https://ned.ipac.caltech.edu/tap/sync";, 0);
-
-  /* Name of default RA Dec columns. */
-  if(p->ra_name==NULL)  p->ra_name="ra";
-  if(p->dec_name==NULL) p->dec_name="dec";
-
-  /* Basic sanity checks. */
-  tap_sanity_checks(p);
+  /* If we should use TAP, do the preparations. */
+  if(p->usetap)
+    {
+      /* Set the URLs, note that this is a simply-linked list, so we need
+         to reverse it in the end (with 'gal_list_str_reverse') to have the
+         same order here. */
+      gal_list_str_add(&p->urls,
+                       "https://ned.ipac.caltech.edu/tap/sync";, 0);
+
+      /* Name of default RA Dec columns. */
+      if(p->ra_name==NULL)  p->ra_name="ra";
+      if(p->dec_name==NULL) p->dec_name="dec";
+
+      /* Basic sanity checks. */
+      tap_sanity_checks(p);
+    }
+  else
+    ned_non_tap(p);
 }
diff --git a/bin/query/query.c b/bin/query/query.c
index b8289b5..cae127d 100644
--- a/bin/query/query.c
+++ b/bin/query/query.c
@@ -284,36 +284,49 @@ query_output_data(struct queryparams *p)
 
 
 void
-query_check_download(struct queryparams *p)
+query_output_finalize(struct queryparams *p)
 {
   size_t len;
-  int status=0;
+  int isxml=0;
   char *logname;
   fitsfile *fptr;
+  int gooddownload=0, status=0;
 
-  /* Open the FITS file and if the status value is still zero, it means
-     everything worked properly. */
-  fits_open_file(&fptr, p->downloadname, READONLY, &status);
-  if(status==0)
+  /* See if its a FITS file or a VOTable. */
+  len=strlen(p->downloadname);
+  if( !strcmp(&p->downloadname[len-4], ".xml") )
+    { isxml=1; gooddownload=1; }
+  else
     {
-      /* Close the FITS file pointer. */
-      fits_close_file(fptr, &status);
+      /* Open the FITS file and if the status value is still zero, it means
+         everything worked properly. */
+      fits_open_file(&fptr, p->downloadname, READONLY, &status);
+      if(status==0)
+        {
+          gooddownload=1;
+          fits_close_file(fptr, &status);
+        }
+    }
 
+  /* If the downloaded file is good, do the preparations. */
+  if(gooddownload)
+    {
       /* Prepare the output dataset. */
       if(p->information)
         {
           if(p->datasetstr)  query_output_meta_dataset(p);
           else               query_output_meta_database(p);
         }
-      else               query_output_data(p);
+      else if(isxml==0)      query_output_data(p);
 
       /* Delete the raw downloaded file if necessary. */
       if(p->keeprawdownload==0) remove(p->downloadname);
     }
+
+  /* If there was an error */
   else
     {
       /* Add a '.log' suffix to the output filename. */
-      len=strlen(p->downloadname);
       logname=gal_pointer_allocate(GAL_TYPE_STRING, len+10, 1,
                                    __func__, "logname");
       sprintf(logname, "%s.log", p->downloadname);
@@ -326,8 +339,8 @@ query_check_download(struct queryparams *p)
             "retrieved! For more, please see '%s'", logname);
     }
 
-  /* Add the query keywords to the first extension (if the output was a
-     FITS file). */
+  /* Add the query keywords to the first extension of the output (if the
+     output was a FITS file). */
   if( p->information==0 && gal_fits_name_is_fits(p->cp.output) )
     {
       gal_fits_key_list_title_add_end(&p->cp.okeys,
@@ -360,12 +373,12 @@ query(struct queryparams *p)
     }
 
   /* Download the requested query. */
-  tap_download(p);
+  if(p->usetap) tap_download(p);
 
   /* Make sure that the result is a readable FITS file, otherwise, abort
      with an error. */
   if(p->dryrun==0)
-    query_check_download(p);
+    query_output_finalize(p);
 
   /* Let the user know that things went well. */
   if(p->dryrun==0 && p->cp.quiet==0)
@@ -379,12 +392,7 @@ query(struct queryparams *p)
           printf("Query's raw downloaded file: %s\n", p->downloadname);
         }
       if(p->information==0)
-        {
-          printf("Query's final output: %s\n", p->cp.output);
-          printf("TIP: use the command below for more on the "
-                 "downloaded table:\n"
-                 "   asttable %s --info\n", p->cp.output);
-        }
+        printf("Query's output written: %s\n", p->cp.output);
     }
 
   /* Clean up. */
diff --git a/bin/query/ui.c b/bin/query/ui.c
index ed37313..204fa8e 100644
--- a/bin/query/ui.c
+++ b/bin/query/ui.c
@@ -266,9 +266,9 @@ ui_read_check_only_options(struct queryparams *p)
 {
   size_t i;
   double *darray;
-  char *basename;
   gal_data_t *tmp;
   int keepinputdir;
+  char *suffix, *rdsuffix, *basename;
 
   /* See if database has been specified. */
   if(p->databasestr==NULL)
@@ -405,26 +405,65 @@ ui_read_check_only_options(struct queryparams *p)
   gal_checkset_writable_remove(p->cp.output, p->cp.keep,
                                p->cp.dontdelete);
 
+  /* Set the suffix of the default names. */
+  if( p->database==QUERY_DATABASE_NED
+      && !strcmp(p->datasetstr, "extinction") )
+    {
+      suffix=".xml";
+      rdsuffix="-raw-download.xml";
+    }
+  else
+    {
+      suffix=".fits";
+      rdsuffix="-raw-download.fits";
+    }
+
+  /* Currently Gnuastro doesn't read or write XML files (VOTable). So if
+     the downloaded file is an XML file but the user hasn't given an XML
+     suffix, abort and inform the user. */
+  if(p->cp.output)
+    {
+      if( !strcmp(suffix,".xml")
+          && strcmp(&p->cp.output[strlen(p->cp.output)-4], ".xml") )
+        error(EXIT_FAILURE, 0, "this dataset's output is a VOTable (with "
+              "an '.xml' suffix). However, Gnuastro doesn't yet support "
+              "VOTable, so it won't do any checks and corrections on "
+              "the downloaded file. Please give an output name with an "
+              "'.xml' suffix to continue");
+    }
+
   /* Set the name for the downloaded and final output name. These are due
      to an internal low-level processing that will be done on the raw
      downloaded file. */
-  if(p->cp.output==NULL)
+  else
     {
-      basename=gal_checkset_malloc_cat(p->databasestr, ".fits");
-      p->cp.output=gal_checkset_make_unique_suffix(basename, ".fits");
+      basename=gal_checkset_malloc_cat(p->databasestr, suffix);
+      p->cp.output=gal_checkset_make_unique_suffix(basename, suffix);
       free(basename);
     }
 
-  /* Make sure the output name doesn't exist (and report an error if
-     '--dontdelete' is called. Just note that for the automatic output, we
-     are basing that on the output, not the input. So we are temporarily
-     activating 'keepinputdir'. */
-  keepinputdir=p->cp.keepinputdir;
-  p->cp.keepinputdir=1;
-  gal_checkset_writable_remove(p->cp.output, 0, p->cp.dontdelete);
-  p->downloadname=gal_checkset_automatic_output(&p->cp, p->cp.output,
-                                                "-raw-download.fits");
-  p->cp.keepinputdir=keepinputdir;
+  /* Currently we don't interally process VOTable (in '.xml' suffix) files,
+     so to keep the next steps un-affected, we'll set Query to not delete
+     the raw download and copy the name of the output into the raw
+     download. */
+  if( !strcmp(suffix, ".xml") )
+    {
+      p->keeprawdownload=1;
+      gal_checkset_allocate_copy(p->cp.output, &p->downloadname);
+    }
+  else
+    {
+      /* Make sure the output name doesn't exist (and report an error if
+         '--dontdelete' is called. Just note that for the automatic output,
+         we are basing that on the output, not the input. So we are
+         temporarily activating 'keepinputdir'. */
+      keepinputdir=p->cp.keepinputdir;
+      p->cp.keepinputdir=1;
+      gal_checkset_writable_remove(p->cp.output, 0, p->cp.dontdelete);
+      p->downloadname=gal_checkset_automatic_output(&p->cp, p->cp.output,
+                                                    rdsuffix);
+      p->cp.keepinputdir=keepinputdir;
+    }
 }
 
 
diff --git a/bin/query/vizier.c b/bin/query/vizier.c
index d7513a0..f2dfcfd 100644
--- a/bin/query/vizier.c
+++ b/bin/query/vizier.c
@@ -143,6 +143,9 @@ vizier_sanity_checks(struct queryparams *p)
           gal_checkset_allocate_copy("II/363/unwise", &p->datasetstr);
         }
     }
+
+  /* Currently we assume VizieR only uses TAP. */
+  p->usetap=1;
 }
 
 
diff --git a/doc/gnuastro.texi b/doc/gnuastro.texi
index 33f4007..4c8b14c 100644
--- a/doc/gnuastro.texi
+++ b/doc/gnuastro.texi
@@ -10536,12 +10536,42 @@ Here is the list of short names for popular datasets 
within Gaia:
 @cindex NED (NASA/IPAC Extragalactic Database)
 The NASA/IPAC Extragalactic Database (NED, @url{http://ned.ipac.caltech.edu}) 
is a fusion database, integrating the information about extra-galactic sources 
from many large sky surveys into a single catalog.
 It covers the full spectrum, from Gamma rays to radio frequencies and is 
updated when new data arrives.
-A query to @code{ned} is submitted to 
@code{https://ned.ipac.caltech.edu/tap/sync}.
+A TAP query to @code{ned} is submitted to 
@code{https://ned.ipac.caltech.edu/tap/sync}.
 
-Currently NED only has its main dataset for TAP access (shown below), more 
datasets will be added for TAP access in the future.
 @itemize
 @item
-@code{objdir --> NEDTAP.objdir}
+@code{objdir --> NEDTAP.objdir}: default TAP-based dataset in NED.
+
+@item
+@cindex VOTable
+@code{extinction}: A command-line interface to the 
@url{https://ned.ipac.caltech.edu/extinction_calculator, NED Extinction 
Calculator}.
+It only takes a central coordinate and returns a VOTable of the calculated 
extinction in many commonly used filters at that point.
+As a result, options like @option{--width} or @option{--radius} are not 
supported.
+However, Gnuastro doesn't yet support the VOTable format.
+Therefore, if you specify an @option{--output} file, it should have an 
@file{.xml} suffix and the downloaded file will not be checked.
+
+Until VOTable support is added to Gnuastro, you can use GREP, AWK and SED to 
convert the VOTable data into a FITS table with a command like below (assuming 
the queried VOTable is called @file{ned-extinction.xml}):
+
+@verbatim
+grep '^<TR><TD>' ned-extinction.xml \
+ | sed -e's|<TR><TD>||' \
+       -e's|</TD></TR>||' \
+       -e's|</TD><TD>|@|g' \
+ | awk 'BEGIN{FS="@"; \
+     print "# Column 1: FILTER [name,str15] Filter name"; \
+     print "# Column 2: CENTRAL [um,f32] Central Wavelength"; \
+     print "# Column 3: EXTINCTION [mag,f32] Galactic Ext."; \
+     print "# Column 4: ADS_REF [ref,str50] ADS reference"} \
+        {printf "%-15s %g %g %s\n", $1, $2, $3, $4}' \
+ | asttable -oned-extinction.fits
+@end verbatim
+
+Once the table is in FITS, you can easily get the extinction for a certain 
filter (for example the @code{SDSS r} filter) like the command below:
+
+@example
+asttable ned-extinction.fits --equal=FILTER,"SDSS r" \
+         -cEXTINCTION
+@end example
 @end itemize
 
 @item vizier