KMOS Pipeline Reference Manual  1.1.5
kmo_sci_red.c
00001 /* $Id: kmo_sci_red.c,v 1.68 2013/06/03 15:00:40 aagudo Exp $
00002  *
00003  * This file is part of the KMOS Pipeline
00004  * Copyright (C) 2002,2003 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., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00019  */
00020 
00021 /*
00022  * $Author: aagudo $
00023  * $Date: 2013/06/03 15:00:40 $
00024  * $Revision: 1.68 $
00025  * $Name: HEAD $
00026  */
00027 
00028 #ifdef HAVE_CONFIG_H
00029 #include <config.h>
00030 #endif
00031 
00032 /*-----------------------------------------------------------------------------
00033  *                              Includes
00034  *----------------------------------------------------------------------------*/
00035 #include <string.h>
00036 #include <math.h>
00037 
00038 #include <cpl.h>
00039 #include "kmclipm_constants.h"
00040 #include "kmclipm_functions.h"
00041 
00042 #include "kmo_debug.h"
00043 #include "kmo_constants.h"
00044 #include "kmo_priv_lcorr.h"
00045 #include "kmo_utils.h"
00046 #include "kmo_error.h"
00047 #include "kmo_dfs.h"
00048 #include "kmo_functions.h"
00049 #include "kmo_priv_arithmetic.h"
00050 #include "kmo_priv_combine.h"
00051 #include "kmo_priv_functions.h"
00052 #include "kmo_priv_reconstruct.h"
00053 #include "kmo_priv_std_star.h"
00054 
00055 /*-----------------------------------------------------------------------------
00056  *              Types
00057  *-----------------------------------------------------------------------------*/
00058 typedef struct {
00059    char *name;
00060    int  count;
00061 } allObjStruct;
00062 
00063 /*-----------------------------------------------------------------------------
00064  *                          Functions prototypes
00065  *----------------------------------------------------------------------------*/
00066 
00067 static int kmo_sci_red_create(cpl_plugin *);
00068 static int kmo_sci_red_exec(cpl_plugin *);
00069 static int kmo_sci_red_destroy(cpl_plugin *);
00070 static int kmo_sci_red(cpl_parameterlist *, cpl_frameset *);
00071 
00072 /*-----------------------------------------------------------------------------
00073  *                          Static variables
00074  *----------------------------------------------------------------------------*/
00075 
00076 static char kmo_sci_red_description[] =
00077 "At least two data frames have to be provided since we need for each IFU poin-\n"
00078 "ting to an object also a sky frame from the same IFU.\n"
00079 "If an OH spectrum is given in the SOF file the lambda axis will be corrected\n"
00080 "using the OH lines as reference.\n"
00081 "Every IFU containing an object will be reconstructed and divided by telluric\n"
00082 "and illumination correction. By default these intermediate cubes are saved to\n"
00083 "disk. Frames just containing skies won’t produce an output here, so the number\n"
00084 "of output frames can be smaller than the number of input frames.\n"
00085 "Then the reconstructed objects with the same object name are combined. These\n"
00086 "outputs are also saved to disk, the number of created files depends on the\n"
00087 "number of reconstructed objects of different name. If the user just wants to\n"
00088 "combine a certain object, the parameters --name or --ifus can be used.\n"
00089 "For exposures taken with the templates KMOS_spec_obs_mapping8 and\n"
00090 "KMOS_spec_obs_mapping24 the recipe behaves a bit different: All active IFUs\n"
00091 "will be combined, regardless of the object names.\n"
00092 "This recipe must be called after the kmo_std_star-recipe.\n"
00093 "\n"
00094 "BASIC PARAMETERS:\n"
00095 "-----------------\n"
00096 "--imethod\n"
00097 "The interpolation method used for reconstruction.\n"
00098 "\n"
00099 "--smethod\n"
00100 "The interpolation method used for shifting.\n"
00101 "\n"
00102 "ADVANCED PARAMETERS\n"
00103 "-------------------\n"
00104 "--flux\n"
00105 "Specify if flux conservation should be applied.\n"
00106 "\n"
00107 "--background\n"
00108 "Specify if background removal should be applied.\n"
00109 "\n"
00110 "--suppress_extension\n"
00111 "If set to TRUE, the arbitrary filename extensions are supressed. If multiple\n"
00112 "products with the same category are produced, they will be numered consecutively\n"
00113 "starting from 0.\n"
00114 "  Advanced reconstruction parameters\n"
00115 "  ----------------------------------\n"
00116 "--neighborhoodRange\n"
00117 "Defines the range to search for neighbors during reconstruction\n"
00118 "\n"
00119 "--b_samples\n"
00120 "The number of samples in spectral direction for the reconstructed cube.\n"
00121 "Ideally this number should be greater than 2048, the detector size.\n"
00122 "\n"
00123 "--b_start\n"
00124 "--b_end\n"
00125 "Used to define manually the start and end wavelength for the reconstructed\n"
00126 "cube. By default the internally defined values are used.\n"
00127 "\n"
00128 "--fast_mode\n"
00129 "If set to TRUE, the reconstructed cubes will be collapsed (using median) and\n"
00130 "only then be shifted and combined.\n"
00131 "\n"
00132 "--pix_scale"
00133 "Change the pixel scale [arcsec]. Default of 0.2\" results into cubes of\n"
00134 "14x14pix, a scale of 0.1\" results into cubes of 28x28pix, etc.\n"
00135 "\n"
00136 "  Advanced combining parameters\n"
00137 "  ----------------------------------\n"
00138 "--edge_nan\n"
00139 "Set borders of two sides of the cubes to NaN before combining them. This minimises\n"
00140 "unwanted border effects when dithering.\n"
00141 "\n"
00142 "--no_combine\n"
00143 "If set to TRUE, the reconstructed cubes will not be combined.\n"
00144 "\n"
00145 "--no_subtract\n"
00146 "If set to TRUE, the found objects and references won't be sky subtracted.\n"
00147 "\n"
00148 "--name\n"
00149 "--ifus\n"
00150 "Since an object can be present only once per exposure and since it can be\n"
00151 "located in different IFUs for the existing exposures, there are two modes to\n"
00152 "identify the objects:\n"
00153 "   * Combine by object names (default)\n"
00154 "   In this case the object name must be provided via the --name parameter. The\n"
00155 "   object name will be searched for in all primary headers of all provided\n"
00156 "   frames in the keyword ESO OCS ARMx NAME.\n"
00157 "\n"
00158 "   * Combine by index (advanced)\n"
00159 "   In this case the --ifus parameter must be provided. The parameter must have\n"
00160 "   the same number of entries as frames are provided, e.g. \"3;1;24\" for 3\n"
00161 "   exposures. The index doesn't reference the extension in the frame but the\n"
00162 "   real index of the IFU as defined in the EXTNAME keyword.\n"
00163 "   (e.g. 'IFU.3.DATA')\n"
00164 "\n"
00165 "--method\n"
00166 "There are following sources to get the shift parameters from:\n"
00167 "   * 'header' (default)\n"
00168 "   The shifts are calculated according to the WCS information stored in the\n"
00169 "   header of every IFU. The output frame will get larger, except the object is\n"
00170 "   at the exact same position for all exposures. The size of the exposures can\n"
00171 "   differ, but the orientation must be the same for all exposures.\n"
00172 "\n"
00173 "   * 'none'\n"
00174 "   The cubes are directly recombined, not shifting at all. The ouput frame\n"
00175 "   will have the same dimensions as the input cubes.\n"
00176 "   If the size differs a warning will be emitted and the cubes will be aligned\n"
00177 "   to the lower left corner. If the orientation differs a warning will be\n"
00178 "   emitted, but the cubes are combined anyway.\n"
00179 "\n"
00180 "   * 'center'\n"
00181 "   The shifts are calculated using a centering algorithm. The cube will be\n"
00182 "   collapsed and a 2D profile will be fitted to it to identify the centre.\n"
00183 "   With the parameter --fmethod the function to fit can be provided. The size\n"
00184 "   of the exposures can differ, but the orientation must be the same for all\n"
00185 "   exposures.\n"
00186 "\n"
00187 "   * 'user'\n"
00188 "   Read the shifts from a user specified file. The path of the file must be\n"
00189 "   provided using the --filename parameter. For every exposure (except the\n"
00190 "   first one) two shift values are expected per line, they have to be separa-\n"
00191 "   ted with simple spaces. The values indicate pixel shifts and are referenced\n"
00192 "   to the first frame. The 1st value is the shift in x-direction to the left,\n"
00193 "   the 2nd the shift in y-direction upwards. The size of the exposures can\n"
00194 "   differ, but the orientation must be the same for all exposures.\n"
00195 "\n"
00196 "--fmethod\n"
00197 "see --method='center'\n"
00198 "The type of function that should be fitted spatially to the collapsed image.\n"
00199 "This fit is used to create a mask to extract the spectrum of the object. Valid\n"
00200 "values are 'gauss' and 'moffat'.\n"
00201 "\n"
00202 "--filename\n"
00203 "see --method='user'\n"
00204 "\n"
00205 "--cpos_rej\n"
00206 "--cneg_rej\n"
00207 "--citer\n"
00208 "see --cmethod='ksigma'\n"
00209 "\n"
00210 "--cmax\n"
00211 "--cmin\n"
00212 "see --cmethod='min_max'\n"
00213 "\n"
00214 "--extrapolate\n"
00215 "By default no extrapolation is applied. This means that the intermediate\n"
00216 "reconstructed cubes will shrink at most one pixel, which is ok for templates\n"
00217 "like KMOS_spec_obs_nodtosky or KMOS_spec_obs_freedither. When the cubes will be\n"
00218 "arranged as a map, a grid is likely to occur between the IFUs. Therefore extra-\n"
00219 "polation during the shifting process can be switched on in order to get IFUs of\n"
00220 "original size. For frames taken with mapping templates, extrapolation is\n"
00221 "switched on automatically.\n"
00222 "\n"
00223 "--xcal_interpolation\n"
00224 "If true interpolate the pixel position in the slitlet (xcal) using the two\n"
00225 "closest rotator angles in the calibration file. Otherwise take the values\n"
00226 "of the closest rotator angle\n"
00227 "\n"
00228 "------------------------------------------------------------------------------\n"
00229 "  Input files:\n"
00230 "\n"
00231 "   DO                    KMOS                                                  \n"
00232 "   category              Type   Explanation                   Required #Frames\n"
00233 "   --------              -----  -----------                   -------- -------\n"
00234 "   SCIENCE               RAW    The science frames                Y      >=1  \n"
00235 "   XCAL                  F2D    x calibration frame               Y       1   \n"
00236 "   YCAL                  F2D    y calibration frame               Y       1   \n"
00237 "   LCAL                  F2D    Wavelength calib. frame           Y       1   \n"
00238 "   MASTER_FLAT           F2D    Master flat                       Y       1   \n"
00239 "   WAVE_BAND             F2L    Table with start-/end-wavelengths Y       1   \n"
00240 "   ILLUM_CORR            F2I    Illumination correction           N      0,1  \n"
00241 "   TELLURIC              F1I    normalised telluric spectrum      N      0,1  \n"
00242 "   OH_SPEC               F1S    Vector holding OH lines           N       1   \n"
00243 "\n"
00244 "  Output files:\n"
00245 "\n"
00246 "   DO                    KMOS\n"
00247 "   category              Type   Explanation\n"
00248 "   --------              -----  -----------\n"
00249 "   SCI_COMBINED          F3I    Combined cubes with noise\n"
00250 "   SCI_RECONSTRUCTED     F3I    Reconstructed cube with noise\n"
00251 "------------------------------------------------------------------------------\n"
00252 "\n";
00253 
00254 /*-----------------------------------------------------------------------------
00255  *                              Functions code
00256  *----------------------------------------------------------------------------*/
00257 
00274 int cpl_plugin_get_info(cpl_pluginlist *list)
00275 {
00276     cpl_recipe *recipe = cpl_calloc(1, sizeof *recipe);
00277     cpl_plugin *plugin = &recipe->interface;
00278 
00279     cpl_plugin_init(plugin,
00280                         CPL_PLUGIN_API,
00281                         KMOS_BINARY_VERSION,
00282                         CPL_PLUGIN_TYPE_RECIPE,
00283                         "kmo_sci_red",
00284                         "Reconstruct and combine data frames dividing "
00285                         "illumination and telluric correction.",
00286                         kmo_sci_red_description,
00287                         "Alex Agudo Berbel",
00288                         "kmos-spark@mpe.mpg.de",
00289                         kmos_get_license(),
00290                         kmo_sci_red_create,
00291                         kmo_sci_red_exec,
00292                         kmo_sci_red_destroy);
00293 
00294     cpl_pluginlist_append(list, plugin);
00295 
00296     return 0;
00297 }
00298 
00306 static int kmo_sci_red_create(cpl_plugin *plugin)
00307 {
00308     cpl_recipe *recipe;
00309     cpl_parameter *p;
00310 
00311     /* Check that the plugin is part of a valid recipe */
00312     if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 
00313         recipe = (cpl_recipe *)plugin;
00314     else
00315         return -1;
00316 
00317     /* Create the parameters list in the cpl_recipe object */
00318     recipe->parameters = cpl_parameterlist_new();
00319 
00320     /* --imethod */
00321     p = cpl_parameter_new_value("kmos.kmo_sci_red.imethod",
00322                                 CPL_TYPE_STRING,
00323                                 "Method to use for interpolation during reconstruction. "
00324                                 "[\"NN\" (nearest neighbour), "
00325                                 "\"lwNN\" (linear weighted nearest neighbor), "
00326                                 "\"swNN\" (square weighted nearest neighbor), "
00327                                 "\"MS\" (Modified Shepard's method)"
00328                                 "\"CS\" (Cubic spline)]",
00329                                 "kmos.kmo_sci_red",
00330                                 "CS");
00331     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "imethod");
00332     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00333     cpl_parameterlist_append(recipe->parameters, p);
00334 
00335     /* --smethod */
00336     p = cpl_parameter_new_value("kmos.kmo_sci_red.smethod",
00337                                 CPL_TYPE_STRING,
00338                                 "Method to use for interpolation during shifting. "
00339                                 "[\"NN\" (nearest neighbour), "
00340                                 "\"CS\" (Cubic spline)]",
00341                                 "kmos.kmo_sci_red",
00342                                 "CS");
00343     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "smethod");
00344     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00345     cpl_parameterlist_append(recipe->parameters, p);
00346 
00347     /* --neighborhoodRange */
00348     p = cpl_parameter_new_value("kmos.kmo_sci_red.neighborhoodRange",
00349                                 CPL_TYPE_DOUBLE,
00350                                 "Defines the range to search for neighbors "
00351                                 "in pixels",
00352                                 "kmos.kmo_sci_red",
00353                                 1.001);
00354     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "neighborhoodRange");
00355     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00356     cpl_parameterlist_append(recipe->parameters, p);
00357 
00358     /* --name */
00359     p = cpl_parameter_new_value("kmos.kmo_sci_red.name",
00360                                 CPL_TYPE_STRING,
00361                                 "Name of the object to combine.",
00362                                 "kmos.kmo_sci_red",
00363                                 "");
00364     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "name");
00365     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00366     cpl_parameterlist_append(recipe->parameters, p);
00367 
00368     /* --ifus */
00369     p = cpl_parameter_new_value("kmos.kmo_sci_red.ifus",
00370                                 CPL_TYPE_STRING,
00371                                 "The indices of the IFUs to combine. "
00372                                 "\"ifu1;ifu2;...\"",
00373                                 "kmos.kmo_sci_red",
00374                                 "");
00375     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "ifus");
00376     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00377     cpl_parameterlist_append(recipe->parameters, p);
00378 
00379     /* --method */
00380     p = cpl_parameter_new_value("kmos.kmo_sci_red.method",
00381                                 CPL_TYPE_STRING,
00382                                 "The shifting method:   "
00383                                 "'none': no shifting, combined directly, "
00384                                 "'header': shift according to WCS (default), "
00385                                 "'center': centering algorithm, "
00386                                 "'user': read shifts from file",
00387                                 "kmos.kmo_sci_red",
00388                                 "header");
00389     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "method");
00390     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00391     cpl_parameterlist_append(recipe->parameters, p);
00392 
00393     /* --fmethod */
00394     p = cpl_parameter_new_value("kmos.kmo_sci_red.fmethod",
00395                                 CPL_TYPE_STRING,
00396                                 "The fitting method (applies only when "
00397                                 "method='center'):   "
00398                                 "'gauss': fit a gauss function to collapsed "
00399                                 "image (default), "
00400                                 "'moffat': fit a moffat function to collapsed"
00401                                 " image",
00402                                 "kmos.kmo_sci_red",
00403                                 "gauss");
00404     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "fmethod");
00405     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00406     cpl_parameterlist_append(recipe->parameters, p);
00407 
00408     /* --filename */
00409     p = cpl_parameter_new_value("kmos.kmo_sci_red.filename",
00410                                 CPL_TYPE_STRING,
00411                                 "The path to the file with the shift vectors."
00412                                 "(Applies only to method='user')",
00413                                 "kmos.kmo_sci_red",
00414                                 "");
00415     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "filename");
00416     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00417     cpl_parameterlist_append(recipe->parameters, p);
00418 
00419     /* --flux */
00420     p = cpl_parameter_new_value("kmos.kmo_sci_red.flux",
00421                                 CPL_TYPE_BOOL,
00422                                 "TRUE: Apply flux conservation. FALSE: otherwise",
00423                                 "kmos.kmo_sci_red",
00424                                 FALSE);
00425     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "flux");
00426     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00427     cpl_parameterlist_append(recipe->parameters, p);
00428 
00429     /* --background */
00430     p = cpl_parameter_new_value("kmos.kmo_sci_red.background",
00431                                 CPL_TYPE_BOOL,
00432                                 "TRUE: Apply background removal. FALSE: otherwise",
00433                                 "kmos.kmo_sci_red",
00434                                 FALSE);
00435     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "background");
00436     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00437     cpl_parameterlist_append(recipe->parameters, p);
00438 
00439     /* --extrapolate */
00440     p = cpl_parameter_new_value("kmos.kmo_sci_red.extrapolate",
00441                                 CPL_TYPE_BOOL,
00442                                 "Applies only to 'method=CS' when doing sub-"
00443                                 "pixel shifts: "
00444                                 "FALSE: shifted IFU will be filled with NaN's "
00445                                 "at the borders,"
00446                                 "TRUE: shifted IFU will be extrapolated at "
00447                                 "the borders",
00448                                 "kmos.kmo_sci_red",
00449                                 FALSE);
00450     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "extrapolate");
00451     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00452     cpl_parameterlist_append(recipe->parameters, p);
00453 
00454     /* --fast_mode */
00455     p = cpl_parameter_new_value("kmos.kmo_sci_red.fast_mode",
00456                                 CPL_TYPE_BOOL,
00457                                 "FALSE: cubes are shifted and combined,"
00458                                 "TRUE: cubes are collapsed and then shifted and combined",
00459                                 "kmos.kmo_sci_red",
00460                                 FALSE);
00461     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "fast_mode");
00462     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00463     cpl_parameterlist_append(recipe->parameters, p);
00464 
00465     /* --edge_nan */
00466     p = cpl_parameter_new_value("kmos.kmo_sci_red.edge_nan",
00467                                 CPL_TYPE_BOOL,
00468                                 "Set borders of cubes to NaN before combining them."
00469                                 "(TRUE (apply) or "
00470                                 "FALSE (don't apply)",
00471                                 "kmos.kmo_sci_red",
00472                                 FALSE);
00473     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "edge_nan");
00474     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00475     cpl_parameterlist_append(recipe->parameters, p);
00476 
00477     /* --no_combine */
00478     p = cpl_parameter_new_value("kmos.kmo_sci_red.no_combine",
00479                                 CPL_TYPE_BOOL,
00480                                 "Don't combine cubes after reconstruction."
00481                                 "(TRUE (apply) or "
00482                                 "FALSE (don't apply)",
00483                                 "kmos.kmo_sci_red",
00484                                 FALSE);
00485     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "no_combine");
00486     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00487     cpl_parameterlist_append(recipe->parameters, p);
00488 
00489     /* --no_subtract */
00490     p = cpl_parameter_new_value("kmos.kmo_sci_red.no_subtract",
00491                                 CPL_TYPE_BOOL,
00492                                 "Don't sky subtract object and references."
00493                                 "(TRUE (apply) or "
00494                                 "FALSE (don't apply)",
00495                                 "kmos.kmo_sci_red",
00496                                 FALSE);
00497     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "no_subtract");
00498     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00499     cpl_parameterlist_append(recipe->parameters, p);
00500 
00501     /* --pix_scale */
00502     p = cpl_parameter_new_value("kmos.kmo_sci_red.pix_scale",
00503                                 CPL_TYPE_DOUBLE,
00504                                 "Change the pixel scale [arcsec]. "
00505                                 "Default of 0.2\" results into cubes of 14x14pix, "
00506                                 "a scale of 0.1\" results into cubes of 28x28pix, "
00507                                 "etc.",
00508                                 "kmos.kmo_sci_red",
00509                                 KMOS_PIX_RESOLUTION);
00510     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "pix_scale");
00511     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00512     cpl_parameterlist_append(recipe->parameters, p);
00513 
00514     /* --xcal_interpolation */
00515     p = cpl_parameter_new_value("kmos.kmo_sci_red.xcal_interpolation",
00516                                 CPL_TYPE_BOOL,
00517                                 "TRUE: Interpolate xcal between rotator angles. FALSE: otherwise",
00518                                 "kmos.kmo_sci_red",
00519                                 TRUE);
00520     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "xcal_interpolation");
00521     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00522     cpl_parameterlist_append(recipe->parameters, p);
00523 
00524     /* --suppress_extension */
00525     p = cpl_parameter_new_value("kmos.kmo_sci_red.suppress_extension",
00526                                 CPL_TYPE_BOOL,
00527                                 "Suppress arbitrary filename extension."
00528                                 "(TRUE (apply) or FALSE (don't apply)",
00529                                 "kmos.kmo_sci_red",
00530                                 FALSE);
00531     cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "suppress_extension");
00532     cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00533     cpl_parameterlist_append(recipe->parameters, p);
00534 
00535     // add parameters for band-definition
00536     kmo_band_pars_create(recipe->parameters,
00537                          "kmos.kmo_sci_red");
00538 
00539     return kmo_combine_pars_create(recipe->parameters,
00540                                    "kmos.kmo_sci_red",
00541                                    DEF_REJ_METHOD,
00542                                    FALSE);
00543 }
00544 
00550 static int kmo_sci_red_exec(cpl_plugin *plugin)
00551 {
00552     cpl_recipe  *recipe;
00553 
00554     /* Get the recipe out of the plugin */
00555     if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 
00556         recipe = (cpl_recipe *)plugin;
00557     else return -1 ;
00558 
00559     return kmo_sci_red(recipe->parameters, recipe->frames);
00560 }
00561 
00567 static int kmo_sci_red_destroy(cpl_plugin *plugin)
00568 {
00569     cpl_recipe *recipe;
00570 
00571     /* Get the recipe out of the plugin */
00572     if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE) 
00573         recipe = (cpl_recipe *)plugin;
00574     else return -1 ;
00575 
00576     cpl_parameterlist_delete(recipe->parameters);
00577     return 0 ;
00578 }
00579 
00580 
00581 
00596 static int kmo_sci_red(cpl_parameterlist *parlist, cpl_frameset *frameset)
00597 {
00598     int                     ret_val                    = 0,
00599                             nr_science_frames          = 0,
00600                             nr_reconstructed_frames    = 0,
00601                             has_illum_corr             = 0,
00602                             has_telluric               = 0,
00603                             *bounds                    = NULL,
00604                             det_nr                     = 0,
00605                             actual_msg_level           = 0,
00606                             print_once                 = FALSE,
00607                             nr_avail_obj_names         = 0,
00608                             found_name                 = FALSE,
00609                             cube_counter_data          = 0,
00610                             cube_counter_noise         = 0,
00611                             citer                      = 0,
00612                             cmin                       = 0,
00613                             cmax                       = 0,
00614                             user_defined_ifu           = 0,
00615                             extrapolate                = 0,
00616                             flux                       = FALSE,
00617                             background                 = FALSE,
00618                             index                      = 0,
00619                             nr_data_alloc              = 0,
00620                             tmp_int                    = 0,
00621                             fast_mode                  = FALSE,
00622                             edge_nan                   = FALSE,
00623                             no_combine                 = FALSE,
00624                             no_subtract                = FALSE,
00625                             xcal_interpolation         = FALSE,
00626                             suppress_extension         = FALSE,
00627                             suppress_index             = 0;
00628     const int               *punused_ifus              = NULL;
00629     double                  neighborhoodRange          = 1.001,
00630                             cpos_rej                   = 0.0,
00631                             cneg_rej                   = 0.0,
00632                             pix_scale                  = 0.0;
00633     char                    *suffix                    = NULL,
00634                             *keyword                   = NULL,
00635                             *extname                   = NULL,
00636                             *fn_suffix                 = NULL,
00637                             *mapping_mode              = NULL,
00638                             **split                    = NULL,
00639                             content[256];
00640     const char              *imethod                   = NULL,
00641                             *smethod                   = NULL,
00642                             *ifus_txt                  = NULL,
00643                             *name                      = NULL,
00644                             *filter_id                 = NULL,
00645                             *tmp_str                   = NULL,
00646                             *filename                  = NULL,
00647                             *fn_out                    = NULL,
00648                             *fn_obj                    = NULL,
00649                             *fn_sky                    = NULL,
00650                             *fn_reconstr               = NULL,
00651                             *comb_method               = NULL,
00652                             *cmethod                   = NULL,
00653                             *fmethod                   = NULL;
00654     cpl_array               **unused_ifus              = NULL;
00655     cpl_frame               *xcal_frame                = NULL,
00656                             *ycal_frame                = NULL,
00657                             *lcal_frame                = NULL,
00658                             *flat_frame                = NULL,
00659                             *illum_frame               = NULL,
00660                             *telluric_frame            = NULL,
00661                             *tmp_frame                 = NULL;
00662     cpl_propertylist        *tmp_header                = NULL,
00663                             *main_header               = NULL,
00664                             **header_data              = NULL,
00665                             **header_noise             = NULL;
00666     cpl_vector              *ifus                      = NULL;
00667     kmclipm_vector          *telluric_data             = NULL,
00668                             *telluric_noise            = NULL;
00669     cpl_image               **lcal                     = NULL,
00670                             *illum_data                = NULL,
00671                             *illum_noise               = NULL,
00672                             *tmpImg                    = NULL;
00673     cpl_imagelist           **cube_data                = NULL,
00674                             **cube_noise               = NULL,
00675                             *combined_data             = NULL,
00676                             *combined_noise            = NULL,
00677                             *tmpCube                   = NULL;
00678     cpl_table               *band_table                = NULL;
00679     cpl_frame               *sky_frame                 = NULL,
00680                             *ref_spectrum_frame        = NULL;;
00681     cpl_bivector            *obj_spectrum              = NULL,
00682                             *ref_spectrum              = NULL;
00683     cpl_vector              *peaks                     = NULL,
00684                             *range                     = NULL;
00685     cpl_polynomial          *lcorr_coeffs              = NULL;
00686     main_fits_desc          desc1,
00687                             desc2;
00688     gridDefinition          gd;
00689     objSkyFrameTableStruct  *obj_sky_struct            = NULL;
00690     allObjStruct            *all_obj                   = NULL;
00691     enum extrapolationType  extrapol_enum              = 0;
00692     enum kmo_frame_type     ft                         = 0;
00693 
00694     KMO_TRY
00695     {
00696 
00697         kmo_init_fits_desc(&desc1);
00698         kmo_init_fits_desc(&desc2);
00699 
00700         //
00701         // check frameset
00702         //
00703         KMO_TRY_ASSURE((parlist != NULL) &&
00704                        (frameset != NULL),
00705                        CPL_ERROR_NULL_INPUT,
00706                        "Not all input data is provided!");
00707 
00708         nr_science_frames = cpl_frameset_count_tags(frameset, SCIENCE);
00709         KMO_TRY_ASSURE(nr_science_frames >= 1,
00710                        CPL_ERROR_ILLEGAL_INPUT,
00711                        "At least one SCIENCE frame is required!");
00712         if (nr_science_frames == 1) {
00713             cpl_msg_warning("", "At least two SCIENCE frames should be provided "
00714                                 "in order to apply sky subtraction!");
00715             cpl_msg_warning("", "All IFUs will be reconstructed regardless if "
00716                                 "they contain object, reference or sky!");
00717         }
00718 
00719         KMO_TRY_ASSURE(cpl_frameset_count_tags(frameset, XCAL) == 1,
00720                        CPL_ERROR_FILE_NOT_FOUND,
00721                        "Exactly one XCAL frame is required!");
00722 
00723         KMO_TRY_ASSURE(cpl_frameset_count_tags(frameset, YCAL) == 1,
00724                        CPL_ERROR_FILE_NOT_FOUND,
00725                        "Exactly one YCAL frame is required!");
00726 
00727         KMO_TRY_ASSURE(cpl_frameset_count_tags(frameset, LCAL) == 1,
00728                        CPL_ERROR_FILE_NOT_FOUND,
00729                        "Exactly one LCAL frame is required!");
00730 
00731         KMO_TRY_ASSURE(cpl_frameset_count_tags(frameset, MASTER_FLAT) == 1,
00732                        CPL_ERROR_FILE_NOT_FOUND,
00733                        "Exactly one MASTER_FLAT frame is required!");
00734 
00735         KMO_TRY_ASSURE(cpl_frameset_count_tags(frameset, WAVE_BAND) == 1,
00736                        CPL_ERROR_FILE_NOT_FOUND,
00737                        "Exactly one WAVE_BAND frame is required!");
00738 
00739         has_illum_corr = cpl_frameset_count_tags(frameset, ILLUM_CORR);
00740         KMO_TRY_ASSURE((has_illum_corr == 0) || (has_illum_corr == 1),
00741                        CPL_ERROR_FILE_NOT_FOUND,
00742                        "At most one ILLUM_CORR frame can be provided!");
00743 
00744         has_telluric = cpl_frameset_count_tags(frameset, TELLURIC);
00745         KMO_TRY_ASSURE((has_telluric == 0) || (has_telluric == 1),
00746                        CPL_ERROR_FILE_NOT_FOUND,
00747                        "At most one TELLURIC frame can be provided!");
00748 
00749         KMO_TRY_ASSURE(kmo_dfs_set_groups(frameset, "kmo_sci_red") == 1,
00750                        CPL_ERROR_ILLEGAL_INPUT,
00751                        "Cannot identify RAW and CALIB frames!");
00752 
00753         KMO_TRY_ASSURE(cpl_frameset_count_tags(frameset, OH_SPEC) == 0 ||
00754                        cpl_frameset_count_tags(frameset, OH_SPEC) == 1,
00755                        CPL_ERROR_ILLEGAL_INPUT,
00756                        "Only a single reference spectrum can be provided!");
00757         //
00758         // get parameters
00759         //
00760         cpl_msg_info("", "--- Parameter setup for kmo_sci_red ------");
00761 
00762         flux = kmo_dfs_get_parameter_bool(parlist,
00763                                           "kmos.kmo_sci_red.flux");
00764 
00765         KMO_TRY_ASSURE((flux == 0) ||
00766                        (flux == 1),
00767                        CPL_ERROR_ILLEGAL_INPUT,
00768                        "flux must be either FALSE or TRUE! %d", flux);
00769 
00770         KMO_TRY_EXIT_IF_ERROR(
00771             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.flux"));
00772 
00773         background = kmo_dfs_get_parameter_bool(parlist,
00774                                           "kmos.kmo_sci_red.background");
00775 
00776         KMO_TRY_ASSURE((background == 0) ||
00777                        (background == 1),
00778                        CPL_ERROR_ILLEGAL_INPUT,
00779                        "background must be either FALSE or TRUE! %d", background);
00780 
00781         KMO_TRY_EXIT_IF_ERROR(
00782             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.background"));
00783 
00784         KMO_TRY_EXIT_IF_NULL(
00785             imethod = kmo_dfs_get_parameter_string(parlist,
00786                                                    "kmos.kmo_sci_red.imethod"));
00787 
00788         KMO_TRY_ASSURE((strcmp(imethod, "NN") == 0) ||
00789                        (strcmp(imethod, "lwNN") == 0) ||
00790                        (strcmp(imethod, "swNN") == 0) ||
00791                        (strcmp(imethod, "MS") == 0) ||
00792                        (strcmp(imethod, "CS") == 0),
00793                        CPL_ERROR_ILLEGAL_INPUT,
00794                        "imethod must be either \"NN\", \"lwNN\", "
00795                        "\"swNN\", \"MS\" or \"CS\"!");
00796 
00797         KMO_TRY_EXIT_IF_ERROR(
00798             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.imethod"));
00799 
00800         KMO_TRY_EXIT_IF_NULL(
00801             smethod = kmo_dfs_get_parameter_string(parlist,
00802                                                    "kmos.kmo_sci_red.smethod"));
00803 
00804         KMO_TRY_ASSURE((strcmp(smethod, "NN") == 0) ||
00805                        (strcmp(smethod, "CS") == 0),
00806                        CPL_ERROR_ILLEGAL_INPUT,
00807                        "smethod must be either \"NN\" or \"CS\"!");
00808 
00809         KMO_TRY_EXIT_IF_ERROR(
00810             kmo_dfs_print_parameter_help(parlist,
00811                                          "kmos.kmo_sci_red.smethod"));
00812 
00813         neighborhoodRange = kmo_dfs_get_parameter_double(parlist,
00814                                           "kmos.kmo_sci_red.neighborhoodRange");
00815         KMO_TRY_CHECK_ERROR_STATE();
00816 
00817         KMO_TRY_ASSURE(neighborhoodRange > 0.0,
00818                        CPL_ERROR_ILLEGAL_INPUT,
00819                        "neighborhoodRange must be greater than 0.0");
00820 
00821         KMO_TRY_EXIT_IF_ERROR(
00822             kmo_dfs_print_parameter_help(parlist,
00823                                          "kmos.kmo_sci_red.neighborhoodRange"));
00824 
00825         KMO_TRY_EXIT_IF_NULL(
00826             comb_method = kmo_dfs_get_parameter_string(parlist,
00827                                            "kmos.kmo_sci_red.method"));
00828 
00829         KMO_TRY_EXIT_IF_NULL(
00830             fmethod = kmo_dfs_get_parameter_string(parlist,
00831                                            "kmos.kmo_sci_red.fmethod"));
00832 
00833         KMO_TRY_ASSURE((strcmp(comb_method, "none") == 0) ||
00834                        (strcmp(comb_method, "header") == 0) ||
00835                        (strcmp(comb_method, "center") == 0) ||
00836                        (strcmp(comb_method, "user") == 0),
00837                        CPL_ERROR_ILLEGAL_INPUT,
00838                        "Following shift methods are available : 'none', "
00839                        "'header', 'center' or 'user'");
00840 
00841         if (strcmp(comb_method, "user") == 0) {
00842             filename = kmo_dfs_get_parameter_string(parlist,
00843                                                    "kmos.kmo_sci_red.filename");
00844             KMO_TRY_CHECK_ERROR_STATE();
00845 
00846             KMO_TRY_ASSURE(strcmp(filename, "") != 0,
00847                            CPL_ERROR_ILLEGAL_INPUT,
00848                            "path of file with shift information must be "
00849                            "provided!");
00850 
00851             KMO_TRY_EXIT_IF_ERROR(
00852                 kmo_dfs_print_parameter_help(parlist,
00853                                              "kmos.kmo_sci_red.filename"));
00854         }
00855 
00856         KMO_TRY_EXIT_IF_ERROR(
00857             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.method"));
00858 
00859         ifus_txt = kmo_dfs_get_parameter_string(parlist,
00860                                                   "kmos.kmo_sci_red.ifus");
00861         KMO_TRY_CHECK_ERROR_STATE();
00862 
00863         name = kmo_dfs_get_parameter_string(parlist, "kmos.kmo_sci_red.name");
00864         KMO_TRY_CHECK_ERROR_STATE();
00865 
00866         if (strcmp(ifus_txt, "") != 0) {
00867             KMO_TRY_ASSURE(strcmp(name, "") == 0,
00868                            CPL_ERROR_ILLEGAL_INPUT,
00869                            "name parameter must be NULL if IFU indices are "
00870                            "provided!");
00871 
00872             KMO_TRY_EXIT_IF_NULL(
00873                 ifus = kmo_identify_values(ifus_txt));
00874 
00875             KMO_TRY_ASSURE(cpl_vector_get_size(ifus) == nr_science_frames,
00876                            CPL_ERROR_ILLEGAL_INPUT,
00877                            "ifus parameter must have the same number of values "
00878                            "than frames provided (for frames just containing "
00879                            "skies insert 0)) (%lld!=%d)",
00880                            cpl_vector_get_size(ifus), nr_science_frames);
00881         }
00882 
00883         if (strcmp(name, "") != 0) {
00884             KMO_TRY_ASSURE(strcmp(ifus_txt, "") == 0,
00885                            CPL_ERROR_ILLEGAL_INPUT,
00886                            "ifus parameter must be NULL if name is provided!");
00887         }
00888 
00889         KMO_TRY_EXIT_IF_ERROR(
00890             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.ifus"));
00891 
00892         KMO_TRY_EXIT_IF_ERROR(
00893             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.name"));
00894 
00895         kmo_band_pars_load(parlist, "kmos.kmo_sci_red");
00896 
00897         extrapolate = kmo_dfs_get_parameter_bool(parlist,
00898                                                 "kmos.kmo_sci_red.extrapolate");
00899         KMO_TRY_CHECK_ERROR_STATE();
00900 
00901         if (strcmp(smethod, "NN") == 0) {
00902             if (extrapolate == TRUE) {
00903                 cpl_msg_warning("", "extrapolation for smethod='NN' not available!");
00904             }
00905             extrapol_enum = NONE_NANS;
00906         } else if (strcmp(smethod, "CS") == 0) {
00907             if (extrapolate == FALSE) {
00908                 extrapol_enum = NONE_NANS;
00909             } else if (extrapolate == TRUE) {
00910                 extrapol_enum = BCS_NATURAL;
00911             } else {
00912                 KMO_TRY_ASSURE(1 == 0,
00913                                CPL_ERROR_ILLEGAL_INPUT,
00914                                "extrapolate must be either FALSE or TRUE!");
00915             }
00916             smethod = "BCS";
00917         } else {
00918             KMO_TRY_ASSURE(1 == 0,
00919                            CPL_ERROR_ILLEGAL_INPUT,
00920                            "method must be either \"CS\" or \"NN\" !");
00921         }
00922         KMO_TRY_CHECK_ERROR_STATE();
00923 
00924         KMO_TRY_EXIT_IF_ERROR(
00925             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.extrapolate"));
00926 
00927         fast_mode = kmo_dfs_get_parameter_bool(parlist,
00928                                                 "kmos.kmo_sci_red.fast_mode");
00929         KMO_TRY_CHECK_ERROR_STATE();
00930         KMO_TRY_ASSURE((fast_mode == TRUE) ||
00931                        (fast_mode == FALSE),
00932                        CPL_ERROR_ILLEGAL_INPUT,
00933                        "fast_mode must be either FALSE or TRUE!");
00934         KMO_TRY_EXIT_IF_ERROR(
00935             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.fast_mode"));
00936 
00937         edge_nan = kmo_dfs_get_parameter_bool(parlist,
00938                                           "kmos.kmo_sci_red.edge_nan");
00939         KMO_TRY_CHECK_ERROR_STATE();
00940         KMO_TRY_EXIT_IF_ERROR(
00941             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.edge_nan"));
00942 
00943         KMO_TRY_ASSURE((edge_nan == TRUE) || (edge_nan == FALSE),
00944                        CPL_ERROR_ILLEGAL_INPUT,
00945                        "edge_nan must be TRUE or FALSE!");
00946 
00947         no_combine = kmo_dfs_get_parameter_bool(parlist,
00948                                                 "kmos.kmo_sci_red.no_combine");
00949         KMO_TRY_CHECK_ERROR_STATE();
00950 
00951         KMO_TRY_EXIT_IF_ERROR(
00952             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.no_combine"));
00953 
00954         KMO_TRY_ASSURE((no_combine == TRUE) || (no_combine == FALSE),
00955                        CPL_ERROR_ILLEGAL_INPUT,
00956                        "no_combine must be TRUE or FALSE!");
00957 
00958         no_subtract = kmo_dfs_get_parameter_bool(parlist,
00959                                                 "kmos.kmo_sci_red.no_subtract");
00960         KMO_TRY_CHECK_ERROR_STATE();
00961 
00962         KMO_TRY_EXIT_IF_ERROR(
00963             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.no_subtract"));
00964 
00965         KMO_TRY_ASSURE((no_subtract == TRUE) || (no_subtract == FALSE),
00966                        CPL_ERROR_ILLEGAL_INPUT,
00967                        "no_subtract must be TRUE or FALSE!");
00968 
00969         pix_scale = kmo_dfs_get_parameter_double(parlist,
00970                                         "kmos.kmo_sci_red.pix_scale");
00971         KMO_TRY_CHECK_ERROR_STATE();
00972         KMO_TRY_EXIT_IF_ERROR(
00973            kmo_dfs_print_parameter_help(parlist,
00974                                        "kmos.kmo_sci_red.pix_scale"));
00975         KMO_TRY_ASSURE((pix_scale >= 0.01) &&
00976                        (pix_scale <= 0.4),
00977                        CPL_ERROR_ILLEGAL_INPUT,
00978                        "pix_scale must be between 0.01 and 0.4 (results in cubes "
00979                        "with 7x7 to 280x280 pixels)!");
00980 
00981         xcal_interpolation = kmo_dfs_get_parameter_bool(parlist,
00982                                            "kmos.kmo_sci_red.xcal_interpolation");
00983         KMO_TRY_CHECK_ERROR_STATE();
00984         KMO_TRY_EXIT_IF_ERROR(
00985             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.xcal_interpolation"));
00986         KMO_TRY_ASSURE((xcal_interpolation == TRUE) ||
00987                        (xcal_interpolation == FALSE),
00988                        CPL_ERROR_ILLEGAL_INPUT,
00989                        "xcal_interpolation must be TRUE or FALSE!");
00990 
00991         suppress_extension = kmo_dfs_get_parameter_bool(parlist,
00992                                           "kmos.kmo_sci_red.suppress_extension");
00993         KMO_TRY_CHECK_ERROR_STATE();
00994         KMO_TRY_EXIT_IF_ERROR(
00995             kmo_dfs_print_parameter_help(parlist, "kmos.kmo_sci_red.suppress_extension"));
00996 
00997         KMO_TRY_ASSURE((suppress_extension == TRUE) || (suppress_extension == FALSE),
00998                        CPL_ERROR_ILLEGAL_INPUT,
00999                        "suppress_extension must be TRUE or FALSE!");
01000 
01001         KMO_TRY_EXIT_IF_ERROR(
01002             kmo_combine_pars_load(parlist,
01003                                   "kmos.kmo_sci_red",
01004                                   &cmethod,
01005                                   &cpos_rej,
01006                                   &cneg_rej,
01007                                   &citer,
01008                                   &cmin,
01009                                   &cmax,
01010                                   FALSE));
01011 
01012         cpl_msg_info("", "-------------------------------------------");
01013 
01014         //
01015         // assure that filters, grating and rotation offsets match for
01016         // XCAL, YCAL, LCAL and for data frame to reconstruct (except DARK
01017         // frames)
01018         //
01019 
01020         // check if filter_id and grating_id match for all detectors
01021         KMO_TRY_EXIT_IF_ERROR(
01022             kmo_check_frameset_setup(frameset, SCIENCE,
01023                                        TRUE, FALSE, TRUE));
01024         KMO_TRY_EXIT_IF_ERROR(
01025             kmo_check_frame_setup(frameset, SCIENCE, YCAL,
01026                                        TRUE, FALSE, TRUE));
01027         KMO_TRY_EXIT_IF_ERROR(
01028             kmo_check_frame_setup(frameset, XCAL, YCAL,
01029                                        TRUE, FALSE, TRUE));
01030         KMO_TRY_EXIT_IF_ERROR(
01031             kmo_check_frame_setup(frameset, XCAL, LCAL,
01032                                        TRUE, FALSE, TRUE));
01033 //        KMO_TRY_EXIT_IF_ERROR(
01034 //            kmo_check_cal_frames_rotangle(frameset, XCAL, YCAL));
01035 //        KMO_TRY_EXIT_IF_ERROR(
01036 //            kmo_check_cal_frames_rotangle(frameset, XCAL, LCAL));
01037         KMO_TRY_EXIT_IF_ERROR(
01038             kmo_check_frame_setup(frameset, XCAL, MASTER_FLAT,
01039                                        TRUE, FALSE, TRUE));
01040 // omit this check because ILLUM_CORR is band-independend
01041 //        if (has_illum_corr) {
01042 //            KMO_TRY_EXIT_IF_ERROR(
01043 //                kmo_check_frame_setup(frameset, XCAL, ILLUM_CORR,
01044 //                                           TRUE, FALSE, TRUE));
01045 //        }
01046 
01047         if (has_telluric) {
01048             KMO_TRY_EXIT_IF_ERROR(
01049                 kmo_check_frame_setup(frameset, XCAL, TELLURIC,
01050                                            TRUE, FALSE, TRUE));
01051         }
01052 
01053         // check descriptors of all frames
01054         KMO_TRY_EXIT_IF_NULL(
01055             xcal_frame = kmo_dfs_get_frame(frameset, XCAL));
01056 
01057         desc1 = kmo_identify_fits_header(cpl_frame_get_filename(xcal_frame));
01058         KMO_TRY_CHECK_ERROR_STATE();
01059 
01060         KMO_TRY_ASSURE((desc1.nr_ext % KMOS_NR_DETECTORS == 0) &&
01061                        (desc1.ex_badpix == FALSE) &&
01062                        (desc1.fits_type == f2d_fits) &&
01063                        (desc1.frame_type == detector_frame),
01064                        CPL_ERROR_ILLEGAL_INPUT,
01065                        "XCAL isn't in the correct format!!!");
01066 
01067         KMO_TRY_EXIT_IF_NULL(
01068             ycal_frame = kmo_dfs_get_frame(frameset, YCAL));
01069         desc2 = kmo_identify_fits_header(cpl_frame_get_filename(ycal_frame));
01070         KMO_TRY_CHECK_ERROR_STATE();
01071 
01072         KMO_TRY_ASSURE((desc1.nr_ext == desc2.nr_ext) &&
01073                        (desc1.ex_badpix == desc2.ex_badpix) &&
01074                        (desc1.fits_type == desc2.fits_type) &&
01075                        (desc1.frame_type == desc2.frame_type),
01076                        CPL_ERROR_ILLEGAL_INPUT,
01077                        "YCAL isn't in the correct format!!!");
01078         kmo_free_fits_desc(&desc2);
01079         kmo_init_fits_desc(&desc2);
01080 
01081         KMO_TRY_EXIT_IF_NULL(
01082             lcal_frame = kmo_dfs_get_frame(frameset, LCAL));
01083         desc2 = kmo_identify_fits_header(cpl_frame_get_filename(lcal_frame));
01084         KMO_TRY_CHECK_ERROR_STATE();
01085 
01086         KMO_TRY_ASSURE((desc2.nr_ext % KMOS_NR_DETECTORS == 0) &&
01087                        (desc1.ex_badpix == desc2.ex_badpix) &&
01088                        (desc1.fits_type == desc2.fits_type) &&
01089                        (desc1.frame_type == desc2.frame_type),
01090                        CPL_ERROR_ILLEGAL_INPUT,
01091                        "LCAL isn't in the correct format!!!");
01092         kmo_free_fits_desc(&desc2);
01093         kmo_init_fits_desc(&desc2);
01094 
01095         KMO_TRY_EXIT_IF_NULL(
01096             flat_frame = kmo_dfs_get_frame(frameset, MASTER_FLAT));
01097         desc2 = kmo_identify_fits_header(cpl_frame_get_filename(flat_frame));
01098         KMO_TRY_CHECK_ERROR_STATE();
01099 
01100         KMO_TRY_ASSURE((desc2.nr_ext % (2*KMOS_NR_DETECTORS) == 0) &&
01101                        (desc1.ex_badpix == desc2.ex_badpix) &&
01102                        (desc1.fits_type == desc2.fits_type) &&
01103                        (desc1.frame_type == desc2.frame_type),
01104                        CPL_ERROR_ILLEGAL_INPUT,
01105                        "MASTER_FLAT isn't in the correct format!!!");
01106         kmo_free_fits_desc(&desc2);
01107         kmo_init_fits_desc(&desc2);
01108 
01109         if (has_illum_corr) {
01110             KMO_TRY_EXIT_IF_NULL(
01111                 illum_frame = kmo_dfs_get_frame(frameset, ILLUM_CORR));
01112             desc2 = kmo_identify_fits_header(
01113                         cpl_frame_get_filename(illum_frame));
01114             KMO_TRY_CHECK_ERROR_STATE();
01115             KMO_TRY_ASSURE(((desc2.nr_ext == 24) || (desc2.nr_ext == 48)) &&
01116                            (desc2.ex_badpix == FALSE) &&
01117                            (desc2.fits_type == f2i_fits) &&
01118                            (desc2.frame_type == ifu_frame),
01119                            CPL_ERROR_ILLEGAL_INPUT,
01120                            "ILLUM_CORR isn't in the correct format!!!");
01121             kmo_free_fits_desc(&desc2);
01122             kmo_init_fits_desc(&desc2);
01123         }
01124 
01125         if (has_telluric) {
01126             KMO_TRY_EXIT_IF_NULL(
01127                 telluric_frame = kmo_dfs_get_frame(frameset, TELLURIC));
01128             desc2 = kmo_identify_fits_header(
01129                         cpl_frame_get_filename(telluric_frame));
01130             KMO_TRY_CHECK_ERROR_STATE();
01131             KMO_TRY_ASSURE(((desc2.nr_ext == 24) || (desc2.nr_ext == 48)) &&
01132                            (desc2.ex_badpix == FALSE) &&
01133                            (desc2.fits_type == f1i_fits) &&
01134                            (desc2.frame_type == ifu_frame),
01135                            CPL_ERROR_ILLEGAL_INPUT,
01136                            "TELLURIC isn't in the correct format!!!");
01137             kmo_free_fits_desc(&desc2);
01138             kmo_init_fits_desc(&desc2);
01139         }
01140 
01141         if (cpl_frameset_count_tags(frameset, OH_SPEC) != 0) {
01142             KMO_TRY_EXIT_IF_NULL(
01143                     ref_spectrum_frame = kmo_dfs_get_frame(frameset, OH_SPEC));
01144         }
01145 
01146         KMO_TRY_EXIT_IF_NULL(
01147             tmp_frame = kmo_dfs_get_frame(frameset, SCIENCE));
01148         while (tmp_frame != NULL ) {
01149             desc2 = kmo_identify_fits_header(cpl_frame_get_filename(tmp_frame));
01150             KMO_TRY_CHECK_ERROR_STATE();
01151             KMO_TRY_ASSURE((desc2.nr_ext == 3) &&
01152                            (desc2.ex_badpix == FALSE) &&
01153                            (desc2.fits_type == raw_fits) &&
01154                            (desc2.frame_type == detector_frame),
01155                            CPL_ERROR_ILLEGAL_INPUT,
01156                            "SCIENCE isn't in the correct format!!!");
01157             kmo_free_fits_desc(&desc2);
01158             kmo_init_fits_desc(&desc2);
01159 
01160             if (mapping_mode == NULL) {
01161                 KMO_TRY_EXIT_IF_NULL(
01162                     tmp_header =
01163                           kmclipm_propertylist_load(
01164                                          cpl_frame_get_filename(tmp_frame), 0));
01165                 if (cpl_propertylist_has(tmp_header, TPL_ID)) {
01166                     KMO_TRY_EXIT_IF_NULL(
01167                         tmp_str = cpl_propertylist_get_string(tmp_header,
01168                                                               TPL_ID));
01169 
01170                     if (strcmp(tmp_str, MAPPING8) == 0)
01171                     {
01172                         mapping_mode = cpl_sprintf("%s", "mapping8");
01173                     }
01174                     if (strcmp(tmp_str, MAPPING24) == 0)
01175                     {
01176                         mapping_mode = cpl_sprintf("%s", "mapping24");
01177                     }
01178                 }
01179                 cpl_propertylist_delete(tmp_header); tmp_header = NULL;
01180             }
01181 
01182             tmp_frame = kmo_dfs_get_frame(frameset, NULL);
01183             KMO_TRY_CHECK_ERROR_STATE();
01184         }
01185 
01186         if (mapping_mode != NULL) {
01187             // we are in mapping mode
01188             if ((ifus != NULL) || (strcmp(name, "") != 0))
01189             {
01190                 cpl_msg_warning("","The SCIENCE frames have been taken in one of the "
01191                                    "mapping modes AND specific IFUs have been "
01192                                    "specified! --> Only processing these!");
01193             } else {
01194                 if (strcmp(smethod, "BCS") == 0) {
01195                     extrapol_enum = BCS_NATURAL;
01196                     cpl_msg_info("","Detected frames taken in mapping mode. "
01197                                     "Changing extrapolation mode to TRUE.");
01198                 }
01199             }
01200             if (fast_mode) {
01201                 cpl_msg_info("", "Creating map in fast_mode.");
01202             }
01203         } else {
01204             if (fast_mode) {
01205                 cpl_msg_info("", "fast_mode has been selected but we aren't in "
01206                              "mapping mode. So your choise for fast_mode is ignored.");
01207             }
01208         }
01209 
01210         KMO_TRY_EXIT_IF_NULL(
01211             suffix = kmo_dfs_get_suffix(xcal_frame, TRUE, TRUE));
01212 
01213         KMO_TRY_EXIT_IF_ERROR(
01214             kmo_check_frame_setup_md5_xycal(frameset));
01215         KMO_TRY_EXIT_IF_ERROR(
01216             kmo_check_frame_setup_md5(frameset));
01217         KMO_TRY_EXIT_IF_ERROR(
01218             kmo_check_frame_setup_sampling(frameset));
01219 
01220         cpl_msg_info("", "Detected instrument setup:   %s", suffix+1);
01221         cpl_msg_info("", "(grating 1, 2 & 3, rotation angle)");
01222 
01223         //
01224         // check which IFUs are active for all frames
01225         //
01226         KMO_TRY_EXIT_IF_NULL(
01227             unused_ifus = kmo_get_unused_ifus(frameset, 1, 1));
01228 
01229 //        KMO_TRY_EXIT_IF_NULL(
01230 //            unused_ifus_after = kmo_duplicate_unused_ifus(unused_ifus_before));
01231 
01232         kmo_print_unused_ifus(unused_ifus, FALSE);
01233 
01234         //
01235         // get bounds, setup grid, setup obj_sky-struct
01236         //
01237 
01238         // get left and right bounds of IFUs
01239         KMO_TRY_EXIT_IF_NULL(
01240             tmp_header =
01241                   kmclipm_propertylist_load(cpl_frame_get_filename(xcal_frame),
01242                                             0));
01243         KMO_TRY_EXIT_IF_NULL(
01244             bounds = kmclipm_extract_bounds(tmp_header));
01245         cpl_propertylist_delete(tmp_header); tmp_header = NULL;
01246 
01247         // setup grid definition, wavelength start and end points will be set
01248         // in the detector loop
01249         KMO_TRY_EXIT_IF_ERROR(
01250             kmclipm_setup_grid(&gd, imethod, neighborhoodRange, pix_scale));
01251 
01252         // get valid STD frames with objects in it and associated sky exposures
01253         KMO_TRY_EXIT_IF_NULL(
01254             obj_sky_struct = kmo_get_obj_sky_frame_table(frameset,
01255                                                          &nr_science_frames,
01256                                                          SCIENCE,
01257                                                          no_subtract));
01258         //
01259         // load lcal-frames
01260         //
01261         KMO_TRY_EXIT_IF_NULL(
01262             lcal = (cpl_image**)
01263                    cpl_calloc(KMOS_NR_DETECTORS, sizeof(cpl_image*)));
01264         for (int i = 0; i < KMOS_NR_DETECTORS; i++) {
01265             KMO_TRY_EXIT_IF_NULL(
01266                 lcal[i] = kmo_dfs_load_image(frameset, LCAL, i+1, FALSE, FALSE, NULL));
01267         }
01268 
01269         //
01270         // allocate intermediate memory
01271         //
01272         KMO_TRY_EXIT_IF_NULL(
01273             all_obj = (allObjStruct*)cpl_calloc(nr_science_frames * KMOS_NR_IFUS, sizeof(allObjStruct)));
01274 
01275         // initialize intermediate memory
01276         for (int i = 0; i < nr_science_frames * KMOS_NR_IFUS; i++) {
01277             all_obj[i].name = NULL;
01278             all_obj[i].count = 0;
01279         }
01280 
01281         nr_data_alloc = KMOS_NR_IFUS;
01282         KMO_TRY_EXIT_IF_NULL(
01283             cube_data =  (cpl_imagelist**)cpl_calloc(nr_data_alloc,
01284                                                      sizeof(cpl_imagelist*)));
01285         KMO_TRY_EXIT_IF_NULL(
01286             cube_noise = (cpl_imagelist**)cpl_calloc(nr_data_alloc,
01287                                                      sizeof(cpl_imagelist*)));
01288         KMO_TRY_EXIT_IF_NULL(
01289             header_data =  (cpl_propertylist**)cpl_calloc(nr_data_alloc,
01290                                                      sizeof(cpl_propertylist*)));
01291         KMO_TRY_EXIT_IF_NULL(
01292             header_noise = (cpl_propertylist**)cpl_calloc(nr_data_alloc,
01293                                                      sizeof(cpl_propertylist*)));
01294 
01295         if (cpl_frameset_count_tags(frameset, SCIENCE) == 1) {
01296             no_combine = TRUE;
01297             cpl_msg_info("", "--no_combine has been set to TRUE since there is only one SCIENCE frame!");
01298         }
01299 
01300         if (no_subtract) {
01301             no_combine = TRUE;
01302             cpl_msg_info("", "--no_combine has been set to TRUE since --no_subtract has been specified by the user!");
01303             cpl_msg_info("", "Combining cubes would combine skies and obejcts which is meaningless.");
01304             cpl_msg_info("", "This can be done manually with the recipe kmo_combine afterwards!");
01305         }
01306 
01307         //
01308         // loop all science frames containing at least one object
01309         //
01310         cpl_msg_info("", "-------------------------------------------");
01311         cpl_msg_info("", "Reconstructing & saving cubes containing objects");
01312         cpl_msg_info("", " ");
01313         for (int sf = 0; sf < nr_science_frames; sf++) {
01314             KMO_TRY_EXIT_IF_NULL(
01315                 fn_obj = cpl_frame_get_filename(
01316                                                obj_sky_struct[sf].objectFrame));
01317 
01318             KMO_TRY_EXIT_IF_NULL(
01319                 main_header = kmclipm_propertylist_load(fn_obj, 0));
01320 
01321             // get ifu and detector number to work at
01322             actual_msg_level = cpl_msg_get_level();
01323             user_defined_ifu = 0;
01324             if ((ifus != NULL) || (strcmp(name, "") != 0)) {
01325                 if (ifus != NULL) {
01326                     // user specified IFUs
01327                     user_defined_ifu = cpl_vector_get(ifus, sf);
01328                     KMO_TRY_CHECK_ERROR_STATE();
01329 
01330                     if (user_defined_ifu < 1) {
01331                         user_defined_ifu = -1;
01332                     }
01333                 } else {
01334                     // user specified an object name
01335                     cpl_msg_set_level(CPL_MSG_OFF);
01336                     user_defined_ifu =
01337                         kmo_get_index_from_ocs_name(obj_sky_struct[sf].objectFrame,
01338                                                     name);
01339                     cpl_msg_set_level(actual_msg_level);
01340                     if (user_defined_ifu == -1) {
01341                         cpl_error_reset();
01342                     }
01343                     KMO_TRY_CHECK_ERROR_STATE();
01344                 }
01345             }
01346 
01347             //
01348             // reconstruct science frame
01349             //
01350             cpl_msg_info("", "   > processing frame: %s", fn_obj);
01351             for (int ifu_nr = 1; ifu_nr <= KMOS_NR_IFUS; ifu_nr++)
01352             {
01353                 det_nr = (ifu_nr - 1)/KMOS_IFUS_PER_DETECTOR + 1;
01354 
01355                 KMO_TRY_ASSURE((det_nr >= 1) &&
01356                                (det_nr <= KMOS_NR_DETECTORS),
01357                                CPL_ERROR_ILLEGAL_INPUT,
01358                                "The provided ifu-numbers are incorrect! They "
01359                                "must be between 1 and %d", KMOS_NR_IFUS);
01360 
01361                 KMO_TRY_EXIT_IF_NULL(
01362                     punused_ifus = cpl_array_get_data_int_const(
01363                                               unused_ifus[det_nr-1]));
01364 
01365                 // get subheader data
01366                 KMO_TRY_EXIT_IF_NULL(
01367                     header_data[ifu_nr-1] = kmclipm_propertylist_load(fn_obj,
01368                                                                       det_nr));
01369                 KMO_TRY_EXIT_IF_NULL(
01370                     extname = kmo_extname_creator(ifu_frame, ifu_nr,
01371                                                   EXT_DATA));
01372                 KMO_TRY_EXIT_IF_ERROR(
01373                     kmclipm_update_property_string(header_data[ifu_nr-1],
01374                                                    EXTNAME, extname,
01375                                                    "FITS extension name"));
01376                 cpl_free(extname); extname = NULL;
01377 
01378                 // Search for keyword ESO OCS ARMi NOTUSED
01379                 // If not present (CPL_ERROR_DATA_NOT_FOUND), do nothing
01380                 KMO_TRY_EXIT_IF_NULL(
01381                     keyword = cpl_sprintf("%s%d%s", IFU_VALID_PREFIX, ifu_nr,
01382                                           IFU_VALID_POSTFIX));
01383                     tmp_str = cpl_propertylist_get_string(main_header, keyword);
01384                 cpl_free(keyword); keyword = NULL;
01385 
01386                 if ((cpl_error_get_code() == CPL_ERROR_DATA_NOT_FOUND) &&
01387                     (bounds[2*(ifu_nr-1)] != -1) &&
01388                     (bounds[2*(ifu_nr-1)+1] != -1) &&
01389                     ((obj_sky_struct[sf].skyframes[ifu_nr-1] != NULL) || (cpl_frameset_count_tags(frameset, SCIENCE) == 1) || no_subtract) &&
01390                     (punused_ifus[(ifu_nr-1) % KMOS_IFUS_PER_DETECTOR] == 0) &&
01391                     ((user_defined_ifu == 0) || (user_defined_ifu == ifu_nr)))
01392                 {
01393                     // IFU is valid
01394                     cpl_error_reset();
01395 
01396                     if ((obj_sky_struct[sf].skyframes[ifu_nr-1] != NO_CORRESPONDING_SKYFRAME) &&
01397                         (cpl_frameset_count_tags(frameset, SCIENCE) != 1) &&
01398                         !no_subtract)
01399                     {
01400                         if (no_subtract) {
01401                             sky_frame = NULL;
01402                             cpl_msg_warning("", "      > Omit sky subtraction on IFU %d", ifu_nr);
01403                         } else {
01404                             sky_frame = obj_sky_struct[sf].skyframes[ifu_nr-1];
01405                             KMO_TRY_EXIT_IF_NULL(
01406                                 fn_sky = cpl_frame_get_filename(sky_frame));
01407                             cpl_msg_info("", "      > IFU %d (with sky in frame: %s)", ifu_nr, fn_sky);
01408                         }
01409                     } else {
01410                         sky_frame = NULL;
01411                         if ((cpl_frameset_count_tags(frameset, SCIENCE) != 1) && (!no_subtract)) {
01412                             cpl_msg_warning("", "      > IFU %d with no corresponding sky frame", ifu_nr);
01413                         }
01414                     }
01415 
01416                     // get filter for this detector and setup grid definition
01417                     // ESO INS FILTi ID
01418                     char *tmp_band_method = getenv("KMO_BAND_METHOD");
01419                     int band_method = 0;
01420                     if (tmp_band_method != NULL) {
01421                         band_method = atoi(tmp_band_method);
01422                     }
01423 
01424                     KMO_TRY_EXIT_IF_NULL(
01425                         keyword = cpl_sprintf("%s%d%s",
01426                                               IFU_FILTID_PREFIX, det_nr,
01427                                               IFU_FILTID_POSTFIX));
01428                     KMO_TRY_EXIT_IF_NULL(
01429                         filter_id = cpl_propertylist_get_string(main_header,
01430                                                                 keyword));
01431                     cpl_free(keyword); keyword = NULL;
01432 
01433                     if (print_once) {
01434                         cpl_msg_set_level(CPL_MSG_WARNING);
01435                     }
01436 
01437                     KMO_TRY_EXIT_IF_NULL(
01438                         band_table = kmo_dfs_load_table(frameset, WAVE_BAND,
01439                                                         1, 0));
01440                     KMO_TRY_EXIT_IF_ERROR(
01441                         kmclipm_setup_grid_band_lcal(&gd, lcal[det_nr-1],
01442                                                      filter_id, band_method,
01443                                                      band_table));
01444                     cpl_table_delete(band_table); band_table = NULL;
01445 
01446                     print_once = TRUE;
01447                     cpl_msg_set_level(actual_msg_level);
01448 
01449                     //
01450                     // calc WCS & update subheader
01451                     //
01452                     KMO_TRY_EXIT_IF_ERROR(
01453                         kmo_calc_wcs_gd(main_header, header_data[ifu_nr-1], ifu_nr, gd));
01454 
01455                     KMO_TRY_EXIT_IF_ERROR(
01456                         kmclipm_update_property_int(header_data[ifu_nr-1],
01457                                                     NAXIS, 3,
01458                                                     "number of data axes"));
01459                     KMO_TRY_EXIT_IF_ERROR(
01460                         kmclipm_update_property_int(header_data[ifu_nr-1],
01461                                                     NAXIS1, gd.x.dim,
01462                                                     "length of data axis 1"));
01463                     KMO_TRY_EXIT_IF_ERROR(
01464                         kmclipm_update_property_int(header_data[ifu_nr-1],
01465                                                     NAXIS2, gd.y.dim,
01466                                                     "length of data axis 2"));
01467                     KMO_TRY_EXIT_IF_ERROR(
01468                         kmclipm_update_property_int(header_data[ifu_nr-1],
01469                                                     NAXIS3, gd.l.dim,
01470                                                     "length of data axis 3"));
01471 
01472                     //
01473                     // reconstruct
01474                     //
01475 
01476                     if (ref_spectrum_frame == NULL) { // no lambda correction using OH lines
01477 
01478                         KMO_TRY_EXIT_IF_ERROR(
01479                                 kmo_reconstruct_sci(ifu_nr,
01480                                         bounds[2*(ifu_nr-1)],
01481                                         bounds[2*(ifu_nr-1)+1],
01482                                         obj_sky_struct[sf].objectFrame,
01483                                         SCIENCE,
01484                                         sky_frame,
01485                                         SCIENCE,
01486                                         flat_frame,
01487                                         xcal_frame,
01488                                         ycal_frame,
01489                                         lcal_frame,
01490                                         NULL,
01491                                         &gd,
01492                                         &cube_data[ifu_nr-1],
01493                                         &cube_noise[ifu_nr-1],
01494                                         flux,
01495                                         background,
01496                                         xcal_interpolation));
01497 
01498                     } else {                    //  lambda correction using OH lines
01499 
01500                         KMO_TRY_EXIT_IF_ERROR(
01501                                 kmo_reconstruct_sci(ifu_nr,
01502                                         bounds[2*(ifu_nr-1)],
01503                                         bounds[2*(ifu_nr-1)+1],
01504                                         obj_sky_struct[sf].objectFrame,
01505                                         SCIENCE,
01506                                         NULL,
01507                                         NULL,
01508                                         flat_frame,
01509                                         xcal_frame,
01510                                         ycal_frame,
01511                                         lcal_frame,
01512                                         NULL,
01513                                         &gd,
01514                                         &cube_data[ifu_nr-1],
01515                                         &cube_noise[ifu_nr-1],
01516                                         FALSE,
01517                                         FALSE,
01518                                         xcal_interpolation));
01519 
01520 
01521                         if (cube_data[ifu_nr-1] != NULL) {
01522                             if (peaks == NULL) {
01523                                 KMO_TRY_EXIT_IF_NULL(
01524                                         range = cpl_vector_new(2));
01525                                 KMO_TRY_EXIT_IF_ERROR(
01526                                         cpl_vector_set(range, 0, gd.l.start));
01527                                 KMO_TRY_EXIT_IF_ERROR(
01528                                         cpl_vector_set(range, 1, gd.l.start + gd.l.delta * gd.l.dim));
01529                                 KMO_TRY_EXIT_IF_NULL(
01530                                         ref_spectrum = kmo_lcorr_read_reference_spectrum(
01531                                                 cpl_frame_get_filename(ref_spectrum_frame),NULL));
01532                                 KMO_TRY_EXIT_IF_NULL(
01533                                         peaks = kmo_lcorr_get_peak_lambdas(ref_spectrum, 0.2, range));
01534 
01535                                 cpl_vector_delete(range);
01536                             }
01537 
01538                             KMO_TRY_EXIT_IF_NULL(
01539                                     obj_spectrum = kmo_lcorr_extract_spectrum(
01540                                             cube_data[ifu_nr-1], header_data[ifu_nr-1], 0.8, NULL));
01541 
01542                             KMO_TRY_EXIT_IF_NULL(
01543                                     lcorr_coeffs = kmo_lcorr_crosscorrelate_spectra(
01544                                             obj_spectrum, ref_spectrum, peaks, filter_id));
01545 
01546                             cpl_bivector_delete(obj_spectrum);
01547 
01548                             const int format_width = 14;
01549                             const int max_coeffs = 6;
01550                             char *coeff_string = NULL;
01551                             char coeff_dump[format_width * max_coeffs + 1];
01552                             cpl_size pows[1];
01553                             coeff_dump[0] = 0;
01554                             for (int ic=0; ic<=cpl_polynomial_get_degree(lcorr_coeffs) && ic < max_coeffs; ic++) {
01555                                 pows[0] = ic;
01556                                 coeff_string = cpl_sprintf(" %*g,",
01557                                         format_width-2, cpl_polynomial_get_coeff(lcorr_coeffs,pows));
01558                                 strncat(coeff_dump, coeff_string, format_width);
01559                                 cpl_free(coeff_string);
01560                             }
01561                             cpl_msg_debug("","Lambda correction coeffs for IFU %d %s",ifu_nr, coeff_dump);
01562 
01563                             cpl_imagelist_delete(cube_data[ifu_nr-1]); cube_data[ifu_nr-1] = NULL;
01564                             if (cube_noise[ifu_nr-1] != NULL) {
01565                                 cpl_imagelist_delete(cube_noise[ifu_nr-1]); cube_noise[ifu_nr-1] = NULL;
01566                             }
01567                             KMO_TRY_EXIT_IF_ERROR(
01568                                     kmo_reconstruct_sci(ifu_nr,
01569                                             bounds[2*(ifu_nr-1)],
01570                                             bounds[2*(ifu_nr-1)+1],
01571                                             obj_sky_struct[sf].objectFrame,
01572                                             SCIENCE,
01573                                             sky_frame,
01574                                             SCIENCE,
01575                                             flat_frame,
01576                                             xcal_frame,
01577                                             ycal_frame,
01578                                             lcal_frame,
01579                                             lcorr_coeffs,
01580                                             &gd,
01581                                             &cube_data[ifu_nr-1],
01582                                             &cube_noise[ifu_nr-1],
01583                                             flux,
01584                                             background,
01585                                             xcal_interpolation));
01586 
01587                             cpl_polynomial_delete(lcorr_coeffs); lcorr_coeffs = NULL;
01588                         }
01589                     }
01590 
01591                     // scale flux according to pixel_scale
01592                     KMO_TRY_EXIT_IF_NULL(
01593                         tmpImg = cpl_imagelist_get(cube_data[ifu_nr-1], 0));
01594                     double scaling = (cpl_image_get_size_x(tmpImg)*cpl_image_get_size_y(tmpImg)) /
01595                                      (KMOS_SLITLET_X*KMOS_SLITLET_Y);
01596                     KMO_TRY_EXIT_IF_ERROR(
01597                         cpl_imagelist_divide_scalar(cube_data[ifu_nr-1], scaling));
01598                     if (cube_noise[ifu_nr-1] != NULL) {
01599                         KMO_TRY_EXIT_IF_ERROR(
01600                             cpl_imagelist_divide_scalar(cube_noise[ifu_nr-1], scaling));
01601                     }
01602 
01603                     //
01604                     // divide cube by telluric correction
01605                     //
01606                     if (has_telluric) {
01607                         telluric_data = kmo_tweak_load_telluric(frameset, ifu_nr, FALSE,
01608                                                                 cube_noise[ifu_nr-1] != NULL);
01609                         KMO_TRY_CHECK_ERROR_STATE();
01610                         if (telluric_data != NULL) {
01611                             // load noise if present
01612                             telluric_noise = kmo_tweak_load_telluric(frameset, ifu_nr, TRUE,
01613                                                                      cube_noise[ifu_nr-1] != NULL);
01614                             KMO_TRY_CHECK_ERROR_STATE();
01615 
01616                             KMO_TRY_EXIT_IF_ERROR(
01617                                 kmo_arithmetic_3D_1D(
01618                                         cube_data[ifu_nr-1], telluric_data,
01619                                         cube_noise[ifu_nr-1], telluric_noise, "/"));
01620 
01621                             kmclipm_vector_delete(telluric_data);
01622                             telluric_data = NULL;
01623                             kmclipm_vector_delete(telluric_noise);
01624                             telluric_noise = NULL;
01625                         }
01626                     }
01627 
01628                     //
01629                     // divide cube by illumination correction
01630                     //
01631                     if (has_illum_corr) {
01632                         cpl_msg_set_level(CPL_MSG_OFF);
01633                         illum_data = kmo_dfs_load_image(frameset, ILLUM_CORR,
01634                                                         ifu_nr, FALSE, FALSE, NULL);
01635                         cpl_msg_set_level(actual_msg_level);
01636                         if (cpl_error_get_code() != CPL_ERROR_NONE) {
01637                             cpl_msg_warning("","No illumination correction for IFU %d available! "
01638                                             "Proceeding anyway.", ifu_nr);
01639                             cpl_error_reset();
01640                         } else {
01641                             cpl_msg_set_level(CPL_MSG_OFF);
01642                             illum_noise = kmo_dfs_load_image(frameset,
01643                                                              ILLUM_CORR,
01644                                                              ifu_nr, TRUE,
01645                                                              FALSE, NULL);
01646                             cpl_msg_set_level(actual_msg_level);
01647                             if (cpl_error_get_code() != CPL_ERROR_NONE) {
01648                                 cpl_msg_warning("","No illumination correction for IFU %d "
01649                                                 "available! Proceeding anyway.", ifu_nr);
01650                                 cpl_image_delete(illum_data); illum_data = NULL;
01651                                 cpl_error_reset();
01652                             }
01653                         }
01654 
01655                         if (illum_data != NULL) {
01656                             KMO_TRY_EXIT_IF_ERROR(
01657                                 kmo_arithmetic_3D_2D(
01658                                             cube_data[ifu_nr-1], illum_data,
01659                                             cube_noise[ifu_nr-1], illum_noise, "/"));
01660                             cpl_image_delete(illum_data); illum_data = NULL;
01661                             cpl_image_delete(illum_noise); illum_noise = NULL;
01662                         }
01663                     }
01664 
01665                     // get object name and store if not already present
01666                     KMO_TRY_EXIT_IF_NULL(
01667                         keyword = cpl_sprintf("%s%d%s", IFU_NAME_PREFIX, ifu_nr,
01668                                               IFU_NAME_POSTFIX));
01669                     KMO_TRY_EXIT_IF_NULL(
01670                         tmp_str = cpl_propertylist_get_string(header_data[ifu_nr-1],
01671                                                               keyword));
01672                     cpl_free(keyword); keyword = NULL;
01673 
01674                     // found keyword, check if it is already in all_obj
01675                     found_name = 0;
01676                     for(int i = 0; i < nr_avail_obj_names; i++) {
01677                         if (strcmp(all_obj[i].name, tmp_str) == 0) {
01678                             found_name = TRUE;
01679                             all_obj[i].count++;
01680                             break;
01681                         }
01682                     }
01683                     if (!found_name) {
01684                         all_obj[nr_avail_obj_names].count++;
01685                         KMO_TRY_EXIT_IF_NULL(
01686                             all_obj[nr_avail_obj_names++].name =
01687                                                     cpl_sprintf("%s", tmp_str));
01688                     }
01689                 } else {
01690                     cpl_error_reset();
01691 
01692                     // IFU is invalid
01693                 }
01694 
01695                 // duplicate subheader data
01696                 KMO_TRY_EXIT_IF_NULL(
01697                     header_noise[ifu_nr-1] =
01698                          cpl_propertylist_duplicate(header_data[ifu_nr-1]));
01699 
01700                 KMO_TRY_EXIT_IF_NULL(
01701                     extname = kmo_extname_creator(ifu_frame, ifu_nr, EXT_NOISE));
01702                 KMO_TRY_EXIT_IF_ERROR(
01703                     kmclipm_update_property_string(header_noise[ifu_nr-1],
01704                                                    EXTNAME, extname,
01705                                                    "FITS extension name"));
01706                 cpl_free(extname); extname = NULL;
01707             } // end for ifu_nr
01708 
01709             //
01710             // count number of reconstructed data- and noise-cubes
01711             //
01712             for (int ifu_nr = 1; ifu_nr <= nr_data_alloc; ifu_nr++) {
01713                 if (cube_data[ifu_nr-1] != NULL) {
01714                     cube_counter_data++;
01715                 }
01716                 if (cube_noise[ifu_nr-1] != NULL) {
01717                     cube_counter_noise++;
01718                 }
01719             }
01720 
01721             //
01722             // save reconstructed cubes of science frame
01723             //
01724             if (cube_counter_data > 0) {
01725                 cpl_msg_info("", "   > saving...");
01726 
01727                 if (!suppress_extension) {
01728                     fn_out = fn_obj;
01729 
01730                     int nr_found = 0;
01731                     // remove any path-elements from filename and use it as
01732                     // suffix
01733                     split = kmo_strsplit(fn_out, "/", &nr_found);
01734 
01735                     fn_suffix = cpl_sprintf("_%s", split[nr_found-1]);
01736                     kmo_strfreev(split);
01737 
01738                     // remove '.fits' at the end if there is any
01739                     char *fff = fn_suffix;
01740                     fff += strlen(fn_suffix)-5;
01741                     if (strcmp(fff, ".fits") == 0) {
01742                         fn_suffix[strlen(fn_suffix)-5] = '\0';
01743                     }
01744                 } else {
01745                     KMO_TRY_EXIT_IF_NULL(
01746                         fn_suffix = cpl_sprintf("_%d", suppress_index++));
01747                 }
01748 
01749                 fn_out = RECONSTRUCTED_CUBE;
01750 
01751                 KMO_TRY_EXIT_IF_ERROR(
01752                     kmo_dfs_save_main_header(frameset, fn_out, fn_suffix,
01753                                              obj_sky_struct[sf].objectFrame,
01754                                              NULL, parlist, cpl_func));
01755 
01756                 for (int ifu_nr = 1; ifu_nr <= KMOS_NR_IFUS; ifu_nr++) {
01757                     KMO_TRY_EXIT_IF_ERROR(
01758                         kmo_dfs_save_cube(cube_data[ifu_nr-1], fn_out,
01759                                           fn_suffix, header_data[ifu_nr-1], 0./0.));
01760 
01761                     if (cube_counter_noise > 0) {
01762                         KMO_TRY_EXIT_IF_ERROR(
01763                             kmo_dfs_save_cube(cube_noise[ifu_nr-1], fn_out,
01764                                               fn_suffix, header_noise[ifu_nr-1],
01765                                               0./0.));
01766                     }
01767 
01768                     cpl_imagelist_delete(cube_data[ifu_nr-1]); cube_data[ifu_nr-1] = NULL;
01769                     cpl_imagelist_delete(cube_noise[ifu_nr-1]); cube_noise[ifu_nr-1] = NULL;
01770                     cpl_propertylist_delete(header_data[ifu_nr-1]); header_data[ifu_nr-1] = NULL;
01771                     cpl_propertylist_delete(header_noise[ifu_nr-1]); header_noise[ifu_nr-1] = NULL;
01772                 } // end for ifu_nr
01773                 cpl_free(fn_suffix); fn_suffix = NULL;
01774             } else {
01775                 cpl_msg_info("", "   > all IFUs invalid, don't save");
01776             } // if (cube_counter_data > 0) {
01777 
01778             cpl_propertylist_delete(main_header); main_header = NULL;
01779         } // end for sf (nr_science_frames)
01780         cpl_free(cube_data);    cube_data = NULL;
01781         cpl_free(cube_noise);    cube_noise = NULL;
01782         cpl_free(header_data);    header_data = NULL;
01783         cpl_free(header_noise);    header_noise = NULL;
01784         cpl_msg_info("", "-------------------------------------------");
01785 
01786         if (lcal != NULL) {
01787             for (int i = 0; i < KMOS_NR_DETECTORS; i++) {
01788                 cpl_image_delete(lcal[i]);
01789             }
01790         }
01791         cpl_free(lcal); lcal = NULL;
01792 
01793         //
01794         // combine
01795         //
01796         suppress_index = 0;
01797         if (!no_combine) {
01798             cpl_msg_info("", "Combining reconstructed objects");
01799             cpl_msg_info("", " ");
01800 
01801             nr_reconstructed_frames = cpl_frameset_count_tags(frameset, RECONSTRUCTED_CUBE);
01802 
01803             if ( (mapping_mode == NULL) || ((mapping_mode != NULL) &&
01804                                             ((ifus != NULL) || (strcmp(name, "") != 0)))
01805                )
01806             {
01807                 // loop all available objects
01808                 for (int i = 0; i < nr_avail_obj_names; i++) {
01809                     cpl_msg_info("", "   > object: %s", all_obj[i].name);
01810                     nr_data_alloc = all_obj[i].count;
01811                     KMO_TRY_EXIT_IF_NULL(
01812                         cube_data =  (cpl_imagelist**)cpl_calloc(nr_data_alloc,
01813                                                            sizeof(cpl_imagelist*)));
01814                     KMO_TRY_EXIT_IF_NULL(
01815                         cube_noise = (cpl_imagelist**)cpl_calloc(nr_data_alloc,
01816                                                            sizeof(cpl_imagelist*)));
01817                     KMO_TRY_EXIT_IF_NULL(
01818                         header_data =  (cpl_propertylist**)cpl_calloc(nr_data_alloc,
01819                                                         sizeof(cpl_propertylist*)));
01820                     KMO_TRY_EXIT_IF_NULL(
01821                         header_noise = (cpl_propertylist**)cpl_calloc(nr_data_alloc,
01822                                                         sizeof(cpl_propertylist*)));
01823 
01824                     // setup cube-list and header-list for kmo_priv_combine()
01825                     cube_counter_data = 0;
01826                     cube_counter_noise = 0;
01827                     KMO_TRY_EXIT_IF_NULL(
01828                         tmp_frame = kmo_dfs_get_frame(frameset, RECONSTRUCTED_CUBE));
01829                     while (tmp_frame != NULL ) {
01830                         KMO_TRY_EXIT_IF_NULL(
01831                             fn_reconstr = cpl_frame_get_filename(tmp_frame));
01832 
01833                         KMO_TRY_EXIT_IF_NULL(
01834                             tmp_header = kmclipm_propertylist_load(fn_reconstr, 0));
01835 
01836                         kmo_free_fits_desc(&desc1);
01837                         kmo_init_fits_desc(&desc1);
01838                         desc1 = kmo_identify_fits_header(fn_reconstr);
01839 
01840                         for (int ifu_nr = 1; ifu_nr <= KMOS_NR_IFUS; ifu_nr++) {
01841                             // check if object-name equals the one in our list
01842                             KMO_TRY_EXIT_IF_NULL(
01843                                 keyword = cpl_sprintf("%s%d%s", IFU_NAME_PREFIX,
01844                                                       ifu_nr, IFU_NAME_POSTFIX));
01845                             KMO_TRY_EXIT_IF_NULL(
01846                                 tmp_str = cpl_propertylist_get_string(tmp_header,
01847                                                                       keyword));
01848                             cpl_free(keyword); keyword = NULL;
01849 
01850                             if (strcmp(all_obj[i].name, tmp_str) == 0) {
01851                                 // found object-IFU with matching name
01852                                 // load data & subheader
01853                                 index = kmo_identify_index(fn_reconstr, ifu_nr, FALSE);
01854                                 KMO_TRY_CHECK_ERROR_STATE();
01855 
01856                                 if (desc1.sub_desc[index-1].valid_data) {
01857                                     KMO_TRY_EXIT_IF_NULL(
01858                                         cube_data[cube_counter_data] =
01859                                             kmclipm_imagelist_load(fn_reconstr,
01860                                                                    CPL_TYPE_FLOAT,
01861                                                                    index));
01862                                     if (edge_nan) {
01863                                         KMO_TRY_EXIT_IF_ERROR(
01864                                             kmo_edge_nan(cube_data[cube_counter_data], ifu_nr));
01865                                     }
01866 
01867                                     KMO_TRY_EXIT_IF_NULL(
01868                                         header_data[cube_counter_data] =
01869                                             kmclipm_propertylist_load(fn_reconstr,
01870                                                                       index));
01871                                     cpl_propertylist_update_string(header_data[cube_counter_data],
01872                                                                    "ESO PRO FRNAME",
01873                                                                    fn_reconstr);
01874                                     cpl_propertylist_update_int(header_data[cube_counter_data],
01875                                                                 "ESO PRO IFUNR",
01876                                                                 index);
01877                                     cube_counter_data++;
01878                                 }
01879 
01880                                 // load noise & subheader (if existing)
01881                                 if (desc1.ex_noise) {
01882                                     index = kmo_identify_index(fn_reconstr, ifu_nr, TRUE);
01883                                     KMO_TRY_CHECK_ERROR_STATE();
01884 
01885                                     if (desc1.sub_desc[index-1].valid_data) {
01886                                         KMO_TRY_EXIT_IF_NULL(
01887                                             cube_noise[cube_counter_noise] =
01888                                                 kmclipm_imagelist_load(fn_reconstr,
01889                                                                        CPL_TYPE_FLOAT,
01890                                                                        index));
01891                                         if (edge_nan) {
01892                                             KMO_TRY_EXIT_IF_ERROR(
01893                                                 kmo_edge_nan(cube_noise[cube_counter_noise], ifu_nr));
01894                                         }
01895                                         KMO_TRY_EXIT_IF_NULL(
01896                                             header_noise[cube_counter_noise] =
01897                                                 kmclipm_propertylist_load(fn_reconstr,
01898                                                                           index));
01899                                         cube_counter_noise++;
01900                                     }
01901                                 }
01902                                 cpl_error_reset();
01903                             } // end if found obj
01904                         } // end for ifu_nr
01905 
01906                         cpl_propertylist_delete(tmp_header); tmp_header = NULL;
01907                         tmp_frame = kmo_dfs_get_frame(frameset, NULL);
01908                         KMO_TRY_CHECK_ERROR_STATE();
01909                     } // end while-loop RECONSTRUCTED_CUBE frames
01910 
01911                     if (cube_counter_data > 1) {
01912                         if (cube_counter_data == cube_counter_noise) {
01913                             KMO_TRY_EXIT_IF_ERROR(
01914                                 kmo_priv_combine(cube_data,
01915                                                  cube_noise,
01916                                                  header_data,
01917                                                  header_noise,
01918                                                  cube_counter_data,
01919                                                  cube_counter_noise,
01920                                                  all_obj[i].name,
01921                                                  "",
01922                                                  comb_method,
01923                                                  smethod,
01924                                                  fmethod,
01925                                                  filename,
01926                                                  cmethod,
01927                                                  cpos_rej,
01928                                                  cneg_rej,
01929                                                  citer,
01930                                                  cmin,
01931                                                  cmax,
01932                                                  extrapol_enum,
01933                                                  flux,
01934                                                  &combined_data,
01935                                                  &combined_noise));
01936                         } else if (cube_counter_noise == 0) {
01937                             // if imethod == "CS"
01938                             KMO_TRY_EXIT_IF_ERROR(
01939                                 kmo_priv_combine(cube_data,
01940                                                  NULL,
01941                                                  header_data,
01942                                                  header_noise,
01943                                                  cube_counter_data,
01944                                                  cube_counter_noise,
01945                                                  all_obj[i].name,
01946                                                  "",
01947                                                  comb_method,
01948                                                  smethod,
01949                                                  fmethod,
01950                                                  filename,
01951                                                  cmethod,
01952                                                  cpos_rej,
01953                                                  cneg_rej,
01954                                                  citer,
01955                                                  cmin,
01956                                                  cmax,
01957                                                  extrapol_enum,
01958                                                  flux,
01959                                                  &combined_data,
01960                                                  &combined_noise));
01961                         } else {
01962                             KMO_TRY_ASSURE(1 == 0,
01963                                            CPL_ERROR_ILLEGAL_INPUT,
01964                                            "The number of cube-data and cube-noise "
01965                                            "isn't the same (%d vs. %d)!",
01966                                            cube_counter_data, cube_counter_noise);
01967                         }
01968                     } else {
01969                         cpl_msg_warning("", "There is only one reconstructed cube with "
01970                                         "this object! Saving it as it is.");
01971                         KMO_TRY_EXIT_IF_NULL(
01972                             combined_data = cpl_imagelist_duplicate(cube_data[0]));
01973 
01974                         if (cube_noise[0] != NULL) {
01975                             KMO_TRY_EXIT_IF_NULL(
01976                                 combined_noise = cpl_imagelist_duplicate(cube_noise[0]));
01977                         }
01978                     } // end if (cube_counter_data > 1)
01979 
01980                     fn_out = COMBINED_CUBE;
01981                     if (!suppress_extension) {
01982                         KMO_TRY_EXIT_IF_NULL(
01983                             fn_suffix = cpl_sprintf("_%s", all_obj[i].name));
01984                     } else {
01985                         KMO_TRY_EXIT_IF_NULL(
01986                             fn_suffix = cpl_sprintf("_%d", suppress_index++));
01987                     }
01988 
01989                     // save combined cube
01990                     KMO_TRY_EXIT_IF_NULL(
01991                         tmp_frame = kmo_dfs_get_frame(frameset, RECONSTRUCTED_CUBE));
01992                     KMO_TRY_EXIT_IF_ERROR(
01993                         kmo_dfs_save_main_header(frameset, fn_out, fn_suffix,
01994                                                  tmp_frame, NULL, parlist, cpl_func));
01995 
01996                     KMO_TRY_EXIT_IF_ERROR(
01997                         kmo_dfs_save_cube(combined_data, fn_out, fn_suffix,
01998                                           header_data[0], 0./0.));
01999 
02000     //                if (combined_noise != NULL) {
02001                         if (header_noise[0] == NULL) {
02002                             KMO_TRY_EXIT_IF_NULL(
02003                                 header_noise[0] =
02004                                      cpl_propertylist_duplicate(header_data[0]));
02005 
02006                             KMO_TRY_EXIT_IF_NULL(
02007                                 tmp_str = cpl_propertylist_get_string(header_data[0],
02008                                                                       EXTNAME));
02009                             KMO_TRY_EXIT_IF_ERROR(
02010                                 kmo_extname_extractor(tmp_str, &ft, &tmp_int, content));
02011                             KMO_TRY_EXIT_IF_NULL(
02012                                 extname = kmo_extname_creator(ifu_frame, tmp_int,
02013                                                               EXT_NOISE));
02014                             KMO_TRY_EXIT_IF_ERROR(
02015                                 kmclipm_update_property_string(header_noise[0],
02016                                                                EXTNAME, extname,
02017                                                                "FITS extension name"));
02018                             cpl_free(extname); extname = NULL;
02019                         }
02020                         KMO_TRY_EXIT_IF_ERROR(
02021                             kmo_dfs_save_cube(combined_noise, fn_out, fn_suffix,
02022                                               header_noise[0], 0./0.));
02023     //                }
02024 
02025                     for (int jj = 0; jj < nr_data_alloc; jj++) {
02026                         cpl_imagelist_delete(cube_data[jj]); cube_data[jj] = NULL;
02027                         cpl_imagelist_delete(cube_noise[jj]); cube_noise[jj] = NULL;
02028                         cpl_propertylist_delete(header_data[jj]); header_data[jj] = NULL;
02029                         cpl_propertylist_delete(header_noise[jj]); header_noise[jj] = NULL;
02030                     }
02031                     cpl_free(cube_data);    cube_data = NULL;
02032                     cpl_free(cube_noise);   cube_noise = NULL;
02033                     cpl_free(header_data);  header_data = NULL;
02034                     cpl_free(header_noise); header_noise = NULL;
02035                     cpl_free(fn_suffix); fn_suffix = NULL;
02036                     cpl_imagelist_delete(combined_data); combined_data = NULL;
02037                     cpl_imagelist_delete(combined_noise); combined_noise = NULL;
02038                 } // for i = nr_avail_obj_names
02039             } else {
02040                 // we are in mapping_mode
02041                 nr_data_alloc = nr_reconstructed_frames*KMOS_NR_IFUS;
02042                 KMO_TRY_EXIT_IF_NULL(
02043                     cube_data = (cpl_imagelist**)cpl_calloc(nr_data_alloc,
02044                                                         sizeof(cpl_imagelist*)));
02045                 KMO_TRY_EXIT_IF_NULL(
02046                     cube_noise = (cpl_imagelist**)cpl_calloc(nr_data_alloc,
02047                                                         sizeof(cpl_imagelist*)));
02048                 KMO_TRY_EXIT_IF_NULL(
02049                     header_data = (cpl_propertylist**)cpl_calloc( nr_data_alloc,
02050                                                         sizeof(cpl_propertylist*)));
02051                 KMO_TRY_EXIT_IF_NULL(
02052                     header_noise = (cpl_propertylist**)cpl_calloc(nr_data_alloc,
02053                                                         sizeof(cpl_propertylist*)));
02054 
02055                 cube_counter_data = 0;
02056                 cube_counter_noise = 0;
02057                 KMO_TRY_EXIT_IF_NULL(
02058                     tmp_frame = kmo_dfs_get_frame(frameset, RECONSTRUCTED_CUBE));
02059                 while (tmp_frame != NULL ) {
02060                     KMO_TRY_EXIT_IF_NULL(
02061                         fn_reconstr = cpl_frame_get_filename(tmp_frame));
02062 
02063                     KMO_TRY_EXIT_IF_NULL(
02064                         tmp_header = kmclipm_propertylist_load(fn_reconstr, 0));
02065 
02066                     kmo_free_fits_desc(&desc1);
02067                     kmo_init_fits_desc(&desc1);
02068                     desc1 = kmo_identify_fits_header(fn_reconstr);
02069                     for (int ifu_nr = 1; ifu_nr <= KMOS_NR_IFUS; ifu_nr++) {
02070                         index = kmo_identify_index(fn_reconstr, ifu_nr, FALSE);
02071                         KMO_TRY_CHECK_ERROR_STATE();
02072 
02073                         if (desc1.sub_desc[index-1].valid_data) {
02074                             KMO_TRY_EXIT_IF_NULL(
02075                                 cube_data[cube_counter_data] =
02076                                     kmclipm_imagelist_load(fn_reconstr, CPL_TYPE_FLOAT,
02077                                                            index));
02078                             if (edge_nan) {
02079                                 KMO_TRY_EXIT_IF_ERROR(
02080                                     kmo_edge_nan(cube_data[cube_counter_data], ifu_nr));
02081                             }
02082 
02083                             if (fast_mode) {
02084                                 KMO_TRY_EXIT_IF_NULL(
02085                                     tmpImg = cpl_imagelist_collapse_median_create(cube_data[cube_counter_data]));
02086                                 KMO_TRY_EXIT_IF_NULL(
02087                                     tmpCube = cpl_imagelist_new());
02088                                 KMO_TRY_EXIT_IF_ERROR(
02089                                     cpl_imagelist_set(tmpCube, tmpImg, 0));
02090                                 cpl_imagelist_delete(cube_data[cube_counter_data]);
02091                                 cube_data[cube_counter_data] = tmpCube;
02092                             }
02093 
02094                             KMO_TRY_EXIT_IF_NULL(
02095                                 header_data[cube_counter_data] =
02096                                     kmclipm_propertylist_load(fn_reconstr, index));
02097                             cpl_propertylist_update_string(header_data[cube_counter_data],
02098                                                     "ESO PRO FRNAME",
02099                                                     fn_reconstr);
02100                             cpl_propertylist_update_int(header_data[cube_counter_data],
02101                                                     "ESO PRO IFUNR",
02102                                                     index);
02103                             cube_counter_data++;
02104                         }
02105 
02106                         // load noise & subheader (if existing)
02107                         if (desc1.ex_noise) {
02108                             index = kmo_identify_index(fn_reconstr, ifu_nr, TRUE);
02109                             KMO_TRY_CHECK_ERROR_STATE();
02110                             if (desc1.sub_desc[index-1].valid_data) {
02111                                 KMO_TRY_EXIT_IF_NULL(
02112                                     cube_noise[cube_counter_noise] =
02113                                         kmclipm_imagelist_load(fn_reconstr, CPL_TYPE_FLOAT,
02114                                                                index));
02115 
02116                                 if (edge_nan) {
02117                                     KMO_TRY_EXIT_IF_ERROR(
02118                                         kmo_edge_nan(cube_noise[cube_counter_noise], ifu_nr));
02119                                 }
02120 
02121                                 if (fast_mode) {
02122                                     KMO_TRY_EXIT_IF_NULL(
02123                                         tmpImg = cpl_imagelist_collapse_median_create(cube_noise[cube_counter_noise]));
02124                                     KMO_TRY_EXIT_IF_NULL(
02125                                         tmpCube = cpl_imagelist_new());
02126                                     KMO_TRY_EXIT_IF_ERROR(
02127                                         cpl_imagelist_set(tmpCube, tmpImg, 0));
02128                                     cpl_imagelist_delete(cube_noise[cube_counter_noise]);
02129                                     cube_noise[cube_counter_noise] = tmpCube;
02130                                 }
02131                                 KMO_TRY_EXIT_IF_NULL(
02132                                     header_noise[cube_counter_noise] =
02133                                         kmclipm_propertylist_load(fn_reconstr, index));
02134                                 cube_counter_noise++;
02135                             }
02136                         }
02137                         cpl_error_reset();
02138                     } // end for ifu_nr
02139 
02140                     cpl_propertylist_delete(tmp_header); tmp_header = NULL;
02141                     tmp_frame = kmo_dfs_get_frame(frameset, NULL);
02142                     KMO_TRY_CHECK_ERROR_STATE();
02143                 } // end while-loop RECONSTRUCTED_CUBE frames
02144 
02145                 if (cube_counter_data > 1) {
02146                     if (cube_counter_data == cube_counter_noise) {
02147                         KMO_TRY_EXIT_IF_ERROR(
02148                             kmo_priv_combine(cube_data,
02149                                              cube_noise,
02150                                              header_data,
02151                                              header_noise,
02152                                              cube_counter_data,
02153                                              cube_counter_noise,
02154                                              mapping_mode,
02155                                              "",
02156                                              comb_method,
02157                                              smethod,
02158                                              fmethod,
02159                                              filename,
02160                                              cmethod,
02161                                              cpos_rej,
02162                                              cneg_rej,
02163                                              citer,
02164                                              cmin,
02165                                              cmax,
02166                                              extrapol_enum,
02167                                              flux,
02168                                              &combined_data,
02169                                              &combined_noise));
02170                     } else if (cube_counter_noise == 0) {
02171                         // if imethod == "CS"
02172                         KMO_TRY_EXIT_IF_ERROR(
02173                             kmo_priv_combine(cube_data,
02174                                              NULL,
02175                                              header_data,
02176                                              header_noise,
02177                                              cube_counter_data,
02178                                              cube_counter_noise,
02179                                              mapping_mode,
02180                                              "",
02181                                              comb_method,
02182                                              smethod,
02183                                              fmethod,
02184                                              filename,
02185                                              cmethod,
02186                                              cpos_rej,
02187                                              cneg_rej,
02188                                              citer,
02189                                              cmin,
02190                                              cmax,
02191                                              extrapol_enum,
02192                                              flux,
02193                                              &combined_data,
02194                                              &combined_noise));
02195                     } else {
02196                         KMO_TRY_ASSURE(1 == 0,
02197                                        CPL_ERROR_ILLEGAL_INPUT,
02198                                        "The number of cube-data and cube-noise "
02199                                        "isn't the same (%d vs. %d)!",
02200                                        cube_counter_data, cube_counter_noise);
02201                     }
02202                 } else {
02203                     cpl_msg_warning("", "There is only one reconstructed cube! "
02204                                         "Saving it as it is.");
02205                     KMO_TRY_EXIT_IF_NULL(
02206                         combined_data = cpl_imagelist_duplicate(cube_data[0]));
02207 
02208                     if (cube_noise[0] != NULL) {
02209                         KMO_TRY_EXIT_IF_NULL(
02210                             combined_noise = cpl_imagelist_duplicate(cube_noise[0]));
02211                     }
02212                 }
02213 
02214                 fn_out = COMBINED_CUBE;
02215                 KMO_TRY_EXIT_IF_NULL(
02216                     fn_suffix = cpl_sprintf("_%s", mapping_mode));
02217 
02218                 // save combined cube
02219                 KMO_TRY_EXIT_IF_NULL(
02220                     tmp_frame = kmo_dfs_get_frame(frameset, RECONSTRUCTED_CUBE));
02221                 KMO_TRY_EXIT_IF_ERROR(
02222                     kmo_dfs_save_main_header(frameset, fn_out, fn_suffix, tmp_frame,
02223                                              NULL, parlist, cpl_func));
02224 
02225                 KMO_TRY_EXIT_IF_ERROR(
02226                     kmo_dfs_save_cube(combined_data, fn_out, fn_suffix,
02227                                       header_data[0], 0./0.));
02228 
02229     //            if (combined_noise != NULL) {
02230                     if (header_noise[0] == NULL) {
02231                         KMO_TRY_EXIT_IF_NULL(
02232                             header_noise[0] =
02233                                  cpl_propertylist_duplicate(header_data[0]));
02234 
02235                         KMO_TRY_EXIT_IF_NULL(
02236                             tmp_str = cpl_propertylist_get_string(header_data[0],
02237                                                                   EXTNAME));
02238                         KMO_TRY_EXIT_IF_ERROR(
02239                             kmo_extname_extractor(tmp_str, &ft, &tmp_int, content));
02240                         KMO_TRY_EXIT_IF_NULL(
02241                             extname = kmo_extname_creator(ifu_frame, tmp_int,
02242                                                           EXT_NOISE));
02243                         KMO_TRY_EXIT_IF_ERROR(
02244                             kmclipm_update_property_string(header_noise[0],
02245                                                            EXTNAME, extname,
02246                                                            "FITS extension name"));
02247                         cpl_free(extname); extname = NULL;
02248                     }
02249                     KMO_TRY_EXIT_IF_ERROR(
02250                         kmo_dfs_save_cube(combined_noise, fn_out, fn_suffix,
02251                                           header_noise[0], 0./0.));
02252     //            }
02253 
02254                 for (int i = 0; i < nr_data_alloc; i++) {
02255                     cpl_imagelist_delete(cube_data[i]); cube_data[i] = NULL;
02256                     cpl_imagelist_delete(cube_noise[i]); cube_noise[i] = NULL;
02257                     cpl_propertylist_delete(header_data[i]); header_data[i] = NULL;
02258                     cpl_propertylist_delete(header_noise[i]); header_noise[i] = NULL;
02259                 }
02260                 cpl_free(cube_data);    cube_data = NULL;
02261                 cpl_free(cube_noise);   cube_noise = NULL;
02262                 cpl_free(header_data);  header_data = NULL;
02263                 cpl_free(header_noise); header_noise = NULL;
02264                 cpl_free(fn_suffix); fn_suffix = NULL;
02265                 cpl_imagelist_delete(combined_data); combined_data = NULL;
02266                 cpl_imagelist_delete(combined_noise); combined_noise = NULL;
02267             } // if mapping_mode
02268         } else {
02269             cpl_msg_info("", "NOT combining reconstructed objects (--no_combine is set)");
02270         } // if (!no_combine)
02271 
02272         cpl_msg_info("", "-------------------------------------------");
02273     }
02274     KMO_CATCH
02275     {
02276         KMO_CATCH_MSG();
02277         ret_val = -1;
02278     }
02279 
02280     if (cube_data != NULL) {
02281         for (int ifu_nr = 1; ifu_nr <= nr_data_alloc; ifu_nr++) {
02282             cpl_imagelist_delete(cube_data[ifu_nr-1]); cube_data[ifu_nr-1] = NULL;
02283         }
02284     }
02285     cpl_free(cube_data);    cube_data = NULL;
02286     if (cube_noise != NULL) {
02287         for (int ifu_nr = 1; ifu_nr <= nr_data_alloc; ifu_nr++) {
02288             cpl_imagelist_delete(cube_noise[ifu_nr-1]); cube_noise[ifu_nr-1] = NULL;
02289         }
02290     }
02291     cpl_free(cube_noise);   cube_noise = NULL;
02292     if (header_data != NULL) {
02293         for (int ifu_nr = 1; ifu_nr <= nr_data_alloc; ifu_nr++) {
02294             cpl_propertylist_delete(header_data[ifu_nr-1]); header_data[ifu_nr-1] = NULL;
02295         }
02296     }
02297     cpl_free(header_data);  header_data = NULL;
02298     if (header_noise != NULL) {
02299         for (int ifu_nr = 1; ifu_nr <= nr_data_alloc; ifu_nr++) {
02300             cpl_propertylist_delete(header_noise[ifu_nr-1]); header_noise[ifu_nr-1] = NULL;
02301         }
02302     }
02303     cpl_free(header_noise); header_noise = NULL;
02304 
02305 
02306     kmo_free_fits_desc(&desc1);
02307     kmo_free_fits_desc(&desc2);
02308     cpl_vector_delete(ifus); ifus = NULL;
02309     cpl_free(mapping_mode); mapping_mode = NULL;
02310     if (unused_ifus != NULL) {
02311         kmo_free_unused_ifus(unused_ifus); unused_ifus = NULL;
02312     }
02313     if (bounds != NULL) {
02314         cpl_free(bounds); bounds = NULL;
02315     }
02316     if (obj_sky_struct != NULL) {
02317         cpl_free(obj_sky_struct); obj_sky_struct = NULL;
02318     }
02319 
02320     // frees for the case of errors
02321     kmclipm_vector_delete(telluric_data); telluric_data = NULL;
02322     kmclipm_vector_delete(telluric_noise); telluric_noise = NULL;
02323     cpl_image_delete(illum_data); illum_data = NULL;
02324     cpl_image_delete(illum_noise); illum_noise = NULL;
02325     cpl_propertylist_delete(tmp_header); tmp_header = NULL;
02326     cpl_table_delete(band_table); band_table = NULL;
02327     cpl_propertylist_delete(main_header); main_header = NULL;
02328     if (lcal != NULL) {
02329         for (int i = 0; i < KMOS_NR_DETECTORS; i++) {
02330             cpl_image_delete(lcal[i]);
02331         }
02332     }
02333     cpl_free(lcal); lcal = NULL;
02334     cpl_free(fn_suffix); fn_suffix = NULL;
02335     cpl_free(suffix); suffix = NULL;
02336 
02337     if (all_obj != NULL) {
02338         for (int i = 0; i < nr_science_frames*KMOS_NR_IFUS; i++) {
02339             cpl_free(all_obj[i].name);
02340         }
02341     }
02342     cpl_free(all_obj); all_obj = NULL;
02343 
02344     return ret_val;
02345 }
02346