fors_science.c

00001 /* $Id: fors_science.c,v 1.26 2011/03/03 14:37:00 cizzo Exp $
00002  *
00003  * This file is part of the FORS Data Reduction Pipeline
00004  * Copyright (C) 2002-2010 European Southern Observatory
00005  *
00006  * This program is free software; you can redistribute it and/or modify
00007  * it under the terms of the GNU General Public License as published by
00008  * the Free Software Foundation; either version 2 of the License, or
00009  * (at your option) any later version.
00010  *
00011  * This program is distributed in the hope that it will be useful,
00012  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00013  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00014  * GNU General Public License for more details.
00015  *
00016  * You should have received a copy of the GNU General Public License
00017  * along with this program; if not, write to the Free Software
00018  * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
00019  */
00020 
00021 /*
00022  * $Author: cizzo $
00023  * $Date: 2011/03/03 14:37:00 $
00024  * $Revision: 1.26 $
00025  * $Name: fors-4_8_6 $
00026  */
00027 
00028 #ifdef HAVE_CONFIG_H
00029 #include <config.h>
00030 #endif
00031 
00032 #include <math.h>
00033 #include <string.h>
00034 #include <cpl.h>
00035 #include <moses.h>
00036 #include <fors_dfs.h>
00037 #include <fors_tools.h>
00038 #include <fors_qc.h>
00039 
00040 static int fors_science_create(cpl_plugin *);
00041 static int fors_science_exec(cpl_plugin *);
00042 static int fors_science_destroy(cpl_plugin *);
00043 static int fors_science(cpl_parameterlist *, cpl_frameset *);
00044 
00045 static char fors_science_description[] =
00046 "This recipe is used to reduce scientific spectra using the extraction\n"
00047 "mask and the products created by the recipe fors_calib. The spectra are\n"
00048 "bias subtracted, flat fielded (if a normalised flat field is specified)\n"
00049 "and remapped eliminating the optical distortions. The wavelength calibration\n"
00050 "can be optionally upgraded using a number of sky lines: if no sky lines\n"
00051 "catalog of wavelengths is specified, an internal one is used instead.\n"
00052 "If the alignment to the sky lines is performed, the input dispersion\n"
00053 "coefficients table is upgraded and saved to disk, and a new CCD wavelengths\n"
00054 "map is created.\n"
00055 "This recipe accepts both FORS1 and FORS2 frames. A grism table (typically\n"
00056 "depending on the instrument mode, and in particular on the grism used)\n"
00057 "may also be specified: this table contains a default recipe parameter\n" 
00058 "setting to control the way spectra are extracted for a specific instrument\n"
00059 "mode, as it is used for automatic run of the pipeline on Paranal and in\n" 
00060 "Garching. If this table is specified, it will modify the default recipe\n" 
00061 "parameter setting, with the exception of those parameters which have been\n" 
00062 "explicitly modifyed on the command line. If a grism table is not specified,\n"
00063 "the input recipe parameters values will always be read from the command\n" 
00064 "line, or from an esorex configuration file if present, or from their\n" 
00065 "generic default values (that are rarely meaningful).\n" 
00066 "In the table below the MXU acronym can be read alternatively as MOS\n"
00067 "and LSS, depending on the instrument mode of the input data. The acronym\n"
00068 "SCI on products should be read STD in case of standard stars observations\n"
00069 "A CURV_COEFF table is not (yet) expected for LSS data.\n"
00070 "Either a scientific or a standard star exposure can be specified in input.\n"
00071 "Only in case of a standard star exposure input, the atmospheric extinction\n"
00072 "table and a table with the physical fluxes of the observed standard star\n"
00073 "must be specified in input, and a spectro-photometric table is created in\n"
00074 "output. This table can then be input again to this recipe, always with an\n"
00075 "atmospheric extinction table, and if a photometric calibration is requested\n"
00076 "then flux calibrated spectra (in units of erg/cm/cm/s/Angstrom) are also\n" 
00077 "written in output.\n\n"
00078 
00079 "Input files:\n\n"
00080 "  DO category:               Type:       Explanation:         Required:\n"
00081 "  SCIENCE_MXU                Raw         Scientific exposure     Y\n"
00082 "  or STANDARD_MXU            Raw         Standard star exposure  Y\n"
00083 "  MASTER_BIAS                Calib       Master bias             Y\n"
00084 "  GRISM_TABLE                Calib       Grism table             .\n"
00085 "  MASTER_SKYLINECAT          Calib       Sky lines catalog       .\n"
00086 "\n"
00087 "  MASTER_NORM_FLAT_MXU       Calib       Normalised flat field   .\n"
00088 "  DISP_COEFF_MXU             Calib       Inverse dispersion      Y\n"
00089 "  CURV_COEFF_MXU             Calib       Spectral curvature      Y\n"
00090 "  SLIT_LOCATION_MXU          Calib       Slits positions table   Y\n"
00091 "\n"
00092 "  or, in case of LSS-like MOS/MXU data,\n"
00093 "\n"
00094 "  MASTER_NORM_FLAT_LONG_MXU  Calib       Normalised flat field   .\n"
00095 "  DISP_COEFF_LONG_MXU        Calib       Inverse dispersion      Y\n"
00096 "\n"
00097 "  In case STANDARD_MXU is specified in input,\n"
00098 "\n"
00099 "  EXTINCT_TABLE              Calib       Atmospheric extinction  Y\n"
00100 "  STD_FLUX_TABLE             Calib       Standard star flux      Y\n"
00101 "\n"
00102 "  In case a photometric calibration is requested for scientific\n"
00103 "  data, the following inputs are mandatory:\n"
00104 "\n"
00105 "  EXTINCT_TABLE              Calib       Atmospheric extinction  Y\n"
00106 "  SPECPHOT_TABLE             Calib       Response curves         Y\n"
00107 "\n"
00108 "  If requested for standard star data, the SPECPHOT_TABLE can be dropped:\n"
00109 "  in this case the correction is applied using the SPECPHOT_TABLE produced\n"
00110 "  in the same run.\n"
00111 "\n"
00112 "Output files:\n"
00113 "\n"
00114 "  DO category:               Data type:  Explanation:\n"
00115 "  REDUCED_SCI_MXU            FITS image  Extracted scientific spectra\n"
00116 "  REDUCED_SKY_SCI_MXU        FITS image  Extracted sky spectra\n"
00117 "  REDUCED_ERROR_SCI_MXU      FITS image  Errors on extracted spectra\n"
00118 "  UNMAPPED_SCI_MXU           FITS image  Sky subtracted scientific spectra\n"
00119 "  MAPPED_SCI_MXU             FITS image  Rectified scientific spectra\n"
00120 "  MAPPED_ALL_SCI_MXU         FITS image  Rectified science spectra with sky\n"
00121 "  MAPPED_SKY_SCI_MXU         FITS image  Rectified sky spectra\n"
00122 "  UNMAPPED_SKY_SCI_MXU       FITS image  Sky on CCD\n"
00123 "  GLOBAL_SKY_SPECTRUM_MXU    FITS table  Global sky spectrum\n"
00124 "  OBJECT_TABLE_SCI_MXU       FITS table  Positions of detected objects\n"
00125 "\n"
00126 "  Only if the sky-alignment of the wavelength solution is requested:\n"
00127 "  SKY_SHIFTS_LONG_SCI_MXU    FITS table  Sky lines offsets (LSS-like data)\n"
00128 "  or SKY_SHIFTS_SLIT_SCI_MXU FITS table  Sky lines offsets (MOS-like data)\n"
00129 "  DISP_COEFF_SCI_MXU         FITS table  Upgraded dispersion coefficients\n"
00130 "  WAVELENGTH_MAP_SCI_MXU     FITS image  Upgraded wavelength map\n"
00131 "\n"
00132 "  Only if a STANDARD_MXU is specified in input:\n"
00133 "  SPECPHOT_TABLE             FITS table  Efficiency and response curves\n"
00134 "\n"
00135 "  Only if a photometric calibration was requested:\n"
00136 "  REDUCED_FLUX_SCI_MXU       FITS image  Flux calibrated scientific spectra\n"
00137 "  REDUCED_FLUX_ERROR_SCI_MXU FITS image  Errors on flux calibrated spectra\n"
00138 "  MAPPED_FLUX_SCI_MXU        FITS image  Flux calibrated slit spectra\n\n";
00139 
00140 #define fors_science_exit(message)            \
00141 {                                             \
00142 if (message) cpl_msg_error(recipe, message);  \
00143 cpl_free(exptime);                            \
00144 cpl_free(instrume);                           \
00145 cpl_image_delete(dummy);                      \
00146 cpl_image_delete(mapped);                     \
00147 cpl_image_delete(mapped_sky);                 \
00148 cpl_image_delete(mapped_cleaned);             \
00149 cpl_image_delete(skylocalmap);                \
00150 cpl_image_delete(skymap);                     \
00151 cpl_image_delete(smapped);                    \
00152 cpl_table_delete(offsets);                    \
00153 cpl_table_delete(photcal);                    \
00154 cpl_table_delete(sky);                        \
00155 cpl_image_delete(bias);                       \
00156 cpl_image_delete(spectra);                    \
00157 cpl_image_delete(coordinate);                 \
00158 cpl_image_delete(norm_flat);                  \
00159 cpl_image_delete(rainbow);                    \
00160 cpl_image_delete(rectified);                  \
00161 cpl_image_delete(wavemap);                    \
00162 cpl_image_delete(wavemaplss);                 \
00163 cpl_propertylist_delete(qclist);              \
00164 cpl_propertylist_delete(header);              \
00165 cpl_propertylist_delete(save_header);         \
00166 cpl_table_delete(grism_table);                \
00167 cpl_table_delete(idscoeff);                   \
00168 cpl_table_delete(maskslits);                  \
00169 cpl_table_delete(overscans);                  \
00170 cpl_table_delete(polytraces);                 \
00171 cpl_table_delete(slits);                      \
00172 cpl_table_delete(wavelengths);                \
00173 cpl_vector_delete(lines);                     \
00174 cpl_msg_indent_less();                        \
00175 return -1;                                    \
00176 }
00177 
00178 
00179 #define fors_science_exit_memcheck(message)   \
00180 {                                             \
00181 if (message) cpl_msg_info(recipe, message);   \
00182 cpl_free(exptime);                            \
00183 cpl_free(instrume);                           \
00184 cpl_image_delete(dummy);                      \
00185 cpl_image_delete(mapped);                     \
00186 cpl_image_delete(mapped_cleaned);             \
00187 cpl_image_delete(mapped_sky);                 \
00188 cpl_image_delete(skylocalmap);                \
00189 cpl_image_delete(skymap);                     \
00190 cpl_image_delete(smapped);                    \
00191 cpl_table_delete(offsets);                    \
00192 cpl_table_delete(photcal);                    \
00193 cpl_table_delete(sky);                        \
00194 cpl_image_delete(bias);                       \
00195 cpl_image_delete(spectra);                    \
00196 cpl_image_delete(coordinate);                 \
00197 cpl_image_delete(norm_flat);                  \
00198 cpl_image_delete(rainbow);                    \
00199 cpl_image_delete(rectified);                  \
00200 cpl_image_delete(wavemap);                    \
00201 cpl_image_delete(wavemaplss);                 \
00202 cpl_propertylist_delete(header);              \
00203 cpl_propertylist_delete(save_header);         \
00204 cpl_table_delete(grism_table);                \
00205 cpl_table_delete(idscoeff);                   \
00206 cpl_table_delete(maskslits);                  \
00207 cpl_table_delete(overscans);                  \
00208 cpl_table_delete(polytraces);                 \
00209 cpl_table_delete(slits);                      \
00210 cpl_table_delete(wavelengths);                \
00211 cpl_vector_delete(lines);                     \
00212 cpl_msg_indent_less();                        \
00213 return 0;                                     \
00214 }
00215 
00216 
00228 int cpl_plugin_get_info(cpl_pluginlist *list)
00229 {
00230     cpl_recipe *recipe = cpl_calloc(1, sizeof *recipe );
00231     cpl_plugin *plugin = &recipe->interface;
00232 
00233     cpl_plugin_init(plugin,
00234                     CPL_PLUGIN_API,
00235                     FORS_BINARY_VERSION,
00236                     CPL_PLUGIN_TYPE_RECIPE,
00237                     "fors_science",
00238                     "Extraction of scientific spectra",
00239                     fors_science_description,
00240                     "Carlo Izzo",
00241                     PACKAGE_BUGREPORT,
00242     "This file is currently part of the FORS Instrument Pipeline\n"
00243     "Copyright (C) 2002-2010 European Southern Observatory\n\n"
00244     "This program is free software; you can redistribute it and/or modify\n"
00245     "it under the terms of the GNU General Public License as published by\n"
00246     "the Free Software Foundation; either version 2 of the License, or\n"
00247     "(at your option) any later version.\n\n"
00248     "This program is distributed in the hope that it will be useful,\n"
00249     "but WITHOUT ANY WARRANTY; without even the implied warranty of\n"
00250     "MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the\n"
00251     "GNU General Public License for more details.\n\n"
00252     "You should have received a copy of the GNU General Public License\n"
00253     "along with this program; if not, write to the Free Software Foundation,\n"
00254     "Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA\n",
00255                     fors_science_create,
00256                     fors_science_exec,
00257                     fors_science_destroy);
00258 
00259     cpl_pluginlist_append(list, plugin);
00260     
00261     return 0;
00262 }
00263 
00264 
00275 static int fors_science_create(cpl_plugin *plugin)
00276 {
00277     cpl_recipe    *recipe;
00278     cpl_parameter *p;
00279 
00280 
00281     /* 
00282      * Check that the plugin is part of a valid recipe 
00283      */
00284 
00285     if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 
00286         recipe = (cpl_recipe *)plugin;
00287     else 
00288         return -1;
00289 
00290     /* 
00291      * Create the parameters list in the cpl_recipe object 
00292      */
00293 
00294     recipe->parameters = cpl_parameterlist_new(); 
00295 
00296 
00297     /*
00298      * Dispersion
00299      */
00300 
00301     p = cpl_parameter_new_value("fors.fors_science.dispersion",
00302                                 CPL_TYPE_DOUBLE,
00303                                 "Resampling step (Angstrom/pixel)",
00304                                 "fors.fors_science",
00305                                 0.0);
00306     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "dispersion");
00307     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00308     cpl_parameterlist_append(recipe->parameters, p);
00309 
00310     /*
00311      * Sky lines alignment
00312      */
00313 
00314     p = cpl_parameter_new_value("fors.fors_science.skyalign",
00315                                 CPL_TYPE_INT,
00316                                 "Polynomial order for sky lines alignment, "
00317                                 "or -1 to avoid alignment",
00318                                 "fors.fors_science",
00319                                 0);
00320     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skyalign");
00321     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00322     cpl_parameterlist_append(recipe->parameters, p);
00323 
00324     /*
00325      * Line catalog table column containing the sky reference wavelengths
00326      */
00327 
00328     p = cpl_parameter_new_value("fors.fors_science.wcolumn",
00329                                 CPL_TYPE_STRING,
00330                                 "Name of sky line catalog table column "
00331                                 "with wavelengths",
00332                                 "fors.fors_science",
00333                                 "WLEN");
00334     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "wcolumn");
00335     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00336     cpl_parameterlist_append(recipe->parameters, p);
00337 
00338     /*
00339      * Start wavelength for spectral extraction
00340      */
00341 
00342     p = cpl_parameter_new_value("fors.fors_science.startwavelength",
00343                                 CPL_TYPE_DOUBLE,
00344                                 "Start wavelength in spectral extraction",
00345                                 "fors.fors_science",
00346                                 0.0);
00347     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "startwavelength");
00348     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00349     cpl_parameterlist_append(recipe->parameters, p);
00350 
00351     /*
00352      * End wavelength for spectral extraction
00353      */
00354 
00355     p = cpl_parameter_new_value("fors.fors_science.endwavelength",
00356                                 CPL_TYPE_DOUBLE,
00357                                 "End wavelength in spectral extraction",
00358                                 "fors.fors_science",
00359                                 0.0);
00360     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "endwavelength");
00361     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00362     cpl_parameterlist_append(recipe->parameters, p);
00363 
00364     /*
00365      * Flux conservation
00366      */
00367 
00368     p = cpl_parameter_new_value("fors.fors_science.flux",
00369                                 CPL_TYPE_BOOL,
00370                                 "Apply flux conservation",
00371                                 "fors.fors_science",
00372                                 TRUE);
00373     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "flux");
00374     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00375     cpl_parameterlist_append(recipe->parameters, p);
00376 
00377     /*
00378      * Apply flat field
00379      */
00380 
00381     p = cpl_parameter_new_value("fors.fors_science.flatfield",
00382                                 CPL_TYPE_BOOL,
00383                                 "Apply flat field",
00384                                 "fors.fors_science",
00385                                 TRUE);
00386     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "flatfield");
00387     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00388     cpl_parameterlist_append(recipe->parameters, p);
00389 
00390     /*
00391      * Global sky subtraction
00392      */
00393 
00394     p = cpl_parameter_new_value("fors.fors_science.skyglobal",
00395                                 CPL_TYPE_BOOL,
00396                                 "Subtract global sky spectrum from CCD",
00397                                 "fors.fors_science",
00398                                 FALSE);
00399     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skyglobal");
00400     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00401     cpl_parameterlist_append(recipe->parameters, p);
00402 
00403     /*
00404      * Local sky subtraction on extracted spectra
00405      */
00406 
00407 /*** New sky subtraction (search NSS)
00408     p = cpl_parameter_new_value("fors.fors_science.skymedian",
00409                                 CPL_TYPE_INT,
00410                                 "Degree of sky fitting polynomial for "
00411                                 "sky subtraction from extracted "
00412                                 "slit spectra (MOS/MXU only, -1 to disable it)",
00413                                 "fors.fors_science",
00414                                 0);
00415     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skymedian");
00416     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00417     cpl_parameterlist_append(recipe->parameters, p);
00418 ***/
00419 
00420     p = cpl_parameter_new_value("fors.fors_science.skymedian",
00421                                 CPL_TYPE_BOOL,
00422                                 "Sky subtraction from extracted slit spectra",
00423                                 "fors.fors_science",
00424                                 FALSE);
00425     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skymedian");
00426     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00427     cpl_parameterlist_append(recipe->parameters, p);
00428 
00429     /*
00430      * Local sky subtraction on CCD spectra
00431      */
00432 
00433     p = cpl_parameter_new_value("fors.fors_science.skylocal",
00434                                 CPL_TYPE_BOOL,
00435                                 "Sky subtraction from CCD slit spectra",
00436                                 "fors.fors_science",
00437                                 TRUE);
00438     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skylocal");
00439     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00440     cpl_parameterlist_append(recipe->parameters, p);
00441 
00442     /*
00443      * Cosmic rays removal
00444      */
00445 
00446     p = cpl_parameter_new_value("fors.fors_science.cosmics",
00447                                 CPL_TYPE_BOOL,
00448                                 "Eliminate cosmic rays hits (only if global "
00449                                 "sky subtraction is also requested)",
00450                                 "fors.fors_science",
00451                                 FALSE);
00452     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "cosmics");
00453     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00454     cpl_parameterlist_append(recipe->parameters, p);
00455 
00456     /*
00457      * Slit margin
00458      */
00459 
00460     p = cpl_parameter_new_value("fors.fors_science.slit_margin",
00461                                 CPL_TYPE_INT,
00462                                 "Number of pixels to exclude at each slit "
00463                                 "in object detection and extraction",
00464                                 "fors.fors_science",
00465                                 3);
00466     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "slit_margin");
00467     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00468     cpl_parameterlist_append(recipe->parameters, p);
00469 
00470     /*
00471      * Extraction radius
00472      */
00473 
00474     p = cpl_parameter_new_value("fors.fors_science.ext_radius",
00475                                 CPL_TYPE_INT,
00476                                 "Maximum extraction radius for detected "
00477                                 "objects (pixel)",
00478                                 "fors.fors_science",
00479                                 12);
00480     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "ext_radius");
00481     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00482     cpl_parameterlist_append(recipe->parameters, p);
00483 
00484     /*
00485      * Contamination radius
00486      */
00487 
00488     p = cpl_parameter_new_value("fors.fors_science.cont_radius",
00489                                 CPL_TYPE_INT,
00490                                 "Minimum distance at which two objects "
00491                                 "of equal luminosity do not contaminate "
00492                                 "each other (pixel)",
00493                                 "fors.fors_science",
00494                                 0);
00495     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "cont_radius");
00496     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00497     cpl_parameterlist_append(recipe->parameters, p);
00498 
00499     /*
00500      * Object extraction method
00501      */
00502 
00503     p = cpl_parameter_new_value("fors.fors_science.ext_mode",
00504                                 CPL_TYPE_INT,
00505                                 "Object extraction method: 0 = aperture, "
00506                                 "1 = Horne optimal extraction",
00507                                 "fors.fors_science",
00508                                 1);
00509     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "ext_mode");
00510     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00511     cpl_parameterlist_append(recipe->parameters, p);
00512 
00513     /*
00514      * Normalise output by exposure time
00515      */
00516 
00517     p = cpl_parameter_new_value("fors.fors_science.time_normalise",
00518                                 CPL_TYPE_BOOL,
00519                                 "Normalise output spectra by the exposure time",
00520                                 "fors.fors_science",
00521                                 TRUE);
00522     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "time_normalise");
00523     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00524     cpl_parameterlist_append(recipe->parameters, p);
00525 
00526     /*
00527      * Order of polynomial modeling the instrument response.
00528      */
00529 
00530     p = cpl_parameter_new_value("fors.fors_science.response",
00531                                 CPL_TYPE_INT,
00532                                 "Order of polynomial modeling the "
00533                                 "instrument response",
00534                                 "fors.fors_science",
00535                                 5);
00536     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "response");
00537     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00538     cpl_parameterlist_append(recipe->parameters, p);
00539 
00540     /*
00541      * Order of polynomial modeling the instrument response.
00542      */
00543 
00544     p = cpl_parameter_new_value("fors.fors_science.photometry",
00545                                 CPL_TYPE_BOOL,
00546                                 "Apply spectrophotometric calibration",
00547                                 "fors.fors_science",
00548                                 FALSE);
00549     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "photometry");
00550     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00551     cpl_parameterlist_append(recipe->parameters, p);
00552 
00553     /*
00554      * Computation of QC1 parameters
00555      */
00556 
00557     p = cpl_parameter_new_value("fors.fors_science.qc",
00558                                 CPL_TYPE_BOOL,
00559                                 "Compute QC1 parameters",
00560                                 "fors.fors_science",
00561                                 TRUE);
00562     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "qc");
00563     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00564     cpl_parameterlist_append(recipe->parameters, p);
00565 
00566     return 0;
00567 }
00568 
00569 
00578 static int fors_science_exec(cpl_plugin *plugin)
00579 {
00580     cpl_recipe *recipe;
00581     
00582     if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 
00583         recipe = (cpl_recipe *)plugin;
00584     else 
00585         return -1;
00586 
00587     return fors_science(recipe->parameters, recipe->frames);
00588 }
00589 
00590 
00599 static int fors_science_destroy(cpl_plugin *plugin)
00600 {
00601     cpl_recipe *recipe;
00602     
00603     if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 
00604         recipe = (cpl_recipe *)plugin;
00605     else 
00606         return -1;
00607 
00608     cpl_parameterlist_delete(recipe->parameters); 
00609 
00610     return 0;
00611 }
00612 
00613 
00623 static int fors_science(cpl_parameterlist *parlist, cpl_frameset *frameset)
00624 {
00625 
00626     const char *recipe = "fors_science";
00627 
00628 
00629     /*
00630      * Input parameters
00631      */
00632 
00633     double      dispersion;
00634     int         skyalign;
00635     const char *wcolumn;
00636     double      startwavelength;
00637     double      endwavelength;
00638     int         flux;
00639     int         flatfield;
00640     int         skyglobal;
00641     int         skylocal;
00642     int         skymedian;
00643     int         cosmics;
00644     int         slit_margin;
00645     int         ext_radius;
00646     int         cont_radius;
00647     int         ext_mode;
00648     int         time_normalise;
00649     int         res_order;
00650     int         photometry;
00651     int         qc;
00652 
00653 
00654     /*
00655      * CPL objects
00656      */
00657 
00658     cpl_imagelist    *all_science;
00659     cpl_image       **images;
00660 
00661     cpl_image        *bias           = NULL;
00662     cpl_image        *norm_flat      = NULL;
00663     cpl_image        *spectra        = NULL;
00664     cpl_image        *rectified      = NULL;
00665     cpl_image        *coordinate     = NULL;
00666     cpl_image        *rainbow        = NULL;
00667     cpl_image        *mapped         = NULL;
00668     cpl_image        *mapped_sky     = NULL;
00669     cpl_image        *mapped_cleaned = NULL;
00670     cpl_image        *smapped        = NULL;
00671     cpl_image        *wavemap        = NULL;
00672     cpl_image        *wavemaplss     = NULL;
00673     cpl_image        *skymap         = NULL;
00674     cpl_image        *skylocalmap    = NULL;
00675     cpl_image        *dummy          = NULL;
00676 
00677     cpl_table        *grism_table    = NULL;
00678     cpl_table        *overscans      = NULL;
00679     cpl_table        *wavelengths    = NULL;
00680     cpl_table        *idscoeff       = NULL;
00681     cpl_table        *slits          = NULL;
00682     cpl_table        *maskslits      = NULL;
00683     cpl_table        *polytraces     = NULL;
00684     cpl_table        *offsets        = NULL;
00685     cpl_table        *sky            = NULL;
00686     cpl_table        *photcal        = NULL;
00687 
00688     cpl_vector       *lines          = NULL;
00689 
00690     cpl_propertylist *header         = NULL;
00691     cpl_propertylist *save_header    = NULL;
00692     cpl_propertylist *qclist         = NULL;
00693 
00694 
00695     /*
00696      * Auxiliary variables
00697      */
00698 
00699     char        version[80];
00700     char       *instrume = NULL;
00701     char       *wheel4;
00702     const char *science_tag;
00703     const char *master_norm_flat_tag;
00704     const char *disp_coeff_tag;
00705     const char *disp_coeff_sky_tag;
00706     const char *wavelength_map_sky_tag;
00707     const char *curv_coeff_tag;
00708     const char *slit_location_tag;
00709     const char *reduced_science_tag;
00710     const char *reduced_flux_science_tag;
00711     const char *reduced_sky_tag;
00712     const char *reduced_error_tag;
00713     const char *reduced_flux_error_tag;
00714     const char *mapped_science_tag;
00715     const char *mapped_flux_science_tag;
00716     const char *unmapped_science_tag;
00717     const char *mapped_science_sky_tag;
00718     const char *mapped_sky_tag;
00719     const char *unmapped_sky_tag;
00720     const char *global_sky_spectrum_tag;
00721     const char *object_table_tag;
00722     const char *skylines_offsets_tag;
00723     const char *specphot_tag;
00724     const char *master_specphot_tag = "MASTER_SPECPHOT_TABLE";
00725     int         mxu, mos, lss;
00726     int         treat_as_lss = 0;
00727     int         have_phot = 0;
00728     int         nslits;
00729     int         nscience;
00730     double     *xpos;
00731     double     *exptime = NULL;
00732     double      alltime;
00733     double      airmass;
00734     double      mxpos;
00735     double      mean_rms;
00736     int         nlines;
00737     int         rebin;
00738     double     *line;
00739     int         nx, ny;
00740     int         ccd_xsize, ccd_ysize;
00741     double      reference;
00742     double      gain;
00743     double      ron;
00744     int         standard;
00745     int         highres;
00746     int         i;
00747     double      wstart;
00748     double      wstep;
00749     int         wcount;
00750 
00751 
00752     snprintf(version, 80, "%s-%s", PACKAGE, PACKAGE_VERSION);
00753 
00754     cpl_msg_set_indentation(2);
00755 
00756     if (dfs_files_dont_exist(frameset))
00757         fors_science_exit(NULL);
00758 
00759 
00760     /* 
00761      * Get configuration parameters
00762      */
00763 
00764     cpl_msg_info(recipe, "Recipe %s configuration parameters:", recipe);
00765     cpl_msg_indent_more();
00766 
00767     if (cpl_frameset_count_tags(frameset, "GRISM_TABLE") > 1)
00768         fors_science_exit("Too many in input: GRISM_TABLE");
00769 
00770     grism_table = dfs_load_table(frameset, "GRISM_TABLE", 1);
00771 
00772     dispersion = dfs_get_parameter_double(parlist, 
00773                     "fors.fors_science.dispersion", grism_table);
00774 
00775     if (dispersion <= 0.0)
00776         fors_science_exit("Invalid resampling step");
00777 
00778     skyalign = dfs_get_parameter_int(parlist, 
00779                     "fors.fors_science.skyalign", NULL);
00780 
00781     if (skyalign > 2)
00782         fors_science_exit("Max polynomial degree for sky alignment is 2");
00783 
00784     wcolumn = dfs_get_parameter_string(parlist, 
00785                     "fors.fors_science.wcolumn", NULL);
00786 
00787     startwavelength = dfs_get_parameter_double(parlist, 
00788                     "fors.fors_science.startwavelength", grism_table);
00789     if (startwavelength < 3000.0 || startwavelength > 13000.0)
00790         fors_science_exit("Invalid wavelength");
00791 
00792     endwavelength = dfs_get_parameter_double(parlist, 
00793                     "fors.fors_science.endwavelength", grism_table);
00794     if (endwavelength < 3000.0 || endwavelength > 13000.0)
00795         fors_science_exit("Invalid wavelength");
00796 
00797     if (endwavelength - startwavelength <= 0.0)
00798         fors_science_exit("Invalid wavelength interval");
00799 
00800     flux = dfs_get_parameter_bool(parlist, "fors.fors_science.flux", NULL);
00801 
00802     flatfield = dfs_get_parameter_bool(parlist, "fors.fors_science.flatfield", 
00803                                        NULL);
00804 
00805     skyglobal = dfs_get_parameter_bool(parlist, "fors.fors_science.skyglobal", 
00806                                        NULL);
00807     skylocal  = dfs_get_parameter_bool(parlist, "fors.fors_science.skylocal", 
00808                                        NULL);
00809     skymedian = dfs_get_parameter_bool(parlist, "fors.fors_science.skymedian", 
00810                                        NULL);
00811 /* NSS
00812     skymedian = dfs_get_parameter_int(parlist, "fors.fors_science.skymedian", 
00813                                        NULL);
00814 */
00815 
00816     if (skylocal && skyglobal)
00817         fors_science_exit("Cannot apply both local and global sky subtraction");
00818 
00819     if (skylocal && skymedian)
00820         fors_science_exit("Cannot apply sky subtraction both on extracted "
00821                           "and non-extracted spectra");
00822 
00823     cosmics = dfs_get_parameter_bool(parlist, 
00824                                      "fors.fors_science.cosmics", NULL);
00825 
00826     if (cosmics)
00827         if (!(skyglobal || skylocal))
00828             fors_science_exit("Cosmic rays correction requires "
00829                               "either skylocal=true or skyglobal=true");
00830 
00831     slit_margin = dfs_get_parameter_int(parlist, 
00832                                         "fors.fors_science.slit_margin",
00833                                         NULL);
00834     if (slit_margin < 0)
00835         fors_science_exit("Value must be zero or positive");
00836 
00837     ext_radius = dfs_get_parameter_int(parlist, "fors.fors_science.ext_radius",
00838                                        NULL);
00839     if (ext_radius < 0)
00840         fors_science_exit("Value must be zero or positive");
00841 
00842     cont_radius = dfs_get_parameter_int(parlist, 
00843                                         "fors.fors_science.cont_radius",
00844                                         NULL);
00845     if (cont_radius < 0)
00846         fors_science_exit("Value must be zero or positive");
00847 
00848     ext_mode = dfs_get_parameter_int(parlist, "fors.fors_science.ext_mode",
00849                                        NULL);
00850     if (ext_mode < 0 || ext_mode > 1)
00851         fors_science_exit("Invalid object extraction mode");
00852 
00853     time_normalise = dfs_get_parameter_bool(parlist, 
00854                              "fors.fors_science.time_normalise", NULL);
00855 
00856     res_order = dfs_get_parameter_int(parlist, "fors.fors_science.response",
00857                                       NULL);
00858     if (res_order < 2 || res_order > 10)
00859         fors_science_exit("Invalid instrument response modeling polynomial");
00860 
00861     photometry = dfs_get_parameter_bool(parlist,
00862                              "fors.fors_science.photometry", NULL);
00863 
00864     qc = dfs_get_parameter_bool(parlist, "fors.fors_science.qc", NULL);
00865 
00866     cpl_table_delete(grism_table); grism_table = NULL;
00867 
00868     if (cpl_error_get_code())
00869         fors_science_exit("Failure getting the configuration parameters");
00870 
00871     
00872     /* 
00873      * Check input set-of-frames
00874      */
00875 
00876     cpl_msg_indent_less();
00877     cpl_msg_info(recipe, "Check input set-of-frames:");
00878     cpl_msg_indent_more();
00879 
00880     if (!dfs_equal_keyword(frameset, "ESO INS GRIS1 ID"))
00881         fors_science_exit("Input frames are not from the same grism");
00882 
00883     if (!dfs_equal_keyword(frameset, "ESO INS FILT1 ID"))
00884         fors_science_exit("Input frames are not from the same filter");
00885 
00886     if (!dfs_equal_keyword(frameset, "ESO DET CHIP1 ID"))
00887         fors_science_exit("Input frames are not from the same chip");
00888 
00889     mxu = cpl_frameset_count_tags(frameset, "SCIENCE_MXU");
00890     mos = cpl_frameset_count_tags(frameset, "SCIENCE_MOS");
00891     lss = cpl_frameset_count_tags(frameset, "SCIENCE_LSS");
00892     standard = 0;
00893 
00894     if (mxu + mos + lss == 0) {
00895         mxu = cpl_frameset_count_tags(frameset, "STANDARD_MXU");
00896         mos = cpl_frameset_count_tags(frameset, "STANDARD_MOS");
00897         lss = cpl_frameset_count_tags(frameset, "STANDARD_LSS");
00898         standard = 1;
00899     }
00900 
00901     if (mxu + mos + lss == 0)
00902         fors_science_exit("Missing input scientific frame");
00903 
00904     nscience = mxu + mos + lss;
00905 
00906     if (mxu && mxu < nscience)
00907         fors_science_exit("Input scientific frames must be of the same type"); 
00908 
00909     if (mos && mos < nscience)
00910         fors_science_exit("Input scientific frames must be of the same type"); 
00911 
00912     if (lss && lss < nscience)
00913         fors_science_exit("Input scientific frames must be of the same type"); 
00914 
00915     if (mxu) {
00916         cpl_msg_info(recipe, "MXU data found");
00917         if (standard) {
00918             science_tag              = "STANDARD_MXU";
00919             reduced_science_tag      = "REDUCED_STD_MXU";
00920             reduced_flux_science_tag = "REDUCED_FLUX_STD_MXU";
00921             unmapped_science_tag     = "UNMAPPED_STD_MXU";
00922             mapped_science_tag       = "MAPPED_STD_MXU";
00923             mapped_flux_science_tag  = "MAPPED_FLUX_STD_MXU";
00924             mapped_science_sky_tag   = "MAPPED_ALL_STD_MXU";
00925             skylines_offsets_tag     = "SKY_SHIFTS_SLIT_STD_MXU";
00926             wavelength_map_sky_tag   = "WAVELENGTH_MAP_STD_MXU";
00927             disp_coeff_sky_tag       = "DISP_COEFF_STD_MXU";
00928             mapped_sky_tag           = "MAPPED_SKY_STD_MXU";
00929             unmapped_sky_tag         = "UNMAPPED_SKY_STD_MXU";
00930             object_table_tag         = "OBJECT_TABLE_STD_MXU";
00931             reduced_sky_tag          = "REDUCED_SKY_STD_MXU";
00932             reduced_error_tag        = "REDUCED_ERROR_STD_MXU";
00933             reduced_flux_error_tag   = "REDUCED_FLUX_ERROR_STD_MXU";
00934             specphot_tag             = "SPECPHOT_TABLE";
00935         }
00936         else {
00937             science_tag              = "SCIENCE_MXU";
00938             reduced_science_tag      = "REDUCED_SCI_MXU";
00939             reduced_flux_science_tag = "REDUCED_FLUX_SCI_MXU";
00940             unmapped_science_tag     = "UNMAPPED_SCI_MXU";
00941             mapped_science_tag       = "MAPPED_SCI_MXU";
00942             mapped_flux_science_tag  = "MAPPED_FLUX_SCI_MXU";
00943             mapped_science_sky_tag   = "MAPPED_ALL_SCI_MXU";
00944             skylines_offsets_tag     = "SKY_SHIFTS_SLIT_SCI_MXU";
00945             wavelength_map_sky_tag   = "WAVELENGTH_MAP_SCI_MXU";
00946             disp_coeff_sky_tag       = "DISP_COEFF_SCI_MXU";
00947             mapped_sky_tag           = "MAPPED_SKY_SCI_MXU";
00948             unmapped_sky_tag         = "UNMAPPED_SKY_SCI_MXU";
00949             object_table_tag         = "OBJECT_TABLE_SCI_MXU";
00950             reduced_sky_tag          = "REDUCED_SKY_SCI_MXU";
00951             reduced_error_tag        = "REDUCED_ERROR_SCI_MXU";
00952             reduced_flux_error_tag   = "REDUCED_FLUX_ERROR_SCI_MXU";
00953             specphot_tag             = "SPECPHOT_TABLE";
00954         }
00955 
00956         master_norm_flat_tag    = "MASTER_NORM_FLAT_MXU";
00957         disp_coeff_tag          = "DISP_COEFF_MXU";
00958         curv_coeff_tag          = "CURV_COEFF_MXU";
00959         slit_location_tag       = "SLIT_LOCATION_MXU";
00960         global_sky_spectrum_tag = "GLOBAL_SKY_SPECTRUM_MXU";
00961 
00962         if (!cpl_frameset_count_tags(frameset, master_norm_flat_tag)) {
00963             master_norm_flat_tag   = "MASTER_NORM_FLAT_LONG_MXU";
00964             disp_coeff_tag         = "DISP_COEFF_LONG_MXU";
00965             slit_location_tag      = "SLIT_LOCATION_LONG_MXU";
00966         }
00967     }
00968 
00969     if (lss) {
00970         cpl_msg_info(recipe, "LSS data found");
00971 
00972         if (cosmics && !skyglobal)
00973             fors_science_exit("Cosmic rays correction for LSS "
00974                               "data requires --skyglobal=true");
00975 
00976         if (standard) {
00977             science_tag              = "STANDARD_LSS";
00978             reduced_science_tag      = "REDUCED_STD_LSS";
00979             reduced_flux_science_tag = "REDUCED_FLUX_STD_LSS";
00980             unmapped_science_tag     = "UNMAPPED_STD_LSS";
00981             mapped_science_tag       = "MAPPED_STD_LSS";
00982             mapped_flux_science_tag  = "MAPPED_FLUX_STD_LSS";
00983             mapped_science_sky_tag   = "MAPPED_ALL_STD_LSS";
00984             skylines_offsets_tag     = "SKY_SHIFTS_LONG_STD_LSS";
00985             wavelength_map_sky_tag   = "WAVELENGTH_MAP_STD_LSS";
00986             disp_coeff_sky_tag       = "DISP_COEFF_STD_LSS";
00987             mapped_sky_tag           = "MAPPED_SKY_STD_LSS";
00988             unmapped_sky_tag         = "UNMAPPED_SKY_STD_LSS";
00989             object_table_tag         = "OBJECT_TABLE_STD_LSS";
00990             reduced_sky_tag          = "REDUCED_SKY_STD_LSS";
00991             reduced_error_tag        = "REDUCED_ERROR_STD_LSS";
00992             reduced_flux_error_tag   = "REDUCED_FLUX_ERROR_STD_LSS";
00993             specphot_tag             = "SPECPHOT_TABLE";
00994         }
00995         else {
00996             science_tag              = "SCIENCE_LSS";
00997             reduced_science_tag      = "REDUCED_SCI_LSS";
00998             reduced_flux_science_tag = "REDUCED_FLUX_SCI_LSS";
00999             unmapped_science_tag     = "UNMAPPED_SCI_LSS";
01000             mapped_science_tag       = "MAPPED_SCI_LSS";
01001             mapped_flux_science_tag  = "MAPPED_FLUX_SCI_LSS";
01002             mapped_science_sky_tag   = "MAPPED_ALL_SCI_LSS";
01003             skylines_offsets_tag     = "SKY_SHIFTS_LONG_SCI_LSS";
01004             wavelength_map_sky_tag   = "WAVELENGTH_MAP_SCI_LSS";
01005             disp_coeff_sky_tag       = "DISP_COEFF_SCI_LSS";
01006             mapped_sky_tag           = "MAPPED_SKY_SCI_LSS";
01007             unmapped_sky_tag         = "UNMAPPED_SKY_SCI_LSS";
01008             object_table_tag         = "OBJECT_TABLE_SCI_LSS";
01009             reduced_sky_tag          = "REDUCED_SKY_SCI_LSS";
01010             reduced_error_tag        = "REDUCED_ERROR_SCI_LSS";
01011             reduced_flux_error_tag   = "REDUCED_FLUX_ERROR_SCI_LSS";
01012             specphot_tag             = "SPECPHOT_TABLE";
01013         }
01014 
01015         master_norm_flat_tag    = "MASTER_NORM_FLAT_LSS";
01016         disp_coeff_tag          = "DISP_COEFF_LSS";
01017         slit_location_tag       = "SLIT_LOCATION_LSS";
01018         global_sky_spectrum_tag = "GLOBAL_SKY_SPECTRUM_LSS";
01019     }
01020 
01021     if (mos) {
01022         cpl_msg_info(recipe, "MOS data found");
01023         if (standard) {
01024             science_tag              = "STANDARD_MOS";
01025             reduced_science_tag      = "REDUCED_STD_MOS";
01026             reduced_flux_science_tag = "REDUCED_FLUX_STD_MOS";
01027             unmapped_science_tag     = "UNMAPPED_STD_MOS";
01028             mapped_science_tag       = "MAPPED_STD_MOS";
01029             mapped_flux_science_tag  = "MAPPED_FLUX_STD_MOS";
01030             mapped_science_sky_tag   = "MAPPED_ALL_STD_MOS";
01031             skylines_offsets_tag     = "SKY_SHIFTS_SLIT_STD_MOS";
01032             wavelength_map_sky_tag   = "WAVELENGTH_MAP_STD_MOS";
01033             disp_coeff_sky_tag       = "DISP_COEFF_STD_MOS";
01034             mapped_sky_tag           = "MAPPED_SKY_STD_MOS";
01035             unmapped_sky_tag         = "UNMAPPED_SKY_STD_MOS";
01036             object_table_tag         = "OBJECT_TABLE_STD_MOS";
01037             reduced_sky_tag          = "REDUCED_SKY_STD_MOS";
01038             reduced_error_tag        = "REDUCED_ERROR_STD_MOS";
01039             reduced_flux_error_tag   = "REDUCED_FLUX_ERROR_STD_MOS";
01040             specphot_tag             = "SPECPHOT_TABLE";
01041         }
01042         else {
01043             science_tag              = "SCIENCE_MOS";
01044             reduced_science_tag      = "REDUCED_SCI_MOS";
01045             reduced_flux_science_tag = "REDUCED_FLUX_SCI_MOS";
01046             unmapped_science_tag     = "UNMAPPED_SCI_MOS";
01047             mapped_science_tag       = "MAPPED_SCI_MOS";
01048             mapped_flux_science_tag  = "MAPPED_FLUX_SCI_MOS";
01049             mapped_science_sky_tag   = "MAPPED_ALL_SCI_MOS";
01050             skylines_offsets_tag     = "SKY_SHIFTS_SLIT_SCI_MOS";
01051             wavelength_map_sky_tag   = "WAVELENGTH_MAP_SCI_MOS";
01052             disp_coeff_sky_tag       = "DISP_COEFF_SCI_MOS";
01053             mapped_sky_tag           = "MAPPED_SKY_SCI_MOS";
01054             unmapped_sky_tag         = "UNMAPPED_SKY_SCI_MOS";
01055             object_table_tag         = "OBJECT_TABLE_SCI_MOS";
01056             reduced_sky_tag          = "REDUCED_SKY_SCI_MOS";
01057             reduced_error_tag        = "REDUCED_ERROR_SCI_MOS";
01058             reduced_flux_error_tag   = "REDUCED_FLUX_ERROR_SCI_MOS";
01059             specphot_tag             = "SPECPHOT_TABLE";
01060         }
01061 
01062         master_norm_flat_tag    = "MASTER_NORM_FLAT_MOS";
01063         disp_coeff_tag          = "DISP_COEFF_MOS";
01064         curv_coeff_tag          = "CURV_COEFF_MOS";
01065         slit_location_tag       = "SLIT_LOCATION_MOS";
01066         global_sky_spectrum_tag = "GLOBAL_SKY_SPECTRUM_MOS";
01067 
01068         if (!cpl_frameset_count_tags(frameset, master_norm_flat_tag)) {
01069             master_norm_flat_tag   = "MASTER_NORM_FLAT_LONG_MOS";
01070             disp_coeff_tag         = "DISP_COEFF_LONG_MOS";
01071             slit_location_tag      = "SLIT_LOCATION_LONG_MOS";
01072         }
01073     }
01074 
01075     if (cpl_frameset_count_tags(frameset, "MASTER_BIAS") == 0)
01076         fors_science_exit("Missing required input: MASTER_BIAS");
01077 
01078     if (cpl_frameset_count_tags(frameset, "MASTER_BIAS") > 1)
01079         fors_science_exit("Too many in input: MASTER_BIAS");
01080 
01081     if (skyalign >= 0)
01082         if (cpl_frameset_count_tags(frameset, "MASTER_SKYLINECAT") > 1)
01083             fors_science_exit("Too many in input: MASTER_SKYLINECAT");
01084 
01085     if (cpl_frameset_count_tags(frameset, disp_coeff_tag) == 0) {
01086         cpl_msg_error(recipe, "Missing required input: %s", disp_coeff_tag);
01087         fors_science_exit(NULL);
01088     }
01089 
01090     if (cpl_frameset_count_tags(frameset, disp_coeff_tag) > 1) {
01091         cpl_msg_error(recipe, "Too many in input: %s", disp_coeff_tag);
01092         fors_science_exit(NULL);
01093     }
01094 
01095     if (cpl_frameset_count_tags(frameset, slit_location_tag) == 0) {
01096         cpl_msg_error(recipe, "Missing required input: %s",
01097                       slit_location_tag);
01098         fors_science_exit(NULL);
01099     }
01100 
01101     if (cpl_frameset_count_tags(frameset, slit_location_tag) > 1) {
01102         cpl_msg_error(recipe, "Too many in input: %s", slit_location_tag);
01103         fors_science_exit(NULL);
01104     }
01105 
01106     if (cpl_frameset_count_tags(frameset, master_norm_flat_tag) > 1) {
01107         if (flatfield) {
01108             cpl_msg_error(recipe, "Too many in input: %s", 
01109                           master_norm_flat_tag);
01110             fors_science_exit(NULL);
01111         }
01112         else {
01113             cpl_msg_warning(recipe, "%s in input are ignored, "
01114                             "since flat field correction was not requested", 
01115                             master_norm_flat_tag);
01116         }
01117     }
01118 
01119     if (cpl_frameset_count_tags(frameset, master_norm_flat_tag) == 1) {
01120         if (!flatfield) {
01121             cpl_msg_warning(recipe, "%s in input is ignored, "
01122                             "since flat field correction was not requested", 
01123                             master_norm_flat_tag);
01124         }
01125     }
01126 
01127     if (cpl_frameset_count_tags(frameset, master_norm_flat_tag) == 0) {
01128         if (flatfield) {
01129             cpl_msg_error(recipe, "Flat field correction was requested, "
01130                           "but no %s are found in input",
01131                           master_norm_flat_tag);
01132             fors_science_exit(NULL);
01133         }
01134     }
01135 
01136     if (standard) {
01137     
01138         if (cpl_frameset_count_tags(frameset, "EXTINCT_TABLE") == 0) {
01139             cpl_msg_warning(recipe, "An EXTINCT_TABLE was not found in input: "
01140                             "instrument response curve will not be produced.");
01141             standard = 0;
01142         }
01143 
01144         if (cpl_frameset_count_tags(frameset, "EXTINCT_TABLE") > 1)
01145             fors_science_exit("Too many in input: EXTINCT_TABLE");
01146     
01147         if (cpl_frameset_count_tags(frameset, "STD_FLUX_TABLE") == 0) {
01148             cpl_msg_warning(recipe, "A STD_FLUX_TABLE was not found in input: "
01149                             "instrument response curve will not be produced.");
01150             standard = 0;
01151         }
01152 
01153         if (cpl_frameset_count_tags(frameset, "STD_FLUX_TABLE") > 1)
01154             fors_science_exit("Too many in input: STD_FLUX_TABLE");
01155 
01156         if (!dfs_equal_keyword(frameset, "ESO OBS TARG NAME")) {
01157             cpl_msg_warning(recipe, "The target name of observation does not "
01158                             "match the standard star catalog: "
01159                             "instrument response curve will not be produced.");
01160             standard = 0;
01161         }
01162     }
01163 
01164     if (photometry) {
01165         if (cpl_frameset_count_tags(frameset, "EXTINCT_TABLE") == 0) {
01166             fors_science_exit("An EXTINCT_TABLE was not found in input: "
01167                               "the requested photometric calibrated spectra "
01168                               "cannot be produced.");
01169         }
01170 
01171         if (cpl_frameset_count_tags(frameset, "EXTINCT_TABLE") > 1)
01172             fors_science_exit("Too many in input: EXTINCT_TABLE");
01173 
01174         have_phot  = cpl_frameset_count_tags(frameset, specphot_tag);
01175         have_phot += cpl_frameset_count_tags(frameset, master_specphot_tag);
01176 
01177         if (!standard) {
01178             if (have_phot == 0) {
01179                 cpl_msg_warning(recipe, 
01180                                 "A SPECPHOT_TABLE was not found in input: "
01181                                 "the requested photometric calibrated "
01182                                 "spectra cannot be produced.");
01183                 photometry = 0;
01184             }
01185 
01186             if (have_phot > 1)
01187                 fors_science_exit("Too many in input: SPECPHOT_TABLE");
01188         }
01189     }
01190 
01191     cpl_msg_indent_less();
01192 
01193 
01194     /*
01195      * Loading input data
01196      */
01197 
01198     exptime = cpl_calloc(nscience, sizeof(double));
01199 
01200     if (nscience > 1) {
01201 
01202         cpl_msg_info(recipe, "Load %d scientific frames and median them...",
01203                      nscience);
01204         cpl_msg_indent_more();
01205 
01206         all_science = cpl_imagelist_new();
01207 
01208         header = dfs_load_header(frameset, science_tag, 0);
01209 
01210         if (header == NULL)
01211             fors_science_exit("Cannot load scientific frame header");
01212 
01213         alltime = exptime[0] = cpl_propertylist_get_double(header, "EXPTIME");
01214 
01215         if (cpl_error_get_code() != CPL_ERROR_NONE)
01216             fors_science_exit("Missing keyword EXPTIME in scientific "
01217                               "frame header");
01218 
01219         if (standard || photometry) {
01220             airmass = fors_get_airmass(header);
01221             if (airmass < 0.0) 
01222                 fors_science_exit("Missing airmass information in "
01223                                   "scientific frame header");
01224         }
01225 
01226         cpl_propertylist_delete(header); header = NULL;
01227 
01228         cpl_msg_info(recipe, "Scientific frame 1 exposure time: %.2f s", 
01229                      exptime[0]);
01230 
01231         for (i = 1; i < nscience; i++) {
01232 
01233             header = dfs_load_header(frameset, NULL, 0);
01234 
01235             if (header == NULL)
01236                 fors_science_exit("Cannot load scientific frame header");
01237     
01238             exptime[i] = cpl_propertylist_get_double(header, "EXPTIME");
01239 
01240             alltime += exptime[i];
01241     
01242             if (cpl_error_get_code() != CPL_ERROR_NONE)
01243                 fors_science_exit("Missing keyword EXPTIME in scientific "
01244                                   "frame header");
01245     
01246             cpl_propertylist_delete(header); header = NULL;
01247 
01248             cpl_msg_info(recipe, "Scientific frame %d exposure time: %.2f s", 
01249                          i+1, exptime[i]);
01250         }
01251 
01252         spectra = dfs_load_image(frameset, science_tag, CPL_TYPE_FLOAT, 0, 0);
01253 
01254         if (spectra == NULL)
01255             fors_science_exit("Cannot load scientific frame");
01256 
01257         cpl_image_divide_scalar(spectra, exptime[0]);
01258         cpl_imagelist_set(all_science, spectra, 0); spectra = NULL;
01259 
01260         for (i = 1; i < nscience; i++) {
01261 
01262             spectra = dfs_load_image(frameset, NULL, CPL_TYPE_FLOAT, 0, 0);
01263 
01264             if (spectra) {
01265                 cpl_image_divide_scalar(spectra, exptime[i]);
01266                 cpl_imagelist_set(all_science, spectra, i); spectra = NULL;
01267             }
01268             else
01269                 fors_science_exit("Cannot load scientific frame");
01270 
01271         }
01272 
01273         spectra = cpl_imagelist_collapse_median_create(all_science);
01274         cpl_image_multiply_scalar(spectra, alltime);
01275 
01276         cpl_imagelist_delete(all_science);
01277     }
01278     else {
01279         cpl_msg_info(recipe, "Load scientific exposure...");
01280         cpl_msg_indent_more();
01281 
01282         header = dfs_load_header(frameset, science_tag, 0);
01283 
01284         if (header == NULL)
01285             fors_science_exit("Cannot load scientific frame header");
01286 
01287         if (standard || photometry) {
01288             airmass = fors_get_airmass(header);
01289             if (airmass < 0.0) 
01290                 fors_science_exit("Missing airmass information in "
01291                                   "scientific frame header");
01292         }
01293 
01294         /*
01295          * Insert here a check on supported filters:
01296          */
01297 
01298         wheel4 = (char *)cpl_propertylist_get_string(header, 
01299                                                      "ESO INS OPTI9 TYPE");
01300         if (cpl_error_get_code() != CPL_ERROR_NONE) {
01301             fors_science_exit("Missing ESO INS OPTI9 TYPE in flat header");
01302         }
01303 
01304         if (strcmp("FILT", wheel4) == 0) {
01305             wheel4 = (char *)cpl_propertylist_get_string(header,
01306                                                          "ESO INS OPTI9 NAME");
01307             cpl_msg_error(recipe, "Unsupported filter: %s", wheel4);
01308             fors_science_exit(NULL);
01309         }
01310 
01311         alltime = exptime[0] = cpl_propertylist_get_double(header, "EXPTIME");
01312 
01313         if (cpl_error_get_code() != CPL_ERROR_NONE)
01314             fors_science_exit("Missing keyword EXPTIME in scientific "
01315                               "frame header");
01316 
01317         cpl_propertylist_delete(header); header = NULL;
01318 
01319         cpl_msg_info(recipe, "Scientific frame exposure time: %.2f s", 
01320                      exptime[0]);
01321 
01322         spectra = dfs_load_image(frameset, science_tag, CPL_TYPE_FLOAT, 0, 0);
01323     }
01324 
01325     if (spectra == NULL)
01326         fors_science_exit("Cannot load scientific frame");
01327 
01328     cpl_free(exptime); exptime = NULL;
01329 
01330     cpl_msg_indent_less();
01331 
01332 
01333     /*
01334      * Get the reference wavelength and the rebin factor along the
01335      * dispersion direction from a scientific exposure
01336      */
01337 
01338     header = dfs_load_header(frameset, science_tag, 0);
01339 
01340     if (header == NULL)
01341         fors_science_exit("Cannot load scientific frame header");
01342 
01343     instrume = (char *)cpl_propertylist_get_string(header, "INSTRUME");
01344     if (instrume == NULL)
01345         fors_science_exit("Missing keyword INSTRUME in scientific header");
01346     instrume = cpl_strdup(instrume);
01347 
01348     if (instrume[4] == '1')
01349         snprintf(version, 80, "%s/%s", "fors1", VERSION);
01350     if (instrume[4] == '2')
01351         snprintf(version, 80, "%s/%s", "fors2", VERSION);
01352 
01353     reference = cpl_propertylist_get_double(header, "ESO INS GRIS1 WLEN");
01354 
01355     if (cpl_error_get_code() != CPL_ERROR_NONE)
01356         fors_science_exit("Missing keyword ESO INS GRIS1 WLEN in scientific "
01357                         "frame header");
01358 
01359     if (reference < 3000.0)   /* Perhaps in nanometers... */
01360         reference *= 10;
01361 
01362     if (reference < 3000.0 || reference > 13000.0) {
01363         cpl_msg_error(recipe, "Invalid central wavelength %.2f read from "
01364                       "keyword ESO INS GRIS1 WLEN in scientific frame header",
01365                       reference);
01366         fors_science_exit(NULL);
01367     }
01368 
01369     cpl_msg_info(recipe, "The central wavelength is: %.2f", reference);
01370 
01371     rebin = cpl_propertylist_get_int(header, "ESO DET WIN1 BINX");
01372 
01373     if (cpl_error_get_code() != CPL_ERROR_NONE)
01374         fors_science_exit("Missing keyword ESO DET WIN1 BINX in scientific "
01375                         "frame header");
01376 
01377     if (rebin != 1) {
01378         dispersion *= rebin;
01379         cpl_msg_warning(recipe, "The rebin factor is %d, and therefore the "
01380                         "resampling step used is %f A/pixel", rebin, 
01381                         dispersion);
01382         ext_radius /= rebin;
01383         cpl_msg_warning(recipe, "The rebin factor is %d, and therefore the "
01384                         "extraction radius used is %d pixel", rebin, 
01385                         ext_radius);
01386     }
01387 
01388     gain = cpl_propertylist_get_double(header, "ESO DET OUT1 CONAD");
01389 
01390     if (cpl_error_get_code() != CPL_ERROR_NONE)
01391         fors_science_exit("Missing keyword ESO DET OUT1 CONAD in scientific "
01392                           "frame header");
01393 
01394     cpl_msg_info(recipe, "The gain factor is: %.2f e-/ADU", gain);
01395 
01396     ron = cpl_propertylist_get_double(header, "ESO DET OUT1 RON");
01397 
01398     if (cpl_error_get_code() != CPL_ERROR_NONE)
01399         fors_science_exit("Missing keyword ESO DET OUT1 RON in scientific "
01400                           "frame header");
01401 
01402     ron /= gain;     /* Convert from electrons to ADU */
01403 
01404     cpl_msg_info(recipe, "The read-out-noise is: %.2f ADU", ron);
01405 
01406     if (mos || mxu) {
01407 /* goto skip; */
01408         if (mos)
01409             maskslits = mos_load_slits_fors_mos(header);
01410         else
01411             maskslits = mos_load_slits_fors_mxu(header);
01412 
01413         /*
01414          * Check if all slits have the same X offset: in such case,
01415          * treat the observation as a long-slit one!
01416          */
01417 
01418         mxpos = cpl_table_get_column_median(maskslits, "xtop");
01419         xpos = cpl_table_get_data_double(maskslits, "xtop");
01420         nslits = cpl_table_get_nrow(maskslits);
01421      
01422         treat_as_lss = 1;
01423         for (i = 0; i < nslits; i++) { 
01424             if (fabs(mxpos-xpos[i]) > 0.01) {
01425                 treat_as_lss = 0;
01426                 break;
01427             }
01428         }
01429 
01430         cpl_table_delete(maskslits); maskslits = NULL;
01431 /* skip: */
01432 
01433         if (treat_as_lss) {
01434             cpl_msg_warning(recipe, "All MOS slits have the same offset: %.2f\n"
01435                             "The LSS data reduction strategy is applied!",
01436                             mxpos);
01437             if (mos) {
01438                 if (standard) {
01439                     skylines_offsets_tag   = "SKY_SHIFTS_LONG_STD_MOS";
01440                 }
01441                 else {
01442                     skylines_offsets_tag   = "SKY_SHIFTS_LONG_SCI_MOS";
01443                 }
01444             }
01445             else {
01446                 if (standard) {
01447                     skylines_offsets_tag   = "SKY_SHIFTS_LONG_STD_MXU";
01448                 }
01449                 else {
01450                     skylines_offsets_tag   = "SKY_SHIFTS_LONG_SCI_MXU";
01451                 }
01452             }
01453         }
01454     }
01455 
01456     if (lss || treat_as_lss) {
01457         if (skylocal) {
01458             if (cosmics)
01459                 fors_science_exit("Cosmic rays correction for LSS or LSS-like "
01460                                   "data requires --skyglobal=true");
01461             skymedian = skylocal;
01462             skylocal = 0;
01463         }
01464     }
01465     else {
01466         if (cpl_frameset_count_tags(frameset, curv_coeff_tag) == 0) {
01467             cpl_msg_error(recipe, "Missing required input: %s", curv_coeff_tag);
01468             fors_science_exit(NULL);
01469         }
01470 
01471         if (cpl_frameset_count_tags(frameset, curv_coeff_tag) > 1) {
01472             cpl_msg_error(recipe, "Too many in input: %s", curv_coeff_tag);
01473             fors_science_exit(NULL);
01474         }
01475     }
01476 
01477     /* Leave the header on for the next step... */
01478 
01479 
01480     /*
01481      * Remove the master bias
01482      */
01483 
01484     cpl_msg_info(recipe, "Remove the master bias...");
01485 
01486     bias = dfs_load_image(frameset, "MASTER_BIAS", CPL_TYPE_FLOAT, 0, 1);
01487 
01488     if (bias == NULL)
01489         fors_science_exit("Cannot load master bias");
01490 
01491     overscans = mos_load_overscans_vimos(header, 1);
01492     cpl_propertylist_delete(header); header = NULL;
01493     dummy = mos_remove_bias(spectra, bias, overscans);
01494     cpl_image_delete(spectra); spectra = dummy; dummy = NULL;
01495     cpl_image_delete(bias); bias = NULL;
01496     cpl_table_delete(overscans); overscans = NULL;
01497 
01498     if (spectra == NULL)
01499         fors_science_exit("Cannot remove bias from scientific frame");
01500 
01501     ccd_xsize = nx = cpl_image_get_size_x(spectra);
01502     ccd_ysize = ny = cpl_image_get_size_y(spectra);
01503 
01504     cpl_msg_indent_less();
01505     cpl_msg_info(recipe, "Load normalised flat field (if present)...");
01506     cpl_msg_indent_more();
01507 
01508     if (flatfield) {
01509 
01510         norm_flat = dfs_load_image(frameset, master_norm_flat_tag, 
01511                                    CPL_TYPE_FLOAT, 0, 1);
01512 
01513         if (norm_flat) {
01514             cpl_msg_info(recipe, "Apply flat field correction...");
01515             if (cpl_image_divide(spectra, norm_flat) != CPL_ERROR_NONE) {
01516                 cpl_msg_error(recipe, "Failure of flat field correction: %s",
01517                               cpl_error_get_message());
01518                 fors_science_exit(NULL);
01519             }
01520             cpl_image_delete(norm_flat); norm_flat = NULL;
01521         }
01522         else {
01523             cpl_msg_error(recipe, "Cannot load input %s for flat field "
01524                           "correction", master_norm_flat_tag);
01525             fors_science_exit(NULL);
01526         }
01527 
01528     }
01529 
01530 
01531     if (skyalign >= 0) {
01532         cpl_msg_indent_less();
01533         cpl_msg_info(recipe, "Load input sky line catalog...");
01534         cpl_msg_indent_more();
01535 
01536         wavelengths = dfs_load_table(frameset, "MASTER_SKYLINECAT", 1);
01537 
01538         if (wavelengths) {
01539 
01540             /*
01541              * Cast the wavelengths into a (double precision) CPL vector
01542              */
01543 
01544             nlines = cpl_table_get_nrow(wavelengths);
01545 
01546             if (nlines == 0)
01547                 fors_science_exit("Empty input sky line catalog");
01548 
01549             if (cpl_table_has_column(wavelengths, wcolumn) != 1) {
01550                 cpl_msg_error(recipe, "Missing column %s in input line "
01551                               "catalog table", wcolumn);
01552                 fors_science_exit(NULL);
01553             }
01554 
01555             line = cpl_malloc(nlines * sizeof(double));
01556     
01557             for (i = 0; i < nlines; i++)
01558                 line[i] = cpl_table_get(wavelengths, wcolumn, i, NULL);
01559 
01560             cpl_table_delete(wavelengths); wavelengths = NULL;
01561 
01562             lines = cpl_vector_wrap(nlines, line);
01563         }
01564         else {
01565             cpl_msg_info(recipe, "No sky line catalog found in input - fine!");
01566         }
01567     }
01568 
01569 
01570     /*
01571      * Load the spectral curvature table, or provide a dummy one
01572      * in case of LSS or LSS-like data (single slit with flat curvature)
01573      */
01574 
01575     if (!(lss || treat_as_lss)) {
01576         polytraces = dfs_load_table(frameset, curv_coeff_tag, 1);
01577         if (polytraces == NULL)
01578             fors_science_exit("Cannot load spectral curvature table");
01579     }
01580 
01581 
01582     /*
01583      * Load the slit location table
01584      */
01585 
01586     slits = dfs_load_table(frameset, slit_location_tag, 1);
01587     if (slits == NULL)
01588         fors_science_exit("Cannot load slits location table");
01589 
01590     if (lss || treat_as_lss) {
01591         int first_row = cpl_table_get_double(slits, "ybottom", 0, NULL);
01592         int last_row = cpl_table_get_double(slits, "ytop", 0, NULL);
01593         int ylow, yhig;
01594 
01595         ylow = first_row + 1;
01596         yhig = last_row + 1;
01597 
01598         dummy = cpl_image_extract(spectra, 1, ylow, nx, yhig);
01599         cpl_image_delete(spectra); spectra = dummy; dummy = NULL;
01600         ny = cpl_image_get_size_y(spectra);
01601     }
01602 
01603 
01604     /*
01605      * Load the wavelength calibration table
01606      */
01607 
01608     idscoeff = dfs_load_table(frameset, disp_coeff_tag, 1);
01609     if (idscoeff == NULL)
01610         fors_science_exit("Cannot load wavelength calibration table");
01611 
01612     cpl_msg_indent_less();
01613     cpl_msg_info(recipe, "Processing scientific spectra...");
01614     cpl_msg_indent_more();
01615 
01616     /*
01617      * This one will also generate the spatial map from the spectral 
01618      * curvature table (in the case of multislit data)
01619      */
01620 
01621     if (lss || treat_as_lss) {
01622         smapped = cpl_image_duplicate(spectra);
01623     }
01624     else {
01625         coordinate = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
01626 
01627         smapped = mos_spatial_calibration(spectra, slits, polytraces, reference,
01628                                           startwavelength, endwavelength,
01629                                           dispersion, flux, coordinate);
01630     }
01631 
01632 
01633     /*
01634      * Generate a rectified wavelength map from the wavelength calibration 
01635      * table
01636      */
01637 
01638     rainbow = mos_map_idscoeff(idscoeff, nx, reference, startwavelength, 
01639                                endwavelength);
01640 
01641     if (dispersion > 1.0)
01642         highres = 0;
01643     else
01644         highres = 1;
01645 
01646     if (skyalign >= 0) {
01647         if (skyalign) {
01648             cpl_msg_info(recipe, "Align wavelength solution to reference "
01649             "skylines applying %d order residual fit...", skyalign);
01650         }
01651         else {
01652             cpl_msg_info(recipe, "Align wavelength solution to reference "
01653             "skylines applying median offset...");
01654         }
01655 
01656         if (lss || treat_as_lss) {
01657             offsets = mos_wavelength_align_lss(smapped, reference, 
01658                                                startwavelength, endwavelength, 
01659                                                idscoeff, lines, highres, 
01660                                                skyalign, rainbow, 4);
01661         }
01662         else {
01663             offsets = mos_wavelength_align(smapped, slits, reference, 
01664                                            startwavelength, endwavelength, 
01665                                            idscoeff, lines, highres, skyalign, 
01666                                            rainbow, 4);
01667         }
01668 
01669         if (offsets) {
01670             if (standard)
01671                 cpl_msg_warning(recipe, "Alignment of the wavelength solution "
01672                                 "to reference sky lines may be unreliable in "
01673                                 "this case!");
01674 
01675             if (dfs_save_table(frameset, offsets, skylines_offsets_tag, NULL, 
01676                                parlist, recipe, version))
01677                 fors_science_exit(NULL);
01678 
01679             cpl_table_delete(offsets); offsets = NULL;
01680         }
01681         else {
01682             if (cpl_error_get_code()) {
01683                 if (cpl_error_get_code() == CPL_ERROR_INCOMPATIBLE_INPUT) {
01684                     cpl_msg_error(recipe, "The IDS coeff table is "
01685                     "incompatible with the input slit position table.");
01686                 }
01687                 cpl_msg_error(cpl_error_get_where(), cpl_error_get_message());
01688                 fors_science_exit(NULL);
01689             }
01690             cpl_msg_warning(recipe, "Alignment of the wavelength solution "
01691                             "to reference sky lines could not be done!");
01692             skyalign = -1;
01693         }
01694 
01695     }
01696 
01697     if (lss || treat_as_lss) {
01698         int first_row = cpl_table_get_double(slits, "ybottom", 0, NULL);
01699         int last_row = cpl_table_get_double(slits, "ytop", 0, NULL);
01700         int ylow, yhig;
01701 
01702         ylow = first_row + 1;
01703         yhig = last_row + 1;
01704 
01705         wavemap = cpl_image_new(ccd_xsize, ccd_ysize, CPL_TYPE_FLOAT);
01706         cpl_image_copy(wavemap, rainbow, 1, ylow);
01707 
01708         wavemaplss = cpl_image_extract(wavemap, 1, ylow, nx, yhig);
01709     }
01710     else {
01711         wavemap = mos_map_wavelengths(coordinate, rainbow, slits, 
01712                                       polytraces, reference, 
01713                                       startwavelength, endwavelength,
01714                                       dispersion);
01715     }
01716 
01717     cpl_image_delete(rainbow); rainbow = NULL;
01718     cpl_image_delete(coordinate); coordinate = NULL;
01719 
01720     /*
01721      * Here the wavelength calibrated slit spectra are created. This frame
01722      * contains sky_science.
01723      */
01724 
01725     mapped_sky = mos_wavelength_calibration(smapped, reference,
01726                                             startwavelength, endwavelength,
01727                                             dispersion, idscoeff, flux);
01728 
01729     cpl_msg_indent_less();
01730     cpl_msg_info(recipe, "Check applied wavelength against skylines...");
01731     cpl_msg_indent_more();
01732 
01733     mean_rms = mos_distortions_rms(mapped_sky, lines, startwavelength,
01734                                    dispersion, 6, highres);
01735 
01736     cpl_vector_delete(lines); lines = NULL;
01737 
01738     cpl_msg_info(recipe, "Mean residual: %f", mean_rms);
01739 
01740     mean_rms = cpl_table_get_column_mean(idscoeff, "error");
01741 
01742     cpl_msg_info(recipe, "Mean model accuracy: %f pixel (%f A)",
01743                  mean_rms, mean_rms * dispersion);
01744 
01745     header = cpl_propertylist_new();
01746     cpl_propertylist_update_double(header, "CRPIX1", 1.0);
01747     cpl_propertylist_update_double(header, "CRPIX2", 1.0);
01748     cpl_propertylist_update_double(header, "CRVAL1", 
01749                                    startwavelength + dispersion/2);
01750     cpl_propertylist_update_double(header, "CRVAL2", 1.0);
01751     /* cpl_propertylist_update_double(header, "CDELT1", dispersion);
01752     cpl_propertylist_update_double(header, "CDELT2", 1.0); */
01753     cpl_propertylist_update_double(header, "CD1_1", dispersion);
01754     cpl_propertylist_update_double(header, "CD1_2", 0.0);
01755     cpl_propertylist_update_double(header, "CD2_1", 0.0);
01756     cpl_propertylist_update_double(header, "CD2_2", 1.0);
01757     cpl_propertylist_update_string(header, "CTYPE1", "LINEAR");
01758     cpl_propertylist_update_string(header, "CTYPE2", "PIXEL");
01759 
01760     if (time_normalise) {
01761         cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
01762         dummy = cpl_image_divide_scalar_create(mapped_sky, alltime);
01763         if (dfs_save_image(frameset, dummy, mapped_science_sky_tag, header, 
01764                            parlist, recipe, version))
01765             fors_science_exit(NULL);
01766         cpl_image_delete(dummy); dummy = NULL;
01767     }
01768     else {
01769         cpl_propertylist_update_string(header, "BUNIT", "ADU");
01770         if (dfs_save_image(frameset, mapped_sky, mapped_science_sky_tag, 
01771                            header, parlist, recipe, version))
01772             fors_science_exit(NULL);
01773     }
01774 
01775 /*    if (skyglobal == 0 && skymedian < 0) {    NSS */
01776     if (skyglobal == 0 && skymedian == 0 && skylocal == 0) {
01777         cpl_image_delete(mapped_sky); mapped_sky = NULL;
01778     }
01779 
01780     if (skyglobal || skylocal) {
01781 
01782         cpl_msg_indent_less();
01783 
01784         if (skyglobal) {
01785             cpl_msg_info(recipe, "Global sky determination...");
01786             cpl_msg_indent_more();
01787             skymap = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
01788             if (lss || treat_as_lss) {
01789                 sky = mos_sky_map_super(spectra, wavemaplss, dispersion, 
01790                                         2.0, 50, skymap);
01791                 cpl_image_delete(wavemaplss); wavemaplss = NULL;
01792             }
01793             else {
01794                 sky = mos_sky_map_super(spectra, wavemap, dispersion, 
01795                                         2.0, 50, skymap);
01796             }
01797             if (sky) {
01798                 cpl_image_subtract(spectra, skymap);
01799             }
01800             else {
01801                 cpl_image_delete(skymap); skymap = NULL;
01802             }
01803         }
01804         else {
01805             cpl_msg_info(recipe, "Local sky determination...");
01806             cpl_msg_indent_more();
01807             skymap = mos_subtract_sky(spectra, slits, polytraces, reference,
01808                            startwavelength, endwavelength, dispersion);
01809         }
01810 
01811         if (skymap) {
01812             if (skyglobal) {
01813                 if (time_normalise)
01814                     cpl_table_divide_scalar(sky, "sky", alltime);
01815                 if (dfs_save_table(frameset, sky, global_sky_spectrum_tag, 
01816                                    NULL, parlist, recipe, version))
01817                     fors_science_exit(NULL);
01818     
01819                 cpl_table_delete(sky); sky = NULL;
01820             }
01821 
01822             save_header = dfs_load_header(frameset, science_tag, 0);
01823 
01824             if (time_normalise)
01825                 cpl_image_divide_scalar(skymap, alltime);
01826             if (dfs_save_image(frameset, skymap, unmapped_sky_tag,
01827                                save_header, parlist, recipe, version))
01828                 fors_science_exit(NULL);
01829 
01830             cpl_image_delete(skymap); skymap = NULL;
01831 
01832             if (dfs_save_image(frameset, spectra, unmapped_science_tag,
01833                                save_header, parlist, recipe, version))
01834                 fors_science_exit(NULL);
01835 
01836             cpl_propertylist_delete(save_header); save_header = NULL;
01837 
01838             if (cosmics) {
01839                 cpl_msg_info(recipe, "Removing cosmic rays...");
01840                 mos_clean_cosmics(spectra, gain, -1., -1.);
01841             }
01842 
01843             /*
01844              * The spatially rectified image, that contained the sky,
01845              * is replaced by a sky-subtracted spatially rectified image:
01846              */
01847 
01848             cpl_image_delete(smapped); smapped = NULL;
01849 
01850             if (lss || treat_as_lss) {
01851                 smapped = cpl_image_duplicate(spectra);
01852             }
01853             else {
01854                 smapped = mos_spatial_calibration(spectra, slits, polytraces, 
01855                                                   reference, startwavelength, 
01856                                                   endwavelength, dispersion, 
01857                                                   flux, NULL);
01858             }
01859         }
01860         else {
01861             cpl_msg_warning(recipe, "Sky subtraction failure");
01862             if (cosmics)
01863                 cpl_msg_warning(recipe, "Cosmic rays removal not performed!");
01864             cosmics = skylocal = skyglobal = 0;
01865         }
01866     }
01867 
01868     cpl_image_delete(spectra); spectra = NULL;
01869     cpl_table_delete(polytraces); polytraces = NULL;
01870 
01871     if (skyalign >= 0) {
01872         save_header = dfs_load_header(frameset, science_tag, 0);
01873         cpl_propertylist_update_string(save_header, "BUNIT", "Angstrom");
01874         if (dfs_save_image(frameset, wavemap, wavelength_map_sky_tag,
01875                            save_header, parlist, recipe, version))
01876             fors_science_exit(NULL);
01877         cpl_propertylist_delete(save_header); save_header = NULL;
01878     }
01879 
01880     cpl_image_delete(wavemap); wavemap = NULL;
01881 
01882     mapped = mos_wavelength_calibration(smapped, reference,
01883                                         startwavelength, endwavelength,
01884                                         dispersion, idscoeff, flux);
01885 
01886     cpl_image_delete(smapped); smapped = NULL;
01887 
01888 /*    if (skymedian >= 0) {    NSS */
01889     if (skymedian) {
01890             cpl_msg_indent_less();
01891             cpl_msg_info(recipe, "Local sky determination...");
01892             cpl_msg_indent_more();
01893        
01894 /*   NSS      skylocalmap = mos_sky_local(mapped, slits, skymedian); */
01895 /*            skylocalmap = mos_sky_local(mapped, slits, 0);        */
01896             skylocalmap = mos_sky_local_old(mapped, slits);       
01897             cpl_image_subtract(mapped, skylocalmap);
01898 /*
01899             if (dfs_save_image(frameset, skylocalmap, mapped_sky_tag, header, 
01900                                parlist, recipe, version))
01901                 fors_science_exit(NULL);
01902 */
01903             cpl_image_delete(skylocalmap); skylocalmap = NULL;
01904     }
01905 
01906 /*    if (skyglobal || skymedian >= 0 || skylocal) {   NSS */
01907     if (skyglobal || skymedian || skylocal) {
01908 
01909         skylocalmap = cpl_image_subtract_create(mapped_sky, mapped);
01910 
01911         cpl_image_delete(mapped_sky); mapped_sky = NULL;
01912 
01913         if (time_normalise) {
01914             cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
01915             dummy = cpl_image_divide_scalar_create(skylocalmap, alltime);
01916             if (dfs_save_image(frameset, dummy, mapped_sky_tag, header,
01917                                parlist, recipe, version))
01918                 fors_science_exit(NULL);
01919             cpl_image_delete(dummy); dummy = NULL;
01920         }
01921         else {
01922             cpl_propertylist_update_string(header, "BUNIT", "ADU");
01923             if (dfs_save_image(frameset, skylocalmap, mapped_sky_tag, header,
01924                                parlist, recipe, version))
01925                 fors_science_exit(NULL);
01926         }
01927 
01928         cpl_msg_indent_less();
01929         cpl_msg_info(recipe, "Object detection...");
01930         cpl_msg_indent_more();
01931 
01932         if (cosmics || nscience > 1) {
01933             dummy = mos_detect_objects(mapped, slits, slit_margin, ext_radius, 
01934                                        cont_radius);
01935         }
01936         else {
01937             mapped_cleaned = cpl_image_duplicate(mapped);
01938             mos_clean_cosmics(mapped_cleaned, gain, -1., -1.);
01939             dummy = mos_detect_objects(mapped_cleaned, slits, slit_margin, 
01940                                        ext_radius, cont_radius);
01941 
01942             cpl_image_delete(mapped_cleaned); mapped_cleaned = NULL;
01943         }
01944 
01945         cpl_image_delete(dummy); dummy = NULL;
01946 
01947         if (dfs_save_table(frameset, slits, object_table_tag, NULL, parlist, 
01948                            recipe, version))
01949             fors_science_exit(NULL);
01950 
01951         cpl_msg_indent_less();
01952         cpl_msg_info(recipe, "Object extraction...");
01953         cpl_msg_indent_more();
01954 
01955         images = mos_extract_objects(mapped, skylocalmap, slits, 
01956                                      ext_mode, ron, gain, 1);
01957 
01958         cpl_image_delete(skylocalmap); skylocalmap = NULL;
01959 
01960         if (images) {
01961             if (standard) {
01962                 cpl_table *ext_table  = NULL;
01963                 cpl_table *flux_table = NULL;
01964 
01965                 ext_table = dfs_load_table(frameset, "EXTINCT_TABLE", 1);
01966                 flux_table = dfs_load_table(frameset, "STD_FLUX_TABLE", 1);
01967 
01968                 photcal = mos_photometric_calibration(images[0], 
01969                                                       startwavelength,
01970                                                       dispersion, gain,
01971                                                       alltime, ext_table,
01972                                                       airmass, flux_table,
01973                                                       res_order);
01974 
01975                 cpl_table_delete(ext_table);
01976                 cpl_table_delete(flux_table);
01977 
01978                 if (photcal) {
01979 
01980                     if (qc) {
01981 
01982                         float *data;
01983                         char  *pipefile = NULL;
01984                         char   keyname[30];
01985 
01986                         qclist = dfs_load_header(frameset, science_tag, 0);
01987 
01988                         if (qclist == NULL)
01989                             fors_science_exit("Cannot reload scientific "
01990                                               "frame header");
01991 
01992                         fors_qc_start_group(qclist, "2.0", instrume);
01993 
01994 
01995                         /*
01996                          * QC1 group header
01997                          */
01998 
01999                         if (fors_qc_write_string("PRO.CATG", specphot_tag,
02000                                                  "Product category", instrume))
02001                             fors_science_exit("Cannot write product category "
02002                                               "to QC log file");
02003 
02004                         if (fors_qc_keyword_to_paf(qclist, "ESO DPR TYPE", NULL,
02005                                                    "DPR type", instrume))
02006                             fors_science_exit("Missing keyword DPR TYPE in "
02007                                               "scientific frame header");
02008 
02009                         if (fors_qc_keyword_to_paf(qclist, "ESO TPL ID", NULL,
02010                                                    "Template", instrume))
02011                             fors_science_exit("Missing keyword TPL ID in "
02012                                               "scientific frame header");
02013 
02014                         if (fors_qc_keyword_to_paf(qclist, 
02015                                                    "ESO INS GRIS1 NAME", NULL,
02016                                                    "Grism name", instrume))
02017                             fors_science_exit("Missing keyword INS GRIS1 NAME "
02018                                               "in scientific frame header");
02019 
02020                         if (fors_qc_keyword_to_paf(qclist, 
02021                                                    "ESO INS GRIS1 ID", NULL,
02022                                                    "Grism identifier", 
02023                                                    instrume))
02024                             fors_science_exit("Missing keyword INS GRIS1 ID "
02025                                               "in scientific frame header");
02026 
02027                         if (cpl_propertylist_has(qclist, "ESO INS FILT1 NAME"))
02028                             fors_qc_keyword_to_paf(qclist, 
02029                                                    "ESO INS FILT1 NAME", NULL,
02030                                                    "Filter name", instrume);
02031 
02032                         if (fors_qc_keyword_to_paf(qclist, 
02033                                                    "ESO INS COLL NAME", NULL,
02034                                                    "Collimator name", 
02035                                                    instrume))
02036                             fors_science_exit("Missing keyword INS COLL NAME "
02037                                               "in scientific frame header");
02038 
02039                         if (fors_qc_keyword_to_paf(qclist, 
02040                                                    "ESO DET CHIP1 ID", NULL,
02041                                                    "Chip identifier", 
02042                                                    instrume))
02043                             fors_science_exit("Missing keyword DET CHIP1 ID "
02044                                               "in scientific frame header");
02045 
02046                         if (fors_qc_keyword_to_paf(qclist, 
02047                                                    "ESO INS MOS10 WID",
02048                                                    "arcsec", "Slit width", 
02049                                                    instrume)) {
02050                             cpl_error_reset();
02051                             if (fors_qc_keyword_to_paf(qclist, 
02052                                                        "ESO INS SLIT WID",
02053                                                        "arcsec", "Slit width", 
02054                                                        instrume)) {
02055                                 if (mos) {
02056                                     fors_science_exit("Missing keyword "
02057                                                   "ESO INS MOS10 WID in "
02058                                                   "scientific frame header");
02059                                 }
02060                                 else {
02061                                     fors_science_exit("Missing keyword "
02062                                                   "ESO INS SLIT WID in "
02063                                                   "scientific frame header");
02064                                 }
02065                             }
02066                         }
02067 
02068 /*
02069                         if (mos) {
02070                             if (fors_qc_keyword_to_paf(qclist, 
02071                                                        "ESO INS MOS10 WID",
02072                                                        "arcsec", "Slit width", 
02073                                                        instrume))
02074                                 fors_science_exit("Missing keyword "
02075                                                   "ESO INS MOS10 WID in "
02076                                                   "scientific frame header");
02077                         }
02078                         else {
02079                             if (fors_qc_keyword_to_paf(qclist, 
02080                                                        "ESO INS SLIT WID",
02081                                                        "arcsec", "Slit width", 
02082                                                        instrume))
02083                                 fors_science_exit("Missing keyword "
02084                                                   "ESO INS SLIT WID in "
02085                                                   "scientific frame header");
02086                         }
02087 */
02088 
02089                         if (fors_qc_keyword_to_paf(qclist, 
02090                                                    "ESO DET WIN1 BINX", NULL,
02091                                                    "Binning factor along X", 
02092                                                    instrume))
02093                             fors_science_exit("Missing keyword ESO "
02094                                               "DET WIN1 BINX "
02095                                               "in scientific frame header");
02096         
02097                         if (fors_qc_keyword_to_paf(qclist, 
02098                                                    "ESO DET WIN1 BINY", NULL,
02099                                                    "Binning factor along Y", 
02100                                                    instrume))
02101                             fors_science_exit("Missing keyword "
02102                                               "ESO DET WIN1 BINY "
02103                                               "in scientific frame header");
02104 
02105                         if (fors_qc_keyword_to_paf(qclist, "ARCFILE", NULL,
02106                                                    "Archive name of input data",
02107                                                    instrume))
02108                             fors_science_exit("Missing keyword ARCFILE in "
02109                                               "scientific frame header");
02110         
02111                         pipefile = dfs_generate_filename_tfits(specphot_tag);
02112                         if (fors_qc_write_string("PIPEFILE", pipefile,
02113                                                  "Pipeline product name", 
02114                                                  instrume))
02115                             fors_science_exit("Cannot write PIPEFILE to "
02116                                               "QC log file");
02117                         cpl_free(pipefile);
02118 
02119 
02120                         /*
02121                          * QC1 parameters
02122                          */
02123 
02124                         wstart = 3700.;
02125                         wstep  = 400.;
02126                         wcount = 15;
02127 
02128                         dummy = cpl_image_new(wcount, 1, CPL_TYPE_FLOAT);
02129                         data = cpl_image_get_data_float(dummy);
02130                         map_table(dummy, wstart, wstep, photcal, 
02131                                   "WAVE", "EFFICIENCY");
02132 
02133                         for (i = 0; i < wcount; i++) {
02134                             sprintf(keyname, "QC.SPEC.EFFICIENCY%d.LAMBDA", 
02135                                     i + 1);
02136                             if (fors_qc_write_qc_double(qclist, 
02137                                                         wstart + wstep * i,
02138                                                         keyname, "Angstrom",
02139                                                         "Wavelength of "
02140                                                         "efficiency evaluation",
02141                                                         instrume)) {
02142                                 fors_science_exit("Cannot write wavelength of "
02143                                                   "efficiency evaluation");
02144                             }
02145 
02146                             sprintf(keyname, "QC.SPEC.EFFICIENCY%d", i + 1);
02147                             if (fors_qc_write_qc_double(qclist,
02148                                                         data[i],
02149                                                         keyname, "e-/photon",
02150                                                         "Efficiency",
02151                                                         instrume)) {
02152                                 fors_science_exit("Cannot write wavelength of "
02153                                                   "efficiency evaluation");
02154                             }
02155                         }
02156 
02157                         cpl_image_delete(dummy); dummy = NULL;
02158 
02159                         fors_qc_end_group();
02160         
02161                     }  /* End of QC1 computation */
02162 
02163                     if (dfs_save_table(frameset, photcal, specphot_tag, qclist,
02164                                        parlist, recipe, version))
02165                         fors_science_exit(NULL);
02166 
02167                     cpl_propertylist_delete(qclist); qclist = NULL;
02168 
02169                     if (have_phot) {
02170                         cpl_table_delete(photcal); photcal = NULL;
02171                     }
02172                 }
02173             }
02174 
02175             if (photometry) {
02176                 cpl_image *calibrated;
02177                 cpl_table *ext_table;
02178 
02179                 ext_table = dfs_load_table(frameset, "EXTINCT_TABLE", 1);
02180 
02181                 if (have_phot) {
02182                     if (cpl_frameset_count_tags(frameset, specphot_tag) == 0) {
02183                         photcal = dfs_load_table(frameset, 
02184                                                  master_specphot_tag, 1);
02185                     }
02186                     else {
02187                         photcal = dfs_load_table(frameset, specphot_tag, 1);
02188                     }
02189                 }
02190 
02191                 calibrated = mos_apply_photometry(images[0], photcal, 
02192                                                   ext_table, startwavelength, 
02193                                                   dispersion, gain, alltime, 
02194                                                   airmass);
02195                 cpl_propertylist_update_string(header, "BUNIT", 
02196                                    "10^(-16) erg/(cm^2 s Angstrom)");
02197 
02198                 if (dfs_save_image(frameset, calibrated,
02199                                    reduced_flux_science_tag, header,
02200                                    parlist, recipe, version)) {
02201                     cpl_image_delete(calibrated);
02202                     fors_science_exit(NULL);
02203                 }
02204 
02205                 cpl_table_delete(ext_table);
02206                 cpl_image_delete(calibrated);
02207             }
02208 
02209             if (time_normalise) {
02210                 cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
02211                 cpl_image_divide_scalar(images[0], alltime);
02212             }
02213             else {
02214                 cpl_propertylist_update_string(header, "BUNIT", "ADU");
02215             }
02216 
02217             if (dfs_save_image(frameset, images[0], reduced_science_tag, header,
02218                                parlist, recipe, version))
02219                 fors_science_exit(NULL);
02220 
02221             if (time_normalise)
02222                 cpl_image_divide_scalar(images[1], alltime);
02223 
02224             if (dfs_save_image(frameset, images[1], reduced_sky_tag, header,
02225                                parlist, recipe, version))
02226                 fors_science_exit(NULL);
02227 
02228             cpl_image_delete(images[1]);
02229 
02230             if (photometry) {
02231                 cpl_image *calibrated;
02232                 cpl_table *ext_table; 
02233  
02234                 ext_table = dfs_load_table(frameset, "EXTINCT_TABLE", 1);
02235 
02236                 calibrated = mos_propagate_photometry_error(images[0], 
02237                                                   images[2], photcal,
02238                                                   ext_table, startwavelength,
02239                                                   dispersion, gain, alltime,
02240                                                   airmass);
02241 
02242                 cpl_propertylist_update_string(header, "BUNIT", 
02243                                    "10^(-16) erg/(cm^2 s Angstrom)");
02244 
02245                 if (dfs_save_image(frameset, calibrated,
02246                                    reduced_flux_error_tag, header,
02247                                    parlist, recipe, version)) {
02248                     cpl_image_delete(calibrated);
02249                     fors_science_exit(NULL);
02250                 }
02251 
02252                 cpl_table_delete(ext_table);
02253                 cpl_image_delete(calibrated);
02254             }
02255 
02256     
02257             if (time_normalise) {
02258                 cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
02259                 cpl_image_divide_scalar(images[2], alltime);
02260             }
02261             else {
02262                 cpl_propertylist_update_string(header, "BUNIT", "ADU");
02263             }
02264 
02265             if (dfs_save_image(frameset, images[2], reduced_error_tag, header,
02266                                parlist, recipe, version))
02267                 fors_science_exit(NULL);
02268 
02269             cpl_image_delete(images[0]);
02270             cpl_image_delete(images[2]);
02271 
02272             cpl_free(images);
02273         }
02274         else {
02275             cpl_msg_warning(recipe, "No objects found: the products "
02276                             "%s, %s, and %s are not created", 
02277                             reduced_science_tag, reduced_sky_tag, 
02278                             reduced_error_tag);
02279         }
02280 
02281     }
02282 
02283     cpl_free(instrume); instrume = NULL;
02284     cpl_table_delete(slits); slits = NULL;
02285 
02286     if (skyalign >= 0) {
02287         if (dfs_save_table(frameset, idscoeff, disp_coeff_sky_tag, NULL, 
02288                            parlist, recipe, version))
02289             fors_science_exit(NULL);
02290     }
02291 
02292     cpl_table_delete(idscoeff); idscoeff = NULL;
02293 
02294 /*    if (skyglobal || skymedian >= 0) {   NSS */
02295     if (skyglobal || skymedian || skylocal) {
02296 
02297         if (photometry) {
02298             cpl_image *calibrated;
02299             cpl_table *ext_table; 
02300 
02301             ext_table = dfs_load_table(frameset, "EXTINCT_TABLE", 1);
02302 
02303             calibrated = mos_apply_photometry(mapped, photcal,
02304                                               ext_table, startwavelength,
02305                                               dispersion, gain, alltime,
02306                                               airmass);
02307 
02308             cpl_propertylist_update_string(header, "BUNIT", 
02309                                "10^(-16) erg/(cm^2 s Angstrom)");
02310 
02311             if (dfs_save_image(frameset, calibrated,
02312                                mapped_flux_science_tag, header,
02313                                parlist, recipe, version)) {
02314                 cpl_image_delete(calibrated);
02315                 fors_science_exit(NULL);
02316             }
02317 
02318             cpl_table_delete(ext_table);
02319             cpl_image_delete(calibrated);
02320         }
02321 
02322         if (time_normalise) {
02323             cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
02324             cpl_image_divide_scalar(mapped, alltime);
02325         }
02326         else {
02327             cpl_propertylist_update_string(header, "BUNIT", "ADU");
02328         }
02329 
02330         if (dfs_save_image(frameset, mapped, mapped_science_tag, header, 
02331                            parlist, recipe, version))
02332             fors_science_exit(NULL);
02333     }
02334 
02335     cpl_table_delete(photcal); photcal = NULL;
02336     cpl_image_delete(mapped); mapped = NULL;
02337     cpl_propertylist_delete(header); header = NULL;
02338 
02339     if (cpl_error_get_code()) {
02340         cpl_msg_error(cpl_error_get_where(), cpl_error_get_message());
02341         fors_science_exit(NULL);
02342     }
02343     else 
02344         return 0;
02345 }

Generated on Fri Mar 4 09:46:01 2011 for FORS Pipeline Reference Manual by  doxygen 1.4.7