CR2RE Pipeline Reference Manual 1.7.0
cr2res_cal_flat.c
1/*
2 * This file is part of the CR2RES Pipeline
3 * Copyright (C) 2002,2003 European Southern Observatory
4 *
5 * This program is free software; you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation; either version 2 of the License, or
8 * (at your option) any later version.
9 *
10 * This program is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this program; if not, write to the Free Software
17 * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02111-1307 USA
18 */
19
20#ifdef HAVE_CONFIG_H
21#include <config.h>
22#endif
23
24/*-----------------------------------------------------------------------------
25 Includes
26 -----------------------------------------------------------------------------*/
27
28#include <string.h>
29#include <cpl.h>
30
31#include "cr2res_utils.h"
32#include "cr2res_calib.h"
33#include "cr2res_pfits.h"
34#include "cr2res_dfs.h"
35#include "cr2res_qc.h"
36#include "cr2res_flat.h"
37#include "cr2res_bpm.h"
38#include "cr2res_trace.h"
39#include "cr2res_extract.h"
40#include "cr2res_io.h"
41
42/*-----------------------------------------------------------------------------
43 Define
44 -----------------------------------------------------------------------------*/
45
46#define RECIPE_STRING "cr2res_cal_flat"
47
48/*-----------------------------------------------------------------------------
49 Plugin registration
50 -----------------------------------------------------------------------------*/
51
52int cpl_plugin_get_info(cpl_pluginlist * list);
53
54/*-----------------------------------------------------------------------------
55 Private function prototypes
56 -----------------------------------------------------------------------------*/
57
58static int cr2res_cal_flat_compare(
59 const cpl_frame * frame1,
60 const cpl_frame * frame2) ;
61static int cr2res_cal_flat_reduce(
62 const cpl_frameset * rawframes,
63 const cpl_frame * tw_frame,
64 const cpl_frame * detlin_frame,
65 const cpl_frame * master_dark_frame,
66 const cpl_frame * bpm_frame,
67 int filter_traces,
68 int subtract_nolight_rows,
69 int cosmics,
70 double bpm_low,
71 double bpm_high,
72 double bpm_linemax,
73 int trace_degree,
74 int trace_min_cluster,
75 int trace_smooth_x,
76 int trace_smooth_y,
77 double trace_threshold,
78 int trace_opening,
79 cr2res_extr_method extr_method,
80 int extract_oversample,
81 int extract_swath_width,
82 int extract_height,
83 double extract_smooth_slit,
84 double extract_smooth_spec,
85 int reduce_det,
86 int reduce_order,
87 int reduce_trace,
88 hdrl_image ** master_flat,
89 cpl_table ** trace_wave,
90 cpl_table ** slit_func,
91 cpl_table ** extract_1d,
92 hdrl_image ** slit_model,
93 cpl_image ** bpm,
94 cpl_propertylist ** ext_plist) ;
95static int cr2res_cal_flat_create(cpl_plugin *);
96static int cr2res_cal_flat_exec(cpl_plugin *);
97static int cr2res_cal_flat_destroy(cpl_plugin *);
98static int cr2res_cal_flat(cpl_frameset *, const cpl_parameterlist *);
99
100/*-----------------------------------------------------------------------------
101 Static variables
102 -----------------------------------------------------------------------------*/
103
104static char cr2res_cal_flat_description[] = "\
105Flat \n\
106 Compute the master flat and perform the traces detection. \n\
107 If a TW file is provided, its traces and slit curvatures are used \n\
108 for the extraction needed to create the master flat. \n\
109 If not provided, the measured traces are used, and the slit is \n\
110 considered vertical. \n\
111 \n\
112 Inputs \n\
113 raw.fits " CR2RES_FLAT_RAW " [1 to n] \n\
114 trace.fits " CR2RES_CAL_FLAT_TW_PROCATG " [0 to 1] \n\
115 or " CR2RES_CAL_FLAT_TW_MERGED_PROCATG " \n\
116 or " CR2RES_UTIL_TRACE_TW_PROCATG " \n\
117 or " CR2RES_UTIL_WAVE_TW_PROCATG " \n\
118 or " CR2RES_CAL_WAVE_TW_PROCATG " \n\
119 or " CR2RES_UTIL_SLIT_CURV_TW_PROCATG " \n\
120 detlin.fits " CR2RES_CAL_DETLIN_COEFFS_PROCATG " [0 to 1] \n\
121 master_dark.fits " CR2RES_CAL_DARK_MASTER_PROCATG " [0 to 1] \n\
122 bpm.fits " CR2RES_CAL_DARK_BPM_PROCATG " [0 to 1] \n\
123 or " CR2RES_CAL_FLAT_BPM_PROCATG " \n\
124 or " CR2RES_CAL_DETLIN_BPM_PROCATG " \n\
125 or " CR2RES_UTIL_BPM_SPLIT_PROCATG " \n\
126 or " CR2RES_UTIL_NORM_BPM_PROCATG " \n\
127 \n\
128 Outputs \n\
129 cr2res_cal_flat_[setting]_[Decker]_bpm.fits "
130 CR2RES_CAL_FLAT_BPM_PROCATG "\n\
131 cr2res_cal_flat_[setting]_[Decker]_blaze.fits "
132 CR2RES_CAL_FLAT_EXTRACT_1D_PROCATG "\n\
133 cr2res_cal_flat_[setting]_[Decker]_slit_model.fits "
134 CR2RES_CAL_FLAT_SLIT_MODEL_PROCATG "\n\
135 cr2res_cal_flat_[setting]_[Decker]_slit_func.fits "
136 CR2RES_CAL_FLAT_SLIT_FUNC_PROCATG "\n\
137 cr2res_cal_flat_[setting]_[Decker]_master_flat.fits "
138 CR2RES_CAL_FLAT_MASTER_PROCATG "\n\
139 cr2res_cal_flat_[setting]_[Decker]_tw.fits "
140 CR2RES_CAL_FLAT_TW_PROCATG "\n\
141 cr2res_cal_flat_[setting]_tw_merged.fits "
142 CR2RES_CAL_FLAT_TW_MERGED_PROCATG "\n\
143 \n\
144 Algorithm \n\
145 group the input frames by different settings \n\
146 loop on groups g: \n\
147 loop on decker positions p: \n\
148 loop on detectors d: \n\
149 cr2res_cal_flat_reduce() computes (master_flat, trace_wave, \n\
150 slit_func,extract_1d,slit_model,bpm)(g,p,d) \n\
151 Save slit_model(g,p) \n\
152 Save extract_1d(g,p) \n\
153 Save master_flat(g,p) \n\
154 Save trace_wave(g,p) \n\
155 Save slit_func(g,p) \n\
156 Save bpm(g,p) \n\
157 Merge the trace_wave(g,p,d) into trace_wave(g,d) \n\
158 Save trace_wave(g) \n\
159 \n\
160 cr2res_cal_flat_reduce() \n\
161 Load the images list \n\
162 Apply the detlin / dark / bpm calibrations \n\
163 Average the images to avg \n\
164 Compute the traces with cr2res_trace(--trace_degree, \n\
165 --trace_min_cluster, --trace_smooth, --trace_opening) \n\
166 on avg \n\
167 For the following step, the computed TW is used if there is no \n\
168 input TW provided, the provided one is used otherwise. \n\
169 loop on the traces t: \n\
170 cr2res_extract_slitdec_curved(--extract_oversample, \n\
171 --extract_swath_width, --extract_height, \n\
172 --extract_smooth_slit, --extract_smooth_spec) \n\
173 -> slit_func(t), extract_1d(t), slit_model(t) \n\
174 Compute the master flat with cr2res_master_flat(avg, slit_model, \n\
175 --bpm_low, --bpm_high, --bpm_lines_ratio) \n\
176 -> master_flat, bpm \n\
177 Merge the bpm with the input bpm \n\
178 Compute QCs \n\
179 store the qc parameters in the returned property list \n\
180 \n\
181 Library functions used: \n\
182 cr2res_io_extract_decker_frameset() \n\
183 cr2res_trace() \n\
184 cr2res_extract_slitdec_curved() \n\
185 cr2res_master_flat() \n\
186 cr2res_qc_flat_lamp_ints() \n\
187 cr2res_qc_flat_mean_level() \n\
188 cr2res_qc_flat_mean_med_flux() \n\
189 cr2res_qc_flat_med_snr() \n\
190 cr2res_qc_overexposed() \n\
191 cr2res_qc_flat_trace_center_y() \n\
192 cr2res_trace_merge() \n\
193 cr2res_io_save_SLIT_MODEL() \n\
194 cr2res_io_save_EXTRACT_1D() \n\
195 cr2res_io_save_MASTER_FLAT() \n\
196 cr2res_io_save_TRACE_WAVE() \n\
197 cr2res_io_save_SLIT_FUNC() \n\
198 cr2res_io_save_BPM() \n\
199";
200
201/*-----------------------------------------------------------------------------
202 Function code
203 -----------------------------------------------------------------------------*/
204
205/*----------------------------------------------------------------------------*/
215/*----------------------------------------------------------------------------*/
216int cpl_plugin_get_info(cpl_pluginlist * list)
217{
218 cpl_recipe * recipe = cpl_calloc(1, sizeof *recipe );
219 cpl_plugin * plugin = &recipe->interface;
220
221 if (cpl_plugin_init(plugin,
222 CPL_PLUGIN_API,
223 CR2RES_BINARY_VERSION,
224 CPL_PLUGIN_TYPE_RECIPE,
225 RECIPE_STRING,
226 "Flat recipe",
227 cr2res_cal_flat_description,
228 CR2RES_PIPELINE_AUTHORS,
229 PACKAGE_BUGREPORT,
230 cr2res_get_license(),
231 cr2res_cal_flat_create,
232 cr2res_cal_flat_exec,
233 cr2res_cal_flat_destroy)) {
234 cpl_msg_error(cpl_func, "Plugin initialization failed");
235 (void)cpl_error_set_where(cpl_func);
236 return 1;
237 }
238
239 if (cpl_pluginlist_append(list, plugin)) {
240 cpl_msg_error(cpl_func, "Error adding plugin to list");
241 (void)cpl_error_set_where(cpl_func);
242 return 1;
243 }
244
245 return 0;
246}
247
248/*----------------------------------------------------------------------------*/
256/*----------------------------------------------------------------------------*/
257static int cr2res_cal_flat_create(cpl_plugin * plugin)
258{
259 cpl_recipe * recipe ;
260 cpl_parameter * p ;
261
262 /* Check that the plugin is part of a valid recipe */
263 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
264 recipe = (cpl_recipe *)plugin;
265 else
266 return -1;
267
268 /* Create the parameters list in the cpl_recipe object */
269 recipe->parameters = cpl_parameterlist_new();
270
271 /* Fill the parameters list */
272 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.subtract_nolight_rows",
273 CPL_TYPE_BOOL, "Subtract the no-light rows",
274 "cr2res.cr2res_cal_flat", FALSE);
275 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "subtract_nolight_rows");
276 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
277 cpl_parameterlist_append(recipe->parameters, p);
278
279 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.cosmics",
280 CPL_TYPE_BOOL, "Find and mark cosmic rays hits as bad",
281 "cr2res.cr2res_cal_flat", FALSE);
282 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "cosmics");
283 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
284 cpl_parameterlist_append(recipe->parameters, p);
285
286 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.bpm_low",
287 CPL_TYPE_DOUBLE, "Low threshold for BPM detection",
288 "cr2res.cr2res_cal_flat", 0.8);
289 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "bpm_low");
290 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
291 cpl_parameterlist_append(recipe->parameters, p);
292
293 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.bpm_high",
294 CPL_TYPE_DOUBLE, "High threshold for BPM detection",
295 "cr2res.cr2res_cal_flat", 1.2);
296 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "bpm_high");
297 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
298 cpl_parameterlist_append(recipe->parameters, p);
299
300 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.bpm_lines_ratio",
301 CPL_TYPE_DOUBLE, "Maximum ratio of bad pixels per line",
302 "cr2res.cr2res_cal_flat", 0.5);
303 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "bpm_lines_ratio");
304 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
305 cpl_parameterlist_append(recipe->parameters, p);
306
307 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.trace_degree",
308 CPL_TYPE_INT, "polynomial degree for the fit to the orders",
309 "cr2res.cr2res_cal_flat", 2);
310 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "trace_degree");
311 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
312 cpl_parameterlist_append(recipe->parameters, p);
313
314 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.trace_min_cluster",
315 CPL_TYPE_INT, "size in pixels of the smallest allowed cluster",
316 "cr2res.cr2res_cal_flat", 100000);
317 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "trace_min_cluster");
318 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
319 cpl_parameterlist_append(recipe->parameters, p);
320
321 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.trace_smooth_x",
322 CPL_TYPE_INT, "Length of the smoothing kernel in x",
323 "cr2res.cr2res_cal_flat", 111);
324 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "trace_smooth_x");
325 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
326 cpl_parameterlist_append(recipe->parameters, p);
327
328 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.trace_smooth_y",
329 CPL_TYPE_INT, "Length of the smoothing kernel in y",
330 "cr2res.cr2res_cal_flat", 401);
331 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "trace_smooth_y");
332 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
333 cpl_parameterlist_append(recipe->parameters, p);
334
335 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.trace_threshold",
336 CPL_TYPE_DOUBLE, "Detection Threshold",
337 "cr2res.cr2res_cal_flat", 1.0);
338 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "trace_threshold");
339 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
340 cpl_parameterlist_append(recipe->parameters, p);
341
342 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.trace_opening",
343 CPL_TYPE_BOOL, "Use a morphological opening to rejoin clusters",
344 "cr2res.cr2res_cal_flat", TRUE);
345 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "trace_opening");
346 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
347 cpl_parameterlist_append(recipe->parameters, p);
348
349 /*
350 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.trace_filter",
351 CPL_TYPE_BOOL, "Only keep the predefined order traces",
352 "cr2res.cr2res_cal_flat", TRUE);
353 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "trace_filter");
354 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
355 cpl_parameterlist_append(recipe->parameters, p);
356 */
357
358 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.extract_method",
359 CPL_TYPE_STRING, "Extraction method (SUM / MEDIAN / TILTSUM / "
360 "OPT_CURV )",
361 "cr2res.cr2res_util_extract", "OPT_CURV");
362 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "extract_method");
363 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
364 cpl_parameterlist_append(recipe->parameters, p);
365
366 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.extract_oversample",
367 CPL_TYPE_INT, "factor by which to oversample the extraction",
368 "cr2res.cr2res_cal_flat", 5);
369 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "extract_oversample");
370 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
371 cpl_parameterlist_append(recipe->parameters, p);
372
373 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.extract_swath_width",
374 CPL_TYPE_INT, "The swath width", "cr2res.cr2res_cal_flat", 800);
375 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "extract_swath_width");
376 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
377 cpl_parameterlist_append(recipe->parameters, p);
378
379 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.extract_height",
380 CPL_TYPE_INT, "Extraction height",
381 "cr2res.cr2res_cal_flat", -1);
382 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "extract_height");
383 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
384 cpl_parameterlist_append(recipe->parameters, p);
385
386 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.extract_smooth_slit",
387 CPL_TYPE_DOUBLE,
388 "Smoothing along the slit",
389 "cr2res.cr2res_cal_flat", 3.0);
390 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "extract_smooth_slit");
391 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
392 cpl_parameterlist_append(recipe->parameters, p);
393
394 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.extract_smooth_spec",
395 CPL_TYPE_DOUBLE,
396 "Smoothing along the spectrum",
397 "cr2res.cr2res_cal_flat", 2.0) ;
398 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "extract_smooth_spec");
399 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
400 cpl_parameterlist_append(recipe->parameters, p);
401
402 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.detector",
403 CPL_TYPE_INT, "Only reduce the specified detector",
404 "cr2res.cr2res_cal_flat", 0);
405 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "detector");
406 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
407 cpl_parameterlist_append(recipe->parameters, p);
408
409 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.order",
410 CPL_TYPE_INT, "Only reduce the specified order",
411 "cr2res.cr2res_cal_flat", -1);
412 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "order");
413 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
414 cpl_parameterlist_append(recipe->parameters, p);
415
416 p = cpl_parameter_new_value("cr2res.cr2res_cal_flat.trace_nb",
417 CPL_TYPE_INT, "Only reduce the specified trace number",
418 "cr2res.cr2res_cal_flat", -1);
419 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "trace_nb");
420 cpl_parameter_disable(p, CPL_PARAMETER_MODE_ENV);
421 cpl_parameterlist_append(recipe->parameters, p);
422
423 return 0;
424}
425
426/*----------------------------------------------------------------------------*/
432/*----------------------------------------------------------------------------*/
433static int cr2res_cal_flat_exec(cpl_plugin * plugin)
434{
435 cpl_recipe *recipe;
436
437 /* Get the recipe out of the plugin */
438 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
439 recipe = (cpl_recipe *)plugin;
440 else return -1;
441
442 return cr2res_cal_flat(recipe->frames, recipe->parameters);
443}
444
445/*----------------------------------------------------------------------------*/
451/*----------------------------------------------------------------------------*/
452static int cr2res_cal_flat_destroy(cpl_plugin * plugin)
453{
454 cpl_recipe *recipe;
455
456 /* Get the recipe out of the plugin */
457 if (cpl_plugin_get_type(plugin) == CPL_PLUGIN_TYPE_RECIPE)
458 recipe = (cpl_recipe *)plugin;
459 else return -1 ;
460
461 cpl_parameterlist_delete(recipe->parameters);
462 return 0 ;
463}
464
465/*----------------------------------------------------------------------------*/
472/*----------------------------------------------------------------------------*/
473static int cr2res_cal_flat(
474 cpl_frameset * frameset,
475 const cpl_parameterlist * parlist)
476{
477 const cpl_parameter * param ;
478 int subtract_nolight_rows, cosmics,
479 trace_degree, trace_min_cluster,
480 trace_opening, trace_filter,
481 extract_oversample, extract_swath_width,
482 extract_height, reduce_det, reduce_order,
483 reduce_trace, trace_smooth_x, trace_smooth_y ;
484 double bpm_low, bpm_high, bpm_lines_ratio,
485 trace_threshold, extract_smooth_slit,
486 extract_smooth_spec ;
487 cr2res_extr_method extr_method;
488 const char * sval ;
489 const cpl_frame * trace_wave_frame ;
490 const cpl_frame * detlin_frame ;
491 const cpl_frame * master_dark_frame ;
492 const cpl_frame * bpm_frame ;
493 cpl_frameset * rawframes ;
494 const char * used_tag ;
495 cpl_frameset * raw_one_setting_decker ;
496 cpl_size * labels ;
497 cpl_size nlabels ;
498 cpl_propertylist * qc_main ;
499
500 hdrl_image * master_flat[CR2RES_NB_DETECTORS] ;
501 cpl_table * trace_wave[CR2RES_NB_DECKER_POSITIONS][CR2RES_NB_DETECTORS] ;
502 cpl_propertylist *
503 ext_plist[CR2RES_NB_DECKER_POSITIONS][CR2RES_NB_DETECTORS] ;
504 cpl_table * trace_wave_merged[CR2RES_NB_DETECTORS];
505 cpl_table * merged ;
506 cpl_table * trace_wave1 ;
507 cpl_table * trace_wave2 ;
508 cpl_table * slit_func[CR2RES_NB_DETECTORS] ;
509 cpl_table * extract_1d[CR2RES_NB_DETECTORS] ;
510 hdrl_image * slit_model[CR2RES_NB_DETECTORS] ;
511 cpl_image * bpm[CR2RES_NB_DETECTORS] ;
512 char * qc_name;
513 char * out_file;
514 int * orders ;
515 cpl_vector * qc_mean,
516 * qc_median,
517 * qc_flux,
518 * qc_rms,
519 * qc_s2n,
520 * qc_nbbad,
521 * qc_centery,
522 * qc_bltot,
523 * qc_blgood,
524 * qc_orderpos,
525 * qc_overexposed ;
526 int l, i, j, det_nr, nb_orders;
527
528 /* Initialise */
529 cr2res_decker decker_values[CR2RES_NB_DECKER_POSITIONS] =
530 {CR2RES_DECKER_NONE, CR2RES_DECKER_1_3, CR2RES_DECKER_2_4} ;
531 char * decker_desc[CR2RES_NB_DECKER_POSITIONS] =
532 {"Open", "Decker1", "Decker2"} ;
533
534 /* RETRIEVE INPUT PARAMETERS */
535 param = cpl_parameterlist_find_const(parlist,
536 "cr2res.cr2res_cal_flat.subtract_nolight_rows");
537 subtract_nolight_rows = cpl_parameter_get_bool(param);
538 param = cpl_parameterlist_find_const(parlist,
539 "cr2res.cr2res_cal_flat.cosmics");
540 cosmics = cpl_parameter_get_bool(param);
541 param = cpl_parameterlist_find_const(parlist,
542 "cr2res.cr2res_cal_flat.bpm_low");
543 bpm_low = cpl_parameter_get_double(param);
544 param = cpl_parameterlist_find_const(parlist,
545 "cr2res.cr2res_cal_flat.bpm_high");
546 bpm_high = cpl_parameter_get_double(param);
547 param = cpl_parameterlist_find_const(parlist,
548 "cr2res.cr2res_cal_flat.bpm_lines_ratio");
549 bpm_lines_ratio = cpl_parameter_get_double(param);
550 param = cpl_parameterlist_find_const(parlist,
551 "cr2res.cr2res_cal_flat.trace_degree");
552 trace_degree = cpl_parameter_get_int(param);
553 param = cpl_parameterlist_find_const(parlist,
554 "cr2res.cr2res_cal_flat.trace_min_cluster");
555 trace_min_cluster = cpl_parameter_get_int(param);
556 param = cpl_parameterlist_find_const(parlist,
557 "cr2res.cr2res_cal_flat.trace_smooth_x");
558 trace_smooth_x = cpl_parameter_get_int(param);
559 param = cpl_parameterlist_find_const(parlist,
560 "cr2res.cr2res_cal_flat.trace_smooth_y");
561 trace_smooth_y = cpl_parameter_get_int(param);
562 param = cpl_parameterlist_find_const(parlist,
563 "cr2res.cr2res_cal_flat.trace_threshold");
564 trace_threshold = cpl_parameter_get_double(param);
565 param = cpl_parameterlist_find_const(parlist,
566 "cr2res.cr2res_cal_flat.trace_opening");
567 trace_opening = cpl_parameter_get_bool(param);
568 /*
569 param = cpl_parameterlist_find_const(parlist,
570 "cr2res.cr2res_cal_flat.trace_filter");
571 trace_filter = cpl_parameter_get_bool(param);
572 */
573 trace_filter = 0 ;
574 param = cpl_parameterlist_find_const(parlist,
575 "cr2res.cr2res_cal_flat.extract_method");
576 sval = cpl_parameter_get_string(param);
577 if (!strcmp(sval, "")) extr_method = CR2RES_EXTR_OPT_CURV;
578 else if (!strcmp(sval, "OPT_CURV")) extr_method = CR2RES_EXTR_OPT_CURV;
579 else if (!strcmp(sval, "SUM")) extr_method = CR2RES_EXTR_SUM;
580 else if (!strcmp(sval, "MEDIAN")) extr_method = CR2RES_EXTR_MEDIAN;
581 else if (!strcmp(sval, "TILTSUM")) extr_method = CR2RES_EXTR_TILTSUM;
582 else {
583 cpl_msg_error(__func__, "Invalid Extraction Method specified");
584 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ;
585 return -1;
586 }
587 param = cpl_parameterlist_find_const(parlist,
588 "cr2res.cr2res_cal_flat.extract_oversample");
589 extract_oversample = cpl_parameter_get_int(param);
590 param = cpl_parameterlist_find_const(parlist,
591 "cr2res.cr2res_cal_flat.extract_swath_width");
592 extract_swath_width = cpl_parameter_get_int(param);
593 param = cpl_parameterlist_find_const(parlist,
594 "cr2res.cr2res_cal_flat.extract_height");
595 extract_height = cpl_parameter_get_int(param);
596 param = cpl_parameterlist_find_const(parlist,
597 "cr2res.cr2res_cal_flat.extract_smooth_slit");
598 extract_smooth_slit = cpl_parameter_get_double(param);
599 param = cpl_parameterlist_find_const(parlist,
600 "cr2res.cr2res_cal_flat.extract_smooth_spec");
601 extract_smooth_spec = cpl_parameter_get_double(param);
602 param = cpl_parameterlist_find_const(parlist,
603 "cr2res.cr2res_cal_flat.detector");
604 reduce_det = cpl_parameter_get_int(param);
605 param = cpl_parameterlist_find_const(parlist,
606 "cr2res.cr2res_cal_flat.order");
607 reduce_order = cpl_parameter_get_int(param);
608 param = cpl_parameterlist_find_const(parlist,
609 "cr2res.cr2res_cal_flat.trace_nb");
610 reduce_trace = cpl_parameter_get_int(param);
611
612 /* Identify the RAW and CALIB frames in the input frameset */
613 if (cr2res_dfs_set_groups(frameset)) {
614 cpl_msg_error(__func__, "Cannot identify RAW and CALIB frames") ;
615 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ;
616 return -1 ;
617 }
618
619 /* Get Calibration frames */
620 trace_wave_frame = cr2res_io_find_TRACE_WAVE(frameset) ;
621 detlin_frame = cpl_frameset_find_const(frameset,
622 CR2RES_CAL_DETLIN_COEFFS_PROCATG);
623 master_dark_frame = cpl_frameset_find_const(frameset,
624 CR2RES_CAL_DARK_MASTER_PROCATG) ;
625 bpm_frame = cr2res_io_find_BPM(frameset) ;
626
627 /* Extract RAW frames */
628 used_tag = CR2RES_FLAT_RAW ;
629 rawframes = cr2res_extract_frameset(frameset, used_tag) ;
630 if (rawframes==NULL || cpl_frameset_get_size(rawframes) <= 0) {
631 cpl_msg_error(__func__, "Cannot find any RAW file") ;
632 cpl_error_set(__func__, CPL_ERROR_DATA_NOT_FOUND) ;
633 return -1 ;
634 }
635
636 /* Label the raw frames with the different settings */
637 if ((labels = cpl_frameset_labelise(rawframes, cr2res_cal_flat_compare,
638 &nlabels)) == NULL) {
639 cpl_msg_error(__func__, "Cannot label input frames") ;
640 cpl_frameset_delete(rawframes) ;
641 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT) ;
642 return -1 ;
643 }
644
645 /* Loop on the settings */
646 for (l = 0; l < (int)nlabels; l++) {
647 cpl_frameset *raw_one_setting;
648 cpl_propertylist *plist;
649 char *setting_id;
650 /* Get the frames for the current setting */
651 raw_one_setting = cpl_frameset_extract(rawframes, labels, (cpl_size)l) ;
652
653 /* Get the current setting */
654 plist = cpl_propertylist_load(cpl_frame_get_filename(
655 cpl_frameset_get_position(raw_one_setting, 0)), 0) ;
656 setting_id = cpl_strdup(cr2res_pfits_get_wlen_id(plist)) ;
657 cr2res_format_setting(setting_id) ;
658 cpl_propertylist_delete(plist) ;
659
660 cpl_msg_info(__func__, "Process SETTING %s", setting_id) ;
661 cpl_msg_indent_more() ;
662
663 /* Loop on the decker positions */
664 for (i=0 ; i<CR2RES_NB_DECKER_POSITIONS ; i++) {
665 /* Initialise */
666 for (det_nr=1 ; det_nr<=CR2RES_NB_DETECTORS ; det_nr++) {
667 trace_wave[i][det_nr-1] = NULL ;
668 ext_plist[i][det_nr-1] = NULL ;
669 }
670
671 /* Get the Frames for the current decker position */
672 raw_one_setting_decker = cr2res_io_extract_decker_frameset(
673 raw_one_setting, used_tag, decker_values[i]) ;
674 if (raw_one_setting_decker == NULL) {
675 cpl_msg_info(__func__, "No files for decker: %s",
676 decker_desc[i]) ;
677 continue ;
678 }
679 cpl_msg_info(__func__, "Reduce %s Frames", decker_desc[i]) ;
680 cpl_msg_indent_more() ;
681
682 /* Loop on the detectors */
683 for (det_nr=1 ; det_nr<=CR2RES_NB_DETECTORS ; det_nr++) {
684 /* Initialise */
685 master_flat[det_nr-1] = NULL ;
686 slit_func[det_nr-1] = NULL ;
687 extract_1d[det_nr-1] = NULL ;
688 slit_model[det_nr-1] = NULL ;
689 bpm[det_nr-1] = NULL ;
690
691 /* Compute only one detector */
692 if (reduce_det != 0 && det_nr != reduce_det) continue ;
693
694 cpl_msg_info(__func__, "Process Detector %d", det_nr) ;
695 cpl_msg_indent_more() ;
696
697 /* Call the reduction function */
698 if (cr2res_cal_flat_reduce(raw_one_setting_decker,
699 trace_wave_frame, detlin_frame, master_dark_frame,
700 bpm_frame, trace_filter, subtract_nolight_rows,
701 cosmics, bpm_low, bpm_high,
702 bpm_lines_ratio, trace_degree, trace_min_cluster,
703 trace_smooth_x, trace_smooth_y, trace_threshold,
704 trace_opening, extr_method, extract_oversample,
705 extract_swath_width, extract_height,
706 extract_smooth_slit, extract_smooth_spec, det_nr,
707 reduce_order, reduce_trace,
708 &(master_flat[det_nr-1]),
709 &(trace_wave[i][det_nr-1]),
710 &(slit_func[det_nr-1]),
711 &(extract_1d[det_nr-1]),
712 &(slit_model[det_nr-1]),
713 &(bpm[det_nr-1]),
714 &(ext_plist[i][det_nr-1])) == -1) {
715 cpl_msg_warning(__func__,
716 "Failed to reduce detector %d of %s Frames",
717 det_nr, decker_desc[i]);
718 cpl_error_reset() ;
719 }
720 cpl_msg_indent_less() ;
721 }
722
723
724 /* Compute QCs for the main header */
725 qc_main = NULL ;
726 if (reduce_det == 0) {
727 qc_mean = cpl_vector_new(CR2RES_NB_DETECTORS) ;
728 qc_median = cpl_vector_new(CR2RES_NB_DETECTORS) ;
729 qc_flux = cpl_vector_new(CR2RES_NB_DETECTORS) ;
730 qc_rms = cpl_vector_new(CR2RES_NB_DETECTORS) ;
731 qc_s2n = cpl_vector_new(CR2RES_NB_DETECTORS) ;
732 qc_nbbad = cpl_vector_new(CR2RES_NB_DETECTORS) ;
733 qc_centery = cpl_vector_new(CR2RES_NB_DETECTORS) ;
734 qc_blgood = cpl_vector_new(CR2RES_NB_DETECTORS) ;
735 qc_bltot = cpl_vector_new(CR2RES_NB_DETECTORS) ;
736 for (det_nr=1 ; det_nr<=CR2RES_NB_DETECTORS ; det_nr++) {
737 cpl_vector_set(qc_mean, det_nr-1,
738 cpl_propertylist_get_double(ext_plist[i][det_nr-1],
739 CR2RES_HEADER_QC_FLAT_MEAN)) ;
740 cpl_vector_set(qc_median, det_nr-1,
741 cpl_propertylist_get_double(ext_plist[i][det_nr-1],
742 CR2RES_HEADER_QC_FLAT_MEDIAN)) ;
743 cpl_vector_set(qc_flux, det_nr-1,
744 cpl_propertylist_get_double(ext_plist[i][det_nr-1],
745 CR2RES_HEADER_QC_FLAT_FLUX)) ;
746 cpl_vector_set(qc_rms, det_nr-1,
747 cpl_propertylist_get_double(ext_plist[i][det_nr-1],
748 CR2RES_HEADER_QC_FLAT_RMS)) ;
749 cpl_vector_set(qc_s2n, det_nr-1,
750 cpl_propertylist_get_double(ext_plist[i][det_nr-1],
751 CR2RES_HEADER_QC_FLAT_S2N)) ;
752 cpl_vector_set(qc_nbbad, det_nr-1,
753 (double)cpl_propertylist_get_int(ext_plist[i][det_nr-1],
754 CR2RES_HEADER_QC_FLAT_NBBAD)) ;
755 cpl_vector_set(qc_centery, det_nr-1,
756 cpl_propertylist_get_double(ext_plist[i][det_nr-1],
757 CR2RES_HEADER_QC_FLAT_CENTERY)) ;
758 cpl_vector_set(qc_blgood, det_nr-1,
759 (double)cpl_propertylist_get_int(ext_plist[i][det_nr-1],
760 CR2RES_HEADER_QC_BLAZE_NGOOD)) ;
761 cpl_vector_set(qc_bltot, det_nr-1,
762 cpl_propertylist_get_double(ext_plist[i][det_nr-1],
763 CR2RES_HEADER_QC_BLAZE_TOT)) ;
764 }
765
766 qc_main = cpl_propertylist_new() ;
767 cpl_propertylist_append_double(qc_main,
768 CR2RES_HEADER_QC_FLAT_MEAN_AVG,
769 cpl_vector_get_mean(qc_mean)) ;
770 cpl_propertylist_append_double(qc_main,
771 CR2RES_HEADER_QC_FLAT_MEAN_RMS,
772 cpl_vector_get_stdev(qc_mean)) ;
773 cpl_propertylist_append_double(qc_main,
774 CR2RES_HEADER_QC_FLAT_MEDIAN_AVG,
775 cpl_vector_get_mean(qc_median)) ;
776 cpl_propertylist_append_double(qc_main,
777 CR2RES_HEADER_QC_FLAT_MEDIAN_RMS,
778 cpl_vector_get_stdev(qc_median)) ;
779 cpl_propertylist_append_double(qc_main,
780 CR2RES_HEADER_QC_FLAT_FLUX_AVG,
781 cpl_vector_get_mean(qc_flux)) ;
782 cpl_propertylist_append_double(qc_main,
783 CR2RES_HEADER_QC_FLAT_FLUX_RMS,
784 cpl_vector_get_stdev(qc_flux)) ;
785 cpl_propertylist_append_double(qc_main,
786 CR2RES_HEADER_QC_FLAT_RMS_AVG,
787 cpl_vector_get_mean(qc_rms)) ;
788 cpl_propertylist_append_double(qc_main,
789 CR2RES_HEADER_QC_FLAT_RMS_RMS,
790 cpl_vector_get_stdev(qc_rms)) ;
791 cpl_propertylist_append_double(qc_main,
792 CR2RES_HEADER_QC_FLAT_S2N_AVG,
793 cpl_vector_get_mean(qc_s2n)) ;
794 cpl_propertylist_append_double(qc_main,
795 CR2RES_HEADER_QC_FLAT_S2N_RMS,
796 cpl_vector_get_stdev(qc_s2n)) ;
797 cpl_propertylist_append_double(qc_main,
798 CR2RES_HEADER_QC_FLAT_NBBAD_AVG,
799 cpl_vector_get_mean(qc_nbbad)) ;
800 cpl_propertylist_append_double(qc_main,
801 CR2RES_HEADER_QC_FLAT_NBBAD_RMS,
802 cpl_vector_get_stdev(qc_nbbad)) ;
803 cpl_propertylist_append_double(qc_main,
804 CR2RES_HEADER_QC_FLAT_CENTERY_AVG,
805 cpl_vector_get_mean(qc_centery)) ;
806 cpl_propertylist_append_double(qc_main,
807 CR2RES_HEADER_QC_FLAT_CENTERY_RMS,
808 cpl_vector_get_stdev(qc_centery)) ;
809 cpl_propertylist_append_double(qc_main,
810 CR2RES_HEADER_QC_BLAZE_NORM,
811 cpl_vector_get_sum(qc_bltot)/cpl_vector_get_sum(qc_blgood)) ;
812 cpl_vector_delete(qc_mean) ;
813 cpl_vector_delete(qc_median) ;
814 cpl_vector_delete(qc_flux) ;
815 cpl_vector_delete(qc_rms) ;
816 cpl_vector_delete(qc_s2n) ;
817 cpl_vector_delete(qc_nbbad) ;
818 cpl_vector_delete(qc_centery) ;
819 cpl_vector_delete(qc_blgood) ;
820 cpl_vector_delete(qc_bltot) ;
821
822 if (trace_wave[i][0] != NULL) {
823 /* Special treatment for ORDERPOSn and OVEREXPOSEDn */
824 /* Get all orders - Use the TW from detector 1 */
825 orders = cr2res_trace_get_order_idx_values(
826 trace_wave[i][0], &nb_orders) ;
827 /* Loop on the orders */
828 for (j=0 ; j<nb_orders ; j++) {
829 qc_orderpos = cpl_vector_new(CR2RES_NB_DETECTORS) ;
830 qc_overexposed = cpl_vector_new(CR2RES_NB_DETECTORS) ;
831 cpl_vector_fill(qc_orderpos, 0.0) ;
832 cpl_vector_fill(qc_overexposed, 0.0) ;
833 for (det_nr = 1; det_nr <= CR2RES_NB_DETECTORS; det_nr++)
834 {
835 // OVEREXPOSEDn
836 qc_name = cpl_sprintf("%s%02d",
837 CR2RES_HEADER_QC_OVEREXPOSED, orders[j]);
838 if (cpl_propertylist_has(ext_plist[i][det_nr - 1], qc_name))
839 {
840 cpl_vector_set(qc_overexposed, det_nr - 1,
841 cpl_propertylist_get_double(
842 ext_plist[i][det_nr - 1], qc_name));
843 }
844 cpl_free(qc_name);
845 // ORDERPOSn
846 qc_name = cpl_sprintf("%s%02d",
847 CR2RES_HEADER_QC_FLAT_ORDERPOS, orders[j]);
848 if (cpl_propertylist_has(ext_plist[i][det_nr - 1], qc_name))
849 {
850 cpl_vector_set(qc_orderpos, det_nr - 1,
851 cpl_propertylist_get_double(
852 ext_plist[i][det_nr - 1], qc_name));
853 }
854 cpl_free(qc_name);
855 }
856
857 // OVEREXPOSEDn.AVG/RMS
858 qc_name = cpl_sprintf("%s%02d AVG",
859 CR2RES_HEADER_QC_OVEREXPOSED, orders[j]) ;
860 cpl_propertylist_append_double(qc_main,
861 qc_name, cpl_vector_get_mean(qc_overexposed)) ;
862 cpl_free(qc_name) ;
863
864 qc_name = cpl_sprintf("%s%02d RMS",
865 CR2RES_HEADER_QC_OVEREXPOSED, orders[j]) ;
866 cpl_propertylist_append_double(qc_main,
867 qc_name, cpl_vector_get_stdev(qc_overexposed)) ;
868 cpl_free(qc_name) ;
869 cpl_vector_delete(qc_overexposed) ;
870
871 // ORDERPOSn.AVG/RMS
872 qc_name = cpl_sprintf("%s%02d AVG",
873 CR2RES_HEADER_QC_FLAT_ORDERPOS, orders[j]) ;
874 cpl_propertylist_append_double(qc_main,
875 qc_name, cpl_vector_get_mean(qc_orderpos)) ;
876 cpl_free(qc_name) ;
877
878 qc_name = cpl_sprintf("%s%02d RMS",
879 CR2RES_HEADER_QC_FLAT_ORDERPOS, orders[j]) ;
880 cpl_propertylist_append_double(qc_main,
881 qc_name, cpl_vector_get_stdev(qc_orderpos)) ;
882 cpl_free(qc_name) ;
883 cpl_vector_delete(qc_orderpos) ;
884 }
885 cpl_free(orders) ;
886 }
887 }
888
889 /* Save Products */
890
891 /* Save only the used RAW ? : raw_one_... instead of 2nd frameset */
892 /* Beware that the calibration PRO RECi CAL will be missing */
893
894 /* MASTER_FLAT */
895 if (nlabels == 1) {
896 out_file = cpl_sprintf("%s_%s_master_flat.fits",
897 RECIPE_STRING, decker_desc[i]) ;
898 } else {
899 out_file = cpl_sprintf("%s_%s_%s_master_flat.fits",
900 RECIPE_STRING, setting_id, decker_desc[i]) ;
901 }
902 cr2res_io_save_MASTER_FLAT(out_file, frameset,
903 frameset, parlist, master_flat, qc_main,
904 ext_plist[i], CR2RES_CAL_FLAT_MASTER_PROCATG,
905 RECIPE_STRING);
906 cpl_free(out_file);
907
908 /* SLIT_MODEL */
909 if (nlabels == 1) {
910 out_file = cpl_sprintf("%s_%s_slit_model.fits",
911 RECIPE_STRING, decker_desc[i]) ;
912 } else {
913 out_file = cpl_sprintf("%s_%s_%s_slit_model.fits",
914 RECIPE_STRING, setting_id, decker_desc[i]) ;
915 }
916 cr2res_io_save_SLIT_MODEL(out_file, frameset,
917 frameset, parlist, slit_model, qc_main,
918 ext_plist[i], CR2RES_CAL_FLAT_SLIT_MODEL_PROCATG,
919 RECIPE_STRING);
920 cpl_free(out_file);
921
922 /* BLAZE */
923 if (nlabels == 1) {
924 out_file = cpl_sprintf("%s_%s_blaze.fits",
925 RECIPE_STRING, decker_desc[i]) ;
926 } else {
927 out_file = cpl_sprintf("%s_%s_%s_blaze.fits",
928 RECIPE_STRING, setting_id, decker_desc[i]) ;
929 }
930 cr2res_io_save_EXTRACT_1D(out_file, frameset,
931 frameset, parlist, extract_1d, qc_main,
932 ext_plist[i], CR2RES_CAL_FLAT_EXTRACT_1D_PROCATG,
933 RECIPE_STRING);
934 cpl_free(out_file);
935
936 /* TRACE_WAVE */
937 if (nlabels == 1) {
938 out_file = cpl_sprintf("%s_%s_tw.fits",
939 RECIPE_STRING, decker_desc[i]) ;
940 } else {
941 out_file = cpl_sprintf("%s_%s_%s_tw.fits",
942 RECIPE_STRING, setting_id, decker_desc[i]) ;
943 }
944 cr2res_io_save_TRACE_WAVE(out_file, frameset,
945 frameset, parlist, trace_wave[i], qc_main,
946 ext_plist[i], CR2RES_CAL_FLAT_TW_PROCATG, RECIPE_STRING);
947 cpl_free(out_file);
948
949 /* SLIT_FUNC */
950 if (nlabels == 1) {
951 out_file = cpl_sprintf("%s_%s_slit_func.fits",
952 RECIPE_STRING, decker_desc[i]) ;
953 } else {
954 out_file = cpl_sprintf("%s_%s_%s_slit_func.fits",
955 RECIPE_STRING, setting_id, decker_desc[i]) ;
956 }
957 cr2res_io_save_SLIT_FUNC(out_file, frameset,
958 frameset, parlist, slit_func, qc_main,
959 ext_plist[i], CR2RES_CAL_FLAT_SLIT_FUNC_PROCATG,
960 RECIPE_STRING);
961 cpl_free(out_file);
962
963 /* BPM */
964 if (nlabels == 1) {
965 out_file = cpl_sprintf("%s_%s_bpm.fits",
966 RECIPE_STRING, decker_desc[i]) ;
967 } else {
968 out_file = cpl_sprintf("%s_%s_%s_bpm.fits",
969 RECIPE_STRING, setting_id, decker_desc[i]) ;
970 }
971 cr2res_io_save_BPM(out_file, frameset,
972 frameset, parlist, bpm, qc_main, ext_plist[i],
973 CR2RES_CAL_FLAT_BPM_PROCATG, RECIPE_STRING) ;
974 cpl_free(out_file);
975
976 /* Free */
977 if (qc_main != NULL) cpl_propertylist_delete(qc_main) ;
978 cpl_frameset_delete(raw_one_setting_decker) ;
979 for (det_nr=1 ; det_nr<=CR2RES_NB_DETECTORS ; det_nr++) {
980 if (master_flat[det_nr-1] != NULL)
981 hdrl_image_delete(master_flat[det_nr-1]) ;
982 if (slit_func[det_nr-1] != NULL)
983 cpl_table_delete(slit_func[det_nr-1]) ;
984 if (extract_1d[det_nr-1] != NULL)
985 cpl_table_delete(extract_1d[det_nr-1]) ;
986 if (slit_model[det_nr-1] != NULL)
987 hdrl_image_delete(slit_model[det_nr-1]) ;
988 if (bpm[det_nr-1] != NULL)
989 cpl_image_delete(bpm[det_nr-1]) ;
990 }
991 cpl_msg_indent_less() ;
992 }
993
994 /* Merge the Decker positions TRACE_WAVE files in a single one */
995 for (det_nr=1 ; det_nr<=CR2RES_NB_DETECTORS ; det_nr++) {
996 /* Initialise */
997 if (trace_wave[0][det_nr-1] != NULL)
998 trace_wave_merged[det_nr-1] =
999 cpl_table_duplicate(trace_wave[0][det_nr-1]) ;
1000 else
1001 trace_wave_merged[det_nr-1] = NULL ;
1002 /* Loop on the other detectors */
1003 for (i=1 ; i<CR2RES_NB_DECKER_POSITIONS ; i++) {
1004 trace_wave1 = trace_wave_merged[det_nr-1] ;
1005 trace_wave2 = trace_wave[i][det_nr-1] ;
1006 if (trace_wave1 == NULL && trace_wave2 == NULL) {
1007 /* Do nothing - go to next iteration */
1008 trace_wave_merged[det_nr-1] = NULL ;
1009 } else if (trace_wave1 == NULL && trace_wave2 != NULL) {
1010 trace_wave_merged[det_nr-1] =
1011 cpl_table_duplicate(trace_wave2) ;
1012 } else if (trace_wave1 != NULL && trace_wave2 == NULL) {
1013 /* Do nothing - go to next iteration */
1014 } else {
1015 merged = cr2res_trace_merge(trace_wave1, trace_wave2) ;
1016 if (merged == NULL) {
1017 cpl_msg_error(__func__, "Failed merging") ;
1018 } else {
1019 cpl_table_delete(trace_wave_merged[det_nr-1]) ;
1020 trace_wave_merged[det_nr-1] = merged ;
1021 }
1022 }
1023 }
1024 }
1025
1026 /* Free the traces */
1027 for (i=0 ; i<CR2RES_NB_DECKER_POSITIONS ; i++)
1028 for (det_nr=1 ; det_nr<=CR2RES_NB_DETECTORS ; det_nr++)
1029 if (trace_wave[i][det_nr-1] != NULL)
1030 cpl_table_delete(trace_wave[i][det_nr-1]) ;
1031
1032
1033 /* Save TRACE_WAVE_MERGED */
1034 if (nlabels == 1) {
1035 out_file = cpl_sprintf("%s_tw_merged.fits",
1036 RECIPE_STRING) ;
1037 } else {
1038 out_file = cpl_sprintf("%s_%s_tw_merged.fits", RECIPE_STRING,
1039 setting_id) ;
1040 }
1041
1042 /* Save only the used RAW ? : raw_one_... instead of 2nd frameset */
1043 /* Beware that the calibration PRO RECi CAL will be missing */
1044
1045 cr2res_io_save_TRACE_WAVE(out_file, frameset, frameset, parlist,
1046 trace_wave_merged, NULL, ext_plist[0],
1047 CR2RES_CAL_FLAT_TW_MERGED_PROCATG, RECIPE_STRING) ;
1048 cpl_free(out_file);
1049 cpl_free(setting_id) ;
1050
1051 /* Free extensions */
1052 for (i=0 ; i<CR2RES_NB_DECKER_POSITIONS ; i++)
1053 for (det_nr=1 ; det_nr<=CR2RES_NB_DETECTORS ; det_nr++)
1054 if (ext_plist[i][det_nr-1] != NULL)
1055 cpl_propertylist_delete(ext_plist[i][det_nr-1]) ;
1056
1057 /* Free the merged traces */
1058 for (det_nr=1 ; det_nr<=CR2RES_NB_DETECTORS ; det_nr++)
1059 if (trace_wave_merged[det_nr-1] != NULL)
1060 cpl_table_delete(trace_wave_merged[det_nr-1]) ;
1061
1062 cpl_frameset_delete(raw_one_setting) ;
1063 cpl_msg_indent_less() ;
1064 }
1065 cpl_free(labels);
1066 cpl_frameset_delete(rawframes) ;
1067
1068 return (int)cpl_error_get_code();
1069}
1070
1071/*----------------------------------------------------------------------------*/
1109/*----------------------------------------------------------------------------*/
1110static int cr2res_cal_flat_reduce(
1111 const cpl_frameset * rawframes,
1112 const cpl_frame * tw_frame,
1113 const cpl_frame * detlin_frame,
1114 const cpl_frame * master_dark_frame,
1115 const cpl_frame * bpm_frame,
1116 int filter_traces,
1117 int subtract_nolight_rows,
1118 int cosmics,
1119 double bpm_low,
1120 double bpm_high,
1121 double bpm_linemax,
1122 int trace_degree,
1123 int trace_min_cluster,
1124 int trace_smooth_x,
1125 int trace_smooth_y,
1126 double trace_threshold,
1127 int trace_opening,
1128 cr2res_extr_method extr_method,
1129 int extract_oversample,
1130 int extract_swath_width,
1131 int extract_height,
1132 double extract_smooth_slit,
1133 double extract_smooth_spec,
1134 int reduce_det,
1135 int reduce_order,
1136 int reduce_trace,
1137 hdrl_image ** master_flat,
1138 cpl_table ** trace_wave,
1139 cpl_table ** slit_func,
1140 cpl_table ** extract_1d,
1141 hdrl_image ** slit_model,
1142 cpl_image ** bpm,
1143 cpl_propertylist ** ext_plist)
1144{
1145 const char * first_file ;
1146 hdrl_image * first_image ;
1147 hdrl_imagelist * imlist ;
1148 hdrl_imagelist * imlist_calibrated ;
1149 hdrl_image * collapsed ;
1150 cpl_image * contrib ;
1151 cpl_propertylist * plist ;
1152 hdrl_image * master_flat_loc ;
1153 cpl_image * my_master_flat ;
1154 cpl_image * bpm_im ;
1155 cpl_image * bpm_flat ;
1156 cpl_mask * bpm_mask ;
1157 cpl_table * computed_traces ;
1158 cpl_table * filtered_traces ;
1159 cpl_table * traces ;
1160 cpl_bivector ** spectrum ;
1161 cpl_vector ** slit_func_vec ;
1162 hdrl_image * model_master;
1163 cpl_table * slit_func_tab ;
1164 cpl_table * extract_tab ;
1165 hdrl_image * model_tmp ;
1166 cpl_vector * dits ;
1167 cpl_vector * ndits ;
1168 cpl_vector * ndsamples ;
1169 int * qc_order_nb ;
1170 double * qc_order_pos ;
1171 int * orders ;
1172 double qc_mean, qc_median, qc_flux, qc_rms, qc_s2n,
1173 qc_trace_centery, dit,
1174 gain, error_factor, bl_tot ;
1175 int i, ext_nr, nb_traces,
1176 nb_orders, qc_nbbad, nbvals, ngood, bl_good ;
1177
1178 /* TODO, make parameters */
1179 int extract_niter = 10;
1180 double extract_kappa = 10;
1181
1182 /* Check Inputs */
1183 if (rawframes == NULL) return -1 ;
1184 if (extr_method != CR2RES_EXTR_OPT_CURV &&
1185 extr_method != CR2RES_EXTR_SUM) {
1186 cpl_msg_error(__func__, "Failed to read the dits") ;
1187 return -1 ;
1188 }
1189
1190 /* Get the Gain */
1191 if (reduce_det == 1) gain = CR2RES_GAIN_CHIP1 ;
1192 else if (reduce_det == 2) gain = CR2RES_GAIN_CHIP2 ;
1193 else if (reduce_det == 3) gain = CR2RES_GAIN_CHIP3 ;
1194 else {
1195 cpl_msg_error(__func__, "Failed to get the Gain value") ;
1196 return -1 ;
1197 }
1198
1199 /* Get the First RAW file */
1200 first_file = cpl_frame_get_filename(
1201 cpl_frameset_get_position_const(rawframes, 0)) ;
1202
1203 /* Get the DIT for the Dark correction */
1204 if ((dits = cr2res_io_read_dits(rawframes)) == NULL) {
1205 cpl_msg_error(__func__, "Failed to read the dits") ;
1206 return -1 ;
1207 }
1208 /*Load the NDITs */
1209 ndits = cr2res_io_read_ndits(rawframes);
1210
1211 /* Load the image list */
1212 imlist = cr2res_io_load_image_list_from_set(rawframes, reduce_det) ;
1213 if (imlist == NULL) {
1214 cpl_msg_error(__func__, "Failed to load the images") ;
1215 cpl_vector_delete(dits) ;
1216 cpl_vector_delete(ndits) ;
1217 return -1 ;
1218 }
1219
1220 /* Set the error factor. */
1221 error_factor = gain * cpl_vector_get(ndits, 0) *
1222 cpl_frameset_get_size(rawframes) ;
1223
1224 for (i=0; i<cpl_vector_get_size(ndits); i++){
1225 if (cpl_vector_get(ndits,i) != cpl_vector_get(ndits, 0))
1226 cpl_msg_warning(__func__, "Raw frames have different NDIT! "
1227 "Error spectrum will likely be scaled incorrectly.");
1228 }
1229
1230 /* Compute traces */
1231 cpl_msg_info(__func__, "Compute the traces") ;
1232 cpl_msg_indent_more() ;
1233 if ((computed_traces = cr2res_trace(
1234 hdrl_image_get_image(hdrl_imagelist_get(imlist, 0)),
1235 trace_smooth_x, trace_smooth_y, trace_threshold,
1236 trace_opening, trace_degree, trace_min_cluster)) == NULL) {
1237 cpl_msg_error(__func__, "Failed compute the traces") ;
1238 cpl_vector_delete(dits) ;
1239 cpl_vector_delete(ndits) ;
1240 hdrl_imagelist_delete(imlist) ;
1241 cpl_msg_indent_less() ;
1242 return -1 ;
1243 }
1244 cpl_msg_indent_less() ;
1245
1246 /* Calibrate the Data */
1247 cpl_msg_info(__func__, "Calibrate the input images") ;
1248 cpl_msg_indent_more() ;
1249 ndsamples = cr2res_io_read_ndsamples(rawframes) ;
1250 imlist_calibrated = cr2res_calib_imagelist(imlist, reduce_det,
1251 0, subtract_nolight_rows, 0, cosmics, NULL,
1252 master_dark_frame, bpm_frame, detlin_frame, dits, ndits,
1253 ndsamples) ;
1254 cpl_vector_delete(ndsamples) ;
1255 if (imlist_calibrated == NULL) {
1256 cpl_msg_error(__func__, "Failed to Calibrate the Data") ;
1257 cpl_vector_delete(dits) ;
1258 cpl_vector_delete(ndits) ;
1259 hdrl_imagelist_delete(imlist) ;
1260 cpl_table_delete(computed_traces) ;
1261 cpl_msg_indent_less() ;
1262 return -1 ;
1263 } else {
1264 /* Replace the calibrated image in the list */
1265 hdrl_imagelist_delete(imlist) ;
1266 imlist = imlist_calibrated ;
1267 }
1268 cpl_vector_delete(dits) ;
1269 cpl_vector_delete(ndits) ;
1270 cpl_msg_indent_less() ;
1271
1272 /* Collapse */
1273 cpl_msg_info(__func__, "Collapse the input images") ;
1274 cpl_msg_indent_more() ;
1275 if (hdrl_imagelist_collapse_mean(imlist, &collapsed, &contrib) !=
1276 CPL_ERROR_NONE) {
1277 cpl_msg_error(__func__, "Failed to Collapse") ;
1278 hdrl_imagelist_delete(imlist) ;
1279 cpl_table_delete(computed_traces) ;
1280 cpl_msg_indent_less() ;
1281 return -1 ;
1282 }
1283 hdrl_imagelist_delete(imlist) ;
1284 cpl_image_delete(contrib) ;
1285 cpl_msg_indent_less() ;
1286
1287 /* Add The remaining Columns to the trace table */
1288 cr2res_trace_add_extra_columns(computed_traces, first_file, reduce_det) ;
1289
1290 /* Filter out traces */
1291 if (filter_traces) {
1292 char *setting_id;
1293 int zp_order;
1294 cpl_msg_info(__func__, "Filter out the traces") ;
1295 cpl_msg_indent_more() ;
1296
1297 /* Get the setting and the zp_order */
1298 plist = cpl_propertylist_load(first_file, 0) ;
1299 setting_id = cpl_strdup(cr2res_pfits_get_wlen_id(plist)) ;
1300 zp_order = cr2res_pfits_get_order_zp(plist) ;
1301 cpl_propertylist_delete(plist) ;
1302 filtered_traces = cr2res_trace_filter(computed_traces,
1303 setting_id, zp_order) ;
1304 cpl_free(setting_id) ;
1305 cpl_table_delete(computed_traces) ;
1306 computed_traces = filtered_traces ;
1307 cpl_msg_indent_less() ;
1308 if (cpl_table_get_nrow(computed_traces) == 0) {
1309 cpl_msg_error(__func__, "All traces are filtered out") ;
1310 hdrl_image_delete(collapsed) ;
1311 cpl_table_delete(computed_traces) ;
1312 cpl_msg_indent_less() ;
1313 return -1 ;
1314 }
1315 }
1316
1317 /* Use the input traces for extraction if they are provided */
1318 if (tw_frame) {
1319 traces = cr2res_io_load_TRACE_WAVE(cpl_frame_get_filename(tw_frame),
1320 reduce_det) ;
1321 if (traces == NULL) {
1322 cpl_msg_error(__func__, "Cannot Load the provided TW file") ;
1323 hdrl_image_delete(collapsed) ;
1324 cpl_table_delete(computed_traces) ;
1325 return -1 ;
1326 }
1327 } else {
1328 traces = cpl_table_duplicate(computed_traces) ;
1329 }
1330 cpl_table_delete(computed_traces) ;
1331
1332 /* Extract */
1333 nb_traces = cpl_table_get_nrow(traces) ;
1334 spectrum = cpl_malloc(nb_traces * sizeof(cpl_bivector *)) ;
1335 slit_func_vec = cpl_malloc(nb_traces * sizeof(cpl_vector *)) ;
1336 model_master = hdrl_image_new(CR2RES_DETECTOR_SIZE, CR2RES_DETECTOR_SIZE) ;
1337 hdrl_image_mul_scalar(model_master, (hdrl_value){0.0, 0.0}) ;
1338
1339 /* Loop over the traces and extract them */
1340 cpl_msg_info(__func__, "Extract the traces") ;
1341 cpl_msg_indent_more() ;
1342 for (i=0 ; i<nb_traces ; i++) {
1343 int order, trace_id;
1344 /* Initialise */
1345 slit_func_vec[i] = NULL ;
1346 spectrum[i] = NULL ;
1347 model_tmp = NULL ;
1348
1349 /* Get Order and trace id */
1350 order = cpl_table_get(traces, CR2RES_COL_ORDER, i, NULL) ;
1351 trace_id = cpl_table_get(traces, CR2RES_COL_TRACENB, i, NULL) ;
1352
1353 /* Check if this order needs to be skipped */
1354 if (reduce_order > -1 && order != reduce_order) continue ;
1355
1356 /* Check if this trace needs to be skipped */
1357 if (reduce_trace > -1 && trace_id != reduce_trace) continue ;
1358
1359 cpl_msg_info(__func__, "Process Order %d/Trace %d", order, trace_id) ;
1360 cpl_msg_indent_more() ;
1361
1362 /* Call the Extraction */
1363 if (extr_method == CR2RES_EXTR_SUM) {
1364 if (cr2res_extract_sum_vert(collapsed, traces, order, trace_id,
1365 extract_height, &(slit_func_vec[i]), &(spectrum[i]),
1366 &model_tmp) != 0) {
1367 cpl_msg_error(__func__, "Cannot (sum-)extract the trace") ;
1368 slit_func_vec[i] = NULL ;
1369 spectrum[i] = NULL ;
1370 model_tmp = NULL ;
1371 cpl_error_reset() ;
1372 cpl_msg_indent_less() ;
1373 continue ;
1374 }
1375 } else if (extr_method == CR2RES_EXTR_MEDIAN) {
1376 if (cr2res_extract_median(collapsed, traces, order, trace_id,
1377 extract_height, &(slit_func_vec[i]), &(spectrum[i]),
1378 &model_tmp) != 0) {
1379 cpl_msg_error(__func__, "Cannot (median-)extract the trace") ;
1380 slit_func_vec[i] = NULL ;
1381 spectrum[i] = NULL ;
1382 model_tmp = NULL ;
1383 cpl_error_reset() ;
1384 cpl_msg_indent_less() ;
1385 continue ;
1386 }
1387 } else if (extr_method == CR2RES_EXTR_TILTSUM) {
1388 if (cr2res_extract_sum_tilt(collapsed, traces, order, trace_id,
1389 extract_height, &(slit_func_vec[i]), &(spectrum[i]),
1390 &model_tmp) != 0) {
1391 cpl_msg_error(__func__, "Cannot (tiltsum-)extract the trace") ;
1392 slit_func_vec[i] = NULL ;
1393 spectrum[i] = NULL ;
1394 model_tmp = NULL ;
1395 cpl_error_reset() ;
1396 cpl_msg_indent_less() ;
1397 continue ;
1398 }
1399 } else if (extr_method == CR2RES_EXTR_OPT_CURV) {
1400 if (cr2res_extract_slitdec_curved(collapsed, traces, NULL,order,
1401 trace_id, extract_height, extract_swath_width,
1402 extract_oversample, 0.0, extract_smooth_slit,
1403 extract_smooth_spec, extract_niter, extract_kappa,
1404 error_factor,
1405 &(slit_func_vec[i]), &(spectrum[i]), &model_tmp) != 0) {
1406 cpl_msg_error(__func__, "Cannot (slitdec-) extract the trace") ;
1407 slit_func_vec[i] = NULL ;
1408 spectrum[i] = NULL ;
1409 model_tmp = NULL ;
1410 cpl_error_reset() ;
1411 cpl_msg_indent_less() ;
1412 continue ;
1413 }
1414 }
1415
1416 /* Update the model global image */
1417 if (model_tmp != NULL) {
1418 hdrl_image_add_image(model_master, model_tmp) ;
1419 hdrl_image_delete(model_tmp) ;
1420 }
1421 cpl_msg_indent_less() ;
1422 }
1423 cpl_msg_indent_less() ;
1424
1425 /* Create the slit_func_tab for the current detector */
1426 slit_func_tab = cr2res_extract_SLITFUNC_create(slit_func_vec, traces) ;
1427
1428 /* Create the extracted_tab for the current detector */
1429 extract_tab = cr2res_extract_EXTRACT1D_create(spectrum, traces) ;
1430
1431 /* Deallocate Vectors */
1432 for (i=0 ; i<nb_traces ; i++) {
1433 if (slit_func_vec[i] != NULL) cpl_vector_delete(slit_func_vec[i]) ;
1434 if (spectrum[i] != NULL) cpl_bivector_delete(spectrum[i]) ;
1435 }
1436 cpl_free(spectrum) ;
1437 cpl_free(slit_func_vec) ;
1438
1439 /* Compute the Master flat */
1440 cpl_msg_info(__func__, "Compute the master flat") ;
1441 cpl_msg_indent_more() ;
1442 if ((master_flat_loc = cr2res_master_flat(collapsed,
1443 model_master, bpm_low, bpm_high, bpm_linemax,
1444 &bpm_flat)) == NULL) {
1445 cpl_msg_error(__func__, "Failed compute the Master Flat") ;
1446 cpl_table_delete(slit_func_tab) ;
1447 cpl_table_delete(traces) ;
1448 cpl_table_delete(extract_tab) ;
1449 hdrl_image_delete(model_master) ;
1450 hdrl_image_delete(collapsed) ;
1451 cpl_msg_indent_less() ;
1452 return -1 ;
1453 }
1454 cpl_msg_indent_less() ;
1455 hdrl_image_delete(collapsed) ;
1456
1457 /* Create BPM image */
1458 bpm_im = NULL ;
1459 if (bpm_frame != NULL) {
1460 if ((bpm_im = cr2res_io_load_BPM(
1461 cpl_frame_get_filename(bpm_frame),
1462 reduce_det, 1)) == NULL) {
1463 cpl_msg_warning(__func__, "Failed to Load the Master BPM") ;
1464 }
1465 }
1466 if (bpm_im == NULL) {
1467 bpm_im = cpl_image_duplicate(bpm_flat) ;
1468 } else {
1469 if (cpl_image_or(bpm_im, NULL, bpm_flat)) {
1470 cpl_msg_error(__func__, "Failed to add the mask to the BPM") ;
1471 cpl_table_delete(traces) ;
1472 cpl_table_delete(slit_func_tab) ;
1473 cpl_table_delete(extract_tab) ;
1474 hdrl_image_delete(model_master) ;
1475 hdrl_image_delete(master_flat_loc) ;
1476 cpl_image_delete(bpm_im) ;
1477 cpl_image_delete(bpm_flat) ;
1478 cpl_msg_indent_less() ;
1479 return -1 ;
1480 }
1481 }
1482 cpl_image_delete(bpm_flat) ;
1483
1484 /* Compute the QC parameters */
1485
1486 /* Load the first RAW image */
1487 first_image = cr2res_io_load_image(first_file, reduce_det) ;
1488
1489 /* Set the BPM */
1490 bpm_mask = cr2res_bpm_extract_mask(bpm_im, CR2RES_BPM_ALL) ;
1491
1492 cpl_image_reject_from_mask(hdrl_image_get_image(first_image), bpm_mask) ;
1493 qc_mean = cpl_image_get_mean(hdrl_image_get_image(first_image)) ;
1494 qc_median = cpl_image_get_median(hdrl_image_get_image(first_image)) ;
1495 plist = cpl_propertylist_load(first_file, 0) ;
1496 dit = cr2res_pfits_get_dit(plist) ;
1497 cpl_propertylist_delete(plist) ;
1498 qc_flux = qc_mean / dit ;
1499
1500 my_master_flat = cpl_image_duplicate(hdrl_image_get_image(master_flat_loc));
1501 cpl_image_reject_from_mask(my_master_flat, bpm_mask) ;
1502 ngood = CR2RES_DETECTOR_SIZE* CR2RES_DETECTOR_SIZE -
1503 cpl_mask_count(bpm_mask) ;
1504 cpl_mask_delete(bpm_mask);
1505 qc_rms = sqrt(cpl_image_get_sqflux(my_master_flat) / ngood) ;
1506 cpl_image_delete(my_master_flat) ;
1507
1508 qc_s2n = cr2res_qc_flat_s2n(extract_tab) ;
1509
1510 cr2res_util_blaze_stat(extract_tab, &bl_good, &bl_tot);
1511
1512 qc_trace_centery = cr2res_qc_flat_trace_center_y(traces) ;
1513 qc_nbbad = cr2res_bpm_count(bpm_im, CR2RES_BPM_FLAT) ;
1514 cr2res_qc_flat_order_positions(traces, &qc_order_nb, &qc_order_pos,&nbvals);
1515
1516 /* Load the extension header for saving */
1517 ext_nr = cr2res_io_get_ext_idx(first_file, reduce_det, 1) ;
1518 plist = cpl_propertylist_load(first_file, ext_nr) ;
1519 if (plist == NULL) {
1520 hdrl_image_delete(first_image) ;
1521 cpl_table_delete(traces) ;
1522 cpl_table_delete(slit_func_tab) ;
1523 cpl_table_delete(extract_tab) ;
1524 hdrl_image_delete(model_master) ;
1525 hdrl_image_delete(master_flat_loc) ;
1526 cpl_image_delete(bpm_im) ;
1527 cpl_msg_error(__func__, "Failed to load the plist") ;
1528 return -1 ;
1529 }
1530
1531 /* QC.OVEREXPOSED */
1532
1533 /* Get all orders */
1534 orders = cr2res_trace_get_order_idx_values(traces, &nb_orders) ;
1535 /* Loop on the orders */
1536 for (i=0 ; i<nb_orders ; i++) {
1537 double qc_overexposed;
1538 qc_overexposed = cr2res_qc_overexposed(
1539 hdrl_image_get_image(first_image), traces, orders[i]) ;
1540 char * qc_name = cpl_sprintf("%s%02d",
1541 CR2RES_HEADER_QC_OVEREXPOSED, orders[i]) ;
1542 cpl_propertylist_append_double(plist, qc_name, qc_overexposed);
1543 cpl_free(qc_name) ;
1544 }
1545 cpl_free(orders) ;
1546 hdrl_image_delete(first_image) ;
1547
1548 /* Store the QC parameters in the plist */
1549 if (qc_order_nb != NULL && qc_order_pos != NULL) {
1550 for (i=0 ; i<nbvals ; i++) {
1551 char * qc_name = cpl_sprintf("%s%02d",
1552 CR2RES_HEADER_QC_FLAT_ORDERPOS, qc_order_nb[i]) ;
1553 cpl_propertylist_append_double(plist, qc_name, qc_order_pos[i]);
1554 cpl_free(qc_name) ;
1555 }
1556 }
1557 if (qc_order_nb != NULL) cpl_free(qc_order_nb) ;
1558 if (qc_order_pos != NULL) cpl_free(qc_order_pos) ;
1559
1560 cpl_propertylist_append_double(plist, CR2RES_HEADER_QC_FLAT_MEAN,
1561 qc_mean) ;
1562 cpl_propertylist_append_double(plist, CR2RES_HEADER_QC_FLAT_MEDIAN,
1563 qc_median);
1564 cpl_propertylist_append_double(plist, CR2RES_HEADER_QC_FLAT_FLUX,
1565 qc_flux) ;
1566 cpl_propertylist_append_double(plist, CR2RES_HEADER_QC_FLAT_RMS,
1567 qc_rms) ;
1568 cpl_propertylist_append_double(plist, CR2RES_HEADER_QC_FLAT_S2N,
1569 qc_s2n) ;
1570 cpl_propertylist_append_double(plist, CR2RES_HEADER_QC_FLAT_CENTERY,
1571 qc_trace_centery) ;
1572 cpl_propertylist_append_int(plist, CR2RES_HEADER_QC_FLAT_NBBAD,
1573 qc_nbbad) ;
1574 cpl_propertylist_append_int(plist, CR2RES_HEADER_QC_BLAZE_NGOOD,
1575 bl_good) ;
1576 cpl_propertylist_append_double(plist, CR2RES_HEADER_QC_BLAZE_TOT,
1577 bl_tot) ;
1578
1579 /* Return the results */
1580 *trace_wave = traces ;
1581 *slit_func = slit_func_tab ;
1582 *extract_1d = extract_tab ;
1583 *slit_model = model_master ;
1584 *master_flat = master_flat_loc ;
1585 *ext_plist = plist ;
1586 *bpm = bpm_im ;
1587 return 0 ;
1588}
1589
1590/*----------------------------------------------------------------------------*/
1597/*----------------------------------------------------------------------------*/
1598static int cr2res_cal_flat_compare(
1599 const cpl_frame * frame1,
1600 const cpl_frame * frame2)
1601{
1602 int comparison ;
1603 cpl_propertylist * plist1 ;
1604 cpl_propertylist * plist2 ;
1605 const char * sval1 ;
1606 const char * sval2 ;
1607
1608 /* Test entries */
1609 if (frame1==NULL || frame2==NULL) return -1 ;
1610
1611 /* Get property lists */
1612 if ((plist1=cpl_propertylist_load(cpl_frame_get_filename(frame1),0))==NULL){
1613 cpl_msg_error(__func__, "getting header from reference frame");
1614 return -1 ;
1615 }
1616 if ((plist2=cpl_propertylist_load(cpl_frame_get_filename(frame2),0))==NULL){
1617 cpl_msg_error(__func__, "getting header from reference frame");
1618 cpl_propertylist_delete(plist1) ;
1619 return -1 ;
1620 }
1621
1622 /* Test status */
1623 if (cpl_error_get_code()) {
1624 cpl_propertylist_delete(plist1) ;
1625 cpl_propertylist_delete(plist2) ;
1626 return -1 ;
1627 }
1628
1629 comparison = 1 ;
1630
1631 /* Compare the SETTING used */
1632 sval1 = cr2res_pfits_get_wlen_id(plist1) ;
1633 sval2 = cr2res_pfits_get_wlen_id(plist2) ;
1634 if (cpl_error_get_code()) {
1635 cpl_msg_error(__func__, "Cannot get the reference wavelength");
1636 cpl_propertylist_delete(plist1) ;
1637 cpl_propertylist_delete(plist2) ;
1638 return -1 ;
1639 }
1640 if (strcmp(sval1, sval2)) comparison = 0 ;
1641
1642 cpl_propertylist_delete(plist1) ;
1643 cpl_propertylist_delete(plist2) ;
1644 return comparison ;
1645}
1646
1647
cr2res_bpm_count
int cr2res_bpm_count(cpl_image *bpm, cr2res_bpm_type type)
Count BPM of a given type.
Definition cr2res_bpm.c:148
cr2res_bpm_extract_mask
cpl_mask * cr2res_bpm_extract_mask(const cpl_image *bpm_ima, cr2res_bpm_type bpm_type)
Extract a mask from a BPM image.
Definition cr2res_bpm.c:247
cr2res_calib_imagelist
hdrl_imagelist * cr2res_calib_imagelist(const hdrl_imagelist *in, int chip, int clean_bad, int subtract_nolight_rows, int subtract_interorder_column, int cosmics_corr, const cpl_frame *flat, const cpl_frame *dark, const cpl_frame *bpm, const cpl_frame *detlin, const cpl_vector *dits, const cpl_vector *ndits, const cpl_vector *ndsamples)
The images calibration routine for a given chip on a list.
Definition cr2res_calib.c:75
cr2res_dfs_set_groups
cpl_error_code cr2res_dfs_set_groups(cpl_frameset *set)
Set the group as RAW or CALIB in a frameset.
Definition cr2res_dfs.c:53
cr2res_extract_EXTRACT1D_create
cpl_table * cr2res_extract_EXTRACT1D_create(cpl_bivector **spectrum, const cpl_table *trace_table)
Create the extract 1D table to be saved.
Definition cr2res_extract.c:1097
cr2res_extract_median
int cr2res_extract_median(const hdrl_image *hdrl_in, const cpl_table *trace_tab, int order, int trace_id, int height, cpl_vector **slit_func, cpl_bivector **spec, hdrl_image **model)
Simple extraction function with the median.
Definition cr2res_extract.c:703
cr2res_extract_sum_tilt
int cr2res_extract_sum_tilt(const hdrl_image *hdrl_in, const cpl_table *trace_tab, int order, int trace_id, int height, cpl_vector **slit_func, cpl_bivector **spec, hdrl_image **model)
Simple extraction function with curvature correction.
Definition cr2res_extract.c:846
cr2res_extract_sum_vert
int cr2res_extract_sum_vert(const hdrl_image *hdrl_in, const cpl_table *trace_tab, int order, int trace_id, int height, cpl_vector **slit_func, cpl_bivector **spec, hdrl_image **model)
Simple extraction function.
Definition cr2res_extract.c:555
cr2res_extract_SLITFUNC_create
cpl_table * cr2res_extract_SLITFUNC_create(cpl_vector **slit_func, const cpl_table *trace_table)
Create the slit functions table to be saved.
Definition cr2res_extract.c:1207
cr2res_extract_slitdec_curved
int cr2res_extract_slitdec_curved(const hdrl_image *img_hdrl, const cpl_table *trace_tab, const cpl_vector *slit_func_vec_in, int order, int trace_id, int height, int swath, int oversample, double pclip, double smooth_slit, double smooth_spec, int niter, double kappa, double error_factor, cpl_vector **slit_func, cpl_bivector **spec, hdrl_image **model)
Extract optimally (slit-decomposition) with polynomial slit.
Definition cr2res_extract.c:1438
cr2res_master_flat
hdrl_image * cr2res_master_flat(const hdrl_image *collapsed, const hdrl_image *model_master, double low, double high, double bad_per_line_limit, cpl_image **bpm)
Compute the Master Flat.
Definition cr2res_flat.c:63
cr2res_io_find_BPM
const cpl_frame * cr2res_io_find_BPM(const cpl_frameset *in)
Get the first CR2RES_BPM_DRSTYPE frame from a frameset.
Definition cr2res_io.c:369
cr2res_io_save_SLIT_MODEL
int cr2res_io_save_SLIT_MODEL(const char *filename, cpl_frameset *allframes, cpl_frameset *inframes, const cpl_parameterlist *parlist, hdrl_image **data, const cpl_propertylist *qc_list, cpl_propertylist **ext_plist, const char *procatg, const char *recipe)
Save a SLIT_MODEL.
Definition cr2res_io.c:1853
cr2res_io_find_TRACE_WAVE
const cpl_frame * cr2res_io_find_TRACE_WAVE(const cpl_frameset *in)
Get the first CR2RES_TW_DRSTYPE frame from a frameset.
Definition cr2res_io.c:231
cr2res_io_read_ndsamples
cpl_vector * cr2res_io_read_ndsamples(const cpl_frameset *in)
Get the NDSAMPLES from a frame set.
Definition cr2res_io.c:492
cr2res_io_save_TRACE_WAVE
int cr2res_io_save_TRACE_WAVE(const char *filename, cpl_frameset *allframes, cpl_frameset *inframes, const cpl_parameterlist *parlist, cpl_table **tables, const cpl_propertylist *qc_list, cpl_propertylist **ext_plist, const char *procatg, const char *recipe)
Save a TRACE_WAVE.
Definition cr2res_io.c:1700
cr2res_io_extract_decker_frameset
cpl_frameset * cr2res_io_extract_decker_frameset(const cpl_frameset *in, const char *tag, cr2res_decker decker)
Extract the frames with the given tag and Decker position.
Definition cr2res_io.c:557
cr2res_io_read_dits
cpl_vector * cr2res_io_read_dits(const cpl_frameset *in)
Get the DITS from a frame set.
Definition cr2res_io.c:432
cr2res_io_save_EXTRACT_1D
int cr2res_io_save_EXTRACT_1D(const char *filename, cpl_frameset *allframes, cpl_frameset *inframes, const cpl_parameterlist *parlist, cpl_table **tables, const cpl_propertylist *qc_list, cpl_propertylist **ext_plist, const char *procatg, const char *recipe)
Save a 1D extracted spectrum.
Definition cr2res_io.c:1762
cr2res_io_load_image_list_from_set
hdrl_imagelist * cr2res_io_load_image_list_from_set(const cpl_frameset *in, int detector)
Load an hdrl image list from an images frameset.
Definition cr2res_io.c:838
cr2res_io_save_MASTER_FLAT
int cr2res_io_save_MASTER_FLAT(const char *filename, cpl_frameset *allframes, cpl_frameset *inframes, const cpl_parameterlist *parlist, hdrl_image **master_flats, const cpl_propertylist *qc_list, cpl_propertylist **ext_plist, const char *procatg, const char *recipe)
Save a MASTER_FLAT.
Definition cr2res_io.c:1669
cr2res_io_get_ext_idx
int cr2res_io_get_ext_idx(const char *filename, int detector, int data)
Get the wished extension number for a detector.
Definition cr2res_io.c:674
cr2res_io_load_image
hdrl_image * cr2res_io_load_image(const char *in, int detector)
Load an hdrl image from a image file.
Definition cr2res_io.c:734
cr2res_io_save_SLIT_FUNC
int cr2res_io_save_SLIT_FUNC(const char *filename, cpl_frameset *allframes, cpl_frameset *inframes, const cpl_parameterlist *parlist, cpl_table **slit_func, const cpl_propertylist *qc_list, cpl_propertylist **ext_plist, const char *procatg, const char *recipe)
Save a SLIT_FUNC.
Definition cr2res_io.c:1822
cr2res_io_save_BPM
int cr2res_io_save_BPM(const char *filename, cpl_frameset *allframes, cpl_frameset *inframes, const cpl_parameterlist *parlist, cpl_image **bpms, const cpl_propertylist *qc_list, cpl_propertylist **ext_plist, const char *procatg, const char *recipe)
Save a BPM.
Definition cr2res_io.c:1585
cr2res_io_load_TRACE_WAVE
cpl_table * cr2res_io_load_TRACE_WAVE(const char *filename, int detector)
Load a table from a TRACE_WAVE.
Definition cr2res_io.c:1139
cr2res_io_load_BPM
cpl_image * cr2res_io_load_BPM(const char *filename, int detector, int data)
Load an image from a BPM.
Definition cr2res_io.c:987
cr2res_io_read_ndits
cpl_vector * cr2res_io_read_ndits(const cpl_frameset *in)
Get the NDITs from a frame set.
Definition cr2res_io.c:462
cr2res_pfits_get_wlen_id
const char * cr2res_pfits_get_wlen_id(const cpl_propertylist *plist)
find out the Setting
Definition cr2res_pfits.c:137
cr2res_pfits_get_dit
double cr2res_pfits_get_dit(const cpl_propertylist *plist)
find out the DIT value
Definition cr2res_pfits.c:199
cr2res_qc_overexposed
double cr2res_qc_overexposed(const cpl_image *ima, const cpl_table *tw, int order_idx)
Computes the Overexposed fraction.
Definition cr2res_qc.c:796
cr2res_qc_flat_s2n
double cr2res_qc_flat_s2n(const cpl_table *extracted)
Computes the S2N on the flat.
Definition cr2res_qc.c:375
cr2res_qc_flat_order_positions
int cr2res_qc_flat_order_positions(const cpl_table *tw, int **order_nb, double **order_pos, int *nbvals)
Compute the central orders positions.
Definition cr2res_qc.c:235
cr2res_qc_flat_trace_center_y
double cr2res_qc_flat_trace_center_y(const cpl_table *trace)
Computes the mean Y coord of the central order.
Definition cr2res_qc.c:325
cr2res_trace_filter
cpl_table * cr2res_trace_filter(const cpl_table *tw, const char *setting, int zp_order)
Only keep the predefined orders of the setting.
Definition cr2res_trace.c:1965
cr2res_trace
cpl_table * cr2res_trace(cpl_image *ima, int smooth_x, int smooth_y, double threshold, int opening, int degree, int min_cluster)
Main function for running all parts of the trace algorithm.
Definition cr2res_trace.c:163
cr2res_trace_add_extra_columns
int cr2res_trace_add_extra_columns(cpl_table *traces, const char *infile, int det_nr)
Add extra columns to the plain trace table.
Definition cr2res_trace.c:972
cr2res_trace_get_order_idx_values
int * cr2res_trace_get_order_idx_values(const cpl_table *trace, int *nb_order_idx_values)
Count and return the different order_idx values from a TW table.
Definition cr2res_trace.c:431
cr2res_trace_merge
cpl_table * cr2res_trace_merge(const cpl_table *trace_wave1, const cpl_table *trace_wave2)
Merge 2 trace_wave tables.
Definition cr2res_trace.c:592
cr2res_format_setting
int cr2res_format_setting(char *setting_id)
Format the setting.
Definition cr2res_utils.c:152
cr2res_util_blaze_stat
double cr2res_util_blaze_stat(const cpl_table *blaze, int *ngood, double *blaze_tot)
Computes the average of the blaze.
Definition cr2res_utils.c:2223
cr2res_extract_frameset
cpl_frameset * cr2res_extract_frameset(const cpl_frameset *in, const char *tag)
Extract the frames with the given tag from a frameset.
Definition cr2res_utils.c:1179
cr2res_get_license
const char * cr2res_get_license(void)
Get the pipeline copyright and license.
Definition cr2res_utils.c:2362
hdrl_image_add_image
cpl_error_code hdrl_image_add_image(hdrl_image *self, const hdrl_image *other)
Add two images, store the result in the first image.
Definition hdrl_image_math.c:59
hdrl_image_mul_scalar
cpl_error_code hdrl_image_mul_scalar(hdrl_image *self, hdrl_value value)
Elementwise multiplication of an image with a scalar.
Definition hdrl_image_math.c:245
hdrl_image_get_image
cpl_image * hdrl_image_get_image(hdrl_image *himg)
get data as cpl image
Definition hdrl_image.c:105
hdrl_image_new
hdrl_image * hdrl_image_new(cpl_size nx, cpl_size ny)
create new zero filled hdrl image
Definition hdrl_image.c:311
hdrl_image_delete
void hdrl_image_delete(hdrl_image *himg)
delete hdrl_image
Definition hdrl_image.c:379
hdrl_imagelist_collapse_mean
cpl_error_code hdrl_imagelist_collapse_mean(const hdrl_imagelist *himlist, hdrl_image **out, cpl_image **contrib)
Mean collapsing of image list.
Definition hdrl_imagelist_basic.c:476
hdrl_imagelist_delete
void hdrl_imagelist_delete(hdrl_imagelist *himlist)
Free all memory used by a hdrl_imagelist object including the images.
Definition hdrl_imagelist_io.c:384
hdrl_imagelist_get
hdrl_image * hdrl_imagelist_get(const hdrl_imagelist *himlist, cpl_size inum)
Get an image from a list of images.
Definition hdrl_imagelist_io.c:210
Generated on Tue Oct 13 2026 12:27:26 for CR2RE Pipeline Reference Manual by doxygen 1.13.2