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