FORS Pipeline Reference Manual  4.12.5
fors_science.c
1 /* $Id: fors_science.c,v 1.42 2013-08-14 15:04:57 cgarcia Exp $
2  *
3  * This file is part of the FORS Data Reduction Pipeline
4  * Copyright (C) 2002-2010 European Southern Observatory
5  *
6  * This program is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 2 of the License, or
9  * (at your option) any later version.
10  *
11  * This program is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with this program; if not, write to the Free Software
18  * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
19  */
20 
21 /*
22  * $Author: cgarcia $
23  * $Date: 2013-08-14 15:04:57 $
24  * $Revision: 1.42 $
25  * $Name: not supported by cvs2svn $
26  */
27 
28 #ifdef HAVE_CONFIG_H
29 #include <config.h>
30 #endif
31 
32 #include <math.h>
33 #include <string.h>
34 #include <cpl.h>
35 #include <moses.h>
36 #include <fors_dfs.h>
37 #include <fors_tools.h>
38 #include <fors_qc.h>
39 #include <fors_utils.h>
40 #include "fors_preprocess.h"
41 #include "fors_subtract_bias.h"
42 
43 static int fors_science_create(cpl_plugin *);
44 static int fors_science_exec(cpl_plugin *);
45 static int fors_science_destroy(cpl_plugin *);
46 static int fors_science(cpl_parameterlist *, cpl_frameset *);
47 
48 
49 static char fors_science_description[] =
50 "This recipe is used to reduce scientific spectra using the extraction\n"
51 "mask and the products created by the recipe fors_calib. The spectra are\n"
52 "bias subtracted, flat fielded (if a normalised flat field is specified)\n"
53 "and remapped eliminating the optical distortions. The wavelength calibration\n"
54 "can be optionally upgraded using a number of sky lines: if no sky lines\n"
55 "catalog of wavelengths is specified, an internal one is used instead.\n"
56 "If the alignment to the sky lines is performed, the input dispersion\n"
57 "coefficients table is upgraded and saved to disk, and a new CCD wavelengths\n"
58 "map is created.\n"
59 "This recipe accepts both FORS1 and FORS2 frames. A grism table (typically\n"
60 "depending on the instrument mode, and in particular on the grism used)\n"
61 "may also be specified: this table contains a default recipe parameter\n"
62 "setting to control the way spectra are extracted for a specific instrument\n"
63 "mode, as it is used for automatic run of the pipeline on Paranal and in\n"
64 "Garching. If this table is specified, it will modify the default recipe\n"
65 "parameter setting, with the exception of those parameters which have been\n"
66 "explicitly modifyed on the command line. If a grism table is not specified,\n"
67 "the input recipe parameters values will always be read from the command\n"
68 "line, or from an esorex configuration file if present, or from their\n"
69 "generic default values (that are rarely meaningful).\n"
70 "In the table below the MXU acronym can be read alternatively as MOS\n"
71 "and LSS, depending on the instrument mode of the input data. The acronym\n"
72 "SCI on products should be read STD in case of standard stars observations\n"
73 "A CURV_COEFF table is not (yet) expected for LSS data.\n"
74 "Either a scientific or a standard star exposure can be specified in input.\n"
75 "Only in case of a standard star exposure input, the atmospheric extinction\n"
76 "table and a table with the physical fluxes of the observed standard star\n"
77 "must be specified in input, and a spectro-photometric table is created in\n"
78 "output. This table can then be input again to this recipe, always with an\n"
79 "atmospheric extinction table, and if a photometric calibration is requested\n"
80 "then flux calibrated spectra (in units of erg/cm/cm/s/Angstrom) are also\n"
81 "written in output.\n\n"
82 
83 "Input files:\n\n"
84 " DO category: Type: Explanation: Required:\n"
85 " SCIENCE_MXU Raw Scientific exposure Y\n"
86 " or STANDARD_MXU Raw Standard star exposure Y\n"
87 " MASTER_BIAS Calib Master bias Y\n"
88 " GRISM_TABLE Calib Grism table .\n"
89 " MASTER_SKYLINECAT Calib Sky lines catalog .\n"
90 "\n"
91 " MASTER_NORM_FLAT_MXU Calib Normalised flat field .\n"
92 " DISP_COEFF_MXU Calib Inverse dispersion Y\n"
93 " CURV_COEFF_MXU Calib Spectral curvature Y\n"
94 " SLIT_LOCATION_MXU Calib Slits positions table Y\n"
95 "\n"
96 " or, in case of LSS-like MOS/MXU data,\n"
97 "\n"
98 " MASTER_NORM_FLAT_LONG_MXU Calib Normalised flat field .\n"
99 " DISP_COEFF_LONG_MXU Calib Inverse dispersion Y\n"
100 "\n"
101 " In case STANDARD_MXU is specified in input,\n"
102 "\n"
103 " EXTINCT_TABLE Calib Atmospheric extinction Y\n"
104 " STD_FLUX_TABLE Calib Standard star flux Y\n"
105 "\n"
106 " The following input files are mandatory if photometric calibrated"
107 " spectra are desired:\n"
108 "\n"
109 " EXTINCT_TABLE Calib Atmospheric extinction Y\n"
110 " SPECPHOT_TABLE Calib Response curves Y\n"
111 "\n"
112 " If requested for standard star data, the SPECPHOT_TABLE can be dropped:\n"
113 " in this case the correction is applied using the SPECPHOT_TABLE produced\n"
114 " in the same run.\n"
115 "\n"
116 "Output files:\n"
117 "\n"
118 " DO category: Data type: Explanation:\n"
119 " REDUCED_SCI_MXU FITS image Extracted scientific spectra\n"
120 " REDUCED_SKY_SCI_MXU FITS image Extracted sky spectra\n"
121 " REDUCED_ERROR_SCI_MXU FITS image Errors on extracted spectra\n"
122 " UNMAPPED_SCI_MXU FITS image Sky subtracted scientific spectra\n"
123 " MAPPED_SCI_MXU FITS image Rectified scientific spectra\n"
124 " MAPPED_ALL_SCI_MXU FITS image Rectified science spectra with sky\n"
125 " MAPPED_SKY_SCI_MXU FITS image Rectified sky spectra\n"
126 " UNMAPPED_SKY_SCI_MXU FITS image Sky on CCD\n"
127 " OBJECT_TABLE_SCI_MXU FITS table Positions of detected objects\n"
128 "\n"
129 " Only if the global sky subtraction is requested:\n"
130 " GLOBAL_SKY_SPECTRUM_MXU FITS table Global sky spectrum\n"
131 "\n"
132 " Only if the sky-alignment of the wavelength solution is requested:\n"
133 " SKY_SHIFTS_LONG_SCI_MXU FITS table Sky lines offsets (LSS-like data)\n"
134 " or SKY_SHIFTS_SLIT_SCI_MXU FITS table Sky lines offsets (MOS-like data)\n"
135 " DISP_COEFF_SCI_MXU FITS table Upgraded dispersion coefficients\n"
136 " WAVELENGTH_MAP_SCI_MXU FITS image Upgraded wavelength map\n"
137 "\n"
138 " Only if a STANDARD_MXU is specified in input:\n"
139 " SPECPHOT_TABLE FITS table Efficiency and response curves\n"
140 "\n"
141 " Only if a photometric calibration was requested:\n"
142 " REDUCED_FLUX_SCI_MXU FITS image Flux calibrated scientific spectra\n"
143 " REDUCED_FLUX_ERROR_SCI_MXU FITS image Errors on flux calibrated spectra\n"
144 " MAPPED_FLUX_SCI_MXU FITS image Flux calibrated slit spectra\n\n";
145 
146 #define fors_science_exit(message) \
147 { \
148 if (message) cpl_msg_error(recipe, message); \
149 cpl_free(exptime); \
150 cpl_free(instrume); \
151 cpl_image_delete(dummy); \
152 cpl_image_delete(mapped); \
153 cpl_image_delete(mapped_sky); \
154 cpl_image_delete(mapped_cleaned); \
155 cpl_image_delete(skylocalmap); \
156 cpl_image_delete(skymap); \
157 cpl_image_delete(smapped); \
158 cpl_table_delete(offsets); \
159 cpl_table_delete(photcal); \
160 cpl_table_delete(sky); \
161 fors_image_delete(&bias); \
162 cpl_image_delete(spectra); \
163 cpl_image_delete(coordinate); \
164 fors_image_delete(&norm_flat); \
165 cpl_image_delete(rainbow); \
166 cpl_image_delete(rectified); \
167 cpl_image_delete(wavemap); \
168 cpl_image_delete(wavemaplss); \
169 cpl_propertylist_delete(qclist); \
170 cpl_propertylist_delete(header); \
171 cpl_propertylist_delete(save_header); \
172 cpl_table_delete(grism_table); \
173 cpl_table_delete(idscoeff); \
174 cpl_table_delete(maskslits); \
175 cpl_table_delete(overscans); \
176 cpl_table_delete(polytraces); \
177 cpl_table_delete(slits); \
178 cpl_table_delete(wavelengths); \
179 cpl_vector_delete(lines); \
180 cpl_msg_indent_less(); \
181 return -1; \
182 }
183 
184 
185 #define fors_science_exit_memcheck(message) \
186 { \
187 if (message) cpl_msg_info(recipe, message); \
188 cpl_free(exptime); \
189 cpl_free(instrume); \
190 cpl_image_delete(dummy); \
191 cpl_image_delete(mapped); \
192 cpl_image_delete(mapped_cleaned); \
193 cpl_image_delete(mapped_sky); \
194 cpl_image_delete(skylocalmap); \
195 cpl_image_delete(skymap); \
196 cpl_image_delete(smapped); \
197 cpl_table_delete(offsets); \
198 cpl_table_delete(photcal); \
199 cpl_table_delete(sky); \
200 fors_image_delete(&bias); \
201 cpl_image_delete(spectra); \
202 cpl_image_delete(coordinate); \
203 fors_image_delete(&norm_flat); \
204 cpl_image_delete(rainbow); \
205 cpl_image_delete(rectified); \
206 cpl_image_delete(wavemap); \
207 cpl_image_delete(wavemaplss); \
208 cpl_propertylist_delete(header); \
209 cpl_propertylist_delete(save_header); \
210 cpl_table_delete(grism_table); \
211 cpl_table_delete(idscoeff); \
212 cpl_table_delete(maskslits); \
213 cpl_table_delete(overscans); \
214 cpl_table_delete(polytraces); \
215 cpl_table_delete(slits); \
216 cpl_table_delete(wavelengths); \
217 cpl_vector_delete(lines); \
218 cpl_msg_indent_less(); \
219 return 0; \
220 }
221 
222 
234 int cpl_plugin_get_info(cpl_pluginlist *list)
235 {
236  cpl_recipe *recipe = cpl_calloc(1, sizeof *recipe );
237  cpl_plugin *plugin = &recipe->interface;
238 
239  cpl_plugin_init(plugin,
240  CPL_PLUGIN_API,
241  FORS_BINARY_VERSION,
242  CPL_PLUGIN_TYPE_RECIPE,
243  "fors_science",
244  "Extraction of scientific spectra",
245  fors_science_description,
246  "Carlo Izzo",
247  PACKAGE_BUGREPORT,
249  fors_science_create,
250  fors_science_exec,
251  fors_science_destroy);
252 
253  cpl_pluginlist_append(list, plugin);
254 
255  return 0;
256 }
257 
258 
269 static int fors_science_create(cpl_plugin *plugin)
270 {
271  cpl_recipe *recipe;
272  cpl_parameter *p;
273 
274 
275  /*
276  * Check that the plugin is part of a valid recipe
277  */
278 
279  if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
280  recipe = (cpl_recipe *)plugin;
281  else
282  return -1;
283 
284  /*
285  * Create the parameters list in the cpl_recipe object
286  */
287 
288  recipe->parameters = cpl_parameterlist_new();
289 
290 
291  /*
292  * Dispersion
293  */
294 
295  p = cpl_parameter_new_value("fors.fors_science.dispersion",
296  CPL_TYPE_DOUBLE,
297  "Resampling step (Angstrom/pixel)",
298  "fors.fors_science",
299  0.0);
300  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "dispersion");
301  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
302  cpl_parameterlist_append(recipe->parameters, p);
303 
304  /*
305  * Sky lines alignment
306  */
307 
308  p = cpl_parameter_new_value("fors.fors_science.skyalign",
309  CPL_TYPE_INT,
310  "Polynomial order for sky lines alignment, "
311  "or -1 to avoid alignment",
312  "fors.fors_science",
313  0);
314  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skyalign");
315  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
316  cpl_parameterlist_append(recipe->parameters, p);
317 
318  /*
319  * Line catalog table column containing the sky reference wavelengths
320  */
321 
322  p = cpl_parameter_new_value("fors.fors_science.wcolumn",
323  CPL_TYPE_STRING,
324  "Name of sky line catalog table column "
325  "with wavelengths",
326  "fors.fors_science",
327  "WLEN");
328  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "wcolumn");
329  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
330  cpl_parameterlist_append(recipe->parameters, p);
331 
332  /*
333  * Start wavelength for spectral extraction
334  */
335 
336  p = cpl_parameter_new_value("fors.fors_science.startwavelength",
337  CPL_TYPE_DOUBLE,
338  "Start wavelength in spectral extraction",
339  "fors.fors_science",
340  0.0);
341  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "startwavelength");
342  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
343  cpl_parameterlist_append(recipe->parameters, p);
344 
345  /*
346  * End wavelength for spectral extraction
347  */
348 
349  p = cpl_parameter_new_value("fors.fors_science.endwavelength",
350  CPL_TYPE_DOUBLE,
351  "End wavelength in spectral extraction",
352  "fors.fors_science",
353  0.0);
354  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "endwavelength");
355  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
356  cpl_parameterlist_append(recipe->parameters, p);
357 
358  /*
359  * Apply flat field
360  */
361 
362  p = cpl_parameter_new_value("fors.fors_science.flatfield",
363  CPL_TYPE_BOOL,
364  "Apply flat field",
365  "fors.fors_science",
366  TRUE);
367  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "flatfield");
368  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
369  cpl_parameterlist_append(recipe->parameters, p);
370 
371  /*
372  * Global sky subtraction
373  */
374 
375  p = cpl_parameter_new_value("fors.fors_science.skyglobal",
376  CPL_TYPE_BOOL,
377  "Subtract global sky spectrum from CCD",
378  "fors.fors_science",
379  FALSE);
380  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skyglobal");
381  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
382  cpl_parameterlist_append(recipe->parameters, p);
383 
384  /*
385  * Local sky subtraction on extracted spectra
386  */
387 
388 /*** New sky subtraction (search NSS)
389  p = cpl_parameter_new_value("fors.fors_science.skymedian",
390  CPL_TYPE_INT,
391  "Degree of sky fitting polynomial for "
392  "sky subtraction from extracted "
393  "slit spectra (MOS/MXU only, -1 to disable it)",
394  "fors.fors_science",
395  0);
396  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skymedian");
397  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
398  cpl_parameterlist_append(recipe->parameters, p);
399 ***/
400 
401  p = cpl_parameter_new_value("fors.fors_science.skymedian",
402  CPL_TYPE_BOOL,
403  "Sky subtraction from extracted slit spectra",
404  "fors.fors_science",
405  FALSE);
406  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skymedian");
407  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
408  cpl_parameterlist_append(recipe->parameters, p);
409 
410  /*
411  * Local sky subtraction on CCD spectra
412  */
413 
414  p = cpl_parameter_new_value("fors.fors_science.skylocal",
415  CPL_TYPE_BOOL,
416  "Sky subtraction from CCD slit spectra",
417  "fors.fors_science",
418  TRUE);
419  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "skylocal");
420  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
421  cpl_parameterlist_append(recipe->parameters, p);
422 
423  /*
424  * Cosmic rays removal
425  */
426 
427  p = cpl_parameter_new_value("fors.fors_science.cosmics",
428  CPL_TYPE_BOOL,
429  "Eliminate cosmic rays hits (only if global "
430  "sky subtraction is also requested)",
431  "fors.fors_science",
432  FALSE);
433  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "cosmics");
434  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
435  cpl_parameterlist_append(recipe->parameters, p);
436 
437  /*
438  * Slit margin
439  */
440 
441  p = cpl_parameter_new_value("fors.fors_science.slit_margin",
442  CPL_TYPE_INT,
443  "Number of pixels to exclude at each slit "
444  "in object detection and extraction",
445  "fors.fors_science",
446  3);
447  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "slit_margin");
448  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
449  cpl_parameterlist_append(recipe->parameters, p);
450 
451  /*
452  * Extraction radius
453  */
454 
455  p = cpl_parameter_new_value("fors.fors_science.ext_radius",
456  CPL_TYPE_INT,
457  "Maximum extraction radius for detected "
458  "objects (pixel)",
459  "fors.fors_science",
460  12);
461  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "ext_radius");
462  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
463  cpl_parameterlist_append(recipe->parameters, p);
464 
465  /*
466  * Contamination radius
467  */
468 
469  p = cpl_parameter_new_value("fors.fors_science.cont_radius",
470  CPL_TYPE_INT,
471  "Minimum distance at which two objects "
472  "of equal luminosity do not contaminate "
473  "each other (pixel)",
474  "fors.fors_science",
475  0);
476  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "cont_radius");
477  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
478  cpl_parameterlist_append(recipe->parameters, p);
479 
480  /*
481  * Object extraction method
482  */
483 
484  p = cpl_parameter_new_value("fors.fors_science.ext_mode",
485  CPL_TYPE_INT,
486  "Object extraction method: 0 = aperture, "
487  "1 = Horne optimal extraction",
488  "fors.fors_science",
489  1);
490  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "ext_mode");
491  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
492  cpl_parameterlist_append(recipe->parameters, p);
493 
494  /*
495  * Order of polynomial modeling the instrument response.
496  */
497 
498  p = cpl_parameter_new_value("fors.fors_science.response",
499  CPL_TYPE_INT,
500  "Order of polynomial modeling the "
501  "instrument response",
502  "fors.fors_science",
503  5);
504  cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "response");
505  cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
506  cpl_parameterlist_append(recipe->parameters, p);
507 
508  return 0;
509 }
510 
511 
520 static int fors_science_exec(cpl_plugin *plugin)
521 {
522  cpl_recipe *recipe;
523 
524  if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
525  recipe = (cpl_recipe *)plugin;
526  else
527  return -1;
528 
529  return fors_science(recipe->parameters, recipe->frames);
530 }
531 
532 
541 static int fors_science_destroy(cpl_plugin *plugin)
542 {
543  cpl_recipe *recipe;
544 
545  if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
546  recipe = (cpl_recipe *)plugin;
547  else
548  return -1;
549 
550  cpl_parameterlist_delete(recipe->parameters);
551 
552  return 0;
553 }
554 
555 
565 static int fors_science(cpl_parameterlist *parlist, cpl_frameset *frameset)
566 {
567 
568  const char *recipe = "fors_science";
569 
570 
571  /*
572  * Input parameters
573  */
574 
575  double dispersion;
576  int skyalign;
577  const char *wcolumn;
578  double startwavelength;
579  double endwavelength;
580  int flatfield;
581  int skyglobal;
582  int skylocal;
583  int skymedian;
584  int cosmics;
585  int slit_margin;
586  int ext_radius;
587  int cont_radius;
588  int ext_mode;
589  int res_order;
590  int photometry;
591 
592 
593  /*
594  * CPL objects
595  */
596 
597  fors_image_list *all_science;
598  cpl_image **images;
599 
600  fors_image *bias = NULL;
601  fors_image *norm_flat = NULL;
602  fors_image *science_ima = NULL;
603  cpl_image *spectra = NULL;
604  cpl_image *rectified = NULL;
605  cpl_image *coordinate = NULL;
606  cpl_image *rainbow = NULL;
607  cpl_image *mapped = NULL;
608  cpl_image *mapped_sky = NULL;
609  cpl_image *mapped_cleaned = NULL;
610  cpl_image *smapped = NULL;
611  cpl_image *wavemap = NULL;
612  cpl_image *wavemaplss = NULL;
613  cpl_image *skymap = NULL;
614  cpl_image *skylocalmap = NULL;
615  cpl_image *dummy = NULL;
616 
617  cpl_table *grism_table = NULL;
618  cpl_table *overscans = NULL;
619  cpl_table *wavelengths = NULL;
620  cpl_table *idscoeff = NULL;
621  cpl_table *slits = NULL;
622  cpl_table *maskslits = NULL;
623  cpl_table *polytraces = NULL;
624  cpl_table *offsets = NULL;
625  cpl_table *sky = NULL;
626  cpl_table *photcal = NULL;
627 
628  cpl_vector *lines = NULL;
629 
630  cpl_propertylist *header = NULL;
631  cpl_propertylist *save_header = NULL;
632  cpl_propertylist *qclist = NULL;
633 
634 
635  /*
636  * Auxiliary variables
637  */
638 
639  char version[80];
640  char *instrume = NULL;
641  char *wheel4;
642  const char *science_tag = NULL;
643  const char *master_norm_flat_tag = NULL;
644  const char *disp_coeff_tag = NULL;
645  const char *disp_coeff_sky_tag = NULL;
646  const char *wavelength_map_sky_tag = NULL;
647  const char *curv_coeff_tag = NULL;
648  const char *slit_location_tag = NULL;
649  const char *reduced_science_tag = NULL;
650  const char *reduced_flux_science_tag = NULL;
651  const char *reduced_sky_tag = NULL;
652  const char *reduced_error_tag = NULL;
653  const char *reduced_flux_error_tag = NULL;
654  const char *mapped_science_tag = NULL;
655  const char *mapped_flux_science_tag = NULL;
656  const char *unmapped_science_tag = NULL;
657  const char *mapped_science_sky_tag = NULL;
658  const char *mapped_sky_tag = NULL;
659  const char *unmapped_sky_tag = NULL;
660  const char *global_sky_spectrum_tag = NULL;
661  const char *object_table_tag = NULL;
662  const char *skylines_offsets_tag = NULL;
663  const char *specphot_tag;
664  const char *master_specphot_tag = "MASTER_SPECPHOT_TABLE";
665  int mxu, mos, lss;
666  int treat_as_lss = 0;
667  int have_phot = 0;
668  int nscience;
669  double *exptime = NULL;
670  double alltime;
671  double airmass = -1;
672  double mxpos;
673  double mean_rms;
674  int nlines;
675  int rebin;
676  double *line;
677  int nx, ny;
678  int ccd_xsize, ccd_ysize;
679  double reference;
680  double gain;
681  double ron;
682  int standard;
683  int highres;
684  int i;
685  double wstart;
686  double wstep;
687  int wcount;
688 
689 
690  snprintf(version, 80, "%s-%s", PACKAGE, PACKAGE_VERSION);
691 
692  cpl_msg_set_indentation(2);
693 
694 
695  /*
696  * Get configuration parameters
697  */
698 
699  cpl_msg_info(recipe, "Recipe %s configuration parameters:", recipe);
700  cpl_msg_indent_more();
701 
702  if (cpl_frameset_count_tags(frameset, "GRISM_TABLE") > 1)
703  fors_science_exit("Too many in input: GRISM_TABLE");
704 
705  grism_table = dfs_load_table(frameset, "GRISM_TABLE", 1);
706 
707  dispersion = dfs_get_parameter_double(parlist,
708  "fors.fors_science.dispersion", grism_table);
709 
710  if (dispersion <= 0.0)
711  fors_science_exit("Invalid resampling step");
712 
713  skyalign = dfs_get_parameter_int(parlist,
714  "fors.fors_science.skyalign", NULL);
715 
716  if (skyalign > 2)
717  fors_science_exit("Max polynomial degree for sky alignment is 2");
718 
719  wcolumn = dfs_get_parameter_string(parlist,
720  "fors.fors_science.wcolumn", NULL);
721 
722  startwavelength = dfs_get_parameter_double(parlist,
723  "fors.fors_science.startwavelength", grism_table);
724  if (startwavelength < 3000.0 || startwavelength > 13000.0)
725  fors_science_exit("Invalid wavelength");
726 
727  endwavelength = dfs_get_parameter_double(parlist,
728  "fors.fors_science.endwavelength", grism_table);
729  if (endwavelength < 3000.0 || endwavelength > 13000.0)
730  fors_science_exit("Invalid wavelength");
731 
732  if (endwavelength - startwavelength <= 0.0)
733  fors_science_exit("Invalid wavelength interval");
734 
735  flatfield = dfs_get_parameter_bool(parlist, "fors.fors_science.flatfield",
736  NULL);
737 
738  skyglobal = dfs_get_parameter_bool(parlist, "fors.fors_science.skyglobal",
739  NULL);
740  skylocal = dfs_get_parameter_bool(parlist, "fors.fors_science.skylocal",
741  NULL);
742  skymedian = dfs_get_parameter_bool(parlist, "fors.fors_science.skymedian",
743  NULL);
744 /* NSS
745  skymedian = dfs_get_parameter_int(parlist, "fors.fors_science.skymedian",
746  NULL);
747 */
748 
749  if (skylocal && skyglobal)
750  fors_science_exit("Cannot apply both local and global sky subtraction");
751 
752  if (skylocal && skymedian)
753  fors_science_exit("Cannot apply sky subtraction both on extracted "
754  "and non-extracted spectra");
755 
756  cosmics = dfs_get_parameter_bool(parlist,
757  "fors.fors_science.cosmics", NULL);
758 
759  if (cosmics)
760  if (!(skyglobal || skylocal))
761  fors_science_exit("Cosmic rays correction requires "
762  "either skylocal=true or skyglobal=true");
763 
764  slit_margin = dfs_get_parameter_int(parlist,
765  "fors.fors_science.slit_margin",
766  NULL);
767  if (slit_margin < 0)
768  fors_science_exit("Value must be zero or positive");
769 
770  ext_radius = dfs_get_parameter_int(parlist, "fors.fors_science.ext_radius",
771  NULL);
772  if (ext_radius < 0)
773  fors_science_exit("Value must be zero or positive");
774 
775  cont_radius = dfs_get_parameter_int(parlist,
776  "fors.fors_science.cont_radius",
777  NULL);
778  if (cont_radius < 0)
779  fors_science_exit("Value must be zero or positive");
780 
781  ext_mode = dfs_get_parameter_int(parlist, "fors.fors_science.ext_mode",
782  NULL);
783  if (ext_mode < 0 || ext_mode > 1)
784  fors_science_exit("Invalid object extraction mode");
785 
786  res_order = dfs_get_parameter_int(parlist, "fors.fors_science.response",
787  NULL);
788  if (res_order < 2 || res_order > 10)
789  fors_science_exit("Invalid instrument response modeling polynomial");
790 
791  cpl_table_delete(grism_table); grism_table = NULL;
792 
793  if (cpl_error_get_code())
794  fors_science_exit("Failure getting the configuration parameters");
795 
796 
797  /*
798  * Check input set-of-frames
799  */
800 
801  cpl_msg_indent_less();
802  cpl_msg_info(recipe, "Check input set-of-frames:");
803  cpl_msg_indent_more();
804 
805  mxu = cpl_frameset_count_tags(frameset, "SCIENCE_MXU");
806  mos = cpl_frameset_count_tags(frameset, "SCIENCE_MOS");
807  lss = cpl_frameset_count_tags(frameset, "SCIENCE_LSS");
808  standard = 0;
809 
810  if (mxu + mos + lss == 0) {
811  mxu = cpl_frameset_count_tags(frameset, "STANDARD_MXU");
812  mos = cpl_frameset_count_tags(frameset, "STANDARD_MOS");
813  lss = cpl_frameset_count_tags(frameset, "STANDARD_LSS");
814  standard = 1;
815  }
816 
817  if (mxu + mos + lss == 0)
818  fors_science_exit("Missing input scientific frame");
819 
820  nscience = mxu + mos + lss;
821 
822  if (mxu && mxu < nscience)
823  fors_science_exit("Input scientific frames must be of the same type");
824 
825  if (mos && mos < nscience)
826  fors_science_exit("Input scientific frames must be of the same type");
827 
828  if (lss && lss < nscience)
829  fors_science_exit("Input scientific frames must be of the same type");
830 
831  if (mxu) {
832  cpl_msg_info(recipe, "MXU data found");
833  if (standard) {
834  science_tag = "STANDARD_MXU";
835  reduced_science_tag = "REDUCED_STD_MXU";
836  reduced_flux_science_tag = "REDUCED_FLUX_STD_MXU";
837  unmapped_science_tag = "UNMAPPED_STD_MXU";
838  mapped_science_tag = "MAPPED_STD_MXU";
839  mapped_flux_science_tag = "MAPPED_FLUX_STD_MXU";
840  mapped_science_sky_tag = "MAPPED_ALL_STD_MXU";
841  skylines_offsets_tag = "SKY_SHIFTS_SLIT_STD_MXU";
842  wavelength_map_sky_tag = "WAVELENGTH_MAP_STD_MXU";
843  disp_coeff_sky_tag = "DISP_COEFF_STD_MXU";
844  mapped_sky_tag = "MAPPED_SKY_STD_MXU";
845  unmapped_sky_tag = "UNMAPPED_SKY_STD_MXU";
846  object_table_tag = "OBJECT_TABLE_STD_MXU";
847  reduced_sky_tag = "REDUCED_SKY_STD_MXU";
848  reduced_error_tag = "REDUCED_ERROR_STD_MXU";
849  reduced_flux_error_tag = "REDUCED_FLUX_ERROR_STD_MXU";
850  specphot_tag = "SPECPHOT_TABLE";
851  }
852  else {
853  science_tag = "SCIENCE_MXU";
854  reduced_science_tag = "REDUCED_SCI_MXU";
855  reduced_flux_science_tag = "REDUCED_FLUX_SCI_MXU";
856  unmapped_science_tag = "UNMAPPED_SCI_MXU";
857  mapped_science_tag = "MAPPED_SCI_MXU";
858  mapped_flux_science_tag = "MAPPED_FLUX_SCI_MXU";
859  mapped_science_sky_tag = "MAPPED_ALL_SCI_MXU";
860  skylines_offsets_tag = "SKY_SHIFTS_SLIT_SCI_MXU";
861  wavelength_map_sky_tag = "WAVELENGTH_MAP_SCI_MXU";
862  disp_coeff_sky_tag = "DISP_COEFF_SCI_MXU";
863  mapped_sky_tag = "MAPPED_SKY_SCI_MXU";
864  unmapped_sky_tag = "UNMAPPED_SKY_SCI_MXU";
865  object_table_tag = "OBJECT_TABLE_SCI_MXU";
866  reduced_sky_tag = "REDUCED_SKY_SCI_MXU";
867  reduced_error_tag = "REDUCED_ERROR_SCI_MXU";
868  reduced_flux_error_tag = "REDUCED_FLUX_ERROR_SCI_MXU";
869  specphot_tag = "SPECPHOT_TABLE";
870  }
871 
872  master_norm_flat_tag = "MASTER_NORM_FLAT_MXU";
873  disp_coeff_tag = "DISP_COEFF_MXU";
874  curv_coeff_tag = "CURV_COEFF_MXU";
875  slit_location_tag = "SLIT_LOCATION_MXU";
876  global_sky_spectrum_tag = "GLOBAL_SKY_SPECTRUM_MXU";
877 
878  if (!cpl_frameset_count_tags(frameset, master_norm_flat_tag)) {
879  master_norm_flat_tag = "MASTER_NORM_FLAT_LONG_MXU";
880  disp_coeff_tag = "DISP_COEFF_LONG_MXU";
881  slit_location_tag = "SLIT_LOCATION_LONG_MXU";
882  }
883  }
884 
885  if (lss) {
886  cpl_msg_info(recipe, "LSS data found");
887 
888  if (cosmics && !skyglobal)
889  fors_science_exit("Cosmic rays correction for LSS "
890  "data requires --skyglobal=true");
891 
892  if (standard) {
893  science_tag = "STANDARD_LSS";
894  reduced_science_tag = "REDUCED_STD_LSS";
895  reduced_flux_science_tag = "REDUCED_FLUX_STD_LSS";
896  unmapped_science_tag = "UNMAPPED_STD_LSS";
897  mapped_science_tag = "MAPPED_STD_LSS";
898  mapped_flux_science_tag = "MAPPED_FLUX_STD_LSS";
899  mapped_science_sky_tag = "MAPPED_ALL_STD_LSS";
900  skylines_offsets_tag = "SKY_SHIFTS_LONG_STD_LSS";
901  wavelength_map_sky_tag = "WAVELENGTH_MAP_STD_LSS";
902  disp_coeff_sky_tag = "DISP_COEFF_STD_LSS";
903  mapped_sky_tag = "MAPPED_SKY_STD_LSS";
904  unmapped_sky_tag = "UNMAPPED_SKY_STD_LSS";
905  object_table_tag = "OBJECT_TABLE_STD_LSS";
906  reduced_sky_tag = "REDUCED_SKY_STD_LSS";
907  reduced_error_tag = "REDUCED_ERROR_STD_LSS";
908  reduced_flux_error_tag = "REDUCED_FLUX_ERROR_STD_LSS";
909  specphot_tag = "SPECPHOT_TABLE";
910  }
911  else {
912  science_tag = "SCIENCE_LSS";
913  reduced_science_tag = "REDUCED_SCI_LSS";
914  reduced_flux_science_tag = "REDUCED_FLUX_SCI_LSS";
915  unmapped_science_tag = "UNMAPPED_SCI_LSS";
916  mapped_science_tag = "MAPPED_SCI_LSS";
917  mapped_flux_science_tag = "MAPPED_FLUX_SCI_LSS";
918  mapped_science_sky_tag = "MAPPED_ALL_SCI_LSS";
919  skylines_offsets_tag = "SKY_SHIFTS_LONG_SCI_LSS";
920  wavelength_map_sky_tag = "WAVELENGTH_MAP_SCI_LSS";
921  disp_coeff_sky_tag = "DISP_COEFF_SCI_LSS";
922  mapped_sky_tag = "MAPPED_SKY_SCI_LSS";
923  unmapped_sky_tag = "UNMAPPED_SKY_SCI_LSS";
924  object_table_tag = "OBJECT_TABLE_SCI_LSS";
925  reduced_sky_tag = "REDUCED_SKY_SCI_LSS";
926  reduced_error_tag = "REDUCED_ERROR_SCI_LSS";
927  reduced_flux_error_tag = "REDUCED_FLUX_ERROR_SCI_LSS";
928  specphot_tag = "SPECPHOT_TABLE";
929  }
930 
931  master_norm_flat_tag = "MASTER_NORM_FLAT_LSS";
932  disp_coeff_tag = "DISP_COEFF_LSS";
933  slit_location_tag = "SLIT_LOCATION_LSS";
934  global_sky_spectrum_tag = "GLOBAL_SKY_SPECTRUM_LSS";
935  }
936 
937  if (mos) {
938  cpl_msg_info(recipe, "MOS data found");
939  if (standard) {
940  science_tag = "STANDARD_MOS";
941  reduced_science_tag = "REDUCED_STD_MOS";
942  reduced_flux_science_tag = "REDUCED_FLUX_STD_MOS";
943  unmapped_science_tag = "UNMAPPED_STD_MOS";
944  mapped_science_tag = "MAPPED_STD_MOS";
945  mapped_flux_science_tag = "MAPPED_FLUX_STD_MOS";
946  mapped_science_sky_tag = "MAPPED_ALL_STD_MOS";
947  skylines_offsets_tag = "SKY_SHIFTS_SLIT_STD_MOS";
948  wavelength_map_sky_tag = "WAVELENGTH_MAP_STD_MOS";
949  disp_coeff_sky_tag = "DISP_COEFF_STD_MOS";
950  mapped_sky_tag = "MAPPED_SKY_STD_MOS";
951  unmapped_sky_tag = "UNMAPPED_SKY_STD_MOS";
952  object_table_tag = "OBJECT_TABLE_STD_MOS";
953  reduced_sky_tag = "REDUCED_SKY_STD_MOS";
954  reduced_error_tag = "REDUCED_ERROR_STD_MOS";
955  reduced_flux_error_tag = "REDUCED_FLUX_ERROR_STD_MOS";
956  specphot_tag = "SPECPHOT_TABLE";
957  }
958  else {
959  science_tag = "SCIENCE_MOS";
960  reduced_science_tag = "REDUCED_SCI_MOS";
961  reduced_flux_science_tag = "REDUCED_FLUX_SCI_MOS";
962  unmapped_science_tag = "UNMAPPED_SCI_MOS";
963  mapped_science_tag = "MAPPED_SCI_MOS";
964  mapped_flux_science_tag = "MAPPED_FLUX_SCI_MOS";
965  mapped_science_sky_tag = "MAPPED_ALL_SCI_MOS";
966  skylines_offsets_tag = "SKY_SHIFTS_SLIT_SCI_MOS";
967  wavelength_map_sky_tag = "WAVELENGTH_MAP_SCI_MOS";
968  disp_coeff_sky_tag = "DISP_COEFF_SCI_MOS";
969  mapped_sky_tag = "MAPPED_SKY_SCI_MOS";
970  unmapped_sky_tag = "UNMAPPED_SKY_SCI_MOS";
971  object_table_tag = "OBJECT_TABLE_SCI_MOS";
972  reduced_sky_tag = "REDUCED_SKY_SCI_MOS";
973  reduced_error_tag = "REDUCED_ERROR_SCI_MOS";
974  reduced_flux_error_tag = "REDUCED_FLUX_ERROR_SCI_MOS";
975  specphot_tag = "SPECPHOT_TABLE";
976  }
977 
978  master_norm_flat_tag = "MASTER_NORM_FLAT_MOS";
979  disp_coeff_tag = "DISP_COEFF_MOS";
980  curv_coeff_tag = "CURV_COEFF_MOS";
981  slit_location_tag = "SLIT_LOCATION_MOS";
982  global_sky_spectrum_tag = "GLOBAL_SKY_SPECTRUM_MOS";
983 
984  if (!cpl_frameset_count_tags(frameset, master_norm_flat_tag)) {
985  master_norm_flat_tag = "MASTER_NORM_FLAT_LONG_MOS";
986  disp_coeff_tag = "DISP_COEFF_LONG_MOS";
987  slit_location_tag = "SLIT_LOCATION_LONG_MOS";
988  }
989  }
990 
991  if (cpl_frameset_count_tags(frameset, "MASTER_BIAS") == 0)
992  fors_science_exit("Missing required input: MASTER_BIAS");
993 
994  if (cpl_frameset_count_tags(frameset, "MASTER_BIAS") > 1)
995  fors_science_exit("Too many in input: MASTER_BIAS");
996 
997  if (skyalign >= 0)
998  if (cpl_frameset_count_tags(frameset, "MASTER_SKYLINECAT") > 1)
999  fors_science_exit("Too many in input: MASTER_SKYLINECAT");
1000 
1001  if (cpl_frameset_count_tags(frameset, disp_coeff_tag) == 0) {
1002  cpl_msg_error(recipe, "Missing required input: %s", disp_coeff_tag);
1003  fors_science_exit(NULL);
1004  }
1005 
1006  if (cpl_frameset_count_tags(frameset, disp_coeff_tag) > 1) {
1007  cpl_msg_error(recipe, "Too many in input: %s", disp_coeff_tag);
1008  fors_science_exit(NULL);
1009  }
1010 
1011  if (cpl_frameset_count_tags(frameset, slit_location_tag) == 0) {
1012  cpl_msg_error(recipe, "Missing required input: %s",
1013  slit_location_tag);
1014  fors_science_exit(NULL);
1015  }
1016 
1017  if (cpl_frameset_count_tags(frameset, slit_location_tag) > 1) {
1018  cpl_msg_error(recipe, "Too many in input: %s", slit_location_tag);
1019  fors_science_exit(NULL);
1020  }
1021 
1022  if (cpl_frameset_count_tags(frameset, master_norm_flat_tag) > 1) {
1023  if (flatfield) {
1024  cpl_msg_error(recipe, "Too many in input: %s",
1025  master_norm_flat_tag);
1026  fors_science_exit(NULL);
1027  }
1028  else {
1029  cpl_msg_warning(recipe, "%s in input are ignored, "
1030  "since flat field correction was not requested",
1031  master_norm_flat_tag);
1032  }
1033  }
1034 
1035  if (cpl_frameset_count_tags(frameset, master_norm_flat_tag) == 1) {
1036  if (!flatfield) {
1037  cpl_msg_warning(recipe, "%s in input is ignored, "
1038  "since flat field correction was not requested",
1039  master_norm_flat_tag);
1040  }
1041  }
1042 
1043  if (cpl_frameset_count_tags(frameset, master_norm_flat_tag) == 0) {
1044  if (flatfield) {
1045  cpl_msg_error(recipe, "Flat field correction was requested, "
1046  "but no %s are found in input",
1047  master_norm_flat_tag);
1048  fors_science_exit(NULL);
1049  }
1050  }
1051 
1052  if (standard) {
1053 
1054  if (cpl_frameset_count_tags(frameset, "EXTINCT_TABLE") == 0) {
1055  cpl_msg_warning(recipe, "An EXTINCT_TABLE was not found in input: "
1056  "instrument response curve will not be produced.");
1057  standard = 0;
1058  }
1059 
1060  if (cpl_frameset_count_tags(frameset, "EXTINCT_TABLE") > 1)
1061  fors_science_exit("Too many in input: EXTINCT_TABLE");
1062 
1063  if (cpl_frameset_count_tags(frameset, "STD_FLUX_TABLE") == 0) {
1064  cpl_msg_warning(recipe, "A STD_FLUX_TABLE was not found in input: "
1065  "instrument response curve will not be produced.");
1066  standard = 0;
1067  }
1068 
1069  if (cpl_frameset_count_tags(frameset, "STD_FLUX_TABLE") > 1)
1070  fors_science_exit("Too many in input: STD_FLUX_TABLE");
1071 
1072  if (!dfs_equal_keyword(frameset, "ESO OBS TARG NAME")) {
1073  cpl_msg_warning(recipe, "The target name of observation does not "
1074  "match the standard star catalog: "
1075  "instrument response curve will not be produced.");
1076  standard = 0;
1077  }
1078  }
1079 
1080  have_phot = cpl_frameset_count_tags(frameset, specphot_tag);
1081  have_phot += cpl_frameset_count_tags(frameset, master_specphot_tag);
1082 
1083  if (!standard) {
1084  if (have_phot == 0) {
1085  cpl_msg_info(recipe,
1086  "A SPECPHOT_TABLE was not found in input: "
1087  "no photometric calibrated "
1088  "spectra will be produced.");
1089  photometry = 0;
1090  }
1091  else
1092  {
1093  cpl_msg_info(recipe,
1094  "A SPECPHOT_TABLE was found in input: "
1095  "photometric calibrated "
1096  "spectra will be produced.");
1097  photometry = 1;
1098  }
1099 
1100  if (photometry) {
1101  if (cpl_frameset_count_tags(frameset, "EXTINCT_TABLE") != 1)
1102  fors_science_exit("One and only one EXTINCT_TABLE is needed "
1103  "to calibrate in photometry");
1104 
1105 
1106  if (have_phot > 1)
1107  fors_science_exit("Too many in input: SPECPHOT_TABLE");
1108  }
1109  }
1110 
1111  if (!dfs_equal_keyword(frameset, "ESO INS GRIS1 ID"))
1112  fors_science_exit("Input frames are not from the same grism");
1113 
1114  if (!dfs_equal_keyword(frameset, "ESO INS FILT1 ID"))
1115  fors_science_exit("Input frames are not from the same filter");
1116 
1117  if (cpl_frameset_count_tags(frameset, "SPECPHOT_TABLE"))
1118  {
1119  cpl_frameset * frameset_nostd = cpl_frameset_duplicate(frameset);
1120  cpl_frameset_erase(frameset_nostd, "SPECPHOT_TABLE");
1121  if (!dfs_equal_keyword(frameset_nostd, "ESO DET CHIP1 ID"))
1122  fors_science_exit("Input frames are not from the same chip."
1123  " This does not apply to specphot table.");
1124  cpl_frameset_delete(frameset_nostd);
1125  }
1126  else
1127  {
1128  if (!dfs_equal_keyword(frameset, "ESO DET CHIP1 ID"))
1129  fors_science_exit("Input frames are not from the same chip");
1130  }
1131 
1132  {
1133  cpl_frameset * frameset_detector = cpl_frameset_duplicate(frameset);
1134  cpl_frameset_erase(frameset_detector, "GRISM_TABLE");
1135  if (!dfs_equal_keyword(frameset_detector, "ESO DET CHIP1 NAME"))
1136  fors_science_exit("Input frames are not from the same chip mosaic");
1137  cpl_frameset_delete(frameset_detector);
1138  }
1139 
1140  cpl_msg_indent_less();
1141 
1142 
1143  /*
1144  * Loading input data
1145  */
1146 
1147  /* Classify frames */
1148  fors_dfs_set_groups(frameset);
1149 
1150  exptime = cpl_calloc(nscience, sizeof(double));
1151 
1152  if (nscience > 1) {
1153 
1154  cpl_msg_info(recipe, "Load %d scientific frames and median them...",
1155  nscience);
1156  cpl_msg_indent_more();
1157 
1158 
1159  header = dfs_load_header(frameset, science_tag, 0);
1160 
1161  if (header == NULL)
1162  fors_science_exit("Cannot load scientific frame header");
1163 
1164  alltime = exptime[0] = cpl_propertylist_get_double(header, "EXPTIME");
1165 
1166  if (cpl_error_get_code() != CPL_ERROR_NONE)
1167  fors_science_exit("Missing keyword EXPTIME in scientific "
1168  "frame header");
1169 
1170  if (standard || photometry) {
1171  airmass = fors_get_airmass(header);
1172  if (airmass < 0.0)
1173  fors_science_exit("Missing airmass information in "
1174  "scientific frame header");
1175  }
1176 
1177  cpl_propertylist_delete(header); header = NULL;
1178 
1179  cpl_msg_info(recipe, "Scientific frame 1 exposure time: %.2f s",
1180  exptime[0]);
1181 
1182  for (i = 1; i < nscience; i++) {
1183 
1184  header = dfs_load_header(frameset, NULL, 0);
1185 
1186  if (header == NULL)
1187  fors_science_exit("Cannot load scientific frame header");
1188 
1189  exptime[i] = cpl_propertylist_get_double(header, "EXPTIME");
1190 
1191  alltime += exptime[i];
1192 
1193  if (cpl_error_get_code() != CPL_ERROR_NONE)
1194  fors_science_exit("Missing keyword EXPTIME in scientific "
1195  "frame header");
1196 
1197  cpl_propertylist_delete(header); header = NULL;
1198 
1199  cpl_msg_info(recipe, "Scientific frame %d exposure time: %.2f s",
1200  i+1, exptime[i]);
1201  }
1202 
1203  all_science = fors_image_list_new();
1204  cpl_frameset * science_frames = fors_frameset_extract(frameset, science_tag);
1205 
1206  for (i = 0; i < nscience; i++)
1207  {
1208  cpl_frame * this_science_frame =
1209  cpl_frameset_get_position(science_frames, i);
1210  fors_image * this_science_ima =
1211  fors_image_load_preprocess(this_science_frame);
1212 
1213  if (this_science_ima) {
1214  fors_image_divide_scalar(this_science_ima, exptime[i], 0);
1215  fors_image_list_insert(all_science, this_science_ima);
1216  }
1217  else
1218  fors_science_exit("Cannot load scientific frame");
1219  }
1220 
1221  hdrl_image *combined;
1222  cpl_image *contrib;
1223  hdrl_imagelist * all_science_hdrl = fors_image_list_to_hdrl(all_science);
1224  hdrl_parameter * combine_par = hdrl_collapse_median_parameter_create();
1225  hdrl_imagelist_collapse(all_science_hdrl, combine_par,
1226  &combined, &contrib);
1227 
1228  science_ima = fors_image_from_hdrl(combined);
1229 
1230  fors_image_multiply_scalar(science_ima, alltime, 0);
1231  spectra = science_ima->data;
1232 
1233  hdrl_imagelist_delete(all_science_hdrl);
1234  fors_image_list_delete(&all_science, fors_image_delete);
1235  hdrl_image_delete(combined);
1236  hdrl_parameter_delete(combine_par);
1237  cpl_image_delete(contrib);
1238  cpl_frameset_delete(science_frames);
1239  }
1240  else {
1241  cpl_msg_info(recipe, "Load scientific exposure...");
1242  cpl_msg_indent_more();
1243 
1244  header = dfs_load_header(frameset, science_tag, 0);
1245 
1246  if (header == NULL)
1247  fors_science_exit("Cannot load scientific frame header");
1248 
1249  if (standard || photometry) {
1250  airmass = fors_get_airmass(header);
1251  if (airmass < 0.0)
1252  fors_science_exit("Missing airmass information in "
1253  "scientific frame header");
1254  }
1255 
1256  /*
1257  * Insert here a check on supported filters:
1258  */
1259 
1260  wheel4 = (char *)cpl_propertylist_get_string(header,
1261  "ESO INS OPTI9 TYPE");
1262  if (cpl_error_get_code() != CPL_ERROR_NONE) {
1263  fors_science_exit("Missing ESO INS OPTI9 TYPE in flat header");
1264  }
1265 
1266  if (strcmp("FILT", wheel4) == 0) {
1267  wheel4 = (char *)cpl_propertylist_get_string(header,
1268  "ESO INS OPTI9 NAME");
1269  cpl_msg_error(recipe, "Unsupported filter: %s", wheel4);
1270  fors_science_exit(NULL);
1271  }
1272 
1273  alltime = exptime[0] = cpl_propertylist_get_double(header, "EXPTIME");
1274 
1275  if (cpl_error_get_code() != CPL_ERROR_NONE)
1276  fors_science_exit("Missing keyword EXPTIME in scientific "
1277  "frame header");
1278 
1279  cpl_propertylist_delete(header); header = NULL;
1280 
1281  cpl_msg_info(recipe, "Scientific frame exposure time: %.2f s",
1282  exptime[0]);
1283 
1284  const cpl_frame * science_frame =
1285  cpl_frameset_find_const(frameset, science_tag);
1286  science_ima = fors_image_load_preprocess(science_frame);
1287  spectra = science_ima->data;
1288  }
1289 
1290  if (spectra == NULL)
1291  fors_science_exit("Cannot load scientific frame");
1292 
1293  cpl_free(exptime); exptime = NULL;
1294 
1295  cpl_msg_indent_less();
1296 
1297 
1298  /*
1299  * Get the reference wavelength and the rebin factor along the
1300  * dispersion direction from a scientific exposure
1301  */
1302 
1303  header = dfs_load_header(frameset, science_tag, 0);
1304 
1305  if (header == NULL)
1306  fors_science_exit("Cannot load scientific frame header");
1307 
1308  instrume = (char *)cpl_propertylist_get_string(header, "INSTRUME");
1309  if (instrume == NULL)
1310  fors_science_exit("Missing keyword INSTRUME in scientific header");
1311  instrume = cpl_strdup(instrume);
1312 
1313  if (instrume[4] == '1')
1314  snprintf(version, 80, "%s/%s", "fors1", VERSION);
1315  if (instrume[4] == '2')
1316  snprintf(version, 80, "%s/%s", "fors2", VERSION);
1317 
1318  reference = cpl_propertylist_get_double(header, "ESO INS GRIS1 WLEN");
1319 
1320  if (cpl_error_get_code() != CPL_ERROR_NONE)
1321  fors_science_exit("Missing keyword ESO INS GRIS1 WLEN in scientific "
1322  "frame header");
1323 
1324  if (reference < 3000.0) /* Perhaps in nanometers... */
1325  reference *= 10;
1326 
1327  if (reference < 3000.0 || reference > 13000.0) {
1328  cpl_msg_error(recipe, "Invalid central wavelength %.2f read from "
1329  "keyword ESO INS GRIS1 WLEN in scientific frame header",
1330  reference);
1331  fors_science_exit(NULL);
1332  }
1333 
1334  cpl_msg_info(recipe, "The central wavelength is: %.2f", reference);
1335 
1336  rebin = cpl_propertylist_get_int(header, "ESO DET WIN1 BINX");
1337 
1338  if (cpl_error_get_code() != CPL_ERROR_NONE)
1339  fors_science_exit("Missing keyword ESO DET WIN1 BINX in scientific "
1340  "frame header");
1341 
1342  if (rebin != 1) {
1343  dispersion *= rebin;
1344  cpl_msg_warning(recipe, "The rebin factor is %d, and therefore the "
1345  "resampling step used is %f A/pixel", rebin,
1346  dispersion);
1347  ext_radius /= rebin;
1348  cpl_msg_warning(recipe, "The rebin factor is %d, and therefore the "
1349  "extraction radius used is %d pixel", rebin,
1350  ext_radius);
1351  }
1352 
1353  gain = cpl_propertylist_get_double(header, "ESO DET OUT1 CONAD");
1354 
1355  if (cpl_error_get_code() != CPL_ERROR_NONE)
1356  fors_science_exit("Missing keyword ESO DET OUT1 CONAD in scientific "
1357  "frame header");
1358 
1359  cpl_msg_info(recipe, "The gain factor is: %.2f e-/ADU", gain);
1360 
1361  ron = cpl_propertylist_get_double(header, "ESO DET OUT1 RON");
1362 
1363  if (cpl_error_get_code() != CPL_ERROR_NONE)
1364  fors_science_exit("Missing keyword ESO DET OUT1 RON in scientific "
1365  "frame header");
1366 
1367  ron /= gain; /* Convert from electrons to ADU */
1368 
1369  cpl_msg_info(recipe, "The read-out-noise is: %.2f ADU", ron);
1370 
1371  if (mos || mxu) {
1372  int nslits_out_det = 0;
1373 /* goto skip; */
1374  if (mos)
1375  maskslits = mos_load_slits_fors_mos(header, &nslits_out_det);
1376  else
1377  maskslits = mos_load_slits_fors_mxu(header);
1378 
1379  /*
1380  * Check if all slits have the same X offset: in such case,
1381  * treat the observation as a long-slit one!
1382  */
1383 
1384  mxpos = cpl_table_get_column_median(maskslits, "xtop");
1385 
1386  treat_as_lss = fors_mos_is_lss_like(maskslits, nslits_out_det);
1387 
1388  cpl_table_delete(maskslits); maskslits = NULL;
1389 /* skip: */
1390 
1391  if (treat_as_lss) {
1392  cpl_msg_warning(recipe, "All MOS slits have the same offset: %.2f\n"
1393  "The LSS data reduction strategy is applied!",
1394  mxpos);
1395  if (mos) {
1396  if (standard) {
1397  skylines_offsets_tag = "SKY_SHIFTS_LONG_STD_MOS";
1398  }
1399  else {
1400  skylines_offsets_tag = "SKY_SHIFTS_LONG_SCI_MOS";
1401  }
1402  }
1403  else {
1404  if (standard) {
1405  skylines_offsets_tag = "SKY_SHIFTS_LONG_STD_MXU";
1406  }
1407  else {
1408  skylines_offsets_tag = "SKY_SHIFTS_LONG_SCI_MXU";
1409  }
1410  }
1411  }
1412  }
1413 
1414  if (lss || treat_as_lss) {
1415  if (skylocal) {
1416  if (cosmics)
1417  fors_science_exit("Cosmic rays correction for LSS or LSS-like "
1418  "data requires --skyglobal=true");
1419  }
1420  }
1421  else {
1422  if (cpl_frameset_count_tags(frameset, curv_coeff_tag) == 0) {
1423  cpl_msg_error(recipe, "Missing required input: %s", curv_coeff_tag);
1424  fors_science_exit(NULL);
1425  }
1426 
1427  if (cpl_frameset_count_tags(frameset, curv_coeff_tag) > 1) {
1428  cpl_msg_error(recipe, "Too many in input: %s", curv_coeff_tag);
1429  fors_science_exit(NULL);
1430  }
1431  }
1432  cpl_propertylist_delete(header); header = NULL;
1433 
1434 
1435  /*
1436  * Remove the master bias
1437  */
1438 
1439  cpl_msg_info(recipe, "Remove the master bias...");
1440 
1441  const cpl_frame * bias_frame =
1442  cpl_frameset_find_const(frameset, "MASTER_BIAS");
1443  bias = fors_image_load(bias_frame);
1444  if (bias == NULL)
1445  fors_science_exit("Cannot load master bias");
1446 
1447  fors_subtract_bias(science_ima, bias);
1448  if (cpl_error_get_code() != CPL_ERROR_NONE)
1449  fors_science_exit("Cannot remove bias from scientific frame");
1450 
1451  fors_image_delete(&bias); bias = NULL;
1452 
1453 
1454  /* Flat-field the science */
1455  cpl_msg_indent_less();
1456  cpl_msg_info(recipe, "Load normalised flat field (if present)...");
1457  cpl_msg_indent_more();
1458 
1459  if (flatfield) {
1460 
1461  const cpl_frame * flat_frame =
1462  cpl_frameset_find_const(frameset, master_norm_flat_tag);
1463  norm_flat = fors_image_load(flat_frame);
1464 
1465  if (norm_flat) {
1466  cpl_msg_info(recipe, "Apply flat field correction...");
1467  fors_image_divide(science_ima, norm_flat);
1468  if (cpl_error_get_code() != CPL_ERROR_NONE) {
1469  cpl_msg_error(recipe, "Failure of flat field correction: %s",
1470  cpl_error_get_message());
1471  fors_science_exit(NULL);
1472  }
1473  fors_image_delete(&norm_flat); norm_flat = NULL;
1474  }
1475  else {
1476  cpl_msg_error(recipe, "Cannot load input %s for flat field "
1477  "correction", master_norm_flat_tag);
1478  fors_science_exit(NULL);
1479  }
1480 
1481  }
1482  spectra = science_ima->data;
1483  ccd_xsize = nx = cpl_image_get_size_x(spectra);
1484  ccd_ysize = ny = cpl_image_get_size_y(spectra);
1485 
1486 
1487  if (skyalign >= 0) {
1488  cpl_msg_indent_less();
1489  cpl_msg_info(recipe, "Load input sky line catalog...");
1490  cpl_msg_indent_more();
1491 
1492  wavelengths = dfs_load_table(frameset, "MASTER_SKYLINECAT", 1);
1493 
1494  if (wavelengths) {
1495 
1496  /*
1497  * Cast the wavelengths into a (double precision) CPL vector
1498  */
1499 
1500  nlines = cpl_table_get_nrow(wavelengths);
1501 
1502  if (nlines == 0)
1503  fors_science_exit("Empty input sky line catalog");
1504 
1505  if (cpl_table_has_column(wavelengths, wcolumn) != 1) {
1506  cpl_msg_error(recipe, "Missing column %s in input line "
1507  "catalog table", wcolumn);
1508  fors_science_exit(NULL);
1509  }
1510 
1511  line = cpl_malloc(nlines * sizeof(double));
1512 
1513  for (i = 0; i < nlines; i++)
1514  line[i] = cpl_table_get(wavelengths, wcolumn, i, NULL);
1515 
1516  cpl_table_delete(wavelengths); wavelengths = NULL;
1517 
1518  lines = cpl_vector_wrap(nlines, line);
1519  }
1520  else {
1521  cpl_msg_info(recipe, "No sky line catalog found in input - fine!");
1522  }
1523  }
1524 
1525 
1526  /*
1527  * Load the slit location table
1528  */
1529 
1530  slits = dfs_load_table(frameset, slit_location_tag, 1);
1531  if (slits == NULL)
1532  fors_science_exit("Cannot load slits location table");
1533 
1534  if (lss || treat_as_lss) {
1535  int first_row = cpl_table_get_double(slits, "ybottom", 0, NULL);
1536  int last_row = cpl_table_get_double(slits, "ytop", 0, NULL);
1537  int ylow, yhig;
1538 
1539  ylow = first_row;
1540  yhig = last_row;
1541 
1542  dummy = cpl_image_extract(spectra, 1, ylow, nx, yhig);
1543  cpl_image_delete(spectra); spectra = dummy; dummy = NULL;
1544  ny = cpl_image_get_size_y(spectra);
1545  }
1546 
1547 
1548  /*
1549  * Load the wavelength calibration table
1550  */
1551 
1552  idscoeff = dfs_load_table(frameset, disp_coeff_tag, 1);
1553  if (idscoeff == NULL)
1554  fors_science_exit("Cannot load wavelength calibration table");
1555 
1556  cpl_msg_indent_less();
1557  cpl_msg_info(recipe, "Processing scientific spectra...");
1558  cpl_msg_indent_more();
1559 
1560  /*
1561  * Load the spectral curvature table, or provide a dummy one
1562  * in case of LSS or LSS-like data (single slit with flat curvature)
1563  */
1564 
1565  if (!(lss || treat_as_lss)) {
1566  polytraces = dfs_load_table(frameset, curv_coeff_tag, 1);
1567  if (polytraces == NULL)
1568  fors_science_exit("Cannot load spectral curvature table");
1569  }
1570  else
1571  {
1572  //Provide a dummy curvature table
1573  int * slit_id = cpl_table_get_data_int(slits, "slit_id");
1574  double * ytop = cpl_table_get_data_double(slits, "ytop");
1575  double * ybottom = cpl_table_get_data_double(slits, "ybottom");
1576  polytraces = cpl_table_new(2);
1577  cpl_table_new_column(polytraces, "slit_id", CPL_TYPE_INT);
1578  cpl_table_new_column(polytraces, "c0", CPL_TYPE_DOUBLE);
1579  cpl_table_set_int(polytraces, "slit_id", 0, slit_id[0]);
1580  cpl_table_set_int(polytraces, "slit_id", 1, slit_id[0]);
1581  cpl_table_set_double(polytraces, "c0", 0, ytop[0]);
1582  cpl_table_set_double(polytraces, "c0", 1, ybottom[0]);
1583  }
1584 
1585  /*
1586  * This one will also generate the spatial map from the spectral
1587  * curvature table (in the case of multislit data)
1588  */
1589 
1590  if (lss || treat_as_lss) {
1591  smapped = cpl_image_duplicate(spectra);
1592  }
1593  else {
1594  coordinate = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
1595 
1596  smapped = mos_spatial_calibration(spectra, slits, polytraces, reference,
1597  startwavelength, endwavelength,
1598  dispersion, 1, coordinate);
1599  }
1600 
1601 
1602  /*
1603  * Generate a rectified wavelength map from the wavelength calibration
1604  * table
1605  */
1606 
1607  rainbow = mos_map_idscoeff(idscoeff, nx, reference, startwavelength,
1608  endwavelength);
1609 
1610  if (dispersion > 1.0)
1611  highres = 0;
1612  else
1613  highres = 1;
1614 
1615  if (skyalign >= 0) {
1616  if (skyalign) {
1617  cpl_msg_info(recipe, "Align wavelength solution to reference "
1618  "skylines applying %d order residual fit...", skyalign);
1619  }
1620  else {
1621  cpl_msg_info(recipe, "Align wavelength solution to reference "
1622  "skylines applying median offset...");
1623  }
1624 
1625  if (lss || treat_as_lss) {
1626  offsets = mos_wavelength_align_lss(smapped, reference,
1627  startwavelength, endwavelength,
1628  idscoeff, lines, highres,
1629  skyalign, rainbow, 4);
1630  }
1631  else {
1632  offsets = mos_wavelength_align(smapped, slits, reference,
1633  startwavelength, endwavelength,
1634  idscoeff, lines, highres, skyalign,
1635  rainbow, 4);
1636  }
1637 
1638  if (offsets) {
1639  if (standard)
1640  cpl_msg_warning(recipe, "Alignment of the wavelength solution "
1641  "to reference sky lines may be unreliable in "
1642  "this case!");
1643 
1644  if (dfs_save_table(frameset, offsets, skylines_offsets_tag, NULL,
1645  parlist, recipe, version))
1646  fors_science_exit(NULL);
1647 
1648  cpl_table_delete(offsets); offsets = NULL;
1649  }
1650  else {
1651  if (cpl_error_get_code()) {
1652  if (cpl_error_get_code() == CPL_ERROR_INCOMPATIBLE_INPUT) {
1653  cpl_msg_error(recipe, "The IDS coeff table is "
1654  "incompatible with the input slit position table.");
1655  }
1656  cpl_msg_error(cpl_func, "Error found in %s: %s",
1657  cpl_error_get_where(), cpl_error_get_message());
1658  fors_science_exit(NULL);
1659  }
1660  cpl_msg_warning(recipe, "Alignment of the wavelength solution "
1661  "to reference sky lines could not be done!");
1662  skyalign = -1;
1663  }
1664 
1665  }
1666 
1667  if (lss || treat_as_lss) {
1668  int first_row = cpl_table_get_double(slits, "ybottom", 0, NULL);
1669  int last_row = cpl_table_get_double(slits, "ytop", 0, NULL);
1670  int ylow, yhig;
1671 
1672  ylow = first_row;
1673  yhig = last_row;
1674 
1675  wavemap = cpl_image_new(ccd_xsize, ccd_ysize, CPL_TYPE_FLOAT);
1676  cpl_image_copy(wavemap, rainbow, 1, ylow);
1677 
1678  wavemaplss = cpl_image_extract(wavemap, 1, ylow, nx, yhig);
1679  }
1680  else {
1681  wavemap = mos_map_wavelengths(coordinate, rainbow, slits,
1682  polytraces, reference,
1683  startwavelength, endwavelength,
1684  dispersion);
1685  }
1686 
1687  cpl_image_delete(rainbow); rainbow = NULL;
1688  cpl_image_delete(coordinate); coordinate = NULL;
1689 
1690  /*
1691  * Here the wavelength calibrated slit spectra are created. This frame
1692  * contains sky_science.
1693  */
1694 
1695  mapped_sky = mos_wavelength_calibration(smapped, reference,
1696  startwavelength, endwavelength,
1697  dispersion, idscoeff, 1);
1698 
1699  cpl_msg_indent_less();
1700  cpl_msg_info(recipe, "Check applied wavelength against skylines...");
1701  cpl_msg_indent_more();
1702 
1703  mean_rms = mos_distortions_rms(mapped_sky, lines, startwavelength,
1704  dispersion, 6, highres);
1705 
1706  cpl_vector_delete(lines); lines = NULL;
1707 
1708  cpl_msg_info(recipe, "Mean residual: %f", mean_rms);
1709 
1710  mean_rms = cpl_table_get_column_mean(idscoeff, "error");
1711 
1712  cpl_msg_info(recipe, "Mean model accuracy: %f pixel (%f A)",
1713  mean_rms, mean_rms * dispersion);
1714 
1715  header = cpl_propertylist_new();
1716  cpl_propertylist_update_double(header, "CRPIX1", 1.0);
1717  cpl_propertylist_update_double(header, "CRPIX2", 1.0);
1718  cpl_propertylist_update_double(header, "CRVAL1",
1719  startwavelength + dispersion/2);
1720  cpl_propertylist_update_double(header, "CRVAL2", 1.0);
1721  /* cpl_propertylist_update_double(header, "CDELT1", dispersion);
1722  cpl_propertylist_update_double(header, "CDELT2", 1.0); */
1723  cpl_propertylist_update_double(header, "CD1_1", dispersion);
1724  cpl_propertylist_update_double(header, "CD1_2", 0.0);
1725  cpl_propertylist_update_double(header, "CD2_1", 0.0);
1726  cpl_propertylist_update_double(header, "CD2_2", 1.0);
1727  cpl_propertylist_update_string(header, "CTYPE1", "LINEAR");
1728  cpl_propertylist_update_string(header, "CTYPE2", "PIXEL");
1729 
1730  /* Perform time normalisation */
1731  cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
1732  dummy = cpl_image_divide_scalar_create(mapped_sky, alltime);
1733  if (dfs_save_image(frameset, dummy, mapped_science_sky_tag, header,
1734  parlist, recipe, version))
1735  fors_science_exit(NULL);
1736  cpl_image_delete(dummy); dummy = NULL;
1737 
1738 /* if (skyglobal == 0 && skymedian < 0) { NSS */
1739  if (skyglobal == 0 && skymedian == 0 && skylocal == 0) {
1740  cpl_image_delete(mapped_sky); mapped_sky = NULL;
1741  }
1742 
1743  if (skyglobal || skylocal) {
1744 
1745  cpl_msg_indent_less();
1746 
1747  if (skyglobal) {
1748  cpl_msg_info(recipe, "Global sky determination...");
1749  cpl_msg_indent_more();
1750  skymap = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
1751  if (lss || treat_as_lss) {
1752  sky = mos_sky_map_super(spectra, wavemaplss, dispersion,
1753  2.0, 50, skymap);
1754  }
1755  else {
1756  sky = mos_sky_map_super(spectra, wavemap, dispersion,
1757  2.0, 50, skymap);
1758  }
1759  if (sky) {
1760  cpl_image_subtract(spectra, skymap);
1761  }
1762  else {
1763  cpl_image_delete(skymap); skymap = NULL;
1764  }
1765  }
1766  else {
1767  cpl_msg_info(recipe, "Local sky determination...");
1768  cpl_msg_indent_more();
1769  skymap = mos_subtract_sky(spectra, slits, polytraces, reference,
1770  startwavelength, endwavelength, dispersion);
1771  }
1772 
1773  if (skymap) {
1774  if (skyglobal) {
1775  /* Perform time normalisation */
1776  cpl_table_divide_scalar(sky, "sky", alltime);
1777 
1778  if (dfs_save_table(frameset, sky, global_sky_spectrum_tag,
1779  NULL, parlist, recipe, version))
1780  fors_science_exit(NULL);
1781 
1782  cpl_table_delete(sky); sky = NULL;
1783  }
1784 
1785  save_header = dfs_load_header(frameset, science_tag, 0);
1786 
1787  /* Perform time normalisation */
1788  cpl_image_divide_scalar(skymap, alltime);
1789 
1790  if (dfs_save_image(frameset, skymap, unmapped_sky_tag,
1791  save_header, parlist, recipe, version))
1792  fors_science_exit(NULL);
1793 
1794  cpl_image_delete(skymap); skymap = NULL;
1795 
1796  if (dfs_save_image(frameset, spectra, unmapped_science_tag,
1797  save_header, parlist, recipe, version))
1798  fors_science_exit(NULL);
1799 
1800  cpl_propertylist_delete(save_header); save_header = NULL;
1801 
1802  if (cosmics) {
1803  cpl_msg_info(recipe, "Removing cosmic rays...");
1804  mos_clean_cosmics(spectra, gain, -1., -1.);
1805  }
1806 
1807  /*
1808  * The spatially rectified image, that contained the sky,
1809  * is replaced by a sky-subtracted spatially rectified image:
1810  */
1811 
1812  cpl_image_delete(smapped); smapped = NULL;
1813 
1814  if (lss || treat_as_lss) {
1815  smapped = cpl_image_duplicate(spectra);
1816  }
1817  else {
1818  smapped = mos_spatial_calibration(spectra, slits, polytraces,
1819  reference, startwavelength,
1820  endwavelength, dispersion,
1821  1, NULL);
1822  }
1823  }
1824  else {
1825  cpl_msg_warning(recipe, "Sky subtraction failure");
1826  if (cosmics)
1827  cpl_msg_warning(recipe, "Cosmic rays removal not performed!");
1828  cosmics = skylocal = skyglobal = 0;
1829  }
1830  }
1831 
1832  cpl_image_delete(spectra); spectra = NULL;
1833  cpl_table_delete(polytraces); polytraces = NULL;
1834  if (lss || treat_as_lss)
1835  cpl_image_delete(wavemaplss); wavemaplss = NULL;
1836 
1837  if (skyalign >= 0) {
1838  save_header = dfs_load_header(frameset, science_tag, 0);
1839  cpl_propertylist_update_string(save_header, "BUNIT", "Angstrom");
1840  if (dfs_save_image(frameset, wavemap, wavelength_map_sky_tag,
1841  save_header, parlist, recipe, version))
1842  fors_science_exit(NULL);
1843  cpl_propertylist_delete(save_header); save_header = NULL;
1844  }
1845 
1846  cpl_image_delete(wavemap); wavemap = NULL;
1847 
1848  mapped = mos_wavelength_calibration(smapped, reference,
1849  startwavelength, endwavelength,
1850  dispersion, idscoeff, 1);
1851 
1852  cpl_image_delete(smapped); smapped = NULL;
1853 
1854 /* if (skymedian >= 0) { NSS */
1855  if (skymedian) {
1856  cpl_msg_indent_less();
1857  cpl_msg_info(recipe, "Local sky determination...");
1858  cpl_msg_indent_more();
1859 
1860 /* NSS skylocalmap = mos_sky_local(mapped, slits, skymedian); */
1861 /* skylocalmap = mos_sky_local(mapped, slits, 0); */
1862  skylocalmap = mos_sky_local_old(mapped, slits);
1863  cpl_image_subtract(mapped, skylocalmap);
1864 /*
1865  if (dfs_save_image(frameset, skylocalmap, mapped_sky_tag, header,
1866  parlist, recipe, version))
1867  fors_science_exit(NULL);
1868 */
1869  cpl_image_delete(skylocalmap); skylocalmap = NULL;
1870  }
1871 
1872  if (photometry) {
1873  if (have_phot && !standard) {
1874  if (cpl_frameset_count_tags(frameset, specphot_tag) == 0) {
1875  photcal = dfs_load_table(frameset,
1876  master_specphot_tag, 1);
1877  }
1878  else {
1879  photcal = dfs_load_table(frameset, specphot_tag, 1);
1880  }
1881  }
1882  }
1883 
1884 
1885 /* if (skyglobal || skymedian >= 0 || skylocal) { NSS */
1886  if (skyglobal || skymedian || skylocal) {
1887 
1888  skylocalmap = cpl_image_subtract_create(mapped_sky, mapped);
1889 
1890  cpl_image_delete(mapped_sky); mapped_sky = NULL;
1891 
1892  /* Perform time normalisation */
1893  cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
1894  dummy = cpl_image_divide_scalar_create(skylocalmap, alltime);
1895  if (dfs_save_image(frameset, dummy, mapped_sky_tag, header,
1896  parlist, recipe, version))
1897  fors_science_exit(NULL);
1898  cpl_image_delete(dummy); dummy = NULL;
1899 
1900  cpl_msg_indent_less();
1901  cpl_msg_info(recipe, "Object detection...");
1902  cpl_msg_indent_more();
1903 
1904  if (cosmics || nscience > 1) {
1905  dummy = mos_detect_objects(mapped, slits, slit_margin, ext_radius,
1906  cont_radius);
1907  }
1908  else {
1909  mapped_cleaned = cpl_image_duplicate(mapped);
1910  mos_clean_cosmics(mapped_cleaned, gain, -1., -1.);
1911  dummy = mos_detect_objects(mapped_cleaned, slits, slit_margin,
1912  ext_radius, cont_radius);
1913 
1914  cpl_image_delete(mapped_cleaned); mapped_cleaned = NULL;
1915  }
1916 
1917  cpl_image_delete(dummy); dummy = NULL;
1918 
1919  if (dfs_save_table(frameset, slits, object_table_tag, NULL, parlist,
1920  recipe, version))
1921  fors_science_exit(NULL);
1922 
1923  cpl_msg_indent_less();
1924  cpl_msg_info(recipe, "Object extraction...");
1925  cpl_msg_indent_more();
1926 
1927  images = mos_extract_objects(mapped, skylocalmap, slits,
1928  ext_mode, ron, gain, 1);
1929 
1930  cpl_image_delete(skylocalmap); skylocalmap = NULL;
1931 
1932  if (images) {
1933  if (standard) {
1934  cpl_table *ext_table = NULL;
1935  cpl_table *flux_table = NULL;
1936 
1937  ext_table = dfs_load_table(frameset, "EXTINCT_TABLE", 1);
1938  flux_table = dfs_load_table(frameset, "STD_FLUX_TABLE", 1);
1939 
1940  photcal = mos_photometric_calibration(images[0],
1941  startwavelength,
1942  dispersion, gain,
1943  alltime, ext_table,
1944  airmass, flux_table,
1945  res_order);
1946 
1947  cpl_table_delete(ext_table);
1948  cpl_table_delete(flux_table);
1949 
1950  if (photcal) {
1951 
1952  float *data;
1953  char *pipefile = NULL;
1954  char keyname[30];
1955 
1956  qclist = dfs_load_header(frameset, science_tag, 0);
1957 
1958  if (qclist == NULL)
1959  fors_science_exit("Cannot reload scientific "
1960  "frame header");
1961 
1962  fors_qc_start_group(qclist, "2.0", instrume);
1963 
1964 
1965  /*
1966  * QC1 group header
1967  */
1968 
1969  if (fors_qc_write_string("PRO.CATG", specphot_tag,
1970  "Product category", instrume))
1971  fors_science_exit("Cannot write product category "
1972  "to QC log file");
1973 
1974  if (fors_qc_keyword_to_paf(qclist, "ESO DPR TYPE", NULL,
1975  "DPR type", instrume))
1976  fors_science_exit("Missing keyword DPR TYPE in "
1977  "scientific frame header");
1978 
1979  if (fors_qc_keyword_to_paf(qclist, "ESO TPL ID", NULL,
1980  "Template", instrume))
1981  fors_science_exit("Missing keyword TPL ID in "
1982  "scientific frame header");
1983 
1984  if (fors_qc_keyword_to_paf(qclist,
1985  "ESO INS GRIS1 NAME", NULL,
1986  "Grism name", instrume))
1987  fors_science_exit("Missing keyword INS GRIS1 NAME "
1988  "in scientific frame header");
1989 
1990  if (fors_qc_keyword_to_paf(qclist,
1991  "ESO INS GRIS1 ID", NULL,
1992  "Grism identifier",
1993  instrume))
1994  fors_science_exit("Missing keyword INS GRIS1 ID "
1995  "in scientific frame header");
1996 
1997  if (cpl_propertylist_has(qclist, "ESO INS FILT1 NAME"))
1998  fors_qc_keyword_to_paf(qclist,
1999  "ESO INS FILT1 NAME", NULL,
2000  "Filter name", instrume);
2001 
2002  if (fors_qc_keyword_to_paf(qclist,
2003  "ESO INS COLL NAME", NULL,
2004  "Collimator name",
2005  instrume))
2006  fors_science_exit("Missing keyword INS COLL NAME "
2007  "in scientific frame header");
2008 
2009  if (fors_qc_keyword_to_paf(qclist,
2010  "ESO DET CHIP1 ID", NULL,
2011  "Chip identifier",
2012  instrume))
2013  fors_science_exit("Missing keyword DET CHIP1 ID "
2014  "in scientific frame header");
2015 
2016  if (fors_qc_keyword_to_paf(qclist,
2017  "ESO INS MOS10 WID",
2018  "arcsec", "Slit width",
2019  instrume)) {
2020  cpl_error_reset();
2021  if (fors_qc_keyword_to_paf(qclist,
2022  "ESO INS SLIT WID",
2023  "arcsec", "Slit width",
2024  instrume)) {
2025  if (mos) {
2026  fors_science_exit("Missing keyword "
2027  "ESO INS MOS10 WID in "
2028  "scientific frame header");
2029  }
2030  else {
2031  fors_science_exit("Missing keyword "
2032  "ESO INS SLIT WID in "
2033  "scientific frame header");
2034  }
2035  }
2036  }
2037 
2038  /*
2039  if (mos) {
2040  if (fors_qc_keyword_to_paf(qclist,
2041  "ESO INS MOS10 WID",
2042  "arcsec", "Slit width",
2043  instrume))
2044  fors_science_exit("Missing keyword "
2045  "ESO INS MOS10 WID in "
2046  "scientific frame header");
2047  }
2048  else {
2049  if (fors_qc_keyword_to_paf(qclist,
2050  "ESO INS SLIT WID",
2051  "arcsec", "Slit width",
2052  instrume))
2053  fors_science_exit("Missing keyword "
2054  "ESO INS SLIT WID in "
2055  "scientific frame header");
2056  }
2057  */
2058 
2059  if (fors_qc_keyword_to_paf(qclist,
2060  "ESO DET WIN1 BINX", NULL,
2061  "Binning factor along X",
2062  instrume))
2063  fors_science_exit("Missing keyword ESO "
2064  "DET WIN1 BINX "
2065  "in scientific frame header");
2066 
2067  if (fors_qc_keyword_to_paf(qclist,
2068  "ESO DET WIN1 BINY", NULL,
2069  "Binning factor along Y",
2070  instrume))
2071  fors_science_exit("Missing keyword "
2072  "ESO DET WIN1 BINY "
2073  "in scientific frame header");
2074 
2075  if (fors_qc_keyword_to_paf(qclist, "ARCFILE", NULL,
2076  "Archive name of input data",
2077  instrume))
2078  fors_science_exit("Missing keyword ARCFILE in "
2079  "scientific frame header");
2080 
2081  pipefile = dfs_generate_filename(specphot_tag);
2082  if (fors_qc_write_string("PIPEFILE", pipefile,
2083  "Pipeline product name",
2084  instrume))
2085  fors_science_exit("Cannot write PIPEFILE to "
2086  "QC log file");
2087  cpl_free(pipefile);
2088 
2089 
2090  /*
2091  * QC1 parameters
2092  */
2093 
2094  wstart = 3700.;
2095  wstep = 400.;
2096  wcount = 15;
2097 
2098  dummy = cpl_image_new(wcount, 1, CPL_TYPE_FLOAT);
2099  data = cpl_image_get_data_float(dummy);
2100  map_table(dummy, wstart, wstep, photcal,
2101  "WAVE", "EFFICIENCY");
2102 
2103  for (i = 0; i < wcount; i++) {
2104  sprintf(keyname, "QC.SPEC.EFFICIENCY%d.LAMBDA",
2105  i + 1);
2106  if (fors_qc_write_qc_double(qclist,
2107  wstart + wstep * i,
2108  keyname, "Angstrom",
2109  "Wavelength of "
2110  "efficiency evaluation",
2111  instrume)) {
2112  fors_science_exit("Cannot write wavelength of "
2113  "efficiency evaluation");
2114  }
2115 
2116  sprintf(keyname, "QC.SPEC.EFFICIENCY%d", i + 1);
2117  if (fors_qc_write_qc_double(qclist,
2118  data[i],
2119  keyname, "e-/photon",
2120  "Efficiency",
2121  instrume)) {
2122  fors_science_exit("Cannot write wavelength of "
2123  "efficiency evaluation");
2124  }
2125  }
2126 
2127  cpl_image_delete(dummy); dummy = NULL;
2128 
2130 
2131  if (dfs_save_table(frameset, photcal, specphot_tag, qclist,
2132  parlist, recipe, version))
2133  fors_science_exit(NULL);
2134 
2135  cpl_propertylist_delete(qclist); qclist = NULL;
2136 
2137  if (have_phot) {
2138  cpl_table_delete(photcal); photcal = NULL;
2139  }
2140  }
2141  }
2142 
2143  if (photometry) {
2144  cpl_image *calibrated;
2145  cpl_table *ext_table;
2146 
2147  ext_table = dfs_load_table(frameset, "EXTINCT_TABLE", 1);
2148 
2149  if (have_phot) {
2150  if (cpl_frameset_count_tags(frameset, specphot_tag) == 0) {
2151  photcal = dfs_load_table(frameset,
2152  master_specphot_tag, 1);
2153  }
2154  else {
2155  photcal = dfs_load_table(frameset, specphot_tag, 1);
2156  }
2157  }
2158 
2159  calibrated = mos_apply_photometry(images[0], photcal,
2160  ext_table, startwavelength,
2161  dispersion, gain, alltime,
2162  airmass);
2163  cpl_propertylist_update_string(header, "BUNIT",
2164  "10^(-16) erg/(cm^2 s Angstrom)");
2165 
2166  if (dfs_save_image(frameset, calibrated,
2167  reduced_flux_science_tag, header,
2168  parlist, recipe, version)) {
2169  cpl_image_delete(calibrated);
2170  fors_science_exit(NULL);
2171  }
2172 
2173  cpl_table_delete(ext_table);
2174  cpl_image_delete(calibrated);
2175  }
2176 
2177  /* Perform time normalisation */
2178  cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
2179  cpl_image_divide_scalar(images[0], alltime);
2180 
2181  if (dfs_save_image(frameset, images[0], reduced_science_tag, header,
2182  parlist, recipe, version))
2183  fors_science_exit(NULL);
2184 
2185  /* Perform time normalisation */
2186  cpl_image_divide_scalar(images[1], alltime);
2187 
2188  if (dfs_save_image(frameset, images[1], reduced_sky_tag, header,
2189  parlist, recipe, version))
2190  fors_science_exit(NULL);
2191 
2192  cpl_image_delete(images[1]);
2193 
2194  if (photometry) {
2195  cpl_image *calibrated;
2196  cpl_table *ext_table;
2197 
2198  ext_table = dfs_load_table(frameset, "EXTINCT_TABLE", 1);
2199 
2200  calibrated = mos_propagate_photometry_error(images[0],
2201  images[2], photcal,
2202  ext_table, startwavelength,
2203  dispersion, gain, alltime,
2204  airmass);
2205 
2206  cpl_propertylist_update_string(header, "BUNIT",
2207  "10^(-16) erg/(cm^2 s Angstrom)");
2208 
2209  if (dfs_save_image(frameset, calibrated,
2210  reduced_flux_error_tag, header,
2211  parlist, recipe, version)) {
2212  cpl_image_delete(calibrated);
2213  fors_science_exit(NULL);
2214  }
2215 
2216  cpl_table_delete(ext_table);
2217  cpl_image_delete(calibrated);
2218  }
2219 
2220 
2221  /* Perform time normalisation */
2222  cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
2223  cpl_image_divide_scalar(images[2], alltime);
2224 
2225  if (dfs_save_image(frameset, images[2], reduced_error_tag, header,
2226  parlist, recipe, version))
2227  fors_science_exit(NULL);
2228 
2229  cpl_image_delete(images[0]);
2230  cpl_image_delete(images[2]);
2231 
2232  cpl_free(images);
2233  }
2234  else {
2235  cpl_msg_warning(recipe, "No objects found: the products "
2236  "%s, %s, and %s are not created",
2237  reduced_science_tag, reduced_sky_tag,
2238  reduced_error_tag);
2239  }
2240 
2241  }
2242 
2243  cpl_free(instrume); instrume = NULL;
2244  cpl_table_delete(slits); slits = NULL;
2245 
2246  if (skyalign >= 0) {
2247  if (dfs_save_table(frameset, idscoeff, disp_coeff_sky_tag, NULL,
2248  parlist, recipe, version))
2249  fors_science_exit(NULL);
2250  }
2251 
2252  cpl_table_delete(idscoeff); idscoeff = NULL;
2253 
2254  if (photometry && photcal) {
2255  cpl_image *calibrated;
2256  cpl_table *ext_table;
2257 
2258 
2259  ext_table = dfs_load_table(frameset, "EXTINCT_TABLE", 1);
2260 
2261  calibrated = mos_apply_photometry(mapped, photcal,
2262  ext_table, startwavelength,
2263  dispersion, gain, alltime,
2264  airmass);
2265 
2266  cpl_propertylist_update_string(header, "BUNIT",
2267  "10^(-16) erg/(cm^2 s Angstrom)");
2268 
2269  if (dfs_save_image(frameset, calibrated,
2270  mapped_flux_science_tag, header,
2271  parlist, recipe, version)) {
2272  cpl_image_delete(calibrated);
2273  fors_science_exit(NULL);
2274  }
2275 
2276  cpl_table_delete(ext_table);
2277  cpl_image_delete(calibrated);
2278  }
2279 
2280  /* Perform time normalisation */
2281  cpl_propertylist_update_string(header, "BUNIT", "ADU/s");
2282  cpl_image_divide_scalar(mapped, alltime);
2283 
2284  if (dfs_save_image(frameset, mapped, mapped_science_tag, header,
2285  parlist, recipe, version))
2286  fors_science_exit(NULL);
2287 
2288  cpl_table_delete(photcal); photcal = NULL;
2289  cpl_image_delete(mapped); mapped = NULL;
2290  cpl_propertylist_delete(header); header = NULL;
2291 
2292  if (cpl_error_get_code()) {
2293  cpl_msg_error(cpl_func, "Error found in %s: %s",
2294  cpl_error_get_where(), cpl_error_get_message());
2295  fors_science_exit(NULL);
2296  }
2297  else
2298  return 0;
2299  return 0;
2300 }
cpl_image * mos_spatial_calibration(cpl_image *spectra, cpl_table *slits, cpl_table *polytraces, double reference, double blue, double red, double dispersion, int flux, cpl_image *calibration)
Spatial remapping of CCD spectra eliminating the spectral curvature.
Definition: moses.c:8264
cpl_table * mos_photometric_calibration(cpl_image *spectra, double startwave, double dispersion, double gain, double exptime, cpl_table *ext_table, double airmass, cpl_table *flux_table, int order)
Produce instrument response curve, with some ancillary information.
Definition: moses.c:17204
void fors_image_multiply_scalar(fors_image *image, double s, double ds)
Multiply by scalar.
Definition: fors_image.c:882
int cpl_plugin_get_info(cpl_pluginlist *list)
Build the list of available plugins, for this module.
Definition: fors_bias.c:62
const char * dfs_get_parameter_string(cpl_parameterlist *parlist, const char *name, const cpl_table *defaults)
Reading a recipe string parameter value.
Definition: fors_dfs.c:587
double mos_distortions_rms(cpl_image *rectified, cpl_vector *lines, double wavestart, double dispersion, int radius, int highres)
Estimate the spectral distortion modeling goodness.
Definition: moses.c:10754
cpl_image * mos_propagate_photometry_error(cpl_image *spectra, cpl_image *errors, cpl_table *response, cpl_table *ext_table, double startwave, double dispersion, double gain, double exptime, double airmass)
Propagate errors from response curve and extracted spectra.
Definition: moses.c:17849
cpl_propertylist * dfs_load_header(cpl_frameset *frameset, const char *category, int ext)
Loading header associated to data of given category.
Definition: fors_dfs.c:951
cpl_table * mos_load_slits_fors_mxu(cpl_propertylist *header)
Create slit location table from FITS header of FORS2-MXU data.
Definition: moses.c:14561
cpl_image ** mos_extract_objects(cpl_image *science, cpl_image *sky, cpl_table *objects, int extraction, double ron, double gain, int ncombined)
Extract detected objects from rectified scientific frame.
Definition: moses.c:14002
cpl_error_code fors_qc_write_qc_double(cpl_propertylist *header, double value, const char *name, const char *unit, const char *comment, const char *instrument)
Write an integer value to the active QC1 PAF object and to a header.
Definition: fors_qc.c:604
cpl_image * mos_map_wavelengths(cpl_image *spatial, cpl_image *calibration, cpl_table *slits, cpl_table *polytraces, double reference, double blue, double red, double dispersion)
Remapping of spatially rectified wavelengths to original CCD pixels.
Definition: moses.c:11087
cpl_error_code fors_qc_keyword_to_paf(cpl_propertylist *header, const char *name, const char *unit, const char *comment, const char *instrument)
Copy a keyword value to the currently active QC1 PAF object.
Definition: fors_qc.c:425
cpl_image * mos_wavelength_calibration(cpl_image *image, double refwave, double firstLambda, double lastLambda, double dispersion, cpl_table *idscoeff, int flux)
Remap at constant wavelength step an image of rectified scientific spectra.
Definition: moses.c:9412
cpl_table * mos_sky_map_super(cpl_image *spectra, cpl_image *wavemap, double dispersion, double factor, int minpoints, cpl_image *skymap)
Create a CCD median sky map.
Definition: moses.c:11849
cpl_table * mos_load_slits_fors_mos(cpl_propertylist *header, int *nslits_out_det)
Create slit location table from FITS header of FORS1/2 MOS data.
Definition: moses.c:14801
void fors_image_delete(fors_image **image)
Deallocate image and set pointer to NULL.
Definition: fors_image.c:162
cpl_table * mos_wavelength_align(cpl_image *image, cpl_table *slits, double refwave, double firstLambda, double lastLambda, cpl_table *idscoeff, cpl_vector *skylines, int highres, int order, cpl_image *calibration, int sradius)
Modify the input wavelength solution to match reference sky lines.
Definition: moses.c:9693
void fors_image_divide_scalar(fors_image *image, double s, double ds)
Divide by scalar.
Definition: fors_image.c:856
cpl_error_code fors_qc_start_group(cpl_propertylist *header, const char *qcdic_version, const char *instrument)
Initiate a new QC1 group.
Definition: fors_qc.c:77
cpl_error_code fors_qc_write_string(const char *name, const char *value, const char *comment, const char *instrument)
Add string parameter to current QC1 group.
Definition: fors_qc.c:235
cpl_image * mos_detect_objects(cpl_image *image, cpl_table *slits, int margin, int maxradius, int conradius)
Detect objects in rectified scientific frame.
Definition: moses.c:13651
cpl_image * mos_sky_local_old(cpl_image *spectra, cpl_table *slits)
Local determination of sky.
Definition: moses.c:12424
void fors_dfs_set_groups(cpl_frameset *set)
Set the group as RAW or CALIB in a frameset.
Definition: fors_dfs.c:257
int dfs_get_parameter_bool(cpl_parameterlist *parlist, const char *name, const cpl_table *defaults)
Reading a recipe boolean parameter value.
Definition: fors_dfs.c:686
cpl_frameset * fors_frameset_extract(const cpl_frameset *frames, const char *tag)
Extract frames with given tag from frameset.
Definition: fors_utils.c:468
int dfs_equal_keyword(cpl_frameset *frameset, const char *keyword)
Saving table data of given category.
Definition: fors_dfs.c:1685
fors_image * fors_image_load(const cpl_frame *frame)
Load image.
Definition: fors_image.c:298
cpl_table * mos_wavelength_align_lss(cpl_image *image, double refwave, double firstLambda, double lastLambda, cpl_table *idscoeff, cpl_vector *skylines, int highres, int order, cpl_image *calibration, int sradius)
Modify the input wavelength solution to match reference sky lines (LSS).
Definition: moses.c:10246
cpl_image * mos_subtract_sky(cpl_image *science, cpl_table *slits, cpl_table *polytraces, double reference, double blue, double red, double dispersion)
Subtract the sky from the scientific CCD exposure.
Definition: moses.c:1969
cpl_image * mos_map_idscoeff(cpl_table *idscoeff, int xsize, double reference, double blue, double red)
Create a wavelengths map from an IDS coefficients table.
Definition: moses.c:10976
int dfs_save_image(cpl_frameset *frameset, const cpl_image *image, const char *category, cpl_propertylist *header, const cpl_parameterlist *parlist, const char *recipename, const char *version)
Saving image data of given category.
Definition: fors_dfs.c:1447
cpl_table * dfs_load_table(cpl_frameset *frameset, const char *category, int ext)
Loading table data of given category.
Definition: fors_dfs.c:901
int dfs_get_parameter_int(cpl_parameterlist *parlist, const char *name, const cpl_table *defaults)
Reading a recipe integer parameter value.
Definition: fors_dfs.c:392
int dfs_save_table(cpl_frameset *frameset, const cpl_table *table, const char *category, cpl_propertylist *header, const cpl_parameterlist *parlist, const char *recipename, const char *version)
Saving table data of given category.
Definition: fors_dfs.c:1574
cpl_image * mos_apply_photometry(cpl_image *spectra, cpl_table *response, cpl_table *ext_table, double startwave, double dispersion, double gain, double exptime, double airmass)
Apply response curve to extracted spectra.
Definition: moses.c:17756
cpl_error_code fors_qc_end_group(void)
Close current QC1 PAF file.
Definition: fors_qc.c:200
void fors_image_divide(fors_image *left, const fors_image *right)
Divide images.
Definition: fors_image.c:733
Definition: list.c:74
double fors_get_airmass(const cpl_propertylist *header)
Compute average airmass.
Definition: fors_tools.c:398
cpl_error_code mos_clean_cosmics(cpl_image *image, float gain, float threshold, float ratio)
Remove cosmic rays from sky-subtracted CCD spectral exposure.
Definition: moses.c:12881
double dfs_get_parameter_double(cpl_parameterlist *parlist, const char *name, const cpl_table *defaults)
Reading a recipe double parameter value.
Definition: fors_dfs.c:489
const char * fors_get_license(void)
Get the pipeline copyright and license.
Definition: fors_utils.c:65