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