00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 #ifdef HAVE_CONFIG_H
00029 #include <config.h>
00030 #endif
00031
00032 #include <math.h>
00033 #include <string.h>
00034 #include <cpl.h>
00035 #include <moses.h>
00036 #include <fors_dfs.h>
00037 #include <fors_qc.h>
00038
00039 static int fors_calib_create(cpl_plugin *);
00040 static int fors_calib_exec(cpl_plugin *);
00041 static int fors_calib_destroy(cpl_plugin *);
00042 static int fors_calib(cpl_parameterlist *, cpl_frameset *);
00043
00044 static char fors_calib_description[] =
00045 "This recipe is used to identify reference lines on LSS, MOS and MXU arc lamp\n"
00046 "exposures, and trace the spectral edges on the corresponding flat field\n"
00047 "exposures. This information is used to determine the spectral extraction\n"
00048 "mask to be applied in the scientific data reduction, performed with the\n"
00049 "recipe fors_science.\n"
00050 "This recipe accepts both FORS1 and FORS2 frames. The input arc lamp and\n"
00051 "flat field exposures are assumed to be obtained quasi-simultaneously, so\n"
00052 "that they would be described by exactly the same instrument distortions.\n"
00053 "A line catalog must be specified, containing the wavelengths of the\n"
00054 "reference arc lamp lines used for the wavelength calibration. A grism\n"
00055 "table (typically depending on the instrument mode, and in particular on\n"
00056 "the grism used) may also be specified: this table contains a default\n"
00057 "recipe parameter setting to control the way spectra are extracted for\n"
00058 "a specific instrument mode, as it is used for automatic run of the\n"
00059 "pipeline on Paranal and in Garching. If this table is specified, it\n"
00060 "will modify the default recipe parameter setting, with the exception of\n"
00061 "those parameters which have been explicitly modifyed on the command line.\n"
00062 "If a grism table is not specified, the input recipe parameters values\n"
00063 "will always be read from the command line, or from an esorex configuration\n"
00064 "file if present, or from their generic default values (that are rarely\n"
00065 "meaningful). Finally a master bias frame must be input to this recipe.\n"
00066 "In the table below the MXU acronym can be read alternatively as MOS\n"
00067 "and LSS, with the exception of CURV_COEFF_LSS, CURV_TRACES_LSS,\n"
00068 "SPATIAL_MAP_LSS, SPECTRA_DETECTION_LSS, and and SLIT_MAP_LSS, which are\n"
00069 "never created. The products SPECTRA_DETECTION_MXU, SLIT_MAP_MXU, and\n"
00070 "DISP_RESIDUALS_MXU, are just created if the --check parameter is set to\n"
00071 "true. The product GLOBAL_DISTORTION_TABLE is just created if more than 12\n"
00072 "separate spectra are found in the CCD.\n\n"
00073 "Input files:\n\n"
00074 " DO category: Type: Explanation: Required:\n"
00075 " SCREEN_FLAT_MXU Raw Flat field exposures Y\n"
00076 " LAMP_MXU Raw Arc lamp exposure Y\n"
00077 " MASTER_BIAS or BIAS Calib Bias frame Y\n"
00078 " MASTER_LINECAT Calib Line catalog Y\n"
00079 " GRISM_TABLE Calib Grism table .\n\n"
00080 "Output files:\n\n"
00081 " DO category: Data type: Explanation:\n"
00082 " MASTER_SCREEN_FLAT_MXU FITS image Combined (sum) flat field\n"
00083 " MASTER_NORM_FLAT_MXU FITS image Normalised flat field\n"
00084 " MAPPED_SCREEN_FLAT_MXU FITS image Wavelength calibrated flat field\n"
00085 " MAPPED_NORM_FLAT_MXU FITS image Wavelength calibrated normalised flat\n"
00086 " REDUCED_LAMP_MXU FITS image Wavelength calibrated arc spectrum\n"
00087 " DISP_COEFF_MXU FITS table Inverse dispersion coefficients\n"
00088 " DISP_RESIDUALS_MXU FITS image Residuals in wavelength calibration\n"
00089 " DISP_RESIDUALS_TABLE_MXU FITS table Residuals in wavelength calibration\n"
00090 " DELTA_IMAGE_MXU FITS image Offset vs linear wavelength calib\n"
00091 " WAVELENGTH_MAP_MXU FITS image Wavelength for each pixel on CCD\n"
00092 " SPECTRA_DETECTION_MXU FITS image Check for preliminary detection\n"
00093 " SLIT_MAP_MXU FITS image Map of central wavelength on CCD\n"
00094 " CURV_TRACES_MXU FITS table Spectral curvature traces\n"
00095 " CURV_COEFF_MXU FITS table Spectral curvature coefficients\n"
00096 " SPATIAL_MAP_MXU FITS image Spatial position along slit on CCD\n"
00097 " SPECTRAL_RESOLUTION_MXU FITS table Resolution at reference arc lines\n"
00098 " SLIT_LOCATION_MXU FITS table Slits on product frames and CCD\n"
00099 " GLOBAL_DISTORTION_TABLE FITS table Global distortions table\n\n";
00100
00101 #define fors_calib_exit(message) \
00102 { \
00103 if (message) cpl_msg_error(recipe, message); \
00104 cpl_free(instrume); \
00105 cpl_free(pipefile); \
00106 cpl_free(fiterror); \
00107 cpl_free(fitlines); \
00108 cpl_image_delete(bias); \
00109 cpl_image_delete(master_bias); \
00110 cpl_image_delete(coordinate); \
00111 cpl_image_delete(checkwave); \
00112 cpl_image_delete(flat); \
00113 cpl_image_delete(master_flat); \
00114 cpl_image_delete(norm_flat); \
00115 cpl_image_delete(mapped_flat); \
00116 cpl_image_delete(mapped_nflat); \
00117 cpl_image_delete(rainbow); \
00118 cpl_image_delete(rectified); \
00119 cpl_image_delete(residual); \
00120 cpl_image_delete(smo_flat); \
00121 cpl_image_delete(spatial); \
00122 cpl_image_delete(spectra); \
00123 cpl_image_delete(wavemap); \
00124 cpl_image_delete(delta); \
00125 cpl_image_delete(rect_flat); \
00126 cpl_image_delete(rect_nflat); \
00127 cpl_image_delete(mapped_flat); \
00128 cpl_image_delete(mapped_nflat); \
00129 cpl_mask_delete(refmask); \
00130 cpl_propertylist_delete(header); \
00131 cpl_propertylist_delete(save_header); \
00132 cpl_propertylist_delete(qclist); \
00133 cpl_table_delete(grism_table); \
00134 cpl_table_delete(idscoeff); \
00135 cpl_table_delete(idscoeff_all); \
00136 cpl_table_delete(restable); \
00137 cpl_table_delete(maskslits); \
00138 cpl_table_delete(overscans); \
00139 cpl_table_delete(traces); \
00140 cpl_table_delete(polytraces); \
00141 cpl_table_delete(slits); \
00142 cpl_table_delete(restab); \
00143 cpl_table_delete(global); \
00144 cpl_table_delete(wavelengths); \
00145 cpl_vector_delete(lines); \
00146 cpl_msg_indent_less(); \
00147 return -1; \
00148 }
00149
00150 #define fors_calib_exit_memcheck(message) \
00151 { \
00152 if (message) cpl_msg_info(recipe, message); \
00153 printf("free instrume (%p)\n", instrume); \
00154 cpl_free(instrume); \
00155 printf("free pipefile (%p)\n", pipefile); \
00156 cpl_free(pipefile); \
00157 printf("free fiterror (%p)\n", fiterror); \
00158 cpl_free(fiterror); \
00159 printf("free fitlines (%p)\n", fitlines); \
00160 cpl_free(fitlines); \
00161 printf("free bias (%p)\n", bias); \
00162 cpl_image_delete(bias); \
00163 printf("free master_bias (%p)\n", master_bias); \
00164 cpl_image_delete(master_bias); \
00165 printf("free coordinate (%p)\n", coordinate); \
00166 cpl_image_delete(coordinate); \
00167 printf("free checkwave (%p)\n", checkwave); \
00168 cpl_image_delete(checkwave); \
00169 printf("free flat (%p)\n", flat); \
00170 cpl_image_delete(flat); \
00171 printf("free master_flat (%p)\n", master_flat); \
00172 cpl_image_delete(master_flat); \
00173 printf("free norm_flat (%p)\n", norm_flat); \
00174 cpl_image_delete(norm_flat); \
00175 printf("free mapped_flat (%p)\n", mapped_flat); \
00176 cpl_image_delete(mapped_flat); \
00177 printf("free mapped_nflat (%p)\n", mapped_nflat); \
00178 cpl_image_delete(mapped_nflat); \
00179 printf("free rainbow (%p)\n", rainbow); \
00180 cpl_image_delete(rainbow); \
00181 printf("free rectified (%p)\n", rectified); \
00182 cpl_image_delete(rectified); \
00183 printf("free residual (%p)\n", residual); \
00184 cpl_image_delete(residual); \
00185 printf("free smo_flat (%p)\n", smo_flat); \
00186 cpl_image_delete(smo_flat); \
00187 printf("free spatial (%p)\n", spatial); \
00188 cpl_image_delete(spatial); \
00189 printf("free spectra (%p)\n", spectra); \
00190 cpl_image_delete(spectra); \
00191 printf("free wavemap (%p)\n", wavemap); \
00192 cpl_image_delete(wavemap); \
00193 printf("free delta (%p)\n", delta); \
00194 cpl_image_delete(delta); \
00195 printf("free rect_flat (%p)\n", rect_flat); \
00196 cpl_image_delete(rect_flat); \
00197 printf("free rect_nflat (%p)\n", rect_nflat); \
00198 cpl_image_delete(rect_nflat); \
00199 printf("free refmask (%p)\n", refmask); \
00200 cpl_mask_delete(refmask); \
00201 printf("free header (%p)\n", header); \
00202 cpl_propertylist_delete(header); \
00203 printf("free save_header (%p)\n", save_header); \
00204 cpl_propertylist_delete(save_header); \
00205 printf("free qclist (%p)\n", qclist); \
00206 cpl_propertylist_delete(qclist); \
00207 printf("free grism_table (%p)\n", grism_table); \
00208 cpl_table_delete(grism_table); \
00209 printf("free idscoeff (%p)\n", idscoeff); \
00210 cpl_table_delete(idscoeff); \
00211 printf("free idscoeff_all (%p)\n", idscoeff_all); \
00212 cpl_table_delete(idscoeff_all); \
00213 printf("free restable (%p)\n", restable); \
00214 cpl_table_delete(restable); \
00215 printf("free maskslits (%p)\n", maskslits); \
00216 cpl_table_delete(maskslits); \
00217 printf("free overscans (%p)\n", overscans); \
00218 cpl_table_delete(overscans); \
00219 printf("free traces (%p)\n", traces); \
00220 cpl_table_delete(traces); \
00221 printf("free polytraces (%p)\n", polytraces); \
00222 cpl_table_delete(polytraces); \
00223 printf("free slits (%p)\n", slits); \
00224 cpl_table_delete(slits); \
00225 printf("free restab (%p)\n", restab); \
00226 cpl_table_delete(restab); \
00227 printf("free global (%p)\n", global); \
00228 cpl_table_delete(global); \
00229 printf("free wavelengths (%p)\n", wavelengths); \
00230 cpl_table_delete(wavelengths); \
00231 printf("free lines (%p)\n", lines); \
00232 cpl_vector_delete(lines); \
00233 cpl_msg_indent_less(); \
00234 return 0; \
00235 }
00236
00237
00249 int cpl_plugin_get_info(cpl_pluginlist *list)
00250 {
00251 cpl_recipe *recipe = cpl_calloc(1, sizeof *recipe );
00252 cpl_plugin *plugin = &recipe->interface;
00253
00254 cpl_plugin_init(plugin,
00255 CPL_PLUGIN_API,
00256 FORS_BINARY_VERSION,
00257 CPL_PLUGIN_TYPE_RECIPE,
00258 "fors_calib",
00259 "Determination of the extraction mask",
00260 fors_calib_description,
00261 "Carlo Izzo",
00262 PACKAGE_BUGREPORT,
00263 "This file is currently part of the FORS Instrument Pipeline\n"
00264 "Copyright (C) 2002-2006 European Southern Observatory\n\n"
00265 "This program is free software; you can redistribute it and/or modify\n"
00266 "it under the terms of the GNU General Public License as published by\n"
00267 "the Free Software Foundation; either version 2 of the License, or\n"
00268 "(at your option) any later version.\n\n"
00269 "This program is distributed in the hope that it will be useful,\n"
00270 "but WITHOUT ANY WARRANTY; without even the implied warranty of\n"
00271 "MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the\n"
00272 "GNU General Public License for more details.\n\n"
00273 "You should have received a copy of the GNU General Public License\n"
00274 "along with this program; if not, write to the Free Software Foundation,\n"
00275 "Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA\n",
00276 fors_calib_create,
00277 fors_calib_exec,
00278 fors_calib_destroy);
00279
00280 cpl_pluginlist_append(list, plugin);
00281
00282 return 0;
00283 }
00284
00285
00296 static int fors_calib_create(cpl_plugin *plugin)
00297 {
00298 cpl_recipe *recipe;
00299 cpl_parameter *p;
00300
00301
00302
00303
00304
00305
00306 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
00307 recipe = (cpl_recipe *)plugin;
00308 else
00309 return -1;
00310
00311
00312
00313
00314
00315 recipe->parameters = cpl_parameterlist_new();
00316
00317
00318
00319
00320
00321
00322 p = cpl_parameter_new_value("fors.fors_calib.dispersion",
00323 CPL_TYPE_DOUBLE,
00324 "Expected spectral dispersion (Angstrom/pixel)",
00325 "fors.fors_calib",
00326 0.0);
00327 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "dispersion");
00328 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00329 cpl_parameterlist_append(recipe->parameters, p);
00330
00331
00332
00333
00334
00335 p = cpl_parameter_new_value("fors.fors_calib.peakdetection",
00336 CPL_TYPE_DOUBLE,
00337 "Initial peak detection threshold (ADU)",
00338 "fors.fors_calib",
00339 0.0);
00340 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "peakdetection");
00341 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00342 cpl_parameterlist_append(recipe->parameters, p);
00343
00344
00345
00346
00347
00348 p = cpl_parameter_new_value("fors.fors_calib.wdegree",
00349 CPL_TYPE_INT,
00350 "Degree of wavelength calibration polynomial",
00351 "fors.fors_calib",
00352 0);
00353 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "wdegree");
00354 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00355 cpl_parameterlist_append(recipe->parameters, p);
00356
00357
00358
00359
00360
00361 p = cpl_parameter_new_value("fors.fors_calib.wradius",
00362 CPL_TYPE_INT,
00363 "Search radius if iterating pattern-matching "
00364 "with first-guess method",
00365 "fors.fors_calib",
00366 4);
00367 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "wradius");
00368 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00369 cpl_parameterlist_append(recipe->parameters, p);
00370
00371
00372
00373
00374
00375 p = cpl_parameter_new_value("fors.fors_calib.wreject",
00376 CPL_TYPE_DOUBLE,
00377 "Rejection threshold in dispersion "
00378 "relation fit (pixel)",
00379 "fors.fors_calib",
00380 0.7);
00381 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "wreject");
00382 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00383 cpl_parameterlist_append(recipe->parameters, p);
00384
00385
00386
00387
00388
00389 p = cpl_parameter_new_value("fors.fors_calib.wmode",
00390 CPL_TYPE_INT,
00391 "Interpolation mode of wavelength solution "
00392 "applicable to LSS-like data (0 = no "
00393 "interpolation, 1 = fill gaps, 2 = global "
00394 "model)",
00395 "fors.fors_calib",
00396 2);
00397 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "wmode");
00398 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00399 cpl_parameterlist_append(recipe->parameters, p);
00400
00401
00402
00403
00404
00405 p = cpl_parameter_new_value("fors.fors_calib.wcolumn",
00406 CPL_TYPE_STRING,
00407 "Name of line catalog table column "
00408 "with wavelengths",
00409 "fors.fors_calib",
00410 "WLEN");
00411 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "wcolumn");
00412 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00413 cpl_parameterlist_append(recipe->parameters, p);
00414
00415
00416
00417
00418
00419 p = cpl_parameter_new_value("fors.fors_calib.cdegree",
00420 CPL_TYPE_INT,
00421 "Degree of spectral curvature polynomial",
00422 "fors.fors_calib",
00423 0);
00424 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "cdegree");
00425 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00426 cpl_parameterlist_append(recipe->parameters, p);
00427
00428
00429
00430
00431
00432 p = cpl_parameter_new_value("fors.fors_calib.cmode",
00433 CPL_TYPE_INT,
00434 "Interpolation mode of curvature solution "
00435 "applicable to MOS-like data (0 = no "
00436 "interpolation, 1 = fill gaps, 2 = global "
00437 "model)",
00438 "fors.fors_calib",
00439 1);
00440 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "cmode");
00441 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00442 cpl_parameterlist_append(recipe->parameters, p);
00443
00444
00445
00446
00447
00448 p = cpl_parameter_new_value("fors.fors_calib.startwavelength",
00449 CPL_TYPE_DOUBLE,
00450 "Start wavelength in spectral extraction",
00451 "fors.fors_calib",
00452 0.0);
00453 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "startwavelength");
00454 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00455 cpl_parameterlist_append(recipe->parameters, p);
00456
00457
00458
00459
00460
00461 p = cpl_parameter_new_value("fors.fors_calib.endwavelength",
00462 CPL_TYPE_DOUBLE,
00463 "End wavelength in spectral extraction",
00464 "fors.fors_calib",
00465 0.0);
00466 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "endwavelength");
00467 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00468 cpl_parameterlist_append(recipe->parameters, p);
00469
00470
00471
00472
00473
00474 p = cpl_parameter_new_value("fors.fors_calib.slit_ident",
00475 CPL_TYPE_BOOL,
00476 "Attempt slit identification for MOS or MXU",
00477 "fors.fors_calib",
00478 TRUE);
00479 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "slit_ident");
00480 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00481 cpl_parameterlist_append(recipe->parameters, p);
00482
00483
00484
00485
00486
00487
00488 p = cpl_parameter_new_value("fors.fors_calib.sdegree",
00489 CPL_TYPE_INT,
00490 "Degree of flat field fitting polynomial "
00491 "along spatial direction (used for LSS "
00492 "data only)",
00493 "fors.fors_calib",
00494 4);
00495 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "sdegree");
00496 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00497 cpl_parameterlist_append(recipe->parameters, p);
00498
00499
00500
00501
00502
00503
00504 p = cpl_parameter_new_value("fors.fors_calib.ddegree",
00505 CPL_TYPE_INT,
00506 "Degree of flat field fitting polynomial "
00507 "along dispersion direction (used for MOS "
00508 "and MXU data only)",
00509 "fors.fors_calib",
00510 -1);
00511 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "ddegree");
00512 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00513 cpl_parameterlist_append(recipe->parameters, p);
00514
00515
00516
00517
00518
00519 p = cpl_parameter_new_value("fors.fors_calib.dradius",
00520 CPL_TYPE_INT,
00521 "Smooth box radius for flat field along "
00522 "dispersion direction",
00523 "fors.fors_calib",
00524 10);
00525 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "dradius");
00526 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00527 cpl_parameterlist_append(recipe->parameters, p);
00528
00529
00530
00531
00532
00533
00534 p = cpl_parameter_new_value("fors.fors_calib.sradius",
00535 CPL_TYPE_INT,
00536 "Smooth box radius for flat field along "
00537 "spatial direction",
00538 "fors.fors_calib",
00539 10);
00540 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "sradius");
00541 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00542 cpl_parameterlist_append(recipe->parameters, p);
00543
00544
00545
00546
00547
00548 p = cpl_parameter_new_value("fors.fors_calib.qc",
00549 CPL_TYPE_BOOL,
00550 "Compute QC1 parameters",
00551 "fors.fors_calib",
00552 TRUE);
00553 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "qc");
00554 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00555 cpl_parameterlist_append(recipe->parameters, p);
00556
00557
00558
00559
00560
00561 p = cpl_parameter_new_value("fors.fors_calib.check",
00562 CPL_TYPE_BOOL,
00563 "Create intermediate products",
00564 "fors.fors_calib",
00565 FALSE);
00566 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "check");
00567 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
00568 cpl_parameterlist_append(recipe->parameters, p);
00569
00570 return 0;
00571 }
00572
00573
00582 static int fors_calib_exec(cpl_plugin *plugin)
00583 {
00584 cpl_recipe *recipe;
00585
00586 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
00587 recipe = (cpl_recipe *)plugin;
00588 else
00589 return -1;
00590
00591 return fors_calib(recipe->parameters, recipe->frames);
00592 }
00593
00594
00603 static int fors_calib_destroy(cpl_plugin *plugin)
00604 {
00605 cpl_recipe *recipe;
00606
00607 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
00608 recipe = (cpl_recipe *)plugin;
00609 else
00610 return -1;
00611
00612 cpl_parameterlist_delete(recipe->parameters);
00613
00614 return 0;
00615 }
00616
00617
00627 static int fors_calib(cpl_parameterlist *parlist, cpl_frameset *frameset)
00628 {
00629
00630 const char *recipe = "fors_calib";
00631
00632
00633
00634
00635
00636
00637 double dispersion;
00638 double peakdetection;
00639 int wdegree;
00640 int wradius;
00641 double wreject;
00642 int wmode;
00643 const char *wcolumn;
00644 int cdegree;
00645 int cmode;
00646 double startwavelength;
00647 double endwavelength;
00648 int slit_ident;
00649 int sdegree;
00650 int ddegree;
00651 int sradius;
00652 int dradius;
00653 int qc;
00654 int check;
00655
00656
00657
00658
00659
00660 cpl_imagelist *biases = NULL;
00661 cpl_image *bias = NULL;
00662 cpl_image *master_bias = NULL;
00663 cpl_image *multi_bias = NULL;
00664 cpl_image *flat = NULL;
00665 cpl_image *master_flat = NULL;
00666 cpl_image *smo_flat = NULL;
00667 cpl_image *norm_flat = NULL;
00668 cpl_image *spectra = NULL;
00669 cpl_image *wavemap = NULL;
00670 cpl_image *delta = NULL;
00671 cpl_image *residual = NULL;
00672 cpl_image *checkwave = NULL;
00673 cpl_image *rectified = NULL;
00674 cpl_image *dummy = NULL;
00675 cpl_image *refimage = NULL;
00676 cpl_image *coordinate = NULL;
00677 cpl_image *rainbow = NULL;
00678 cpl_image *spatial = NULL;
00679 cpl_image *rect_flat = NULL;
00680 cpl_image *rect_nflat = NULL;
00681 cpl_image *mapped_flat = NULL;
00682 cpl_image *mapped_nflat = NULL;
00683
00684 cpl_mask *refmask = NULL;
00685
00686 cpl_table *grism_table = NULL;
00687 cpl_table *overscans = NULL;
00688 cpl_table *wavelengths = NULL;
00689 cpl_table *idscoeff = NULL;
00690 cpl_table *idscoeff_all = NULL;
00691 cpl_table *restable = NULL;
00692 cpl_table *slits = NULL;
00693 cpl_table *positions = NULL;
00694 cpl_table *maskslits = NULL;
00695 cpl_table *traces = NULL;
00696 cpl_table *polytraces = NULL;
00697 cpl_table *restab = NULL;
00698 cpl_table *global = NULL;
00699
00700 cpl_vector *lines = NULL;
00701
00702 cpl_propertylist *header = NULL;
00703 cpl_propertylist *save_header = NULL;
00704 cpl_propertylist *qclist = NULL;
00705
00706
00707
00708
00709
00710 char version[80];
00711 const char *arc_tag;
00712 const char *flat_tag;
00713 const char *master_screen_flat_tag;
00714 const char *master_norm_flat_tag;
00715 const char *reduced_lamp_tag;
00716 const char *disp_residuals_tag;
00717 const char *disp_coeff_tag;
00718 const char *wavelength_map_tag;
00719 const char *spectra_detection_tag;
00720 const char *spectral_resolution_tag;
00721 const char *slit_map_tag;
00722 const char *curv_traces_tag;
00723 const char *curv_coeff_tag;
00724 const char *spatial_map_tag;
00725 const char *slit_location_tag;
00726 const char *global_distortion_tag = "GLOBAL_DISTORTION_TABLE";
00727 const char *disp_residuals_table_tag;
00728 const char *delta_image_tag;
00729 const char *mapped_screen_flat_tag;
00730 const char *mapped_norm_flat_tag;
00731 const char *keyname;
00732 int mxu, mos, lss;
00733 int treat_as_lss = 0;
00734 int nslits;
00735 float *data;
00736 double *xpos;
00737 double mxpos;
00738 double mean_rms;
00739 double mean_rms_err;
00740 double alltime;
00741 int nflats;
00742 int nbias;
00743 int nlines;
00744 int rebin;
00745 double *line;
00746 double *fiterror = NULL;
00747 int *fitlines = NULL;
00748 int nx, ny;
00749 double reference;
00750 double gain;
00751 int compute_central_wave;
00752 int ccd_xsize, ccd_ysize;
00753 int i;
00754
00755 char *instrume = NULL;
00756 char *pipefile = NULL;
00757 char *wheel4;
00758
00759
00760 snprintf(version, 80, "%s-%s", PACKAGE, PACKAGE_VERSION);
00761
00762 cpl_msg_set_indentation(2);
00763
00764 if (dfs_files_dont_exist(frameset))
00765 fors_calib_exit(NULL);
00766
00767
00768
00769
00770
00771 cpl_msg_info(recipe, "Recipe %s configuration parameters:", recipe);
00772 cpl_msg_indent_more();
00773
00774 if (cpl_frameset_count_tags(frameset, "GRISM_TABLE") > 1)
00775 fors_calib_exit("Too many in input: GRISM_TABLE");
00776
00777 grism_table = dfs_load_table(frameset, "GRISM_TABLE", 1);
00778
00779 dispersion = dfs_get_parameter_double(parlist,
00780 "fors.fors_calib.dispersion", grism_table);
00781
00782 if (dispersion <= 0.0)
00783 fors_calib_exit("Invalid spectral dispersion value");
00784
00785 peakdetection = dfs_get_parameter_double(parlist,
00786 "fors.fors_calib.peakdetection", grism_table);
00787 if (peakdetection <= 0.0)
00788 fors_calib_exit("Invalid peak detection level");
00789
00790 wdegree = dfs_get_parameter_int(parlist,
00791 "fors.fors_calib.wdegree", grism_table);
00792
00793 if (wdegree < 1)
00794 fors_calib_exit("Invalid polynomial degree");
00795
00796 if (wdegree > 5)
00797 fors_calib_exit("Max allowed polynomial degree is 5");
00798
00799 wradius = dfs_get_parameter_int(parlist, "fors.fors_calib.wradius", NULL);
00800
00801 if (wradius < 0)
00802 fors_calib_exit("Invalid search radius");
00803
00804 wreject = dfs_get_parameter_double(parlist, "fors.fors_calib.wreject", NULL);
00805
00806 if (wreject <= 0.0)
00807 fors_calib_exit("Invalid rejection threshold");
00808
00809 wmode = dfs_get_parameter_int(parlist, "fors.fors_calib.wmode", NULL);
00810
00811 if (wmode < 0 || wmode > 2)
00812 fors_calib_exit("Invalid wavelength solution interpolation mode");
00813
00814 wcolumn = dfs_get_parameter_string(parlist, "fors.fors_calib.wcolumn",
00815 NULL);
00816
00817 cdegree = dfs_get_parameter_int(parlist, "fors.fors_calib.cdegree",
00818 grism_table);
00819
00820 if (cdegree < 1)
00821 fors_calib_exit("Invalid polynomial degree");
00822
00823 if (cdegree > 5)
00824 fors_calib_exit("Max allowed polynomial degree is 5");
00825
00826 cmode = dfs_get_parameter_int(parlist, "fors.fors_calib.cmode", NULL);
00827
00828 if (cmode < 0 || cmode > 2)
00829 fors_calib_exit("Invalid curvature solution interpolation mode");
00830
00831 startwavelength = dfs_get_parameter_double(parlist,
00832 "fors.fors_calib.startwavelength", grism_table);
00833 if (startwavelength > 1.0)
00834 if (startwavelength < 3000.0 || startwavelength > 13000.0)
00835 fors_calib_exit("Invalid wavelength");
00836
00837 endwavelength = dfs_get_parameter_double(parlist,
00838 "fors.fors_calib.endwavelength", grism_table);
00839 if (endwavelength > 1.0) {
00840 if (endwavelength < 3000.0 || endwavelength > 13000.0)
00841 fors_calib_exit("Invalid wavelength");
00842 if (startwavelength < 1.0)
00843 fors_calib_exit("Invalid wavelength interval");
00844 }
00845
00846 if (startwavelength > 1.0)
00847 if (endwavelength - startwavelength <= 0.0)
00848 fors_calib_exit("Invalid wavelength interval");
00849
00850 slit_ident = dfs_get_parameter_bool(parlist,
00851 "fors.fors_calib.slit_ident", NULL);
00852
00853 sdegree = dfs_get_parameter_int(parlist, "fors.fors_calib.sdegree", NULL);
00854 ddegree = dfs_get_parameter_int(parlist, "fors.fors_calib.ddegree", NULL);
00855 sradius = dfs_get_parameter_int(parlist, "fors.fors_calib.sradius", NULL);
00856 dradius = dfs_get_parameter_int(parlist, "fors.fors_calib.dradius", NULL);
00857
00858 if (sradius < 1 || dradius < 1)
00859 fors_calib_exit("Invalid smoothing box radius");
00860
00861 qc = dfs_get_parameter_bool(parlist, "fors.fors_calib.qc", NULL);
00862
00863 check = dfs_get_parameter_bool(parlist, "fors.fors_calib.check", NULL);
00864
00865 cpl_table_delete(grism_table); grism_table = NULL;
00866
00867 if (cpl_error_get_code())
00868 fors_calib_exit("Failure getting the configuration parameters");
00869
00870
00871
00872
00873
00874
00875 cpl_msg_indent_less();
00876 cpl_msg_info(recipe, "Check input set-of-frames:");
00877 cpl_msg_indent_more();
00878
00879 if (!dfs_equal_keyword(frameset, "ESO INS GRIS1 ID"))
00880 fors_calib_exit("Input frames are not from the same grism");
00881
00882 if (!dfs_equal_keyword(frameset, "ESO INS FILT1 ID"))
00883 fors_calib_exit("Input frames are not from the same filter");
00884
00885 if (!dfs_equal_keyword(frameset, "ESO DET CHIP1 ID"))
00886 fors_calib_exit("Input frames are not from the same chip");
00887
00888 mxu = cpl_frameset_count_tags(frameset, "LAMP_MXU");
00889 mos = cpl_frameset_count_tags(frameset, "LAMP_MOS");
00890 lss = cpl_frameset_count_tags(frameset, "LAMP_LSS");
00891
00892 if (mxu + mos + lss == 0)
00893 fors_calib_exit("Missing input arc lamp frame");
00894
00895 if (mxu + mos + lss > 1)
00896 fors_calib_exit("Just one input arc lamp frame is allowed");
00897
00898 if (mxu) {
00899 cpl_msg_info(recipe, "MXU data found");
00900 arc_tag = "LAMP_MXU";
00901 flat_tag = "SCREEN_FLAT_MXU";
00902 master_screen_flat_tag = "MASTER_SCREEN_FLAT_MXU";
00903 master_norm_flat_tag = "MASTER_NORM_FLAT_MXU";
00904 reduced_lamp_tag = "REDUCED_LAMP_MXU";
00905 disp_residuals_tag = "DISP_RESIDUALS_MXU";
00906 disp_coeff_tag = "DISP_COEFF_MXU";
00907 wavelength_map_tag = "WAVELENGTH_MAP_MXU";
00908 spectra_detection_tag = "SPECTRA_DETECTION_MXU";
00909 spectral_resolution_tag = "SPECTRAL_RESOLUTION_MXU";
00910 slit_map_tag = "SLIT_MAP_MXU";
00911 curv_traces_tag = "CURV_TRACES_MXU";
00912 curv_coeff_tag = "CURV_COEFF_MXU";
00913 spatial_map_tag = "SPATIAL_MAP_MXU";
00914 slit_location_tag = "SLIT_LOCATION_MXU";
00915 disp_residuals_table_tag = "DISP_RESIDUALS_TABLE_MXU";
00916 delta_image_tag = "DELTA_IMAGE_MXU";
00917 mapped_screen_flat_tag = "MAPPED_SCREEN_FLAT_MXU";
00918 mapped_norm_flat_tag = "MAPPED_NORM_FLAT_MXU";
00919 }
00920
00921 if (lss) {
00922 cpl_msg_info(recipe, "LSS data found");
00923 arc_tag = "LAMP_LSS";
00924 flat_tag = "SCREEN_FLAT_LSS";
00925 master_screen_flat_tag = "MASTER_SCREEN_FLAT_LSS";
00926 master_norm_flat_tag = "MASTER_NORM_FLAT_LSS";
00927 reduced_lamp_tag = "REDUCED_LAMP_LSS";
00928 spectral_resolution_tag = "SPECTRAL_RESOLUTION_LSS";
00929 disp_residuals_tag = "DISP_RESIDUALS_LSS";
00930 disp_coeff_tag = "DISP_COEFF_LSS";
00931 slit_location_tag = "SLIT_LOCATION_LSS";
00932 wavelength_map_tag = "WAVELENGTH_MAP_LSS";
00933 slit_map_tag = "SLIT_MAP_LSS";
00934 disp_residuals_table_tag = "DISP_RESIDUALS_TABLE_LSS";
00935 delta_image_tag = "DELTA_IMAGE_LSS";
00936 mapped_screen_flat_tag = "MAPPED_SCREEN_FLAT_LSS";
00937 mapped_norm_flat_tag = "MAPPED_NORM_FLAT_LSS";
00938 }
00939
00940 if (mos) {
00941 cpl_msg_info(recipe, "MOS data found");
00942 arc_tag = "LAMP_MOS";
00943 flat_tag = "SCREEN_FLAT_MOS";
00944 master_screen_flat_tag = "MASTER_SCREEN_FLAT_MOS";
00945 master_norm_flat_tag = "MASTER_NORM_FLAT_MOS";
00946 reduced_lamp_tag = "REDUCED_LAMP_MOS";
00947 disp_residuals_tag = "DISP_RESIDUALS_MOS";
00948 disp_coeff_tag = "DISP_COEFF_MOS";
00949 wavelength_map_tag = "WAVELENGTH_MAP_MOS";
00950 spectra_detection_tag = "SPECTRA_DETECTION_MOS";
00951 spectral_resolution_tag = "SPECTRAL_RESOLUTION_MOS";
00952 slit_map_tag = "SLIT_MAP_MOS";
00953 curv_traces_tag = "CURV_TRACES_MOS";
00954 curv_coeff_tag = "CURV_COEFF_MOS";
00955 spatial_map_tag = "SPATIAL_MAP_MOS";
00956 slit_location_tag = "SLIT_LOCATION_MOS";
00957 disp_residuals_table_tag = "DISP_RESIDUALS_TABLE_MOS";
00958 delta_image_tag = "DELTA_IMAGE_MOS";
00959 mapped_screen_flat_tag = "MAPPED_SCREEN_FLAT_MOS";
00960 mapped_norm_flat_tag = "MAPPED_NORM_FLAT_MOS";
00961 }
00962
00963 nbias = 0;
00964 if (cpl_frameset_count_tags(frameset, "MASTER_BIAS") == 0) {
00965 if (cpl_frameset_count_tags(frameset, "BIAS") == 0)
00966 fors_calib_exit("Missing required input: MASTER_BIAS or BIAS");
00967 nbias = cpl_frameset_count_tags(frameset, "BIAS");
00968 }
00969
00970 if (cpl_frameset_count_tags(frameset, "MASTER_BIAS") > 1)
00971 fors_calib_exit("Too many in input: MASTER_BIAS");
00972
00973 if (cpl_frameset_count_tags(frameset, "MASTER_LINECAT") == 0)
00974 fors_calib_exit("Missing required input: MASTER_LINECAT");
00975
00976 if (cpl_frameset_count_tags(frameset, "MASTER_LINECAT") > 1)
00977 fors_calib_exit("Too many in input: MASTER_LINECAT");
00978
00979 nflats = cpl_frameset_count_tags(frameset, flat_tag);
00980
00981 if (nflats < 1) {
00982 cpl_msg_error(recipe, "Missing required input: %s", flat_tag);
00983 fors_calib_exit(NULL);
00984 }
00985
00986 cpl_msg_indent_less();
00987
00988 if (nflats > 1)
00989 cpl_msg_info(recipe, "Load %d flat field frames and sum them...",
00990 nflats);
00991 else
00992 cpl_msg_info(recipe, "Load flat field exposure...");
00993
00994 cpl_msg_indent_more();
00995
00996 header = dfs_load_header(frameset, flat_tag, 0);
00997
00998 if (header == NULL)
00999 fors_calib_exit("Cannot load flat field frame header");
01000
01001
01002
01003
01004
01005 wheel4 = (char *)cpl_propertylist_get_string(header, "ESO INS OPTI9 TYPE");
01006 if (cpl_error_get_code() != CPL_ERROR_NONE) {
01007 fors_calib_exit("Missing keyword ESO INS OPTI9 TYPE in flat header");
01008 }
01009
01010 if (strcmp("FILT", wheel4) == 0) {
01011 wheel4 = (char *)cpl_propertylist_get_string(header,
01012 "ESO INS OPTI9 NAME");
01013 cpl_msg_error(recipe, "Unsupported filter: %s", wheel4);
01014 fors_calib_exit(NULL);
01015 }
01016
01017 alltime = cpl_propertylist_get_double(header, "EXPTIME");
01018
01019 if (cpl_error_get_code() != CPL_ERROR_NONE)
01020 fors_calib_exit("Missing keyword EXPTIME in flat field frame header");
01021
01022 cpl_propertylist_delete(header);
01023
01024 for (i = 1; i < nflats; i++) {
01025
01026 header = dfs_load_header(frameset, NULL, 0);
01027
01028 if (header == NULL)
01029 fors_calib_exit("Cannot load flat field frame header");
01030
01031 alltime += cpl_propertylist_get_double(header, "EXPTIME");
01032
01033 if (cpl_error_get_code() != CPL_ERROR_NONE)
01034 fors_calib_exit("Missing keyword EXPTIME in flat field "
01035 "frame header");
01036
01037 cpl_propertylist_delete(header);
01038
01039 }
01040
01041 master_flat = dfs_load_image(frameset, flat_tag, CPL_TYPE_FLOAT, 0, 0);
01042
01043 if (master_flat == NULL)
01044 fors_calib_exit("Cannot load flat field");
01045
01046 for (i = 1; i < nflats; i++) {
01047 flat = dfs_load_image(frameset, NULL, CPL_TYPE_FLOAT, 0, 0);
01048 if (flat) {
01049 cpl_image_add(master_flat, flat);
01050 cpl_image_delete(flat); flat = NULL;
01051 }
01052 else
01053 fors_calib_exit("Cannot load flat field");
01054 }
01055
01056
01057
01058
01059
01060
01061
01062
01063
01064
01065
01066
01067 header = dfs_load_header(frameset, arc_tag, 0);
01068
01069 if (header == NULL)
01070 fors_calib_exit("Cannot load arc lamp header");
01071
01072 instrume = (char *)cpl_propertylist_get_string(header, "INSTRUME");
01073 if (instrume == NULL)
01074 fors_calib_exit("Missing keyword INSTRUME in arc lamp header");
01075 instrume = cpl_strdup(instrume);
01076
01077 if (instrume[4] == '1')
01078 snprintf(version, 80, "%s/%s", "fors1", VERSION);
01079 if (instrume[4] == '2')
01080 snprintf(version, 80, "%s/%s", "fors2", VERSION);
01081
01082 reference = cpl_propertylist_get_double(header, "ESO INS GRIS1 WLEN");
01083
01084 if (cpl_error_get_code() != CPL_ERROR_NONE)
01085 fors_calib_exit("Missing keyword ESO INS GRIS1 WLEN in arc lamp "
01086 "frame header");
01087
01088 if (reference < 3000.0)
01089 reference *= 10;
01090
01091 if (reference < 3000.0 || reference > 13000.0) {
01092 cpl_msg_error(recipe, "Invalid central wavelength %.2f read from "
01093 "keyword ESO INS GRIS1 WLEN in arc lamp frame header",
01094 reference);
01095 fors_calib_exit(NULL);
01096 }
01097
01098 cpl_msg_info(recipe, "The central wavelength is: %.2f", reference);
01099
01100 rebin = cpl_propertylist_get_int(header, "ESO DET WIN1 BINX");
01101
01102 if (cpl_error_get_code() != CPL_ERROR_NONE)
01103 fors_calib_exit("Missing keyword ESO DET WIN1 BINX in arc lamp "
01104 "frame header");
01105
01106 if (rebin != 1) {
01107 dispersion *= rebin;
01108 cpl_msg_warning(recipe, "The rebin factor is %d, and therefore the "
01109 "working dispersion used is %f A/pixel", rebin,
01110 dispersion);
01111 }
01112
01113 gain = cpl_propertylist_get_double(header, "ESO DET OUT1 CONAD");
01114
01115 if (cpl_error_get_code() != CPL_ERROR_NONE)
01116 fors_calib_exit("Missing keyword ESO DET OUT1 CONAD in arc lamp "
01117 "frame header");
01118
01119 cpl_msg_info(recipe, "The gain factor is: %.2f e-/ADU", gain);
01120
01121 if (mos || mxu) {
01122
01123 cpl_msg_info(recipe, "Produce mask slit position table...");
01124 if (mos)
01125 maskslits = mos_load_slits_fors_mos(header);
01126 else
01127 maskslits = mos_load_slits_fors_mxu(header);
01128
01129
01130
01131
01132
01133
01134 mxpos = cpl_table_get_column_median(maskslits, "xtop");
01135 xpos = cpl_table_get_data_double(maskslits, "xtop");
01136 nslits = cpl_table_get_nrow(maskslits);
01137
01138 treat_as_lss = 1;
01139 for (i = 0; i < nslits; i++) {
01140 if (fabs(mxpos-xpos[i]) > 0.01) {
01141 treat_as_lss = 0;
01142 break;
01143 }
01144 }
01145
01146
01147 if (treat_as_lss) {
01148 cpl_msg_warning(recipe, "All MOS slits have the same offset: %.2f\n"
01149 "The LSS data reduction strategy is applied!",
01150 mxpos);
01151 cpl_table_delete(maskslits); maskslits = NULL;
01152 if (mos) {
01153 master_screen_flat_tag = "MASTER_SCREEN_FLAT_LONG_MOS";
01154 master_norm_flat_tag = "MASTER_NORM_FLAT_LONG_MOS";
01155 disp_residuals_table_tag = "DISP_RESIDUALS_TABLE_LONG_MOS";
01156 delta_image_tag = "DELTA_IMAGE_LONG_MOS";
01157 spectral_resolution_tag = "SPECTRAL_RESOLUTION_LONG_MOS";
01158 reduced_lamp_tag = "REDUCED_LAMP_LONG_MOS";
01159 disp_coeff_tag = "DISP_COEFF_LONG_MOS";
01160 wavelength_map_tag = "WAVELENGTH_MAP_LONG_MOS";
01161 slit_location_tag = "SLIT_LOCATION_LONG_MOS";
01162 mapped_screen_flat_tag = "MAPPED_SCREEN_FLAT_LONG_MOS";
01163 mapped_norm_flat_tag = "MAPPED_NORM_FLAT_LONG_MOS";
01164 }
01165 else {
01166 master_screen_flat_tag = "MASTER_SCREEN_FLAT_LONG_MXU";
01167 master_norm_flat_tag = "MASTER_NORM_FLAT_LONG_MXU";
01168 disp_residuals_table_tag = "DISP_RESIDUALS_TABLE_LONG_MXU";
01169 delta_image_tag = "DELTA_IMAGE_LONG_MXU";
01170 spectral_resolution_tag = "SPECTRAL_RESOLUTION_LONG_MXU";
01171 reduced_lamp_tag = "REDUCED_LAMP_LONG_MXU";
01172 disp_coeff_tag = "DISP_COEFF_LONG_MXU";
01173 wavelength_map_tag = "WAVELENGTH_MAP_LONG_MXU";
01174 slit_location_tag = "SLIT_LOCATION_LONG_MXU";
01175 mapped_screen_flat_tag = "MAPPED_SCREEN_FLAT_LONG_MXU";
01176 mapped_norm_flat_tag = "MAPPED_NORM_FLAT_LONG_MXU";
01177 }
01178 }
01179 }
01180
01181 if (slit_ident == 0) {
01182 cpl_table_delete(maskslits); maskslits = NULL;
01183 }
01184
01185
01186
01187
01188
01189
01190
01191
01192
01193 if (nbias) {
01194
01195
01196
01197
01198
01199 cpl_msg_info(recipe, "Generate the master from input raw biases...");
01200
01201 if (nbias > 1) {
01202
01203 biases = cpl_imagelist_new();
01204
01205 bias = dfs_load_image(frameset, "BIAS", CPL_TYPE_FLOAT, 0, 0);
01206
01207 if (bias == NULL)
01208 fors_calib_exit("Cannot load bias frame");
01209
01210 cpl_imagelist_set(biases, bias, 0);
01211
01212 for (i = 1; i < nbias; i++) {
01213 bias = dfs_load_image(frameset, NULL, CPL_TYPE_FLOAT, 0, 0);
01214 if (bias)
01215 cpl_imagelist_set(biases, bias, i);
01216 else
01217 fors_calib_exit("Cannot load bias frame");
01218 }
01219
01220 master_bias = cpl_imagelist_collapse_median_create(biases);
01221
01222 cpl_imagelist_delete(biases);
01223 }
01224 else {
01225 master_bias = dfs_load_image(frameset, "BIAS",
01226 CPL_TYPE_FLOAT, 0, 1);
01227 if (master_bias == NULL)
01228 fors_calib_exit("Cannot load bias");
01229 }
01230
01231 }
01232 else {
01233 master_bias = dfs_load_image(frameset, "MASTER_BIAS",
01234 CPL_TYPE_FLOAT, 0, 1);
01235 if (master_bias == NULL)
01236 fors_calib_exit("Cannot load master bias");
01237 }
01238
01239 cpl_msg_info(recipe, "Remove the master bias...");
01240
01241 overscans = mos_load_overscans_vimos(header, 1);
01242 cpl_propertylist_delete(header); header = NULL;
01243
01244 if (nbias) {
01245 int xlow = cpl_table_get_int(overscans, "xlow", 0, NULL);
01246 int ylow = cpl_table_get_int(overscans, "ylow", 0, NULL);
01247 int xhig = cpl_table_get_int(overscans, "xhig", 0, NULL);
01248 int yhig = cpl_table_get_int(overscans, "yhig", 0, NULL);
01249 dummy = cpl_image_extract(master_bias, xlow+1, ylow+1, xhig, yhig);
01250 cpl_image_delete(master_bias); master_bias = dummy;
01251
01252 if (dfs_save_image(frameset, master_bias, "MASTER_BIAS",
01253 NULL, parlist, recipe, version))
01254 fors_calib_exit(NULL);
01255 }
01256
01257 if (nflats > 1) {
01258 multi_bias = cpl_image_multiply_scalar_create(master_bias, nflats);
01259 dummy = mos_remove_bias(master_flat, multi_bias, overscans);
01260 cpl_image_delete(multi_bias);
01261 }
01262 else {
01263 dummy = mos_remove_bias(master_flat, master_bias, overscans);
01264 }
01265 cpl_image_delete(master_flat);
01266 master_flat = dummy;
01267
01268 if (master_flat == NULL)
01269 fors_calib_exit("Cannot remove bias from flat field");
01270
01271 cpl_msg_indent_less();
01272 cpl_msg_info(recipe, "Load arc lamp exposure...");
01273 cpl_msg_indent_more();
01274
01275 spectra = dfs_load_image(frameset, arc_tag, CPL_TYPE_FLOAT, 0, 0);
01276
01277 if (spectra == NULL)
01278 fors_calib_exit("Cannot load arc lamp exposure");
01279
01280 cpl_msg_info(recipe, "Remove the master bias...");
01281
01282 dummy = mos_remove_bias(spectra, master_bias, overscans);
01283 cpl_table_delete(overscans); overscans = NULL;
01284 cpl_image_delete(master_bias); master_bias = NULL;
01285 cpl_image_delete(spectra); spectra = dummy;
01286
01287 if (spectra == NULL)
01288 fors_calib_exit("Cannot remove bias from arc lamp exposure");
01289
01290 cpl_msg_indent_less();
01291 cpl_msg_info(recipe, "Load input line catalog...");
01292 cpl_msg_indent_more();
01293
01294 wavelengths = dfs_load_table(frameset, "MASTER_LINECAT", 1);
01295
01296 if (wavelengths == NULL)
01297 fors_calib_exit("Cannot load line catalog");
01298
01299
01300
01301
01302
01303
01304 nlines = cpl_table_get_nrow(wavelengths);
01305
01306 if (nlines == 0)
01307 fors_calib_exit("Empty input line catalog");
01308
01309 if (cpl_table_has_column(wavelengths, wcolumn) != 1) {
01310 cpl_msg_error(recipe, "Missing column %s in input line catalog table",
01311 wcolumn);
01312 fors_calib_exit(NULL);
01313 }
01314
01315 line = cpl_malloc(nlines * sizeof(double));
01316
01317 for (i = 0; i < nlines; i++)
01318 line[i] = cpl_table_get(wavelengths, wcolumn, i, NULL);
01319
01320 cpl_table_delete(wavelengths); wavelengths = NULL;
01321
01322 lines = cpl_vector_wrap(nlines, line);
01323
01324
01325 if (lss || treat_as_lss) {
01326
01327 int first_row, last_row;
01328 int ylow, yhig;
01329
01330
01331
01332
01333
01334
01335
01336
01337
01338
01339
01340
01341
01342
01343
01344 cpl_msg_indent_less();
01345 cpl_msg_info(recipe, "Perform flat field normalisation...");
01346 cpl_msg_indent_more();
01347
01348 norm_flat = cpl_image_duplicate(master_flat);
01349
01350 smo_flat = mos_normalise_longflat(norm_flat, sradius, dradius,
01351 sdegree);
01352
01353 cpl_image_delete(smo_flat); smo_flat = NULL;
01354
01355
01356
01357 if (1) {
01358 save_header = dfs_load_header(frameset, flat_tag, 0);
01359 cpl_propertylist_update_int(save_header,
01360 "ESO PRO DATANCOM", nflats);
01361
01362 if (dfs_save_image(frameset, master_flat, master_screen_flat_tag,
01363 save_header, parlist, recipe, version))
01364 fors_calib_exit(NULL);
01365
01366 }
01367
01368 if (dfs_save_image(frameset, norm_flat, master_norm_flat_tag,
01369 save_header, parlist, recipe, version))
01370 fors_calib_exit(NULL);
01371
01372 cpl_propertylist_delete(save_header); save_header = NULL;
01373
01374
01375
01376
01377
01378
01379
01380
01381
01382 cpl_msg_indent_less();
01383 cpl_msg_info(recipe, "Perform wavelength calibration...");
01384 cpl_msg_indent_more();
01385
01386 nx = cpl_image_get_size_x(spectra);
01387 ny = cpl_image_get_size_y(spectra);
01388
01389 wavemap = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
01390 idscoeff_all = cpl_table_new(ny);
01391
01392 if (mos_saturation_process(spectra))
01393 fors_calib_exit("Cannot process saturation");
01394
01395 if (mos_subtract_background(spectra))
01396 fors_calib_exit("Cannot subtract the background");
01397
01398 rectified = mos_wavelength_calibration_raw(spectra, lines, dispersion,
01399 peakdetection, wradius,
01400 wdegree, wreject, reference,
01401 &startwavelength,
01402 &endwavelength, NULL,
01403 NULL, idscoeff_all, wavemap,
01404 NULL, NULL, NULL);
01405
01406 if (rectified == NULL)
01407 fors_calib_exit("Wavelength calibration failure.");
01408
01409 cpl_image_delete(rectified); rectified = NULL;
01410
01411 first_row = 0;
01412 while (!cpl_table_is_valid(idscoeff_all, "c0", first_row))
01413 first_row++;
01414
01415 last_row = ny - 1;
01416 while (!cpl_table_is_valid(idscoeff_all, "c0", last_row))
01417 last_row--;
01418
01419 ylow = first_row + 1;
01420 yhig = last_row + 1;
01421
01422 if (ylow >= yhig) {
01423 cpl_error_reset();
01424 fors_calib_exit("No spectra could be detected.");
01425 }
01426
01427 cpl_msg_info(recipe,
01428 "Spectral pattern was detected on %d out of %d CCD rows",
01429 yhig - ylow, ny);
01430
01431 dummy = cpl_image_extract(spectra, 1, ylow, nx, yhig);
01432 cpl_image_delete(spectra); spectra = dummy;
01433
01434 ccd_ysize = ny;
01435 ny = cpl_image_get_size_y(spectra);
01436
01437 if (check)
01438 residual = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
01439
01440 fiterror = cpl_calloc(ny, sizeof(double));
01441 fitlines = cpl_calloc(ny, sizeof(int));
01442 idscoeff = cpl_table_new(ny);
01443 restable = cpl_table_new(nlines);
01444
01445 if (mos_saturation_process(spectra))
01446 fors_calib_exit("Cannot process saturation");
01447
01448 if (mos_subtract_background(spectra))
01449 fors_calib_exit("Cannot subtract the background");
01450
01451 rectified = mos_wavelength_calibration_raw(spectra, lines, dispersion,
01452 peakdetection, wradius,
01453 wdegree, wreject, reference,
01454 &startwavelength,
01455 &endwavelength, fitlines,
01456 fiterror, idscoeff, NULL,
01457 residual, restable, NULL);
01458
01459 if (rectified == NULL)
01460 fors_calib_exit("Wavelength calibration failure.");
01461
01462
01463
01464
01465
01466 slits = cpl_table_new(1);
01467 cpl_table_new_column(slits, "slit_id", CPL_TYPE_INT);
01468 cpl_table_new_column(slits, "xtop", CPL_TYPE_DOUBLE);
01469 cpl_table_new_column(slits, "ytop", CPL_TYPE_DOUBLE);
01470 cpl_table_new_column(slits, "xbottom", CPL_TYPE_DOUBLE);
01471 cpl_table_new_column(slits, "ybottom", CPL_TYPE_DOUBLE);
01472 cpl_table_new_column(slits, "position", CPL_TYPE_INT);
01473 cpl_table_new_column(slits, "length", CPL_TYPE_INT);
01474 cpl_table_set_column_unit(slits, "xtop", "pixel");
01475 cpl_table_set_column_unit(slits, "ytop", "pixel");
01476 cpl_table_set_column_unit(slits, "xbottom", "pixel");
01477 cpl_table_set_column_unit(slits, "ybottom", "pixel");
01478 cpl_table_set_column_unit(slits, "position", "pixel");
01479 cpl_table_set_column_unit(slits, "length", "pixel");
01480 cpl_table_set_int(slits, "slit_id", 0, 0);
01481 cpl_table_set_double(slits, "xtop", 0, 0);
01482 cpl_table_set_double(slits, "ytop", 0, last_row);
01483 cpl_table_set_double(slits, "xbottom", 0, 0);
01484 cpl_table_set_double(slits, "ybottom", 0, first_row);
01485 cpl_table_set_int(slits, "position", 0, 0);
01486 cpl_table_set_int(slits, "length", 0, ny);
01487
01488 if (dfs_save_table(frameset, slits, slit_location_tag, NULL,
01489 parlist, recipe, version))
01490 fors_calib_exit(NULL);
01491
01492 cpl_table_delete(slits); slits = NULL;
01493
01494 if (dfs_save_table(frameset, restable, disp_residuals_table_tag, NULL,
01495 parlist, recipe, version))
01496 fors_calib_exit(NULL);
01497
01498 cpl_table_delete(restable); restable = NULL;
01499
01500 if (wmode) {
01501 cpl_image_delete(rectified); rectified = NULL;
01502 cpl_image_delete(wavemap); wavemap = NULL;
01503 mos_interpolate_wavecalib(idscoeff, wavemap, wmode);
01504 mos_interpolate_wavecalib(idscoeff_all, wavemap, wmode);
01505 wavemap = mos_map_idscoeff(idscoeff_all, nx, reference,
01506 startwavelength, endwavelength);
01507 rectified = mos_wavelength_calibration(spectra, reference,
01508 startwavelength,
01509 endwavelength, dispersion,
01510 idscoeff, 0);
01511 }
01512
01513 cpl_table_delete(idscoeff_all); idscoeff_all = NULL;
01514
01515 cpl_table_wrap_double(idscoeff, fiterror, "error"); fiterror = NULL;
01516 cpl_table_set_column_unit(idscoeff, "error", "pixel");
01517 cpl_table_wrap_int(idscoeff, fitlines, "nlines"); fitlines = NULL;
01518
01519 for (i = 0; i < ny; i++)
01520 if (!cpl_table_is_valid(idscoeff, "c0", i))
01521 cpl_table_set_invalid(idscoeff, "error", i);
01522
01523 delta = mos_map_pixel(idscoeff, reference, startwavelength,
01524 endwavelength, dispersion, 2);
01525
01526
01527 dummy = cpl_image_extract(master_flat, 1, ylow, nx, yhig);
01528 cpl_image_delete(master_flat); master_flat = dummy;
01529
01530 mapped_flat = mos_wavelength_calibration(master_flat, reference,
01531 startwavelength, endwavelength,
01532 dispersion, idscoeff, 0);
01533
01534 cpl_image_delete(master_flat); master_flat = NULL;
01535
01536 dummy = cpl_image_extract(norm_flat, 1, ylow, nx, yhig);
01537 cpl_image_delete(norm_flat); norm_flat = dummy;
01538
01539 mapped_nflat = mos_wavelength_calibration(norm_flat, reference,
01540 startwavelength, endwavelength,
01541 dispersion, idscoeff, 0);
01542
01543 cpl_image_delete(norm_flat); norm_flat = NULL;
01544
01545 header = cpl_propertylist_new();
01546 cpl_propertylist_update_double(header, "CRPIX1", 1.0);
01547 cpl_propertylist_update_double(header, "CRPIX2", 1.0);
01548 cpl_propertylist_update_double(header, "CRVAL1",
01549 startwavelength + dispersion/2);
01550 cpl_propertylist_update_double(header, "CRVAL2", 1.0);
01551
01552
01553 cpl_propertylist_update_double(header, "CD1_1", dispersion);
01554 cpl_propertylist_update_double(header, "CD1_2", 0.0);
01555 cpl_propertylist_update_double(header, "CD2_1", 0.0);
01556 cpl_propertylist_update_double(header, "CD2_2", 1.0);
01557 cpl_propertylist_update_string(header, "CTYPE1", "LINEAR");
01558 cpl_propertylist_update_string(header, "CTYPE2", "PIXEL");
01559
01560 if (dfs_save_image(frameset, delta, delta_image_tag,
01561 header, parlist, recipe, version))
01562 fors_calib_exit(NULL);
01563
01564 cpl_image_delete(delta); delta = NULL;
01565
01566 if (dfs_save_image(frameset, mapped_flat, mapped_screen_flat_tag,
01567 header, parlist, recipe, version))
01568 fors_calib_exit(NULL);
01569
01570 cpl_image_delete(mapped_flat); mapped_flat = NULL;
01571
01572 if (dfs_save_image(frameset, mapped_nflat, mapped_norm_flat_tag,
01573 header, parlist, recipe, version))
01574 fors_calib_exit(NULL);
01575
01576 cpl_image_delete(mapped_nflat); mapped_nflat = NULL;
01577
01578 cpl_propertylist_delete(header); header = NULL;
01579
01580 cpl_msg_info(recipe, "Valid solutions found: %d out of %d rows",
01581 ny - cpl_table_count_invalid(idscoeff, "c0"), ny);
01582
01583 cpl_image_delete(spectra); spectra = NULL;
01584
01585 mean_rms = mos_distortions_rms(rectified, lines, startwavelength,
01586 dispersion, 6, 0);
01587
01588 cpl_msg_info(recipe, "Mean residual: %f pixel", mean_rms);
01589
01590 mean_rms = cpl_table_get_column_mean(idscoeff, "error");
01591 mean_rms_err = cpl_table_get_column_stdev(idscoeff, "error");
01592
01593 cpl_msg_info(recipe, "Mean model accuracy: %f pixel (%f A)",
01594 mean_rms, mean_rms * dispersion);
01595
01596 restab = mos_resolution_table(rectified, startwavelength, dispersion,
01597 60000, lines);
01598
01599 if (restab) {
01600 cpl_msg_info(recipe, "Mean spectral resolution: %.2f",
01601 cpl_table_get_column_mean(restab, "resolution"));
01602 cpl_msg_info(recipe,
01603 "Mean reference lines FWHM: %.2f +/- %.2f pixel",
01604 cpl_table_get_column_mean(restab, "fwhm") / dispersion,
01605 cpl_table_get_column_mean(restab, "fwhm_rms") / dispersion);
01606
01607 if (qc) {
01608
01609 header = dfs_load_header(frameset, arc_tag, 0);
01610
01611 if (header == NULL)
01612 fors_calib_exit("Cannot reload arc lamp header");
01613
01614 qclist = cpl_propertylist_new();
01615
01616 fors_qc_start_group(qclist, "2.0", instrume);
01617
01618
01619
01620
01621
01622
01623 if (fors_qc_write_string("PRO.CATG", spectral_resolution_tag,
01624 "Product category", instrume))
01625 fors_calib_exit("Cannot write product category to "
01626 "QC log file");
01627
01628 if (fors_qc_keyword_to_paf(header, "ESO DPR TYPE", NULL,
01629 "DPR type", instrume))
01630 fors_calib_exit("Missing keyword DPR TYPE in arc "
01631 "lamp header");
01632
01633 if (fors_qc_keyword_to_paf(header, "ESO TPL ID", NULL,
01634 "Template", instrume))
01635 fors_calib_exit("Missing keyword TPL ID in arc "
01636 "lamp header");
01637
01638 if (fors_qc_keyword_to_paf(header, "ESO INS GRIS1 NAME", NULL,
01639 "Grism name", instrume))
01640 fors_calib_exit("Missing keyword INS GRIS1 NAME in arc "
01641 "lamp header");
01642
01643 if (fors_qc_keyword_to_paf(header, "ESO INS GRIS1 ID", NULL,
01644 "Grism identifier", instrume))
01645 fors_calib_exit("Missing keyword INS GRIS1 ID in arc "
01646 "lamp header");
01647
01648 if (cpl_propertylist_has(header, "ESO INS FILT1 NAME"))
01649 fors_qc_keyword_to_paf(header, "ESO INS FILT1 NAME", NULL,
01650 "Filter name", instrume);
01651
01652 if (fors_qc_keyword_to_paf(header, "ESO INS COLL NAME", NULL,
01653 "Collimator name", instrume))
01654 fors_calib_exit("Missing keyword INS COLL NAME in arc "
01655 "lamp header");
01656
01657 if (fors_qc_keyword_to_paf(header, "ESO DET CHIP1 ID", NULL,
01658 "Chip identifier", instrume))
01659 fors_calib_exit("Missing keyword DET CHIP1 ID in arc "
01660 "lamp header");
01661
01662 if (mos) {
01663 if (fors_qc_keyword_to_paf(header, "ESO INS MOS10 WID",
01664 "arcsec", "Slit width", instrume))
01665 fors_calib_exit("Missing keyword ESO INS MOS1 WID in "
01666 "arc lamp header");
01667 }
01668 else {
01669 if (fors_qc_keyword_to_paf(header, "ESO INS SLIT WID",
01670 "arcsec", "Slit width", instrume))
01671 fors_calib_exit("Missing keyword ESO INS SLIT WID in "
01672 "arc lamp header");
01673 }
01674
01675 if (fors_qc_keyword_to_paf(header, "ESO DET WIN1 BINX", NULL,
01676 "Binning factor along X", instrume))
01677 fors_calib_exit("Missing keyword ESO DET WIN1 BINX "
01678 "in frame header");
01679
01680 if (fors_qc_keyword_to_paf(header, "ESO DET WIN1 BINY", NULL,
01681 "Binning factor along Y", instrume))
01682 fors_calib_exit("Missing keyword ESO DET WIN1 BINY "
01683 "in frame header");
01684
01685 if (fors_qc_keyword_to_paf(header, "ARCFILE", NULL,
01686 "Archive name of input data",
01687 instrume))
01688 fors_calib_exit("Missing keyword ARCFILE in arc "
01689 "lamp header");
01690
01691 cpl_propertylist_delete(header); header = NULL;
01692
01693 pipefile = dfs_generate_filename_tfits(spectral_resolution_tag);
01694 if (fors_qc_write_string("PIPEFILE", pipefile,
01695 "Pipeline product name", instrume))
01696 fors_calib_exit("Cannot write PIPEFILE to QC log file");
01697 cpl_free(pipefile); pipefile = NULL;
01698
01699
01700
01701
01702
01703
01704 if (mos)
01705 keyname = "QC.MOS.RESOLUTION";
01706 else
01707 keyname = "QC.LSS.RESOLUTION";
01708
01709 if (fors_qc_write_qc_double(qclist,
01710 cpl_table_get_column_mean(restab,
01711 "resolution"),
01712 keyname, NULL,
01713 "Mean spectral resolution",
01714 instrume)) {
01715 fors_calib_exit("Cannot write mean spectral resolution to "
01716 "QC log file");
01717 }
01718
01719 if (mos)
01720 keyname = "QC.MOS.RESOLUTION.RMS";
01721 else
01722 keyname = "QC.LSS.RESOLUTION.RMS";
01723
01724 if (fors_qc_write_qc_double(qclist,
01725 cpl_table_get_column_stdev(restab,
01726 "resolution"),
01727 keyname, NULL,
01728 "Scatter of spectral resolution",
01729 instrume)) {
01730 fors_calib_exit("Cannot write spectral resolution scatter "
01731 "to QC log file");
01732 }
01733
01734 if (mos)
01735 keyname = "QC.MOS.RESOLUTION.NWAVE";
01736 else
01737 keyname = "QC.LSS.RESOLUTION.NWAVE";
01738
01739 if (fors_qc_write_qc_int(qclist, cpl_table_get_nrow(restab) -
01740 cpl_table_count_invalid(restab,
01741 "resolution"),
01742 keyname, NULL,
01743 "Number of examined wavelengths "
01744 "for resolution computation",
01745 instrume)) {
01746 fors_calib_exit("Cannot write number of lines used in "
01747 "spectral resolution computation "
01748 "to QC log file");
01749 }
01750
01751 if (mos)
01752 keyname = "QC.MOS.RESOLUTION.MEANRMS";
01753 else
01754 keyname = "QC.LSS.RESOLUTION.MEANRMS";
01755
01756 if (fors_qc_write_qc_double(qclist,
01757 cpl_table_get_column_mean(restab,
01758 "resolution_rms"),
01759 keyname, NULL,
01760 "Mean error on spectral "
01761 "resolution computation",
01762 instrume)) {
01763 fors_calib_exit("Cannot write mean error in "
01764 "spectral resolution computation "
01765 "to QC log file");
01766 }
01767
01768 if (mos)
01769 keyname = "QC.MOS.RESOLUTION.NLINES";
01770 else
01771 keyname = "QC.LSS.RESOLUTION.NLINES";
01772
01773 if (fors_qc_write_qc_int(qclist,
01774 cpl_table_get_column_mean(restab, "nlines") *
01775 cpl_table_get_nrow(restab),
01776 keyname, NULL,
01777 "Number of lines for spectral "
01778 "resolution computation",
01779 instrume)) {
01780 fors_calib_exit("Cannot write number of examined "
01781 "wavelengths in spectral resolution computation "
01782 "to QC log file");
01783 }
01784
01785 fors_qc_end_group();
01786
01787 }
01788
01789 if (dfs_save_table(frameset, restab, spectral_resolution_tag,
01790 qclist, parlist, recipe, version))
01791 fors_calib_exit(NULL);
01792
01793 cpl_table_delete(restab); restab = NULL;
01794 cpl_propertylist_delete(qclist); qclist = NULL;
01795
01796 }
01797 else
01798 fors_calib_exit("Cannot compute the spectral resolution table");
01799
01800 cpl_vector_delete(lines); lines = NULL;
01801
01802
01803
01804
01805
01806
01807 header = cpl_propertylist_new();
01808 cpl_propertylist_update_double(header, "CRPIX1", 1.0);
01809 cpl_propertylist_update_double(header, "CRPIX2", 1.0);
01810 cpl_propertylist_update_double(header, "CRVAL1",
01811 startwavelength + dispersion/2);
01812 cpl_propertylist_update_double(header, "CRVAL2", 1.0);
01813
01814
01815 cpl_propertylist_update_double(header, "CD1_1", dispersion);
01816 cpl_propertylist_update_double(header, "CD1_2", 0.0);
01817 cpl_propertylist_update_double(header, "CD2_1", 0.0);
01818 cpl_propertylist_update_double(header, "CD2_2", 1.0);
01819 cpl_propertylist_update_string(header, "CTYPE1", "LINEAR");
01820 cpl_propertylist_update_string(header, "CTYPE2", "PIXEL");
01821 cpl_propertylist_update_int(header, "ESO PRO DATANCOM", 1);
01822
01823 if (dfs_save_image(frameset, rectified, reduced_lamp_tag, header,
01824 parlist, recipe, version))
01825 fors_calib_exit(NULL);
01826
01827 cpl_image_delete(rectified); rectified = NULL;
01828 cpl_propertylist_delete(header); header = NULL;
01829
01830 if (dfs_save_table(frameset, idscoeff, disp_coeff_tag, NULL,
01831 parlist, recipe, version))
01832 fors_calib_exit(NULL);
01833
01834 cpl_table_delete(idscoeff); idscoeff = NULL;
01835
01836 header = dfs_load_header(frameset, arc_tag, 0);
01837
01838 if (header == NULL)
01839 fors_calib_exit("Cannot reload arc lamp header");
01840
01841 if (qc) {
01842
01843 compute_central_wave = 0;
01844 if (lss) {
01845
01846
01847
01848
01849 compute_central_wave = 1;
01850 }
01851 else {
01852 if (fabs(mxpos) < 0.05)
01853 compute_central_wave = 1;
01854 }
01855
01856 fors_qc_start_group(header, "2.0", instrume);
01857
01858
01859
01860
01861
01862 if (fors_qc_write_string("PRO.CATG", wavelength_map_tag,
01863 "Product category", instrume))
01864 fors_calib_exit("Cannot write product category to "
01865 "QC log file");
01866
01867 if (fors_qc_keyword_to_paf(header, "ESO DPR TYPE", NULL,
01868 "DPR type", instrume))
01869 fors_calib_exit("Missing keyword DPR TYPE in arc "
01870 "lamp header");
01871
01872 if (fors_qc_keyword_to_paf(header, "ESO TPL ID", NULL,
01873 "Template", instrume))
01874 fors_calib_exit("Missing keyword TPL ID in arc "
01875 "lamp header");
01876
01877 if (fors_qc_keyword_to_paf(header, "ESO INS GRIS1 NAME", NULL,
01878 "Grism name", instrume))
01879 fors_calib_exit("Missing keyword INS GRIS1 NAME in arc "
01880 "lamp header");
01881
01882 if (fors_qc_keyword_to_paf(header, "ESO INS GRIS1 ID", NULL,
01883 "Grism identifier", instrume))
01884 fors_calib_exit("Missing keyword INS GRIS1 ID in arc "
01885 "lamp header");
01886
01887 if (cpl_propertylist_has(header, "ESO INS FILT1 NAME"))
01888 fors_qc_keyword_to_paf(header, "ESO INS FILT1 NAME", NULL,
01889 "Filter name", instrume);
01890
01891 if (fors_qc_keyword_to_paf(header, "ESO INS COLL NAME", NULL,
01892 "Collimator name", instrume))
01893 fors_calib_exit("Missing keyword INS COLL NAME in arc "
01894 "lamp header");
01895
01896 if (fors_qc_keyword_to_paf(header, "ESO DET CHIP1 ID", NULL,
01897 "Chip identifier", instrume))
01898 fors_calib_exit("Missing keyword DET CHIP1 ID in arc "
01899 "lamp header");
01900
01901 if (mos) {
01902 if (fors_qc_keyword_to_paf(header, "ESO INS MOS10 WID",
01903 "arcsec", "Slit width", instrume))
01904 fors_calib_exit("Missing keyword ESO INS MOS1 WID in "
01905 "arc lamp header");
01906 }
01907 else {
01908 if (fors_qc_keyword_to_paf(header, "ESO INS SLIT WID",
01909 "arcsec", "Slit width", instrume))
01910 fors_calib_exit("Missing keyword ESO INS SLIT WID in "
01911 "arc lamp header");
01912 }
01913
01914 if (fors_qc_keyword_to_paf(header, "ESO DET WIN1 BINX", NULL,
01915 "Binning factor along X", instrume))
01916 fors_calib_exit("Missing keyword ESO DET WIN1 BINX "
01917 "in frame header");
01918
01919 if (fors_qc_keyword_to_paf(header, "ESO DET WIN1 BINY", NULL,
01920 "Binning factor along Y", instrume))
01921 fors_calib_exit("Missing keyword ESO DET WIN1 BINY "
01922 "in frame header");
01923
01924 if (fors_qc_keyword_to_paf(header, "ARCFILE", NULL,
01925 "Archive name of input data",
01926 instrume))
01927 fors_calib_exit("Missing keyword ARCFILE in arc "
01928 "lamp header");
01929
01930 pipefile = dfs_generate_filename(wavelength_map_tag);
01931 if (fors_qc_write_string("PIPEFILE", pipefile,
01932 "Pipeline product name", instrume))
01933 fors_calib_exit("Cannot write PIPEFILE to QC log file");
01934 cpl_free(pipefile); pipefile = NULL;
01935
01936
01937
01938
01939
01940
01941 if (fors_qc_write_qc_double(header,
01942 mean_rms,
01943 "QC.WAVE.ACCURACY",
01944 "pixel",
01945 "Mean accuracy of wavecalib model",
01946 instrume)) {
01947 fors_calib_exit("Cannot write mean wavelength calibration "
01948 "accuracy to QC log file");
01949 }
01950
01951 if (fors_qc_write_qc_double(header,
01952 mean_rms_err,
01953 "QC.WAVE.ACCURACY.ERROR",
01954 "pixel",
01955 "Error on accuracy of wavecalib model",
01956 instrume)) {
01957 fors_calib_exit("Cannot write error on wavelength calibration "
01958 "accuracy to QC log file");
01959 }
01960
01961 if (compute_central_wave) {
01962
01963 data = cpl_image_get_data(wavemap);
01964
01965 if (lss) {
01966 if (fors_qc_write_qc_double(header,
01967 data[nx/2 + ccd_ysize*nx/2],
01968 "QC.LSS.CENTRAL.WAVELENGTH",
01969 "Angstrom",
01970 "Wavelength at CCD center",
01971 instrume)) {
01972 fors_calib_exit("Cannot write central wavelength to QC "
01973 "log file");
01974 }
01975 }
01976 else {
01977 if (fors_qc_write_qc_double(header,
01978 data[nx/2 + ccd_ysize*nx/2],
01979 "QC.MOS.CENTRAL.WAVELENGTH",
01980 "Angstrom",
01981 "Wavelength at CCD center",
01982 instrume)) {
01983 fors_calib_exit("Cannot write central wavelength to QC "
01984 "log file");
01985 }
01986 }
01987 }
01988
01989 fors_qc_end_group();
01990
01991 }
01992
01993 if (dfs_save_image(frameset, wavemap, wavelength_map_tag, header,
01994 parlist, recipe, version))
01995 fors_calib_exit(NULL);
01996
01997 cpl_image_delete(wavemap); wavemap = NULL;
01998
01999 cpl_propertylist_erase_regexp(header, "^ESO QC ", 0);
02000
02001 if (check) {
02002 cpl_propertylist_update_double(header, "CRPIX2", 1.0);
02003 cpl_propertylist_update_double(header, "CRVAL2", 1.0);
02004
02005 cpl_propertylist_update_double(header, "CD1_1", 1.0);
02006 cpl_propertylist_update_double(header, "CD1_2", 0.0);
02007 cpl_propertylist_update_double(header, "CD2_1", 0.0);
02008 cpl_propertylist_update_double(header, "CD2_2", 1.0);
02009 cpl_propertylist_update_string(header, "CTYPE1", "LINEAR");
02010 cpl_propertylist_update_string(header, "CTYPE2", "PIXEL");
02011
02012 if (dfs_save_image(frameset, residual, disp_residuals_tag, header,
02013 parlist, recipe, version))
02014 fors_calib_exit(NULL);
02015
02016 cpl_image_delete(residual); residual = NULL;
02017 }
02018
02019 cpl_propertylist_delete(header); header = NULL;
02020 cpl_free(instrume); instrume = NULL;
02021
02022 return 0;
02023
02024 }
02025
02026
02027
02028
02029
02030
02031
02032
02033
02034
02035 cpl_msg_indent_less();
02036 cpl_msg_info(recipe, "Detecting spectra on CCD...");
02037 cpl_msg_indent_more();
02038
02039 ccd_xsize = nx = cpl_image_get_size_x(spectra);
02040 ccd_ysize = ny = cpl_image_get_size_y(spectra);
02041
02042 refmask = cpl_mask_new(nx, ny);
02043
02044 if (mos_saturation_process(spectra))
02045 fors_calib_exit("Cannot process saturation");
02046
02047 if (mos_subtract_background(spectra))
02048 fors_calib_exit("Cannot subtract the background");
02049
02050 checkwave = mos_wavelength_calibration_raw(spectra, lines, dispersion,
02051 peakdetection, wradius,
02052 wdegree, wreject, reference,
02053 &startwavelength, &endwavelength,
02054 NULL, NULL, NULL, NULL, NULL,
02055 NULL, refmask);
02056
02057 if (checkwave == NULL)
02058 fors_calib_exit("Wavelength calibration failure.");
02059
02060
02061
02062
02063
02064 header = cpl_propertylist_new();
02065 cpl_propertylist_update_double(header, "CRPIX1", 1.0);
02066 cpl_propertylist_update_double(header, "CRPIX2", 1.0);
02067 cpl_propertylist_update_double(header, "CRVAL1",
02068 startwavelength + dispersion/2);
02069 cpl_propertylist_update_double(header, "CRVAL2", 1.0);
02070
02071
02072 cpl_propertylist_update_double(header, "CD1_1", dispersion);
02073 cpl_propertylist_update_double(header, "CD1_2", 0.0);
02074 cpl_propertylist_update_double(header, "CD2_1", 0.0);
02075 cpl_propertylist_update_double(header, "CD2_2", 1.0);
02076 cpl_propertylist_update_string(header, "CTYPE1", "LINEAR");
02077 cpl_propertylist_update_string(header, "CTYPE2", "PIXEL");
02078
02079 if (check) {
02080 if (dfs_save_image(frameset, checkwave, spectra_detection_tag, header,
02081 parlist, recipe, version))
02082 fors_calib_exit(NULL);
02083 }
02084
02085 cpl_image_delete(checkwave); checkwave = NULL;
02086 cpl_propertylist_delete(header); header = NULL;
02087
02088 cpl_msg_info(recipe, "Locate slits at reference wavelength on CCD...");
02089 slits = mos_locate_spectra(refmask);
02090
02091 if (!slits) {
02092 cpl_msg_error(cpl_error_get_where(), cpl_error_get_message());
02093 fors_calib_exit("No slits could be detected!");
02094 }
02095
02096 refimage = cpl_image_new_from_mask(refmask);
02097 cpl_mask_delete(refmask); refmask = NULL;
02098
02099 if (check) {
02100 save_header = dfs_load_header(frameset, arc_tag, 0);
02101 if (dfs_save_image(frameset, refimage, slit_map_tag, save_header,
02102 parlist, recipe, version))
02103 fors_calib_exit(NULL);
02104 cpl_propertylist_delete(save_header); save_header = NULL;
02105 }
02106
02107 cpl_image_delete(refimage); refimage = NULL;
02108
02109 if (slit_ident) {
02110
02111
02112
02113
02114
02115
02116
02117
02118
02119
02120
02121
02122
02123
02124
02125 cpl_msg_indent_less();
02126 cpl_msg_info(recipe, "Attempt slit identification (optional)...");
02127 cpl_msg_indent_more();
02128
02129 positions = mos_identify_slits(slits, maskslits, NULL);
02130
02131 if (positions) {
02132 cpl_table_delete(slits);
02133 slits = positions;
02134
02135
02136
02137
02138
02139 cpl_table_and_selected_double(slits,
02140 "ybottom", CPL_GREATER_THAN, ny-1);
02141 cpl_table_or_selected_double(slits,
02142 "ytop", CPL_LESS_THAN, 0);
02143 cpl_table_erase_selected(slits);
02144
02145 nslits = cpl_table_get_nrow(slits);
02146
02147 if (nslits == 0)
02148 fors_calib_exit("No slits found on the CCD");
02149
02150 cpl_msg_info(recipe, "%d slits are entirely or partially "
02151 "contained in CCD", nslits);
02152
02153 }
02154 else {
02155 slit_ident = 0;
02156 cpl_msg_info(recipe, "Global distortion model cannot be computed");
02157 if (cpl_error_get_code() != CPL_ERROR_NONE) {
02158 fors_calib_exit(NULL);
02159 }
02160 }
02161 }
02162
02163
02164
02165
02166
02167
02168 cpl_msg_indent_less();
02169 cpl_msg_info(recipe, "Determining spectral curvature...");
02170 cpl_msg_indent_more();
02171
02172 cpl_msg_info(recipe, "Tracing master flat field spectra edges...");
02173 traces = mos_trace_flat(master_flat, slits, reference,
02174 startwavelength, endwavelength, dispersion);
02175
02176 if (!traces)
02177 fors_calib_exit("Tracing failure");
02178
02179 cpl_msg_info(recipe, "Fitting flat field spectra edges...");
02180 polytraces = mos_poly_trace(slits, traces, cdegree);
02181
02182 if (!polytraces)
02183 fors_calib_exit("Trace fitting failure");
02184
02185 if (cmode) {
02186 cpl_msg_info(recipe, "Computing global spectral curvature model...");
02187 mos_global_trace(slits, polytraces, cmode);
02188 }
02189
02190 if (dfs_save_table(frameset, traces, curv_traces_tag, NULL, parlist,
02191 recipe, version))
02192 fors_calib_exit(NULL);
02193
02194 cpl_table_delete(traces); traces = NULL;
02195
02196 coordinate = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
02197 spatial = mos_spatial_calibration(spectra, slits, polytraces, reference,
02198 startwavelength, endwavelength,
02199 dispersion, 0, coordinate);
02200
02201 if (!slit_ident) {
02202 cpl_image_delete(spectra); spectra = NULL;
02203 }
02204
02205
02206
02207
02208
02209
02210
02211 cpl_msg_indent_less();
02212 cpl_msg_info(recipe, "Perform flat field normalisation...");
02213 cpl_msg_indent_more();
02214
02215 norm_flat = cpl_image_duplicate(master_flat);
02216
02217 smo_flat = mos_normalise_flat(norm_flat, coordinate, slits, polytraces,
02218 reference, startwavelength, endwavelength,
02219 dispersion, dradius, ddegree);
02220
02221 cpl_image_delete(smo_flat); smo_flat = NULL;
02222
02223
02224 save_header = dfs_load_header(frameset, flat_tag, 0);
02225 cpl_propertylist_update_int(save_header, "ESO PRO DATANCOM", nflats);
02226
02227 rect_flat = mos_spatial_calibration(master_flat, slits, polytraces,
02228 reference, startwavelength,
02229 endwavelength, dispersion, 0, NULL);
02230 rect_nflat = mos_spatial_calibration(norm_flat, slits, polytraces,
02231 reference, startwavelength,
02232 endwavelength, dispersion, 0, NULL);
02233
02234
02235 if (dfs_save_image(frameset, master_flat, master_screen_flat_tag,
02236 save_header, parlist, recipe, version))
02237 fors_calib_exit(NULL);
02238
02239 cpl_image_delete(master_flat); master_flat = NULL;
02240
02241 if (dfs_save_image(frameset, norm_flat, master_norm_flat_tag,
02242 save_header, parlist, recipe, version))
02243 fors_calib_exit(NULL);
02244
02245 cpl_image_delete(norm_flat); norm_flat = NULL;
02246 cpl_propertylist_delete(save_header); save_header = NULL;
02247
02248
02249
02250
02251
02252
02253
02254 cpl_msg_indent_less();
02255 cpl_msg_info(recipe, "Perform final wavelength calibration...");
02256 cpl_msg_indent_more();
02257
02258 nx = cpl_image_get_size_x(spatial);
02259 ny = cpl_image_get_size_y(spatial);
02260
02261 idscoeff = cpl_table_new(ny);
02262 restable = cpl_table_new(nlines);
02263 rainbow = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
02264 if (check)
02265 residual = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
02266 fiterror = cpl_calloc(ny, sizeof(double));
02267 fitlines = cpl_calloc(ny, sizeof(int));
02268
02269 rectified = mos_wavelength_calibration_final(spatial, slits, lines,
02270 dispersion, peakdetection,
02271 wradius, wdegree, wreject,
02272 reference, &startwavelength,
02273 &endwavelength, fitlines,
02274 fiterror, idscoeff, rainbow,
02275 residual, restable);
02276
02277
02278
02279
02280
02281 cpl_image_delete(spatial); spatial = NULL;
02282
02283 if (rectified == NULL)
02284 fors_calib_exit("Wavelength calibration failure.");
02285
02286 if (dfs_save_table(frameset, restable, disp_residuals_table_tag, NULL,
02287 parlist, recipe, version))
02288 fors_calib_exit(NULL);
02289
02290 cpl_table_delete(restable); restable = NULL;
02291
02292 cpl_table_wrap_double(idscoeff, fiterror, "error"); fiterror = NULL;
02293 cpl_table_set_column_unit(idscoeff, "error", "pixel");
02294 cpl_table_wrap_int(idscoeff, fitlines, "nlines"); fitlines = NULL;
02295
02296 for (i = 0; i < ny; i++)
02297 if (!cpl_table_is_valid(idscoeff, "c0", i))
02298 cpl_table_set_invalid(idscoeff, "error", i);
02299
02300 delta = mos_map_pixel(idscoeff, reference, startwavelength,
02301 endwavelength, dispersion, 2);
02302
02303 header = cpl_propertylist_new();
02304 cpl_propertylist_update_double(header, "CRPIX1", 1.0);
02305 cpl_propertylist_update_double(header, "CRPIX2", 1.0);
02306 cpl_propertylist_update_double(header, "CRVAL1",
02307 startwavelength + dispersion/2);
02308 cpl_propertylist_update_double(header, "CRVAL2", 1.0);
02309
02310
02311 cpl_propertylist_update_double(header, "CD1_1", dispersion);
02312 cpl_propertylist_update_double(header, "CD1_2", 0.0);
02313 cpl_propertylist_update_double(header, "CD2_1", 0.0);
02314 cpl_propertylist_update_double(header, "CD2_2", 1.0);
02315 cpl_propertylist_update_string(header, "CTYPE1", "LINEAR");
02316 cpl_propertylist_update_string(header, "CTYPE2", "PIXEL");
02317
02318 if (dfs_save_image(frameset, delta, delta_image_tag,
02319 header, parlist, recipe, version))
02320 fors_calib_exit(NULL);
02321
02322 cpl_image_delete(delta); delta = NULL;
02323 cpl_propertylist_delete(header); header = NULL;
02324
02325 mean_rms = mos_distortions_rms(rectified, lines, startwavelength,
02326 dispersion, 6, 0);
02327
02328 cpl_msg_info(recipe, "Mean residual: %f pixel", mean_rms);
02329
02330 mean_rms = cpl_table_get_column_mean(idscoeff, "error");
02331 mean_rms_err = cpl_table_get_column_stdev(idscoeff, "error");
02332
02333 cpl_msg_info(recipe, "Mean model accuracy: %f pixel (%f A)",
02334 mean_rms, mean_rms * dispersion);
02335
02336 restab = mos_resolution_table(rectified, startwavelength, dispersion,
02337 60000, lines);
02338
02339 if (restab) {
02340 cpl_msg_info(recipe, "Mean spectral resolution: %.2f",
02341 cpl_table_get_column_mean(restab, "resolution"));
02342 cpl_msg_info(recipe, "Mean reference lines FWHM: %.2f +/- %.2f pixel",
02343 cpl_table_get_column_mean(restab, "fwhm") / dispersion,
02344 cpl_table_get_column_mean(restab, "fwhm_rms") / dispersion);
02345
02346 if (qc) {
02347
02348 header = dfs_load_header(frameset, arc_tag, 0);
02349
02350 if (header == NULL)
02351 fors_calib_exit("Cannot reload arc lamp header");
02352
02353 qclist = cpl_propertylist_new();
02354
02355 fors_qc_start_group(qclist, "2.0", instrume);
02356
02357
02358
02359
02360
02361
02362 if (fors_qc_write_string("PRO.CATG", spectral_resolution_tag,
02363 "Product category", instrume))
02364 fors_calib_exit("Cannot write product category to "
02365 "QC log file");
02366
02367 if (fors_qc_keyword_to_paf(header, "ESO DPR TYPE", NULL,
02368 "DPR type", instrume))
02369 fors_calib_exit("Missing keyword DPR TYPE in arc "
02370 "lamp header");
02371
02372 if (fors_qc_keyword_to_paf(header, "ESO TPL ID", NULL,
02373 "Template", instrume))
02374 fors_calib_exit("Missing keyword TPL ID in arc "
02375 "lamp header");
02376
02377 if (fors_qc_keyword_to_paf(header, "ESO INS GRIS1 NAME", NULL,
02378 "Grism name", instrume))
02379 fors_calib_exit("Missing keyword INS GRIS1 NAME in arc "
02380 "lamp header");
02381
02382 if (fors_qc_keyword_to_paf(header, "ESO INS GRIS1 ID", NULL,
02383 "Grism identifier", instrume))
02384 fors_calib_exit("Missing keyword INS GRIS1 ID in arc "
02385 "lamp header");
02386
02387 if (cpl_propertylist_has(header, "ESO INS FILT1 NAME"))
02388 fors_qc_keyword_to_paf(header, "ESO INS FILT1 NAME", NULL,
02389 "Filter name", instrume);
02390
02391 if (fors_qc_keyword_to_paf(header, "ESO INS COLL NAME", NULL,
02392 "Collimator name", instrume))
02393 fors_calib_exit("Missing keyword INS COLL NAME in arc "
02394 "lamp header");
02395
02396 if (fors_qc_keyword_to_paf(header, "ESO DET CHIP1 ID", NULL,
02397 "Chip identifier", instrume))
02398 fors_calib_exit("Missing keyword DET CHIP1 ID in arc "
02399 "lamp header");
02400
02401 if (fors_qc_keyword_to_paf(header, "ARCFILE", NULL,
02402 "Archive name of input data",
02403 instrume))
02404 fors_calib_exit("Missing keyword ARCFILE in arc "
02405 "lamp header");
02406
02407 cpl_propertylist_delete(header); header = NULL;
02408
02409 pipefile = dfs_generate_filename_tfits(spectral_resolution_tag);
02410 if (fors_qc_write_string("PIPEFILE", pipefile,
02411 "Pipeline product name", instrume))
02412 fors_calib_exit("Cannot write PIPEFILE to QC log file");
02413 cpl_free(pipefile); pipefile = NULL;
02414
02415
02416
02417
02418
02419
02420 if (mos)
02421 keyname = "QC.MOS.RESOLUTION";
02422 else
02423 keyname = "QC.MXU.RESOLUTION";
02424
02425 if (fors_qc_write_qc_double(qclist,
02426 cpl_table_get_column_mean(restab,
02427 "resolution"),
02428 keyname,
02429 "Angstrom",
02430 "Mean spectral resolution",
02431 instrume)) {
02432 fors_calib_exit("Cannot write mean spectral resolution to QC "
02433 "log file");
02434 }
02435
02436 if (mos)
02437 keyname = "QC.MOS.RESOLUTION.RMS";
02438 else
02439 keyname = "QC.MXU.RESOLUTION.RMS";
02440
02441 if (fors_qc_write_qc_double(qclist,
02442 cpl_table_get_column_stdev(restab,
02443 "resolution"),
02444 keyname,
02445 "Angstrom",
02446 "Scatter of spectral resolution",
02447 instrume)) {
02448 fors_calib_exit("Cannot write spectral resolution scatter "
02449 "to QC log file");
02450 }
02451
02452 if (mos)
02453 keyname = "QC.MOS.RESOLUTION.NWAVE";
02454 else
02455 keyname = "QC.MXU.RESOLUTION.NWAVE";
02456
02457 if (fors_qc_write_qc_int(qclist, cpl_table_get_nrow(restab) -
02458 cpl_table_count_invalid(restab,
02459 "resolution"),
02460 keyname,
02461 NULL,
02462 "Number of examined wavelengths "
02463 "for resolution computation",
02464 instrume)) {
02465 fors_calib_exit("Cannot write number of lines used in "
02466 "spectral resolution computation "
02467 "to QC log file");
02468 }
02469
02470 if (mos)
02471 keyname = "QC.MOS.RESOLUTION.MEANRMS";
02472 else
02473 keyname = "QC.MXU.RESOLUTION.MEANRMS";
02474
02475 if (fors_qc_write_qc_double(qclist,
02476 cpl_table_get_column_mean(restab,
02477 "resolution_rms"),
02478 keyname, NULL,
02479 "Mean error on spectral "
02480 "resolution computation",
02481 instrume)) {
02482 fors_calib_exit("Cannot write mean error in "
02483 "spectral resolution computation "
02484 "to QC log file");
02485 }
02486
02487 if (mos)
02488 keyname = "QC.MOS.RESOLUTION.NLINES";
02489 else
02490 keyname = "QC.MXU.RESOLUTION.NLINES";
02491
02492 if (fors_qc_write_qc_int(qclist,
02493 cpl_table_get_column_mean(restab, "nlines") *
02494 cpl_table_get_nrow(restab),
02495 keyname, NULL,
02496 "Number of lines for spectral "
02497 "resolution computation",
02498 instrume)) {
02499 fors_calib_exit("Cannot write number of examined "
02500 "wavelengths in spectral resolution computation "
02501 "to QC log file");
02502 }
02503
02504 fors_qc_end_group();
02505
02506 }
02507
02508 if (dfs_save_table(frameset, restab, spectral_resolution_tag, qclist,
02509 parlist, recipe, version))
02510 fors_calib_exit(NULL);
02511
02512 cpl_table_delete(restab); restab = NULL;
02513 cpl_propertylist_delete(qclist); qclist = NULL;
02514
02515 }
02516 else
02517 fors_calib_exit("Cannot compute the spectral resolution table");
02518
02519 cpl_vector_delete(lines); lines = NULL;
02520
02521 if (dfs_save_table(frameset, idscoeff, disp_coeff_tag, NULL,
02522 parlist, recipe, version))
02523 fors_calib_exit(NULL);
02524
02525
02526 mapped_flat = mos_wavelength_calibration(rect_flat, reference,
02527 startwavelength, endwavelength,
02528 dispersion, idscoeff, 0);
02529
02530 mapped_nflat = mos_wavelength_calibration(rect_nflat, reference,
02531 startwavelength, endwavelength,
02532 dispersion, idscoeff, 0);
02533
02534 cpl_image_delete(rect_flat); rect_flat = NULL;
02535 cpl_image_delete(rect_nflat); rect_nflat = NULL;
02536
02537
02538
02539
02540
02541
02542 if (slit_ident) {
02543
02544 cpl_msg_info(recipe, "Computing global distortions model");
02545 global = mos_global_distortion(slits, maskslits, idscoeff,
02546 polytraces, reference);
02547
02548 if (global && 0) {
02549 cpl_table *stest;
02550 cpl_table *ctest;
02551 cpl_table *dtest;
02552 cpl_image *itest;
02553
02554 stest = mos_build_slit_location(global, maskslits, ccd_ysize);
02555
02556 ctest = mos_build_curv_coeff(global, maskslits, stest);
02557 if (dfs_save_table(frameset, ctest, "CURVS", NULL,
02558 parlist, recipe, version))
02559 fors_calib_exit(NULL);
02560
02561 itest = mos_spatial_calibration(spectra, stest, ctest,
02562 reference, startwavelength,
02563 endwavelength, dispersion,
02564 0, NULL);
02565 cpl_table_delete(ctest); ctest = NULL;
02566 cpl_image_delete(itest); itest = NULL;
02567 if (dfs_save_table(frameset, stest, "SLITS", NULL,
02568 parlist, recipe, version))
02569 fors_calib_exit(NULL);
02570
02571 dtest = mos_build_disp_coeff(global, stest);
02572 if (dfs_save_table(frameset, dtest, "DISPS", NULL,
02573 parlist, recipe, version))
02574 fors_calib_exit(NULL);
02575
02576 cpl_table_delete(dtest); dtest = NULL;
02577 cpl_table_delete(stest); stest = NULL;
02578 }
02579
02580 if (global) {
02581 if (dfs_save_table(frameset, global, global_distortion_tag, NULL,
02582 parlist, recipe, version))
02583 fors_calib_exit(NULL);
02584 cpl_table_delete(global); global = NULL;
02585 }
02586
02587 cpl_image_delete(spectra); spectra = NULL;
02588 cpl_table_delete(maskslits); maskslits = NULL;
02589 }
02590
02591 cpl_table_delete(idscoeff); idscoeff = NULL;
02592
02593 header = cpl_propertylist_new();
02594 cpl_propertylist_update_double(header, "CRPIX1", 1.0);
02595 cpl_propertylist_update_double(header, "CRPIX2", 1.0);
02596 cpl_propertylist_update_double(header, "CRVAL1",
02597 startwavelength + dispersion/2);
02598 cpl_propertylist_update_double(header, "CRVAL2", 1.0);
02599
02600
02601 cpl_propertylist_update_double(header, "CD1_1", dispersion);
02602 cpl_propertylist_update_double(header, "CD1_2", 0.0);
02603 cpl_propertylist_update_double(header, "CD2_1", 0.0);
02604 cpl_propertylist_update_double(header, "CD2_2", 1.0);
02605 cpl_propertylist_update_string(header, "CTYPE1", "LINEAR");
02606 cpl_propertylist_update_string(header, "CTYPE2", "PIXEL");
02607 cpl_propertylist_update_int(header, "ESO PRO DATANCOM", 1);
02608
02609 if (dfs_save_image(frameset, rectified, reduced_lamp_tag, header,
02610 parlist, recipe, version))
02611 fors_calib_exit(NULL);
02612
02613 cpl_image_delete(rectified); rectified = NULL;
02614
02615 cpl_propertylist_update_int(header, "ESO PRO DATANCOM", nflats);
02616
02617 if (dfs_save_image(frameset, mapped_flat, mapped_screen_flat_tag, header,
02618 parlist, recipe, version))
02619 fors_calib_exit(NULL);
02620
02621 cpl_image_delete(mapped_flat); mapped_flat = NULL;
02622
02623 if (dfs_save_image(frameset, mapped_nflat, mapped_norm_flat_tag, header,
02624 parlist, recipe, version))
02625 fors_calib_exit(NULL);
02626
02627 cpl_image_delete(mapped_nflat); mapped_nflat = NULL;
02628
02629 cpl_propertylist_delete(header); header = NULL;
02630
02631 if (check) {
02632 save_header = dfs_load_header(frameset, arc_tag, 0);
02633
02634 cpl_propertylist_update_double(save_header, "CRPIX2", 1.0);
02635 cpl_propertylist_update_double(save_header, "CRVAL2", 1.0);
02636
02637 cpl_propertylist_update_double(save_header, "CD1_1", 1.0);
02638 cpl_propertylist_update_double(save_header, "CD1_2", 0.0);
02639 cpl_propertylist_update_double(save_header, "CD2_1", 0.0);
02640 cpl_propertylist_update_double(save_header, "CD2_2", 1.0);
02641 cpl_propertylist_update_string(save_header, "CTYPE1", "LINEAR");
02642 cpl_propertylist_update_string(save_header, "CTYPE2", "PIXEL");
02643
02644 if (dfs_save_image(frameset, residual, disp_residuals_tag, save_header,
02645 parlist, recipe, version))
02646 fors_calib_exit(NULL);
02647
02648 cpl_image_delete(residual); residual = NULL;
02649 cpl_propertylist_delete(save_header); save_header = NULL;
02650 }
02651
02652 wavemap = mos_map_wavelengths(coordinate, rainbow, slits, polytraces,
02653 reference, startwavelength, endwavelength,
02654 dispersion);
02655
02656 cpl_image_delete(rainbow); rainbow = NULL;
02657
02658 save_header = dfs_load_header(frameset, arc_tag, 0);
02659
02660 if (qc) {
02661 fors_qc_start_group(save_header, "2.0", instrume);
02662
02663
02664
02665
02666
02667 if (fors_qc_write_string("PRO.CATG", wavelength_map_tag,
02668 "Product category", instrume))
02669 fors_calib_exit("Cannot write product category to "
02670 "QC log file");
02671
02672 if (fors_qc_keyword_to_paf(save_header, "ESO DPR TYPE", NULL,
02673 "DPR type", instrume))
02674 fors_calib_exit("Missing keyword DPR TYPE in arc "
02675 "lamp header");
02676
02677 if (fors_qc_keyword_to_paf(save_header, "ESO TPL ID", NULL,
02678 "Template", instrume))
02679 fors_calib_exit("Missing keyword TPL ID in arc "
02680 "lamp header");
02681
02682 if (fors_qc_keyword_to_paf(save_header, "ESO INS GRIS1 NAME", NULL,
02683 "Grism name", instrume))
02684 fors_calib_exit("Missing keyword INS GRIS1 NAME in arc "
02685 "lamp header");
02686
02687 if (fors_qc_keyword_to_paf(save_header, "ESO INS GRIS1 ID", NULL,
02688 "Grism identifier", instrume))
02689 fors_calib_exit("Missing keyword INS GRIS1 ID in arc "
02690 "lamp header");
02691
02692 if (cpl_propertylist_has(save_header, "ESO INS FILT1 NAME"))
02693 fors_qc_keyword_to_paf(save_header, "ESO INS FILT1 NAME", NULL,
02694 "Filter name", instrume);
02695
02696 if (fors_qc_keyword_to_paf(save_header, "ESO INS COLL NAME", NULL,
02697 "Collimator name", instrume))
02698 fors_calib_exit("Missing keyword INS COLL NAME in arc "
02699 "lamp header");
02700
02701 if (fors_qc_keyword_to_paf(save_header, "ESO DET CHIP1 ID", NULL,
02702 "Chip identifier", instrume))
02703 fors_calib_exit("Missing keyword DET CHIP1 ID in arc "
02704 "lamp header");
02705
02706 if (fors_qc_keyword_to_paf(save_header, "ESO DET WIN1 BINX", NULL,
02707 "Binning factor along X", instrume))
02708 fors_calib_exit("Missing keyword ESO DET WIN1 BINX "
02709 "in frame header");
02710
02711 if (fors_qc_keyword_to_paf(save_header, "ESO DET WIN1 BINY", NULL,
02712 "Binning factor along Y", instrume))
02713 fors_calib_exit("Missing keyword ESO DET WIN1 BINY "
02714 "in frame header");
02715
02716 if (fors_qc_keyword_to_paf(save_header, "ARCFILE", NULL,
02717 "Archive name of input data",
02718 instrume))
02719 fors_calib_exit("Missing keyword ARCFILE in arc "
02720 "lamp header");
02721
02722 pipefile = dfs_generate_filename(wavelength_map_tag);
02723 if (fors_qc_write_string("PIPEFILE", pipefile,
02724 "Pipeline product name", instrume))
02725 fors_calib_exit("Cannot write PIPEFILE to QC log file");
02726 cpl_free(pipefile); pipefile = NULL;
02727
02728
02729
02730
02731
02732
02733 if (fors_qc_write_qc_double(save_header,
02734 mean_rms,
02735 "QC.WAVE.ACCURACY",
02736 "pixel",
02737 "Mean accuracy of wavecalib model",
02738 instrume)) {
02739 fors_calib_exit("Cannot write mean wavelength calibration "
02740 "accuracy to QC log file");
02741 }
02742
02743
02744 if (fors_qc_write_qc_double(save_header,
02745 mean_rms_err,
02746 "QC.WAVE.ACCURACY.ERROR",
02747 "pixel",
02748 "Error on accuracy of wavecalib model",
02749 instrume)) {
02750 fors_calib_exit("Cannot write error on wavelength calibration "
02751 "accuracy to QC log file");
02752 }
02753
02754 fors_qc_end_group();
02755
02756 }
02757
02758 if (dfs_save_image(frameset, wavemap, wavelength_map_tag, save_header,
02759 parlist, recipe, version))
02760 fors_calib_exit(NULL);
02761
02762 cpl_image_delete(wavemap); wavemap = NULL;
02763
02764 cpl_propertylist_erase_regexp(save_header, "^ESO QC ", 0);
02765
02766 if (dfs_save_image(frameset, coordinate, spatial_map_tag, save_header,
02767 parlist, recipe, version))
02768 fors_calib_exit(NULL);
02769
02770 cpl_image_delete(coordinate); coordinate = NULL;
02771 cpl_propertylist_delete(save_header); save_header = NULL;
02772
02773 header = NULL;
02774
02775 if (qc) {
02776
02777 double maxpos, maxneg, maxcurve, maxslope;
02778
02779 header = dfs_load_header(frameset, arc_tag, 0);
02780
02781 fors_qc_start_group(header, "2.0", instrume);
02782
02783
02784
02785
02786
02787 if (fors_qc_write_string("PRO.CATG", curv_coeff_tag,
02788 "Product category", instrume))
02789 fors_calib_exit("Cannot write product category to "
02790 "QC log file");
02791
02792 if (fors_qc_keyword_to_paf(header, "ESO DPR TYPE", NULL,
02793 "DPR type", instrume))
02794 fors_calib_exit("Missing keyword DPR TYPE in arc "
02795 "lamp header");
02796
02797 if (fors_qc_keyword_to_paf(header, "ESO TPL ID", NULL,
02798 "Template", instrume))
02799 fors_calib_exit("Missing keyword TPL ID in arc "
02800 "lamp header");
02801 if (fors_qc_keyword_to_paf(header, "ARCFILE", NULL,
02802 "Archive name of input data",
02803 instrume))
02804 fors_calib_exit("Missing keyword ARCFILE in arc "
02805 "lamp header");
02806
02807 pipefile = dfs_generate_filename(curv_coeff_tag);
02808 if (fors_qc_write_string("PIPEFILE", pipefile,
02809 "Pipeline product name", instrume))
02810 fors_calib_exit("Cannot write PIPEFILE to QC log file");
02811 cpl_free(pipefile); pipefile = NULL;
02812
02813
02814
02815
02816
02817 maxpos = fabs(cpl_table_get_column_max(polytraces, "c2"));
02818 maxneg = fabs(cpl_table_get_column_min(polytraces, "c2"));
02819 maxcurve = maxpos > maxneg ? maxpos : maxneg;
02820 if (fors_qc_write_qc_double(header,
02821 maxcurve,
02822 "QC.TRACE.MAX.CURVATURE",
02823 "Y pixel / X pixel ^2",
02824 "Max observed curvature in spectral tracing",
02825 instrume)) {
02826 fors_calib_exit("Cannot write max observed curvature in spectral "
02827 "tracing to QC log file");
02828 }
02829
02830 maxpos = fabs(cpl_table_get_column_max(polytraces, "c1"));
02831 maxneg = fabs(cpl_table_get_column_min(polytraces, "c1"));
02832 maxslope = maxpos > maxneg ? maxpos : maxneg;
02833 if (fors_qc_write_qc_double(header,
02834 maxslope,
02835 "QC.TRACE.MAX.SLOPE",
02836 "Y pixel / X pixel",
02837 "Max observed slope in spectral tracing",
02838 instrume)) {
02839 fors_calib_exit("Cannot write max observed slope in spectral "
02840 "tracing to QC log file");
02841 }
02842
02843 fors_qc_end_group();
02844 }
02845
02846 cpl_free(instrume); instrume = NULL;
02847
02848 if (dfs_save_table(frameset, polytraces, curv_coeff_tag, header,
02849 parlist, recipe, version))
02850 fors_calib_exit(NULL);
02851
02852 cpl_propertylist_delete(header); header = NULL;
02853 cpl_table_delete(polytraces); polytraces = NULL;
02854
02855 if (dfs_save_table(frameset, slits, slit_location_tag, NULL,
02856 parlist, recipe, version))
02857 fors_calib_exit(NULL);
02858
02859 cpl_table_delete(slits); slits = NULL;
02860
02861 if (cpl_error_get_code()) {
02862 cpl_msg_error(cpl_error_get_where(), cpl_error_get_message());
02863 fors_calib_exit(NULL);
02864 }
02865
02866 return 0;
02867 }