GIRAFFE Pipeline Reference Manual
Overview Modules Structures Source Files
GIRAFFE
Credits
giscience.c
1/*
2 * This file is part of the GIRAFFE Pipeline
3 * Copyright (C) 2002-2019 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 02110-1301 USA
18 */
19
20#ifdef HAVE_CONFIG_H
21# include <config.h>
22#endif
23
24#include <assert.h>
25#include <math.h>
26
27#include <cxslist.h>
28#include <cxmessages.h>
29
30#include <cpl_recipe.h>
31#include <cpl_plugininfo.h>
32#include <cpl_parameterlist.h>
33#include <cpl_frameset.h>
34#include <cpl_msg.h>
35#include <cpl_errorstate.h>
36
37#include <irplib_sdp_spectrum.h>
38
39#include "gialias.h"
40#include "giimage.h"
41#include "giframe.h"
42#include "gifibers.h"
43#include "gifiberutils.h"
44#include "gislitgeometry.h"
45#include "gipsfdata.h"
46#include "gibias.h"
47#include "gidark.h"
48#include "giextract.h"
49#include "giflat.h"
50#include "gitransmission.h"
51#include "girebinning.h"
52#include "gisgcalibration.h"
53#include "giastrometry.h"
54#include "gifov.h"
55#include "gimessages.h"
56#include "gierror.h"
57#include "giutils.h"
58
59const char *sciqcpar[] = {
60 GIALIAS_WLSTART, GIALIAS_WLEND,
61 GIALIAS_WLSTEP, GIALIAS_QCAIR, GIALIAS_QCFWHM,
62 GIALIAS_QCNFIB, GIALIAS_QCNFIBSCI, GIALIAS_QCNFIBSKY,
63 GIALIAS_QCMEANRED, GIALIAS_QCNSATSCI, GIALIAS_QCSNR,
64 GIALIAS_QCMAG, GIALIAS_QCDLTTEMP, GIALIAS_QCDLTTIME,
65 GIALIAS_QCBRIGHTFLG
66};
67
68static cxint giscience(cpl_parameterlist*, cpl_frameset*);
69
70static cxint _giraffe_make_sdp_spectra(const cxchar* flux_filename,
71 const cxchar* err_filename,
72 cxint nassoc_keys,
73 cpl_frameset* allframes,
74 const cpl_parameterlist* parlist,
75 const cxchar* recipe_id);
76
77const cxdouble saturation = 60000.;
78
79/*-----------------------------------------------------------------------------
80 Utility functions
81 -----------------------------------------------------------------------------*/
82
83 static void replace_spaces_with_underscores(char *str)
84 {
85
86 for (int i = 0; str[i] != '\0'; i++) {
87 if (str[i] == ' ') {
88 str[i] = '_'; // Replace space with an underscore
89 }
90 }
91 }
92
93/*
94 * Create the recipe instance, i.e. setup the parameter list for this
95 * recipe and make it available to the application using the interface.
96 */
97
98static cxint
99giscience_create(cpl_plugin* plugin)
100{
101
102 cpl_recipe* recipe = (cpl_recipe*)plugin;
103
104 cpl_parameter* p = NULL;
105
106
107 giraffe_error_init();
108
109
110 /*
111 * We have to provide the option we accept to the application. We
112 * need to setup our parameter list and hook it into the recipe
113 * interface.
114 */
115
116 recipe->parameters = cpl_parameterlist_new();
117 cx_assert(recipe->parameters != NULL);
118
119
120 /*
121 * Fill the parameter list.
122 */
123
124 /* Bias removal */
125
126 giraffe_bias_config_add(recipe->parameters);
127
128 /* Dark subtraction */
129
130 /* TBD */
131
132 /* Spectrum extraction */
133
134 giraffe_extract_config_add(recipe->parameters);
135
136 /* Flat fielding and relative fiber transmission correction */
137
138 giraffe_flat_config_add(recipe->parameters);
139
140 /* Spectrum rebinning */
141
142 giraffe_rebin_config_add(recipe->parameters);
143
144 /* Simultaneous wavelength calibration correction */
145
146 p = cpl_parameter_new_value("giraffe.siwc.apply",
147 CPL_TYPE_BOOL,
148 "Enable simultaneous wavelength calibration "
149 "correction.",
150 "giraffe.siwc",
151 TRUE);
152
153 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "siwc-apply");
154 cpl_parameterlist_append(recipe->parameters, p);
155
156 giraffe_sgcalibration_config_add(recipe->parameters);
157
158 /* Image reconstruction (IFU and Argus only) */
159
160 giraffe_fov_config_add(recipe->parameters);
161
162 /* Science Data Product format generation parameters: */
163
164 p = cpl_parameter_new_value("giraffe.sdp.format.generate",
165 CPL_TYPE_BOOL,
166 "TRUE if additional files should be generated"
167 " in Science Data Product (SDP) format.",
168 "giraffe.sdp",
169 FALSE);
170
171 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI, "generate-SDP-format");
172 cpl_parameterlist_append(recipe->parameters, p);
173
174 p = cpl_parameter_new_value("giraffe.sdp.nassoc.keys",
175 CPL_TYPE_INT,
176 "Sets the number of dummy (empty) ASSONi,"
177 " ASSOCi and ASSOMi keywords to create.",
178 "giraffe.sdp",
179 (int)0);
180
181 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
182 "dummy-association-keys");
183 cpl_parameterlist_append(recipe->parameters, p);
184
185 return 0;
186
187}
188
189
190/*
191 * Execute the plugin instance given by the interface.
192 */
193
194static cxint
195giscience_exec(cpl_plugin* plugin)
196{
197 cxint result;
198 cpl_errorstate prev_state;
199
200 cpl_recipe* recipe = (cpl_recipe*)plugin;
201
202
203 cx_assert(recipe->parameters != NULL);
204 cx_assert(recipe->frames != NULL);
205
206 prev_state = cpl_errorstate_get();
207 result = giscience(recipe->parameters, recipe->frames);
208 if (result != 0) {
209 cpl_errorstate_dump(prev_state, CPL_FALSE, cpl_errorstate_dump_one);
210 }
211 return result;
212}
213
214
215static cxint
216giscience_destroy(cpl_plugin* plugin)
217{
218
219 cpl_recipe* recipe = (cpl_recipe*)plugin;
220
221
222 /*
223 * We just destroy what was created during the plugin initialization
224 * phase, i.e. the parameter list. The frame set is managed by the
225 * application which called us, so we must not touch it,
226 */
227
228 cpl_parameterlist_delete(recipe->parameters);
229
230 giraffe_error_clear();
231
232 return 0;
233
234}
235
236
237/*
238 * The actual recipe starts here.
239 */
240
241static cxint
242giscience(cpl_parameterlist* config, cpl_frameset* set)
243{
244
245 const cxchar* const _id = "giscience";
246
247
248 const cxchar* filename = NULL;
249
250 cxbool siwc = FALSE;
251 cxbool calsim = FALSE;
252
253 cxbool gensdp = FALSE;
254 cxint nassoc_keys = 0;
255 cxchar* flux_filename = NULL;
256 cxchar* err_filename = NULL;
257
258 cxint status = 0;
259
260 cxlong i;
261 cxlong nscience = 0;
262
263 cxdouble exptime = 0.;
264
265 cx_slist* slist = NULL;
266
267 cpl_propertylist* properties = NULL;
268
269 cpl_matrix* biasareas = NULL;
270
271 cpl_frame* science_frame = NULL;
272 cpl_frame* mbias_frame = NULL;
273 cpl_frame* mdark_frame = NULL;
274 cpl_frame* bpixel_frame = NULL;
275 cpl_frame* slight_frame = NULL;
276 cpl_frame* locy_frame = NULL;
277 cpl_frame* locw_frame = NULL;
278 cpl_frame* psfdata_frame = NULL;
279 cpl_frame* grating_frame = NULL;
280 cpl_frame* linemask_frame = NULL;
281 cpl_frame* slit_frame = NULL;
282 cpl_frame* wcal_frame = NULL;
283 cpl_frame* rscience_frame = NULL;
284 cpl_frame* sext_frame = NULL;
285 cpl_frame* rbin_frame = NULL;
286
287 cpl_parameter* p = NULL;
288
289 GiImage* mbias = NULL;
290 GiImage* mdark = NULL;
291 GiImage* bpixel = NULL;
292 GiImage* slight = NULL;
293 GiImage* sscience = NULL;
294 GiImage* rscience = NULL;
295
296 GiTable* fibers = NULL;
297 GiTable* slitgeometry = NULL;
298 GiTable* grating = NULL;
299 GiTable* wcalcoeff = NULL;
300
301 GiLocalization* localization = NULL;
302 GiExtraction* extraction = NULL;
303 GiRebinning* rebinning = NULL;
304
305 GiBiasConfig* bias_config = NULL;
306 GiExtractConfig* extract_config = NULL;
307 GiFlatConfig* flat_config = NULL;
308 GiRebinConfig* rebin_config = NULL;
309
310 GiInstrumentMode mode;
311
312 GiRecipeInfo info = {(cxchar*)_id, 1, NULL, config};
313
314 GiGroupInfo groups[] = {
315 {GIFRAME_SCIENCE, CPL_FRAME_GROUP_RAW},
316 {GIFRAME_BADPIXEL_MAP, CPL_FRAME_GROUP_CALIB},
317 {GIFRAME_BIAS_MASTER, CPL_FRAME_GROUP_CALIB},
318 {GIFRAME_DARK_MASTER, CPL_FRAME_GROUP_CALIB},
319 {GIFRAME_FIBER_FLAT_EXTSPECTRA, CPL_FRAME_GROUP_CALIB},
320 {GIFRAME_FIBER_FLAT_EXTERRORS, CPL_FRAME_GROUP_CALIB},
321 {GIFRAME_SCATTERED_LIGHT_MODEL, CPL_FRAME_GROUP_CALIB},
322 {GIFRAME_LOCALIZATION_CENTROID, CPL_FRAME_GROUP_CALIB},
323 {GIFRAME_LOCALIZATION_WIDTH, CPL_FRAME_GROUP_CALIB},
324 {GIFRAME_PSF_CENTROID, CPL_FRAME_GROUP_CALIB},
325 {GIFRAME_PSF_WIDTH, CPL_FRAME_GROUP_CALIB},
326 {GIFRAME_PSF_DATA, CPL_FRAME_GROUP_CALIB},
327 {GIFRAME_WAVELENGTH_SOLUTION, CPL_FRAME_GROUP_CALIB},
328 {GIFRAME_LINE_MASK, CPL_FRAME_GROUP_CALIB},
329 {GIFRAME_SLITSETUP, CPL_FRAME_GROUP_CALIB},
330 {GIFRAME_SLITMASTER, CPL_FRAME_GROUP_CALIB},
331 {GIFRAME_GRATING, CPL_FRAME_GROUP_CALIB},
332 {NULL, CPL_FRAME_GROUP_NONE}
333 };
334
335
336
337 if (!config) {
338 cpl_msg_error(_id, "Invalid parameter list! Aborting ...");
339 return 1;
340 }
341
342 if (!set) {
343 cpl_msg_error(_id, "Invalid frame set! Aborting ...");
344 return 1;
345 }
346
347 status = giraffe_frameset_set_groups(set, groups);
348
349 if (status != 0) {
350 cpl_msg_error(_id, "Setting frame group information failed!");
351 return 1;
352 }
353
354
355 /*
356 * Verify the frame set contents
357 */
358
359 nscience = cpl_frameset_count_tags(set, GIFRAME_SCIENCE);
360
361 if (nscience < 1) {
362 cpl_msg_error(_id, "Too few (%ld) raw frames (%s) present in "
363 "frame set! Aborting ...", nscience, GIFRAME_SCIENCE);
364 return 1;
365 }
366
367 locy_frame = cpl_frameset_find(set, GIFRAME_PSF_CENTROID);
368
369 if (locy_frame == NULL) {
370
371 locy_frame = cpl_frameset_find(set, GIFRAME_LOCALIZATION_CENTROID);
372
373 if (locy_frame == NULL) {
374 cpl_msg_info(_id, "No master localization (centroid position) "
375 "present in frame set. Aborting ...");
376 return 1;
377 }
378
379 }
380
381 locw_frame = cpl_frameset_find(set, GIFRAME_PSF_WIDTH);
382
383 if (locw_frame == NULL) {
384
385 locw_frame = cpl_frameset_find(set, GIFRAME_LOCALIZATION_WIDTH);
386
387 if (locw_frame == NULL) {
388 cpl_msg_info(_id, "No master localization (spectrum width) "
389 "present in frame set. Aborting ...");
390 return 1;
391 }
392
393 }
394
395 grating_frame = cpl_frameset_find(set, GIFRAME_GRATING);
396
397 if (!grating_frame) {
398 cpl_msg_error(_id, "No grating data present in frame set. "
399 "Aborting ...");
400 return 1;
401 }
402
403 slit_frame = giraffe_get_slitgeometry(set);
404
405 if (!slit_frame) {
406 cpl_msg_error(_id, "No slit geometry present in frame set. "
407 "Aborting ...");
408 return 1;
409 }
410
411 wcal_frame = cpl_frameset_find(set, GIFRAME_WAVELENGTH_SOLUTION);
412
413 if (!wcal_frame) {
414 cpl_msg_error(_id, "No dispersion solution present in frame set. "
415 "Aborting ...");
416 return 1;
417 }
418
419 linemask_frame = cpl_frameset_find(set, GIFRAME_LINE_MASK);
420
421 if (!linemask_frame) {
422 cpl_msg_warning(_id, "No reference line mask present in frame set.");
423 }
424
425 bpixel_frame = cpl_frameset_find(set, GIFRAME_BADPIXEL_MAP);
426
427 if (!bpixel_frame) {
428 cpl_msg_info(_id, "No bad pixel map present in frame set.");
429 }
430
431 mbias_frame = cpl_frameset_find(set, GIFRAME_BIAS_MASTER);
432
433 if (!mbias_frame) {
434 cpl_msg_info(_id, "No master bias present in frame set.");
435 }
436
437 mdark_frame = cpl_frameset_find(set, GIFRAME_DARK_MASTER);
438
439 if (!mdark_frame) {
440 cpl_msg_info(_id, "No master dark present in frame set.");
441 }
442
443 slight_frame = cpl_frameset_find(set, GIFRAME_SCATTERED_LIGHT_MODEL);
444
445 if (!slight_frame) {
446 cpl_msg_info(_id, "No scattered light model present in frame set.");
447 }
448
449 psfdata_frame = cpl_frameset_find(set, GIFRAME_PSF_DATA);
450
451 if (!psfdata_frame) {
452 cpl_msg_info(_id, "No PSF profile parameters present in frame set.");
453 }
454
455
456 /*
457 * Load raw images
458 */
459
460 slist = cx_slist_new();
461
462 science_frame = cpl_frameset_find(set, GIFRAME_SCIENCE);
463
464 for (i = 0; i < nscience; i++) {
465
466 filename = cpl_frame_get_filename(science_frame);
467
468 GiImage* raw = giraffe_image_new(CPL_TYPE_DOUBLE);
469
470
471 status = giraffe_image_load(raw, filename, 0);
472
473 if (status) {
474 cpl_msg_error(_id, "Cannot load raw science frame from '%s'. "
475 "Aborting ...", filename);
476
477 cx_slist_destroy(slist, (cx_free_func) giraffe_image_delete);
478
479 return 1;
480 }
481
482 cx_slist_push_back(slist, raw);
483
484 science_frame = cpl_frameset_find(set, NULL);
485
486 }
487
488 nscience = (cxint)cx_slist_size(slist);
489 sscience = cx_slist_pop_front(slist);
490
491 properties = giraffe_image_get_properties(sscience);
492 cx_assert(properties != NULL);
493
494 cpl_errorstate tempes = cpl_errorstate_get();
495 int nsaturated = -1;
496 cpl_image *sciim = giraffe_image_get(sscience);
497 if (sciim != NULL) {
498 double *scipix = cpl_image_get_data_double(sciim);
499 if (scipix != NULL) {
500 const size_t nxy = (size_t)(cpl_image_get_size_x(sciim) *
501 cpl_image_get_size_y(sciim));
502 size_t scii;
503 nsaturated = 0;
504 for (scii = 0; scii < nxy; scii++) {
505 if (scipix[scii] > saturation)
506 nsaturated++;
507 }
508 }
509 }
510 cpl_errorstate_set(tempes);
511
512
513 cpl_propertylist_update_int(properties, GIALIAS_QCNSATSCI, nsaturated);
514
515 if (nscience > 1) {
516
517 /*
518 * Create a stacked science image from the list of raw images.
519 * Each raw image is disposed when it is no longer needed.
520 */
521
522 cpl_msg_info(_id, "Averaging science frames ...");
523
524 exptime = cpl_propertylist_get_double(properties, GIALIAS_EXPTIME);
525
526 for (i = 1; i < nscience; i++) {
527
528 cpl_propertylist* _properties;
529
530 GiImage* science = cx_slist_pop_front(slist);
531
532
533 cpl_image_add(giraffe_image_get(sscience),
534 giraffe_image_get(science));
535
536 _properties = giraffe_image_get_properties(science);
537 cx_assert(_properties != NULL);
538
539 exptime += cpl_propertylist_get_double(_properties, GIALIAS_EXPTIME);
540
541 giraffe_image_delete(science);
542
543 }
544
545 cpl_image_divide_scalar(giraffe_image_get(sscience), nscience);
546 }
547
548 cx_assert(cx_slist_empty(slist));
549 cx_slist_delete(slist);
550 slist = NULL;
551
552
553 if (nscience > 1) {
554
555 /*
556 * Update stacked science image properties
557 */
558
559 cpl_msg_info(_id, "Updating stacked science image properties ...");
560
561 cpl_propertylist_set_double(properties, GIALIAS_EXPTIME,
562 exptime / nscience);
563
564 cpl_propertylist_append_double(properties, GIALIAS_EXPTTOT, exptime);
565 cpl_propertylist_set_comment(properties, GIALIAS_EXPTTOT,
566 "Total exposure time of all frames "
567 "combined");
568
569 cpl_propertylist_erase(properties, GIALIAS_TPLEXPNO);
570
571 }
572
573 cpl_propertylist_append_int(properties, GIALIAS_DATANCOM, nscience);
574 cpl_propertylist_set_comment(properties, GIALIAS_DATANCOM,
575 "Number of combined frames");
576
577
578
579 /*
580 * Prepare for bias subtraction
581 */
582
583 bias_config = giraffe_bias_config_create(config);
584
585 /*
586 * Setup user defined areas to use for the bias computation
587 */
588
589 if (bias_config->method == GIBIAS_METHOD_MASTER ||
590 bias_config->method == GIBIAS_METHOD_ZMASTER) {
591
592 if (!mbias_frame) {
593 cpl_msg_error(_id, "Missing master bias frame! Selected bias "
594 "removal method requires a master bias frame!");
595
596 giraffe_bias_config_destroy(bias_config);
597 giraffe_image_delete(sscience);
598
599 return 1;
600 }
601 else {
602 filename = cpl_frame_get_filename(mbias_frame);
603
604
605 mbias = giraffe_image_new(CPL_TYPE_DOUBLE);
606 status = giraffe_image_load(mbias, filename, 0);
607
608 if (status) {
609 cpl_msg_error(_id, "Cannot load master bias from '%s'. "
610 "Aborting ...", filename);
611
612 giraffe_bias_config_destroy(bias_config);
613 giraffe_image_delete(sscience);
614
615 return 1;
616 }
617 }
618 }
619
620
621 /*
622 * Load bad pixel map if it is present in the frame set.
623 */
624
625 if (bpixel_frame) {
626
627 filename = cpl_frame_get_filename(bpixel_frame);
628
629
630 bpixel = giraffe_image_new(CPL_TYPE_INT);
631 status = giraffe_image_load(bpixel, filename, 0);
632
633 if (status) {
634 cpl_msg_error(_id, "Cannot load bad pixel map from '%s'. "
635 "Aborting ...", filename);
636
637 giraffe_image_delete(bpixel);
638 bpixel = NULL;
639
640 if (mbias != NULL) {
641 giraffe_image_delete(mbias);
642 mbias = NULL;
643 }
644
645 giraffe_bias_config_destroy(bias_config);
646 bias_config = NULL;
647
648 giraffe_image_delete(sscience);
649 sscience = NULL;
650
651 return 1;
652 }
653
654 }
655
656
657 /*
658 * Compute and remove the bias from the stacked flat field frame.
659 */
660
661 rscience = giraffe_image_new(CPL_TYPE_DOUBLE);
662
663 status = giraffe_bias_remove(rscience, sscience, mbias, bpixel, biasareas,
664 bias_config);
665
666 giraffe_image_delete(sscience);
667
668 if (mbias) {
669 giraffe_image_delete(mbias);
670 mbias = NULL;
671 }
672
673 giraffe_bias_config_destroy(bias_config);
674
675 if (status) {
676 cpl_msg_error(_id, "Bias removal failed. Aborting ...");
677
678 giraffe_image_delete(rscience);
679 rscience = NULL;
680
681 if (bpixel != NULL) {
682 giraffe_image_delete(bpixel);
683 bpixel = NULL;
684 }
685
686 return 1;
687 }
688
689
690 /*
691 * Load master dark if it is present in the frame set and correct
692 * the master flat field for the dark current.
693 */
694
695 if (mdark_frame) {
696
697 GiDarkConfig dark_config = {GIDARK_METHOD_ZMASTER, 0.};
698
699
700 cpl_msg_info(_id, "Correcting for dark current ...");
701
702 filename = cpl_frame_get_filename(mdark_frame);
703
704 mdark = giraffe_image_new(CPL_TYPE_DOUBLE);
705 status = giraffe_image_load(mdark, filename, 0);
706
707 if (status != 0) {
708 cpl_msg_error(_id, "Cannot load master dark from '%s'. "
709 "Aborting ...", filename);
710
711 giraffe_image_delete(rscience);
712 rscience = NULL;
713
714 if (bpixel != NULL) {
715 giraffe_image_delete(bpixel);
716 bpixel = NULL;
717 }
718
719 return 1;
720 }
721
722 status = giraffe_subtract_dark(rscience, mdark, bpixel, NULL,
723 &dark_config);
724
725 if (status != 0) {
726 cpl_msg_error(_id, "Dark subtraction failed! Aborting ...");
727
728 giraffe_image_delete(mdark);
729 mdark = NULL;
730
731 giraffe_image_delete(rscience);
732 rscience = NULL;
733
734 if (bpixel != NULL) {
735 giraffe_image_delete(bpixel);
736 bpixel = NULL;
737 }
738
739 return 1;
740 }
741
742 giraffe_image_delete(mdark);
743 mdark = NULL;
744
745 }
746
747
748 /*
749 * Update the reduced science properties, save the reduced science frame
750 * and register it as product.
751 */
752
753 cpl_msg_info(_id, "Writing pre-processed science image ...");
754
755 giraffe_image_add_info(rscience, &info, set);
756
757 rscience_frame = giraffe_frame_create_image(rscience,
758 GIFRAME_SCIENCE_REDUCED,
759 CPL_FRAME_LEVEL_INTERMEDIATE,
760 TRUE, TRUE);
761
762 if (rscience_frame == NULL) {
763 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
764
765 giraffe_image_delete(rscience);
766
767 return 1;
768 }
769
770 cpl_frameset_insert(set, rscience_frame);
771
772
773 /*
774 * Determine fiber setup
775 */
776
777 science_frame = cpl_frameset_find(set, GIFRAME_SCIENCE);
778
779 cpl_msg_info(_id, "Building fiber setup for frame '%s'.",
780 cpl_frame_get_filename(science_frame));
781
782 fibers = giraffe_fibers_setup(science_frame, locy_frame);
783
784 if (!fibers) {
785 cpl_msg_error(_id, "Cannot create fiber setup for frame '%s'! "
786 "Aborting ...", cpl_frame_get_filename(science_frame));
787
788 if (bpixel) {
789 giraffe_image_delete(bpixel);
790 bpixel = NULL;
791 }
792
793 giraffe_image_delete(rscience);
794 rscience = NULL;
795
796 return 1;
797 }
798
799 cpl_msg_info(_id, "Fiber reference setup taken from localization "
800 "frame '%s'.", cpl_frame_get_filename(locy_frame));
801
802
803 /*
804 * Load fiber localization
805 */
806
807 localization = giraffe_localization_new();
808
809 filename = cpl_frame_get_filename(locy_frame);
810
811 localization->locy = giraffe_image_new(CPL_TYPE_DOUBLE);
812 status = giraffe_image_load(localization->locy, filename, 0);
813
814 if (status) {
815 cpl_msg_error(_id, "Cannot load localization (centroid "
816 "position) frame from '%s'. Aborting ...",
817 filename);
818
819 giraffe_localization_destroy(localization);
820
821 if (bpixel) {
822 giraffe_image_delete(bpixel);
823 bpixel = NULL;
824 }
825
826 giraffe_table_delete(fibers);
827 giraffe_image_delete(rscience);
828
829 return 1;
830 }
831
832
833 filename = cpl_frame_get_filename(locw_frame);
834
835 localization->locw = giraffe_image_new(CPL_TYPE_DOUBLE);
836 status = giraffe_image_load(localization->locw, filename, 0);
837
838 if (status) {
839 cpl_msg_error(_id, "Cannot load localization (spectrum width) "
840 "frame from '%s'. Aborting ...", filename);
841
842 giraffe_localization_destroy(localization);
843
844 if (bpixel) {
845 giraffe_image_delete(bpixel);
846 bpixel = NULL;
847 }
848
849 giraffe_table_delete(fibers);
850 giraffe_image_delete(rscience);
851
852 return 1;
853 }
854
855
856 /*
857 * Spectrum extraction
858 */
859
860 if (slight_frame) {
861
862 filename = cpl_frame_get_filename(slight_frame);
863
864
865 slight = giraffe_image_new(CPL_TYPE_DOUBLE);
866 status = giraffe_image_load(slight, filename, 0);
867
868 if (status) {
869 cpl_msg_error(_id, "Cannot load scattered light model from '%s'. "
870 "Aborting ...", filename);
871
872 giraffe_image_delete(slight);
873
874 giraffe_localization_destroy(localization);
875
876 if (bpixel) {
877 giraffe_image_delete(bpixel);
878 bpixel = NULL;
879 }
880
881 giraffe_table_delete(fibers);
882 giraffe_image_delete(rscience);
883
884 return 1;
885
886 }
887
888 }
889
890
891 extract_config = giraffe_extract_config_create(config);
892
893 if ((extract_config->emethod == GIEXTRACT_OPTIMAL) ||
894 (extract_config->emethod == GIEXTRACT_HORNE)) {
895
896 if (psfdata_frame == NULL) {
897
898 const cxchar* emethod = "Optimal";
899
900 if (extract_config->emethod == GIEXTRACT_HORNE) {
901 emethod = "Horne";
902 }
903
904 cpl_msg_error(_id, "%s spectrum extraction requires PSF "
905 "profile data. Aborting ...", emethod);
906
907 giraffe_extract_config_destroy(extract_config);
908 extract_config = NULL;
909
910 if (slight != NULL) {
911 giraffe_image_delete(slight);
912 slight = NULL;
913 }
914
915 giraffe_localization_destroy(localization);
916 localization = NULL;
917
918 if (bpixel) {
919 giraffe_image_delete(bpixel);
920 bpixel = NULL;
921 }
922
923 giraffe_table_delete(fibers);
924 fibers = NULL;
925
926 giraffe_image_delete(rscience);
927 rscience = NULL;
928
929 return 1;
930
931 }
932 else {
933
934 filename = cpl_frame_get_filename(psfdata_frame);
935
936 localization->psf = giraffe_psfdata_new();
937 status = giraffe_psfdata_load(localization->psf, filename);
938
939 if (status) {
940 cpl_msg_error(_id, "Cannot load PSF profile data frame from "
941 "'%s'. Aborting ...", filename);
942
943 giraffe_extract_config_destroy(extract_config);
944 extract_config = NULL;
945
946 if (slight != NULL) {
947 giraffe_image_delete(slight);
948 slight = NULL;
949 }
950
951 giraffe_localization_destroy(localization);
952 localization = NULL;
953
954 if (bpixel) {
955 giraffe_image_delete(bpixel);
956 bpixel = NULL;
957 }
958
959 giraffe_table_delete(fibers);
960 fibers = NULL;
961
962 giraffe_image_delete(rscience);
963 rscience = NULL;
964
965 return 1;
966
967 }
968
969 }
970
971 }
972
973
974 extraction = giraffe_extraction_new();
975
976 status = giraffe_extract_spectra(extraction, rscience, fibers,
977 localization, bpixel, slight,
978 extract_config);
979
980 if (status) {
981 cpl_msg_error(_id, "Spectrum extraction failed! Aborting ...");
982
983 giraffe_extraction_destroy(extraction);
984 giraffe_extract_config_destroy(extract_config);
985
986 giraffe_image_delete(slight);
987
988 giraffe_localization_destroy(localization);
989
990 if (bpixel) {
991 giraffe_image_delete(bpixel);
992 bpixel = NULL;
993 }
994
995 giraffe_table_delete(fibers);
996 giraffe_image_delete(rscience);
997
998 return 1;
999 }
1000
1001 giraffe_image_delete(slight);
1002 slight = NULL;
1003
1004 if (bpixel) {
1005 giraffe_image_delete(bpixel);
1006 bpixel = NULL;
1007 }
1008
1009 giraffe_image_delete(rscience);
1010 rscience = NULL;
1011
1012 giraffe_extract_config_destroy(extract_config);
1013
1014
1015 /*
1016 * Apply flat field and apply the relative fiber transmission correction.
1017 */
1018
1019 flat_config = giraffe_flat_config_create(config);
1020
1021 if (flat_config->load == TRUE) {
1022
1023 cpl_frame* flat_frame = NULL;
1024
1025 GiImage* flat = NULL;
1026
1027
1028 flat_frame = cpl_frameset_find(set, GIFRAME_FIBER_FLAT_EXTSPECTRA);
1029
1030 if (flat_frame == NULL) {
1031 cpl_msg_error(_id, "Missing flat field spectra frame!");
1032
1033 giraffe_flat_config_destroy(flat_config);
1034
1035 giraffe_extraction_destroy(extraction);
1036 giraffe_localization_destroy(localization);
1037
1038 giraffe_table_delete(wcalcoeff);
1039
1040 giraffe_table_delete(grating);
1041 giraffe_table_delete(fibers);
1042
1043 return 1;
1044 }
1045
1046 filename = cpl_frame_get_filename(flat_frame);
1047
1048 flat = giraffe_image_new(CPL_TYPE_DOUBLE);
1049 status = giraffe_image_load(flat, filename, 0);
1050
1051 if (status) {
1052 cpl_msg_error(_id, "Cannot load flat field spectra from '%s'. "
1053 "Aborting ...", filename);
1054
1055 giraffe_image_delete(flat);
1056
1057 giraffe_flat_config_destroy(flat_config);
1058
1059 giraffe_extraction_destroy(extraction);
1060 giraffe_localization_destroy(localization);
1061
1062 giraffe_table_delete(wcalcoeff);
1063
1064 giraffe_table_delete(grating);
1065 giraffe_table_delete(fibers);
1066
1067 return 1;
1068 }
1069
1070 if (flat_config->apply == TRUE) {
1071
1072 GiImage* errors = NULL;
1073
1074
1075 flat_frame = cpl_frameset_find(set, GIFRAME_FIBER_FLAT_EXTERRORS);
1076
1077 if (flat_frame == NULL) {
1078 cpl_msg_warning(_id, "Missing flat field spectra errors "
1079 "frame!");
1080 }
1081 else {
1082
1083 filename = cpl_frame_get_filename(flat_frame);
1084
1085 errors = giraffe_image_new(CPL_TYPE_DOUBLE);
1086 status = giraffe_image_load(errors, filename, 0);
1087
1088 if (status) {
1089 cpl_msg_error(_id, "Cannot load flat field spectra "
1090 "errors from '%s'. Aborting ...",
1091 filename);
1092
1093 giraffe_image_delete(errors);
1094 giraffe_image_delete(flat);
1095
1096 giraffe_flat_config_destroy(flat_config);
1097
1098 giraffe_extraction_destroy(extraction);
1099 giraffe_localization_destroy(localization);
1100
1101 giraffe_table_delete(wcalcoeff);
1102
1103 giraffe_table_delete(grating);
1104 giraffe_table_delete(fibers);
1105
1106 return 1;
1107 }
1108
1109 }
1110
1111 cpl_msg_info(_id, "Applying flat field correction ...");
1112
1113 status = giraffe_flat_apply(extraction, fibers, flat, errors,
1114 flat_config);
1115
1116 if (status) {
1117 cpl_msg_error(_id, "Flat field correction failed! "
1118 "Aborting ...");
1119
1120 giraffe_image_delete(errors);
1121 giraffe_image_delete(flat);
1122
1123 giraffe_flat_config_destroy(flat_config);
1124
1125 giraffe_extraction_destroy(extraction);
1126 giraffe_localization_destroy(localization);
1127
1128 giraffe_table_delete(wcalcoeff);
1129
1130 giraffe_table_delete(grating);
1131 giraffe_table_delete(fibers);
1132
1133 return 1;
1134 }
1135
1136 giraffe_image_delete(errors);
1137 errors = NULL;
1138
1139 }
1140
1141 if (flat_config->transmission == TRUE) {
1142
1143 const cxchar* _filename = cpl_frame_get_filename(flat_frame);
1144
1145 GiTable* _fibers = NULL;
1146
1147
1148 cpl_msg_info(_id, "Loading fiber setup for frame '%s'.",
1149 _filename);
1150
1151 _fibers = giraffe_fiberlist_load(_filename, 1, "FIBER_SETUP");
1152
1153 if (!_fibers) {
1154 cpl_msg_error(_id, "Cannot create fiber setup for "
1155 "frame '%s'! Aborting ...", _filename);
1156
1157 giraffe_image_delete(flat);
1158
1159 giraffe_flat_config_destroy(flat_config);
1160
1161 giraffe_extraction_destroy(extraction);
1162 giraffe_localization_destroy(localization);
1163
1164 giraffe_table_delete(wcalcoeff);
1165
1166 giraffe_table_delete(grating);
1167 giraffe_table_delete(fibers);
1168
1169 return 1;
1170 }
1171
1172 cpl_msg_info(_id, "Applying relative fiber transmission "
1173 "correction");
1174
1175 status = giraffe_transmission_setup(fibers, _fibers);
1176 giraffe_table_delete(_fibers);
1177
1178 if (status == 0) {
1179 status = giraffe_transmission_apply(extraction, fibers);
1180 }
1181
1182 if (status) {
1183
1184 cpl_msg_error(_id, "Relative transmission correction failed! "
1185 "Aborting ...");
1186
1187 giraffe_image_delete(flat);
1188
1189 giraffe_flat_config_destroy(flat_config);
1190
1191 giraffe_extraction_destroy(extraction);
1192 giraffe_localization_destroy(localization);
1193
1194 giraffe_table_delete(wcalcoeff);
1195
1196 giraffe_table_delete(grating);
1197 giraffe_table_delete(fibers);
1198
1199 return 1;
1200
1201 }
1202
1203 }
1204
1205 giraffe_image_delete(flat);
1206
1207 }
1208
1209 giraffe_flat_config_destroy(flat_config);
1210
1211
1212 /*
1213 * Save the spectrum extraction results and register them as
1214 * products.
1215 */
1216
1217 cpl_msg_info(_id, "Writing extracted spectra ...");
1218
1219 /* Extracted spectra */
1220
1221 giraffe_image_add_info(extraction->spectra, &info, set);
1222
1223 cpl_propertylist* extprop = giraffe_image_get_properties(extraction->spectra);
1224
1225 giraffe_qc_update_sci_props(extprop, fibers);
1226
1227
1228 sext_frame = giraffe_frame_create_image(extraction->spectra,
1229 GIFRAME_SCIENCE_EXTSPECTRA,
1230 CPL_FRAME_LEVEL_FINAL,
1231 TRUE, TRUE);
1232
1233 if (sext_frame == NULL) {
1234 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
1235
1236 giraffe_extraction_destroy(extraction);
1237 giraffe_localization_destroy(localization);
1238
1239 giraffe_table_delete(wcalcoeff);
1240
1241 giraffe_table_delete(grating);
1242 giraffe_table_delete(fibers);
1243
1244 return 1;
1245 }
1246
1247 status = giraffe_fiberlist_attach(sext_frame, fibers);
1248
1249 if (status) {
1250 cpl_msg_error(_id, "Cannot attach fiber setup to local file '%s'! "
1251 "Aborting ...", cpl_frame_get_filename(sext_frame));
1252
1253 cpl_frame_delete(sext_frame);
1254
1255 giraffe_extraction_destroy(extraction);
1256 giraffe_localization_destroy(localization);
1257
1258 giraffe_table_delete(wcalcoeff);
1259
1260 giraffe_table_delete(grating);
1261 giraffe_table_delete(fibers);
1262
1263 return 1;
1264 }
1265
1266 cpl_frameset_insert(set, sext_frame);
1267
1268 /* Extracted spectra errors */
1269
1270
1271
1272 cpl_propertylist* exterrprop = giraffe_image_get_properties(extraction->error);
1273
1274 giraffe_qc_update_sci_props(exterrprop, fibers);
1275
1276
1277 giraffe_image_add_info(extraction->error, &info, set);
1278
1279 sext_frame = giraffe_frame_create_image(extraction->error,
1280 GIFRAME_SCIENCE_EXTERRORS,
1281 CPL_FRAME_LEVEL_FINAL,
1282 TRUE, TRUE);
1283
1284 if (sext_frame == NULL) {
1285 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
1286
1287 giraffe_extraction_destroy(extraction);
1288 giraffe_localization_destroy(localization);
1289
1290 giraffe_table_delete(wcalcoeff);
1291
1292 giraffe_table_delete(grating);
1293 giraffe_table_delete(fibers);
1294
1295 return 1;
1296 }
1297
1298 status = giraffe_fiberlist_attach(sext_frame, fibers);
1299
1300 if (status) {
1301 cpl_msg_error(_id, "Cannot attach fiber setup to local file '%s'! "
1302 "Aborting ...", cpl_frame_get_filename(sext_frame));
1303
1304 cpl_frame_delete(sext_frame);
1305
1306 giraffe_extraction_destroy(extraction);
1307 giraffe_localization_destroy(localization);
1308
1309 giraffe_table_delete(wcalcoeff);
1310
1311 giraffe_table_delete(grating);
1312 giraffe_table_delete(fibers);
1313
1314 return 1;
1315 }
1316
1317 cpl_frameset_insert(set, sext_frame);
1318
1319 /* Extracted spectra pixels */
1320
1321 if (extraction->npixels != NULL) {
1322
1323 cpl_propertylist* extpixprop = giraffe_image_get_properties(extraction->npixels);
1324
1325 giraffe_qc_update_sci_props(extpixprop, fibers);
1326
1327 giraffe_image_add_info(extraction->npixels, &info, set);
1328
1329 sext_frame = giraffe_frame_create_image(extraction->npixels,
1330 GIFRAME_SCIENCE_EXTPIXELS,
1331 CPL_FRAME_LEVEL_FINAL,
1332 TRUE, TRUE);
1333
1334 if (sext_frame == NULL) {
1335 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
1336
1337 giraffe_extraction_destroy(extraction);
1338 giraffe_localization_destroy(localization);
1339
1340 giraffe_table_delete(wcalcoeff);
1341
1342 giraffe_table_delete(grating);
1343 giraffe_table_delete(fibers);
1344
1345 return 1;
1346 }
1347
1348 status = giraffe_fiberlist_attach(sext_frame, fibers);
1349
1350 if (status) {
1351 cpl_msg_error(_id, "Cannot attach fiber setup to local file '%s'! "
1352 "Aborting ...", cpl_frame_get_filename(sext_frame));
1353
1354 cpl_frame_delete(sext_frame);
1355
1356 giraffe_extraction_destroy(extraction);
1357 giraffe_localization_destroy(localization);
1358
1359 giraffe_table_delete(wcalcoeff);
1360
1361 giraffe_table_delete(grating);
1362 giraffe_table_delete(fibers);
1363
1364 return 1;
1365 }
1366
1367 cpl_frameset_insert(set, sext_frame);
1368
1369 }
1370
1371 /* Extracted spectra centroids */
1372
1373 giraffe_image_add_info(extraction->centroid, &info, set);
1374
1375 sext_frame = giraffe_frame_create_image(extraction->centroid,
1376 GIFRAME_SCIENCE_EXTTRACE,
1377 CPL_FRAME_LEVEL_FINAL,
1378 TRUE, TRUE);
1379
1380 if (sext_frame == NULL) {
1381 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
1382
1383 giraffe_extraction_destroy(extraction);
1384 giraffe_localization_destroy(localization);
1385
1386 giraffe_table_delete(wcalcoeff);
1387
1388 giraffe_table_delete(grating);
1389 giraffe_table_delete(fibers);
1390
1391 return 1;
1392 }
1393
1394 status = giraffe_fiberlist_attach(sext_frame, fibers);
1395
1396 if (status) {
1397 cpl_msg_error(_id, "Cannot attach fiber setup to local file '%s'! "
1398 "Aborting ...", cpl_frame_get_filename(sext_frame));
1399
1400 cpl_frame_delete(sext_frame);
1401
1402 giraffe_extraction_destroy(extraction);
1403 giraffe_localization_destroy(localization);
1404
1405 giraffe_table_delete(wcalcoeff);
1406
1407 giraffe_table_delete(grating);
1408 giraffe_table_delete(fibers);
1409
1410 return 1;
1411 }
1412
1413 cpl_frameset_insert(set, sext_frame);
1414
1415 /* Extraction model spectra */
1416
1417 if (extraction->model != NULL) {
1418
1419 giraffe_image_add_info(extraction->model, &info, set);
1420
1421 sext_frame = giraffe_frame_create_image(extraction->model,
1422 GIFRAME_SCIENCE_EXTMODEL,
1423 CPL_FRAME_LEVEL_FINAL,
1424 TRUE, TRUE);
1425
1426 if (sext_frame == NULL) {
1427 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
1428
1429 giraffe_extraction_destroy(extraction);
1430 giraffe_localization_destroy(localization);
1431
1432 giraffe_table_delete(wcalcoeff);
1433
1434 giraffe_table_delete(grating);
1435 giraffe_table_delete(fibers);
1436
1437 return 1;
1438 }
1439
1440 status = giraffe_fiberlist_attach(sext_frame, fibers);
1441
1442 if (status != 0) {
1443 cpl_msg_error(_id, "Cannot attach fiber setup to local file '%s'! "
1444 "Aborting ...", cpl_frame_get_filename(sext_frame));
1445
1446 cpl_frame_delete(sext_frame);
1447
1448 giraffe_extraction_destroy(extraction);
1449 giraffe_localization_destroy(localization);
1450
1451 giraffe_table_delete(wcalcoeff);
1452
1453 giraffe_table_delete(grating);
1454 giraffe_table_delete(fibers);
1455
1456 return 1;
1457 }
1458
1459 cpl_frameset_insert(set, sext_frame);
1460
1461 }
1462
1463
1464 /*
1465 * Load dispersion solution
1466 */
1467
1468
1469 filename = (cxchar *)cpl_frame_get_filename(wcal_frame);
1470
1471 tempes = cpl_errorstate_get();
1472
1473 const char qcdltcopy[] = "^(ESO INS TEMP53 VAL|MJD-OBS)$";
1474 cpl_propertylist * qcdltlist = cpl_propertylist_load_regexp(filename, 0,
1475 qcdltcopy, 0);
1476 cpl_errorstate_set(tempes);
1477
1478 wcalcoeff = giraffe_table_new();
1479 status = giraffe_table_load(wcalcoeff, filename, 1, NULL);
1480
1481 if (status) {
1482 cpl_msg_error(_id, "Cannot load dispersion solution from "
1483 "'%s'. Aborting ...", filename);
1484
1485 giraffe_extraction_destroy(extraction);
1486 giraffe_localization_destroy(localization);
1487
1488 giraffe_table_delete(wcalcoeff);
1489
1490 giraffe_table_delete(grating);
1491 giraffe_table_delete(fibers);
1492
1493 return 1;
1494 }
1495
1496
1497 /*
1498 * Load grating data
1499 */
1500
1501 filename = (cxchar *)cpl_frame_get_filename(grating_frame);
1502
1503
1504 grating = giraffe_table_new();
1505 status = giraffe_table_load(grating, filename, 1, NULL);
1506
1507 if (status) {
1508 cpl_msg_error(_id, "Cannot load grating data from '%s'. "
1509 "Aborting ...", filename);
1510
1511 giraffe_extraction_destroy(extraction);
1512 giraffe_localization_destroy(localization);
1513
1514 giraffe_table_delete(wcalcoeff);
1515
1516 giraffe_table_delete(grating);
1517 giraffe_table_delete(fibers);
1518
1519 return 1;
1520 }
1521
1522
1523 /*
1524 * Load slit geometry data
1525 */
1526
1527
1528 filename = (cxchar *)cpl_frame_get_filename(slit_frame);
1529
1530 slitgeometry = giraffe_slitgeometry_load(fibers, filename, 1, NULL);
1531
1532 if (slitgeometry == NULL) {
1533 cpl_msg_error(_id, "Cannot load slit geometry data from '%s'. "
1534 "Aborting ...", filename);
1535
1536 giraffe_table_delete(wcalcoeff);
1537
1538 giraffe_extraction_destroy(extraction);
1539 giraffe_localization_destroy(localization);
1540
1541 giraffe_table_delete(wcalcoeff);
1542
1543 giraffe_table_delete(grating);
1544 giraffe_table_delete(fibers);
1545
1546 return 1;
1547 }
1548 else {
1549
1550 /*
1551 * Check whether the contains the positions for all fibers
1552 * provided by the fiber setup. If this is not the case
1553 * this is an error.
1554 */
1555
1556 if (giraffe_fiberlist_compare(slitgeometry, fibers) != 1) {
1557 cpl_msg_error(_id, "Slit geometry data from '%s' is not "
1558 "applicable for current fiber setup! "
1559 "Aborting ...", filename);
1560
1561 giraffe_table_delete(slitgeometry);
1562 giraffe_table_delete(wcalcoeff);
1563
1564 giraffe_extraction_destroy(extraction);
1565 giraffe_localization_destroy(localization);
1566
1567 giraffe_table_delete(wcalcoeff);
1568
1569 giraffe_table_delete(grating);
1570 giraffe_table_delete(fibers);
1571
1572 return 1;
1573 }
1574
1575 }
1576
1577
1578
1579 /*
1580 * Spectrum rebinning
1581 */
1582
1583 cpl_msg_info(_id, "Spectrum rebinning");
1584
1585 rebin_config = giraffe_rebin_config_create(config);
1586
1587 rebinning = giraffe_rebinning_new();
1588
1589 status = giraffe_rebin_spectra(rebinning, extraction, fibers,
1590 localization, grating, slitgeometry,
1591 wcalcoeff, rebin_config);
1592
1593 if (status) {
1594 cpl_msg_error(_id, "Rebinning of science spectra failed! Aborting...");
1595
1596 giraffe_rebinning_destroy(rebinning);
1597
1598 giraffe_extraction_destroy(extraction);
1599 giraffe_localization_destroy(localization);
1600
1601 giraffe_table_delete(wcalcoeff);
1602
1603 giraffe_table_delete(slitgeometry);
1604 giraffe_table_delete(grating);
1605 giraffe_table_delete(fibers);
1606
1607 giraffe_rebin_config_destroy(rebin_config);
1608
1609 return 1;
1610
1611 }
1612
1613
1614 /*
1615 * Optionally compute and apply spectral shifts from the simultaneous
1616 * calibration fibers. This is only done if the simultaneous calibration
1617 * fibers were used.
1618 */
1619
1620 p = cpl_parameterlist_find(config, "giraffe.siwc.apply");
1621 cx_assert(p != NULL);
1622
1623 siwc = cpl_parameter_get_bool(p);
1624 p = NULL;
1625
1626 properties = giraffe_image_get_properties(rebinning->spectra);
1627 cx_assert(properties != NULL);
1628
1629
1630 if (cpl_propertylist_has(properties, GIALIAS_STSCTAL) == TRUE) {
1631 calsim = cpl_propertylist_get_bool(properties, GIALIAS_STSCTAL);
1632 }
1633
1634
1635 if ((siwc == TRUE) && (calsim == TRUE) && (linemask_frame != NULL)) {
1636
1637 GiTable* linemask = giraffe_table_new();
1638
1639 GiSGCalConfig* siwc_config = NULL;
1640
1641
1642 siwc_config = giraffe_sgcalibration_config_create(config);
1643
1644 if (siwc_config == NULL) {
1645
1646 giraffe_table_delete(linemask);
1647 linemask = NULL;
1648
1649 giraffe_rebinning_destroy(rebinning);
1650
1651 giraffe_extraction_destroy(extraction);
1652 giraffe_localization_destroy(localization);
1653
1654 giraffe_table_delete(wcalcoeff);
1655
1656 giraffe_table_delete(slitgeometry);
1657 giraffe_table_delete(grating);
1658 giraffe_table_delete(fibers);
1659
1660 giraffe_rebin_config_destroy(rebin_config);
1661
1662 return 1;
1663
1664 }
1665
1666 filename = cpl_frame_get_filename(linemask_frame);
1667
1668 status = giraffe_table_load(linemask, filename, 1, NULL);
1669
1670 if (status) {
1671 cpl_msg_error(_id, "Cannot load line reference mask from '%s'. "
1672 "Aborting ...", filename);
1673
1674 giraffe_sgcalibration_config_destroy(siwc_config);
1675 siwc_config = NULL;
1676
1677 giraffe_table_delete(linemask);
1678 linemask = NULL;
1679
1680 giraffe_rebinning_destroy(rebinning);
1681
1682 giraffe_extraction_destroy(extraction);
1683 giraffe_localization_destroy(localization);
1684
1685 giraffe_table_delete(wcalcoeff);
1686
1687 giraffe_table_delete(slitgeometry);
1688 giraffe_table_delete(grating);
1689 giraffe_table_delete(fibers);
1690
1691 giraffe_rebin_config_destroy(rebin_config);
1692
1693 return 1;
1694
1695 }
1696
1697
1698 status = giraffe_compute_offsets(fibers, rebinning, grating,
1699 linemask, siwc_config);
1700
1701 if (status != 0) {
1702 cpl_msg_error(_id, "Applying simultaneous wavelength "
1703 "calibration correction failed! Aborting...");
1704
1705 giraffe_sgcalibration_config_destroy(siwc_config);
1706 siwc_config = NULL;
1707
1708 giraffe_table_delete(linemask);
1709 linemask = NULL;
1710
1711 giraffe_rebinning_destroy(rebinning);
1712
1713 giraffe_extraction_destroy(extraction);
1714 giraffe_localization_destroy(localization);
1715
1716 giraffe_table_delete(wcalcoeff);
1717
1718 giraffe_table_delete(slitgeometry);
1719 giraffe_table_delete(grating);
1720 giraffe_table_delete(fibers);
1721
1722 giraffe_rebin_config_destroy(rebin_config);
1723
1724 return 1;
1725
1726 }
1727
1728 giraffe_sgcalibration_config_destroy(siwc_config);
1729 siwc_config = NULL;
1730
1731 giraffe_table_delete(linemask);
1732 linemask = NULL;
1733
1734 giraffe_rebinning_destroy(rebinning);
1735 rebinning = giraffe_rebinning_new();
1736
1737 status = giraffe_rebin_spectra(rebinning, extraction, fibers,
1738 localization, grating, slitgeometry,
1739 wcalcoeff, rebin_config);
1740
1741 if (status) {
1742 cpl_msg_error(_id, "Rebinning of science spectra failed! "
1743 "Aborting...");
1744
1745 giraffe_rebinning_destroy(rebinning);
1746
1747 giraffe_extraction_destroy(extraction);
1748 giraffe_localization_destroy(localization);
1749
1750 giraffe_table_delete(wcalcoeff);
1751
1752 giraffe_table_delete(slitgeometry);
1753 giraffe_table_delete(grating);
1754 giraffe_table_delete(fibers);
1755
1756 giraffe_rebin_config_destroy(rebin_config);
1757
1758 return 1;
1759
1760 }
1761
1762 }
1763
1764 giraffe_extraction_destroy(extraction);
1765 extraction = NULL;
1766
1767 giraffe_localization_destroy(localization);
1768 localization = NULL;
1769
1770 giraffe_rebin_config_destroy(rebin_config);
1771 rebin_config = NULL;
1772
1773
1774 /*
1775 * Compute barycentric correction for each object spectrum (fiber)
1776 */
1777
1778 status = giraffe_add_rvcorrection(fibers, rebinning->spectra);
1779
1780 switch (status) {
1781 case 0:
1782 {
1783 break;
1784 }
1785
1786 case 1:
1787 {
1788 cpl_msg_warning(_id, "Missing observation time properties! "
1789 "Barycentric correction computation "
1790 "skipped!");
1791 status = 0;
1792 break;
1793 }
1794 case 2:
1795 {
1796 cpl_msg_warning(_id, "Missing telescope location properties! "
1797 "Barycentric correction computation "
1798 "skipped!");
1799 status = 0;
1800 break;
1801 }
1802 case 3:
1803 {
1804 cpl_msg_warning(_id, "Object positions are not available "
1805 "Barycentric correction computation "
1806 "skipped!");
1807 status = 0;
1808 break;
1809 }
1810 default:
1811 {
1812 cpl_msg_error(_id, "Barycentric correction computation "
1813 "failed! Aborting...");
1814
1815 giraffe_rebinning_destroy(rebinning);
1816
1817 giraffe_table_delete(wcalcoeff);
1818
1819 giraffe_table_delete(slitgeometry);
1820 giraffe_table_delete(grating);
1821 giraffe_table_delete(fibers);
1822
1823 return 1;
1824 break;
1825 }
1826
1827 }
1828
1829
1830 /*
1831 * Save and register the results of the spectrum rebinning.
1832 */
1833
1834 /* Rebinned spectra */
1835
1836 /*
1837 * Add the QC parameters to the rebinning spectra
1838 */
1839
1840 tempes = cpl_errorstate_get();
1841
1842 cpl_propertylist* rbprop = giraffe_image_get_properties(rebinning->spectra);
1843
1844
1845 giraffe_propertylist_copy(rbprop, GIALIAS_WLSTART, rbprop,
1846 GIALIAS_BINWLMIN);
1847 giraffe_propertylist_copy(rbprop, GIALIAS_WLEND, rbprop,
1848 GIALIAS_BINWLMAX);
1849 giraffe_propertylist_copy(rbprop, GIALIAS_WLSTEP, rbprop,
1850 GIALIAS_BINSTEP);
1851
1852 double dlttemp = cpl_propertylist_get_double(rbprop, "ESO INS TEMP53 VAL");
1853 dlttemp -= cpl_propertylist_get_double(qcdltlist, "ESO INS TEMP53 VAL");
1854
1855 double dltdate = cpl_propertylist_get_double(rbprop, "MJD-OBS");
1856 dltdate -= cpl_propertylist_get_double(qcdltlist, "MJD-OBS");
1857
1858 cpl_propertylist_delete(qcdltlist);
1859
1860 cpl_propertylist_append_double(rbprop, GIALIAS_QCDLTTEMP, dlttemp);
1861 cpl_propertylist_append_double(rbprop, GIALIAS_QCDLTTIME, dltdate);
1862
1863
1864 cpl_errorstate_set(tempes);
1865
1866 giraffe_image_add_info(rebinning->spectra, &info, set);
1867
1868 rbin_frame = giraffe_frame_create_image(rebinning->spectra,
1869 GIFRAME_SCIENCE_RBNSPECTRA,
1870 CPL_FRAME_LEVEL_FINAL,
1871 TRUE, TRUE);
1872
1873 if (rbin_frame == NULL) {
1874 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
1875
1876 giraffe_rebinning_destroy(rebinning);
1877
1878 giraffe_table_delete(wcalcoeff);
1879
1880 giraffe_table_delete(slitgeometry);
1881 giraffe_table_delete(grating);
1882 giraffe_table_delete(fibers);
1883
1884 return 1;
1885 }
1886
1887 status = giraffe_fiberlist_attach(rbin_frame, fibers);
1888
1889 if (status) {
1890 cpl_msg_error(_id, "Cannot attach fiber setup to local "
1891 "file '%s'! Aborting ...",
1892 cpl_frame_get_filename(rbin_frame));
1893
1894 giraffe_rebinning_destroy(rebinning);
1895 giraffe_table_delete(wcalcoeff);
1896
1897 giraffe_table_delete(slitgeometry);
1898 giraffe_table_delete(grating);
1899 giraffe_table_delete(fibers);
1900
1901 cpl_frame_delete(rbin_frame);
1902
1903 return 1;
1904 }
1905
1906 cpl_frameset_insert(set, rbin_frame);
1907 flux_filename = cpl_strdup(cpl_frame_get_filename(rbin_frame));
1908
1909 /* Rebinned spectra errors */
1910
1911 tempes = cpl_errorstate_get();
1912
1913 cpl_propertylist* rberrprop = giraffe_image_get_properties(rebinning->errors);
1914
1915
1916 giraffe_propertylist_copy(rberrprop, GIALIAS_WLSTART, rberrprop,
1917 GIALIAS_BINWLMIN);
1918 giraffe_propertylist_copy(rberrprop, GIALIAS_WLEND, rberrprop,
1919 GIALIAS_BINWLMAX);
1920 giraffe_propertylist_copy(rberrprop, GIALIAS_WLSTEP, rberrprop,
1921 GIALIAS_BINSTEP);
1922
1923 cpl_propertylist_append_double(rberrprop, GIALIAS_QCDLTTEMP, dlttemp);
1924 cpl_propertylist_append_double(rberrprop, GIALIAS_QCDLTTIME, dltdate);
1925
1926
1927 giraffe_image_add_info(rebinning->errors, &info, set);
1928
1929 cpl_errorstate_set(tempes);
1930
1931 rbin_frame = giraffe_frame_create_image(rebinning->errors,
1932 GIFRAME_SCIENCE_RBNERRORS,
1933 CPL_FRAME_LEVEL_FINAL,
1934 TRUE, TRUE);
1935
1936 if (rbin_frame == NULL) {
1937 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
1938
1939 giraffe_rebinning_destroy(rebinning);
1940
1941 giraffe_table_delete(wcalcoeff);
1942
1943 giraffe_table_delete(slitgeometry);
1944 giraffe_table_delete(grating);
1945 giraffe_table_delete(fibers);
1946 cpl_free(flux_filename);
1947
1948 return 1;
1949 }
1950
1951 status = giraffe_fiberlist_attach(rbin_frame, fibers);
1952
1953 if (status) {
1954 cpl_msg_error(_id, "Cannot attach fiber setup to local "
1955 "file '%s'! Aborting ...",
1956 cpl_frame_get_filename(rbin_frame));
1957
1958 giraffe_rebinning_destroy(rebinning);
1959
1960 giraffe_table_delete(wcalcoeff);
1961
1962 giraffe_table_delete(slitgeometry);
1963 giraffe_table_delete(grating);
1964 giraffe_table_delete(fibers);
1965
1966 cpl_frame_delete(rbin_frame);
1967 cpl_free(flux_filename);
1968
1969 return 1;
1970 }
1971
1972 cpl_frameset_insert(set, rbin_frame);
1973 err_filename = cpl_strdup(cpl_frame_get_filename(rbin_frame));
1974
1975 properties = giraffe_image_get_properties(rebinning->spectra);
1976 mode = giraffe_get_mode(properties);
1977
1978 /*
1979 * Optionally generate spectra in Science Data Product (SDP) format.
1980 */
1981 p = cpl_parameterlist_find(config, "giraffe.sdp.format.generate");
1982 cx_assert(p != NULL);
1983 gensdp = cpl_parameter_get_bool(p);
1984 p = cpl_parameterlist_find(config, "giraffe.sdp.nassoc.keys");
1985 cx_assert(p != NULL);
1986 nassoc_keys = cpl_parameter_get_int(p);
1987 p = NULL;
1988 if (gensdp) {
1989 if (mode == GIMODE_MEDUSA) {
1990 status = _giraffe_make_sdp_spectra(flux_filename, err_filename,
1991 nassoc_keys, set, config, _id);
1992 if (status) {
1993 cpl_msg_error(_id, "Failed to generate spectra in Science Data"
1994 " Product format.");
1995
1996 giraffe_rebinning_destroy(rebinning);
1997
1998 giraffe_table_delete(wcalcoeff);
1999
2000 giraffe_table_delete(slitgeometry);
2001 giraffe_table_delete(grating);
2002 giraffe_table_delete(fibers);
2003
2004 cpl_free(flux_filename);
2005 cpl_free(err_filename);
2006
2007 return 1;
2008 }
2009 } else {
2010 cpl_msg_warning(_id, "Requested to generate SDP 1D spectra, but"
2011 " this is currently only supported for the MEDUSA"
2012 " mode. Skipping SDP generation.");
2013 }
2014 }
2015 cpl_free(flux_filename);
2016 cpl_free(err_filename);
2017
2018
2019 /*
2020 * Optional image and data cube construction (only for IFU and Argus)
2021 */
2022
2023 if (mode == GIMODE_IFU || mode == GIMODE_ARGUS) {
2024
2025 cpl_frame* rimg_frame = NULL;
2026
2027 GiFieldOfView* fov = NULL;
2028
2029 GiFieldOfViewConfig* fov_config = NULL;
2030
2031 GiFieldOfViewCubeFormat cube_format = GIFOV_FORMAT_ESO3D;
2032
2033
2034 fov_config = giraffe_fov_config_create(config);
2035
2036 cube_format = fov_config->format;
2037
2038
2039 cpl_msg_info(_id, "Reconstructing image and data cube from rebinned "
2040 "spectra ...");
2041
2042 fov = giraffe_fov_new();
2043
2044 status = giraffe_fov_build(fov, rebinning, fibers, wcalcoeff, grating,
2045 slitgeometry, fov_config);
2046
2047 if (status) {
2048
2049 if (status == -2) {
2050 cpl_msg_warning(_id, "No reconstructed image was built. "
2051 "Fiber list has no fiber position "
2052 "information.");
2053 }
2054 else {
2055 cpl_msg_error(_id, "Image reconstruction failed! Aborting...");
2056
2057 giraffe_fov_delete(fov);
2058 giraffe_rebinning_destroy(rebinning);
2059
2060 giraffe_table_delete(wcalcoeff);
2061
2062 giraffe_table_delete(slitgeometry);
2063 giraffe_table_delete(grating);
2064 giraffe_table_delete(fibers);
2065
2066 giraffe_fov_config_destroy(fov_config);
2067
2068 return 1;
2069 }
2070
2071 }
2072
2073 giraffe_fov_config_destroy(fov_config);
2074
2075
2076 /*
2077 * Save and register the results of the image reconstruction.
2078 */
2079
2080 /* Reconstructed image */
2081
2082 giraffe_image_add_info(fov->fov.spectra, &info, set);
2083
2084 rimg_frame = giraffe_frame_create_image(fov->fov.spectra,
2085 GIFRAME_SCIENCE_RCSPECTRA,
2086 CPL_FRAME_LEVEL_FINAL,
2087 TRUE, TRUE);
2088
2089 if (rimg_frame == NULL) {
2090 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
2091
2092 giraffe_fov_delete(fov);
2093 giraffe_rebinning_destroy(rebinning);
2094
2095 giraffe_table_delete(wcalcoeff);
2096
2097 giraffe_table_delete(slitgeometry);
2098 giraffe_table_delete(grating);
2099 giraffe_table_delete(fibers);
2100
2101 return 1;
2102 }
2103
2104 cpl_frameset_insert(set, rimg_frame);
2105
2106
2107 /* Reconstructed image errors */
2108
2109 giraffe_image_add_info(fov->fov.errors, &info, set);
2110
2111 rimg_frame = giraffe_frame_create_image(fov->fov.errors,
2112 GIFRAME_SCIENCE_RCERRORS,
2113 CPL_FRAME_LEVEL_FINAL,
2114 TRUE, TRUE);
2115
2116 if (rimg_frame == NULL) {
2117 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
2118
2119 giraffe_fov_delete(fov);
2120 giraffe_rebinning_destroy(rebinning);
2121
2122 giraffe_table_delete(wcalcoeff);
2123
2124 giraffe_table_delete(slitgeometry);
2125 giraffe_table_delete(grating);
2126 giraffe_table_delete(fibers);
2127
2128 return 1;
2129 }
2130
2131 cpl_frameset_insert(set, rimg_frame);
2132
2133
2134 /* Save data cubes according to format selection */
2135
2136 if (cube_format == GIFOV_FORMAT_SINGLE) {
2137
2138 /* Spectrum cube */
2139
2140 if (fov->cubes.spectra != NULL) {
2141
2142 cxint component = 0;
2143
2144 GiFrameCreator creator = (GiFrameCreator) giraffe_fov_save_cubes;
2145
2146
2147 properties = giraffe_image_get_properties(rebinning->spectra);
2148 properties = cpl_propertylist_duplicate(properties);
2149
2150 giraffe_add_frameset_info(properties, set, info.sequence);
2151
2152 rimg_frame = giraffe_frame_create(GIFRAME_SCIENCE_CUBE_SPECTRA,
2153 CPL_FRAME_LEVEL_FINAL,
2154 properties,
2155 fov,
2156 &component,
2157 creator);
2158
2159 cpl_propertylist_delete(properties);
2160 properties = NULL;
2161
2162 if (rimg_frame == NULL) {
2163 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
2164
2165 giraffe_fov_delete(fov);
2166 fov = NULL;
2167
2168 giraffe_rebinning_destroy(rebinning);
2169 rebinning = NULL;
2170
2171 giraffe_table_delete(wcalcoeff);
2172 wcalcoeff = NULL;
2173
2174 giraffe_table_delete(slitgeometry);
2175 slitgeometry = NULL;
2176
2177 giraffe_table_delete(grating);
2178 grating = NULL;
2179
2180 giraffe_table_delete(fibers);
2181 fibers = NULL;
2182
2183 return 1;
2184 }
2185
2186 status = giraffe_fiberlist_attach(rimg_frame, fibers);
2187
2188 if (status != 0) {
2189 cpl_msg_error(_id, "Cannot attach fiber setup to local "
2190 "file '%s'! Aborting ...",
2191 cpl_frame_get_filename(rimg_frame));
2192
2193 cpl_frame_delete(rimg_frame);
2194
2195 giraffe_fov_delete(fov);
2196 fov = NULL;
2197
2198 giraffe_rebinning_destroy(rebinning);
2199 rebinning = NULL;
2200
2201 giraffe_table_delete(wcalcoeff);
2202 wcalcoeff = NULL;
2203
2204 giraffe_table_delete(slitgeometry);
2205 slitgeometry = NULL;
2206
2207 giraffe_table_delete(grating);
2208 grating = NULL;
2209
2210 giraffe_table_delete(fibers);
2211 fibers = NULL;
2212
2213 return 1;
2214 }
2215
2216 cpl_frameset_insert(set, rimg_frame);
2217
2218 }
2219
2220 /* Error cube */
2221
2222 if (fov->cubes.errors != NULL) {
2223
2224 cxint component = 1;
2225
2226 GiFrameCreator creator = (GiFrameCreator) giraffe_fov_save_cubes;
2227
2228
2229 properties = giraffe_image_get_properties(rebinning->errors);
2230 properties = cpl_propertylist_duplicate(properties);
2231
2232 giraffe_add_frameset_info(properties, set, info.sequence);
2233
2234 rimg_frame = giraffe_frame_create(GIFRAME_SCIENCE_CUBE_ERRORS,
2235 CPL_FRAME_LEVEL_FINAL,
2236 properties,
2237 fov,
2238 &component,
2239 creator);
2240
2241 cpl_propertylist_delete(properties);
2242 properties = NULL;
2243
2244 if (rimg_frame == NULL) {
2245 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
2246
2247 giraffe_fov_delete(fov);
2248 fov = NULL;
2249
2250 giraffe_rebinning_destroy(rebinning);
2251 rebinning = NULL;
2252
2253 giraffe_table_delete(wcalcoeff);
2254 wcalcoeff = NULL;
2255
2256 giraffe_table_delete(slitgeometry);
2257 slitgeometry = NULL;
2258
2259 giraffe_table_delete(grating);
2260 grating = NULL;
2261
2262 giraffe_table_delete(fibers);
2263 fibers = NULL;
2264
2265 return 1;
2266 }
2267
2268 status = giraffe_fiberlist_attach(rimg_frame, fibers);
2269
2270 if (status != 0) {
2271 cpl_msg_error(_id, "Cannot attach fiber setup to local "
2272 "file '%s'! Aborting ...",
2273 cpl_frame_get_filename(rimg_frame));
2274
2275 cpl_frame_delete(rimg_frame);
2276
2277 giraffe_fov_delete(fov);
2278 fov = NULL;
2279
2280 giraffe_rebinning_destroy(rebinning);
2281 rebinning = NULL;
2282
2283 giraffe_table_delete(wcalcoeff);
2284 wcalcoeff = NULL;
2285
2286 giraffe_table_delete(slitgeometry);
2287 slitgeometry = NULL;
2288
2289 giraffe_table_delete(grating);
2290 grating = NULL;
2291
2292 giraffe_table_delete(fibers);
2293 fibers = NULL;
2294
2295 return 1;
2296 }
2297
2298 cpl_frameset_insert(set, rimg_frame);
2299 }
2300
2301 }
2302 else {
2303
2304 /* Data Cube (ESO 3D format) */
2305
2306 GiFrameCreator creator = (GiFrameCreator) giraffe_fov_save_cubes_eso3d;
2307
2308 properties = giraffe_image_get_properties(rebinning->spectra);
2309 properties = cpl_propertylist_duplicate(properties);
2310
2311 giraffe_add_frameset_info(properties, set, info.sequence);
2312
2313 rimg_frame = giraffe_frame_create(GIFRAME_SCIENCE_CUBE,
2314 CPL_FRAME_LEVEL_FINAL,
2315 properties,
2316 fov,
2317 NULL,
2318 creator);
2319
2320 cpl_propertylist_delete(properties);
2321 properties = NULL;
2322
2323 if (rimg_frame == NULL) {
2324 cpl_msg_error(_id, "Cannot create local file! Aborting ...");
2325
2326 giraffe_fov_delete(fov);
2327 fov = NULL;
2328
2329 giraffe_rebinning_destroy(rebinning);
2330 rebinning = NULL;
2331
2332 giraffe_table_delete(wcalcoeff);
2333 wcalcoeff = NULL;
2334
2335 giraffe_table_delete(slitgeometry);
2336 slitgeometry = NULL;
2337
2338 giraffe_table_delete(grating);
2339 grating = NULL;
2340
2341 giraffe_table_delete(fibers);
2342 fibers = NULL;
2343
2344 return 1;
2345 }
2346
2347 status = giraffe_fiberlist_attach(rimg_frame, fibers);
2348
2349 if (status != 0) {
2350 cpl_msg_error(_id, "Cannot attach fiber setup to local "
2351 "file '%s'! Aborting ...",
2352 cpl_frame_get_filename(rimg_frame));
2353
2354 cpl_frame_delete(rimg_frame);
2355
2356 giraffe_fov_delete(fov);
2357 fov = NULL;
2358
2359 giraffe_rebinning_destroy(rebinning);
2360 rebinning = NULL;
2361
2362 giraffe_table_delete(wcalcoeff);
2363 wcalcoeff = NULL;
2364
2365 giraffe_table_delete(slitgeometry);
2366 slitgeometry = NULL;
2367
2368 giraffe_table_delete(grating);
2369 grating = NULL;
2370
2371 giraffe_table_delete(fibers);
2372 fibers = NULL;
2373
2374 return 1;
2375 }
2376
2377 cpl_frameset_insert(set, rimg_frame);
2378
2379 }
2380
2381 giraffe_fov_delete(fov);
2382 fov = NULL;
2383
2384 }
2385
2386
2387 /*
2388 * Cleanup
2389 */
2390
2391 giraffe_table_delete(wcalcoeff);
2392
2393 giraffe_table_delete(slitgeometry);
2394 giraffe_table_delete(grating);
2395 giraffe_table_delete(fibers);
2396
2397 giraffe_rebinning_destroy(rebinning);
2398
2399 return 0;
2400
2401}
2402
2403
2404/*
2405 * Build table of contents, i.e. the list of available plugins, for
2406 * this module. This function is exported.
2407 */
2408
2409int
2410cpl_plugin_get_info(cpl_pluginlist* list)
2411{
2412
2413 cpl_recipe* recipe = cx_calloc(1, sizeof *recipe);
2414 cpl_plugin* plugin = &recipe->interface;
2415
2416
2417 cpl_plugin_init(plugin,
2418 CPL_PLUGIN_API,
2419 GIRAFFE_BINARY_VERSION,
2420 CPL_PLUGIN_TYPE_RECIPE,
2421 "giscience",
2422 "Process a science observation.",
2423 "For detailed information please refer to the "
2424 "GIRAFFE pipeline user manual.\nIt is available at "
2425 "http://www.eso.org/pipelines.",
2426 "Giraffe Pipeline",
2427 PACKAGE_BUGREPORT,
2428 giraffe_get_license(),
2429 giscience_create,
2430 giscience_exec,
2431 giscience_destroy);
2432
2433 cpl_pluginlist_append(list, plugin);
2434
2435 return 0;
2436
2437}
2438
2439
2440/* Structure of a lookup table entry for the LUT used to get values for
2441 * SPEC_RES. */
2442typedef struct _giraffe_lut_entry {
2443 const char* expmode;
2444 double specres;
2445 double lamrms; /* given in (px), to convert to (nm) one needs:
2446 lamrms * spec_length / 4100 */
2447 double spec_length;
2448 int lamnlin;
2449} giraffe_lut_entry;
2450
2451
2452#ifndef NDEBUG
2453
2458static cpl_boolean _giraffe_lut_is_sorted(const giraffe_lut_entry* lut,
2459 size_t nentries)
2460{
2461 size_t i;
2462 if (nentries < 2) return CPL_TRUE;
2463 for (i = 0; i < nentries - 1; ++i) {
2464 if (strcmp(lut[i].expmode, lut[i+1].expmode) >= 0) {
2465 return CPL_FALSE;
2466 }
2467 }
2468 return CPL_TRUE;
2469}
2470
2471#endif /* NDEBUG */
2472
2473
2481static const giraffe_lut_entry* _giraffe_find_lut_entry(const char* expmode)
2482{
2483 static giraffe_lut_entry lut[] = {
2484 /* INS_EXP_MODE lambda LAMRMS spec_length LAMNLIN
2485 * R = ------------- (px) (nm) line count
2486 * delta(lambda)
2487 */
2488 {"H379.", 25000, 0.2250, 20.0, 60},
2489 {"H379.0", 25000, 0.2250, 20.0, 60},
2490 {"H395.8", 22000, 0.1913, 21.1, 53},
2491 {"H412.4", 30000, 0.1326, 22.3, 48},
2492 {"H429.7", 23000, 0.1471, 23.5, 67},
2493 {"H447.1", 20000, 0.0906, 24.9, 63},
2494 {"H447.1A", 20000, 0.0906, 24.9, 63},
2495 {"H447.1B", 31000, 0.1297, 23.5, 58},
2496 {"H465.6", 23000, 0.1573, 22.4, 57},
2497 {"H484.5", 19000, 0.1175, 23.5, 57},
2498 {"H484.5A", 19000, 0.1175, 23.5, 57},
2499 {"H484.5B", 33000, 0.0907, 24.6, 61},
2500 {"H504.8", 22000, 0.1726, 25.2, 56},
2501 {"H525.8", 17000, 0.1397, 25.8, 49},
2502 {"H525.8A", 17000, 0.1397, 25.8, 49},
2503 {"H525.8B", 29000, 0.2014, 26.4, 44},
2504 {"H548.8", 20000, 0.2389, 26.9, 51},
2505 {"H572.8", 27000, 0.1844, 27.5, 49},
2506 {"H599.3", 18000, 0.1683, 28.3, 57},
2507 {"H627.3", 24000, 0.1268, 29.0, 50},
2508 {"H651.5", 17000, 0.1185, 30.0, 42},
2509 {"H651.5A", 17000, 0.1185, 30.0, 42},
2510 {"H651.5B", 35000, 0.1253, 31.5, 37},
2511 {"H665.", 17000, 0.2076, 33.0, 45},
2512 {"H665.0", 17000, 0.2076, 33.0, 45},
2513 {"H679.7", 19000, 0.1621, 34.5, 51},
2514 {"H710.5", 25000, 0.1495, 36.0, 48},
2515 {"H737.", 16000, 0.1456, 37.5, 47},
2516 {"H737.0", 16000, 0.1456, 37.5, 47},
2517 {"H737.0A", 16000, 0.1456, 37.5, 47},
2518 {"H737.0B", 35000, 0.1147, 39.5, 40},
2519 {"H769.1", 19000, 0.2811, 41.8, 35},
2520 {"H805.3", 14000, 0.2662, 42.4, 39},
2521 {"H805.3A", 14000, 0.2662, 42.4, 39},
2522 {"H805.3B", 25000, 0.2342, 42.9, 28},
2523 {"H836.6", 16000, 0.2032, 43.3, 21},
2524 {"H836.6A", 16000, 0.2032, 43.3, 21},
2525 {"H836.6B", 34000, 0.1073, 43.7, 14},
2526 {"H875.7", 18000, 0.2026, 44.3, 29},
2527 {"H920.5", 12000, 0.1568, 44.9, 32},
2528 {"H920.5A", 12000, 0.1568, 44.9, 32},
2529 {"H920.5B", 24000, 0.2531, 45.9, 28},
2530 {"L385.7", 7500, 0.3358, 58.0, 43},
2531 {"L427.2", 6300, 0.2152, 61.1, 62},
2532 {"L479.7", 7500, 0.1554, 71.0, 61},
2533 {"L543.1", 5800, 0.2065, 82.1, 58},
2534 {"L614.2", 6600, 0.1803, 79.0, 51},
2535 {"L682.2", 8100, 0.1843, 74.0, 50},
2536 {"L773.4", 5400, 0.1617, 94.0, 44},
2537 {"L881.7", 6600, 0.1614, 119.0, 29}
2538 };
2539 static const size_t nentries = sizeof(lut) / sizeof(giraffe_lut_entry);
2540 int low = 0; /* Bottom of search region. */
2541 int high = (int)nentries - 1; /* Top of search region. */
2542
2543 assert(_giraffe_lut_is_sorted(lut, nentries));
2544 assert(expmode != NULL);
2545
2546 /* Perform a binary search for the entry. */
2547 do {
2548 int mid = (low + high) >> 1; /* Find mid point of search range. */
2549 int result = strcmp(expmode, lut[mid].expmode);
2550 if (result == 0) {
2551 return &lut[mid];
2552 } else if (result < 0) {
2553 high = mid - 1;
2554 } else {
2555 low = mid + 1;
2556 }
2557 } while (high >= low);
2558 return NULL;
2559}
2560
2567static double _giraffe_lookup_specres(const char* expmode)
2568{
2569 const giraffe_lut_entry* entry = _giraffe_find_lut_entry(expmode);
2570 if (entry == NULL) return NAN;
2571 return entry->specres;
2572}
2573
2580static double _giraffe_lookup_lamrms(const char* expmode)
2581{
2582 const giraffe_lut_entry* entry = _giraffe_find_lut_entry(expmode);
2583 if (entry == NULL) return NAN;
2584 if (isnan(entry->lamrms) || isnan(entry->spec_length)) return NAN;
2585 return entry->lamrms * entry->spec_length / 4100.;
2586}
2587
2594static double _giraffe_lookup_lamnlin(const char* expmode)
2595{
2596 const giraffe_lut_entry* entry = _giraffe_find_lut_entry(expmode);
2597 if (entry == NULL) return -1;
2598 return entry->lamnlin;
2599}
2600
2601
2612static cpl_type _giraffe_calc_wave_type(double crval2, double crpix2,
2613 double cdelt2, cpl_size naxis2)
2614{
2615 static const double errfrac = 0.02;
2616 static const double single_precision_digits = 7;
2617 double lo = (1.0 - crpix2) * cdelt2;
2618 double hi = ((double)naxis2 - crpix2) * cdelt2;
2619 double maxwave = crval2 + (hi > lo ? hi : lo);
2620 double binfrac = (maxwave != 0.0) ? fabs(cdelt2 / maxwave) : 0.0;
2621 if (binfrac * errfrac < pow(10, -single_precision_digits)) {
2622 return CPL_TYPE_DOUBLE;
2623 } else {
2624 return CPL_TYPE_FLOAT;
2625 }
2626}
2627
2628
2638static cpl_boolean _giraffe_ancillary_data_available(const char* filename,
2639 const GiTable* fibertable)
2640{
2641 cpl_table* tbl = giraffe_table_get(fibertable);
2642 cpl_size i;
2643 const char** spectypes = NULL;
2644
2645 assert(filename != NULL);
2646
2647 cpl_error_ensure(tbl != NULL, cpl_error_get_code(), return CPL_FALSE,
2648 "The fiber table is not available for '%s'.", filename);
2649
2650 spectypes = cpl_table_get_data_string_const(tbl, GIALIAS_COLUMN_TYPE);
2651 cpl_error_ensure(spectypes != NULL, cpl_error_get_code(), return CPL_FALSE,
2652 "Could not fetch the '%s' column from the fiber setup"
2653 " table in '%s'.", GIALIAS_COLUMN_TYPE, filename);
2654 /*
2655 * Step through the table and check for fiber types other than Medusa
2656 * fibers. Also the simultaneous calibration fibers are not considered
2657 * as ancillary spectra and are ignored.
2658 */
2659 for (i = 0; i < cpl_table_get_nrow(tbl); ++i) {
2660 if ((spectypes[i][0] != '\0') && (strcmp(spectypes[i], "M") != 0)) {
2661 return CPL_TRUE;
2662 }
2663 }
2664 return CPL_FALSE;
2665}
2666
2667
2687static cpl_error_code _giraffe_make_ancillary_file(cpl_frameset* allframes,
2688 const char* outputfile,
2689 const char* infilename,
2690 const GiImage* fluximage,
2691 const GiTable* fibertable)
2692{
2693 cpl_error_code error = CPL_ERROR_NONE;
2694 cxint retcode;
2695 cpl_frame* frame = NULL;
2696 cpl_image* srcimg = giraffe_image_get(fluximage);
2697 GiImage* image = NULL;
2698 GiTable* table = giraffe_table_duplicate(fibertable);
2699 cpl_image* img = NULL;
2700 cpl_table* tbl = giraffe_table_get(table);
2701 cpl_image* subimg = NULL;
2702 cpl_propertylist* comments = NULL;
2703 cpl_size ny, i;
2704 int* indices = NULL;
2705
2706 assert(allframes != NULL);
2707 assert(outputfile != NULL);
2708 assert(infilename != NULL);
2709
2710 cpl_error_ensure(srcimg != NULL && table != NULL && tbl != NULL,
2711 cpl_error_get_code(), goto cleanup,
2712 "The image or table are not available for '%s'.",
2713 infilename);
2714
2715 /* Setup a new frame for the output file. */
2716 frame = cpl_frame_new();
2717 error |= cpl_frame_set_filename(frame, outputfile);
2718 error |= cpl_frame_set_tag(frame, GIALIAS_ASSO_PROCATG_VALUE);
2719 error |= cpl_frame_set_type(frame, CPL_FRAME_TYPE_IMAGE);
2720 error |= cpl_frame_set_group(frame, CPL_FRAME_GROUP_PRODUCT);
2721 error |= cpl_frame_set_level(frame, CPL_FRAME_LEVEL_FINAL);
2722 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2723 "Failed to setup a new output frame for '%s'.",
2724 outputfile);
2725
2726 /* First step to filter out the science entries in the image and fiber setup
2727 * table is to go through the table entries, mark the entries to remove, and
2728 * delete the ones that are marked. i.e. unselect everything, select all
2729 * non-science spectra, ignoring the simultaneous calibration spectra,
2730 * invert the selection and remove everything that remains marked as
2731 * selected. */
2732 error |= cpl_table_unselect_all(tbl);
2733 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2734 "Failed to unselect all entries in fiber setup table from"
2735 " '%s'.", infilename);
2736 cpl_table_or_selected_string(tbl, GIALIAS_COLUMN_TYPE, CPL_NOT_EQUAL_TO,
2737 "M");
2738 cpl_table_and_selected_int(tbl, GIALIAS_COLUMN_RP, CPL_NOT_EQUAL_TO, -1);
2739 cpl_table_not_selected(tbl);
2740 error |= cpl_table_erase_selected(tbl);
2741 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2742 "Failed to erase selected entries in fiber setup table"
2743 " from '%s'.", infilename);
2744
2745 /* Create a new output image which is wide enough to store the data with the
2746 * stripped out spectra removed. */
2747 ny = cpl_image_get_size_y(srcimg);
2748 image = giraffe_image_create(cpl_image_get_type(srcimg),
2749 cpl_table_get_nrow(tbl), ny);
2750 img = giraffe_image_get(image);
2751 retcode = giraffe_image_set_properties(
2752 image, giraffe_image_get_properties(fluximage));
2753 cpl_error_ensure(image != NULL && img != NULL && retcode == 0,
2754 cpl_error_get_code(), goto cleanup,
2755 "Failed to create image for output file '%s'.",
2756 outputfile);
2757
2758 error |= cpl_propertylist_update_string(giraffe_image_get_properties(image),
2759 GIALIAS_PROCATG,
2760 GIALIAS_ASSO_PROCATG_VALUE);
2761 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2762 "Could not update keyword '%s' for output file '%s'.",
2763 GIALIAS_PROCATG, outputfile);
2764
2765
2766 error |= cpl_propertylist_update_string(giraffe_image_get_properties(image),
2767 GIALIAS_PRODCATG,
2768 GIALIAS_PRODCATG_MOSSKY);
2769 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2770 "Could not update keyword '%s' for output file '%s'.",
2771 GIALIAS_PRODCATG, outputfile);
2772
2773 /* Strip out extra comments which will be added by CPL anyway. */
2774 if (cpl_propertylist_has(giraffe_image_get_properties(image), "COMMENT")) {
2775 cpl_size ic;
2776 const char* comments_to_remove[] = {
2777 " FITS (Flexible Image Transport System) format is defined in"
2778 " 'Astronomy",
2779 " and Astrophysics', volume 376, page 359; bibcode: 2001A&A..."
2780 "376..359H",
2781 NULL
2782 };
2783 cpl_propertylist* props = giraffe_image_get_properties(image);
2784 comments = cpl_propertylist_new();
2785 error |= cpl_propertylist_copy_property_regexp(comments, props,
2786 "^COMMENT$", 0);
2787 cpl_propertylist_erase_regexp(props, "^COMMENT$", 0);
2788 for (ic = 0; ic < cpl_propertylist_get_size(comments); ++ic) {
2789 const char** cmnt_str;
2790 cpl_property* p = cpl_propertylist_get(comments, ic);
2791 for (cmnt_str = comments_to_remove; *cmnt_str != NULL; ++cmnt_str) {
2792 if (strcmp(cpl_property_get_string(p), *cmnt_str) == 0) {
2793 goto dont_add_comment;
2794 }
2795 }
2796 /* Add back comments that should not be removed. */
2797 error |= cpl_propertylist_append_property(props, p);
2798 dont_add_comment:
2799 /* Land up here if the comment was found in the comments_to_remove
2800 * list of strings. */
2801 ;
2802 }
2803 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2804 "Failed to cleanup comments in primary HDU for '%s'.",
2805 outputfile);
2806 cpl_propertylist_delete(comments);
2807 comments = NULL;
2808 }
2809
2810 /* We now have to relabel the index numbers in the fiber setup table and
2811 * copy corresponding image columns to the new image. */
2812 indices = cpl_table_get_data_int(tbl, GIALIAS_COLUMN_INDEX);
2813 cpl_error_ensure(indices != NULL, cpl_error_get_code(), goto cleanup,
2814 "Could not fetch the '%s' column from the fiber setup"
2815 " table in '%s'.", GIALIAS_COLUMN_INDEX, infilename);
2816 for (i = 0; i < cpl_table_get_nrow(tbl); ++i) {
2817 cpl_size oldindex = indices[i];
2818 cpl_size newindex = i+1;
2819 indices[i] = newindex;
2820 subimg = cpl_image_extract(srcimg, oldindex, 1, oldindex, ny);
2821 cpl_error_ensure(subimg != NULL, cpl_error_get_code(), goto cleanup,
2822 "Could not extract sub image from '%s' at column %"
2823 CPL_SIZE_FORMAT".", infilename, oldindex);
2824 error |= cpl_image_copy(img, subimg, newindex, 1);
2825 cpl_image_delete(subimg);
2826 subimg = NULL;
2827 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2828 "Could write sub image from '%s' at column %"
2829 CPL_SIZE_FORMAT" to new image in '%s' at column %"
2830 CPL_SIZE_FORMAT".", infilename, oldindex, outputfile,
2831 newindex);
2832 }
2833
2834 /* Now write the actual FITS file. */
2835 retcode = giraffe_image_save(image, outputfile);
2836 cpl_error_ensure(retcode == 0, cpl_error_get_code(), goto cleanup,
2837 "Failed to write image to file '%s'.", outputfile);
2838 retcode = giraffe_fiberlist_attach(frame, table);
2839 cpl_error_ensure(retcode == 0, cpl_error_get_code(), goto cleanup,
2840 "Failed to attach the fiber setup table to file '%s'.",
2841 outputfile);
2842
2843 /* Add the new frame to the output frame set for the ancillary file.
2844 * Note: this should be the last step so that we do not delete this frame
2845 * anymore once in the frame set. */
2846 error |= cpl_frameset_insert(allframes, frame);
2847 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2848 "Could not add a new frame to the frame list for '%s'.",
2849 outputfile);
2850
2851 giraffe_image_delete(image);
2852 giraffe_table_delete(table);
2853 return CPL_ERROR_NONE;
2854
2855cleanup:
2856 /* Error handling. Note: NULL pointer checks done by delete functions. */
2857 cpl_image_delete(subimg);
2858 giraffe_image_delete(image);
2859 giraffe_table_delete(table);
2860 cpl_frame_delete(frame);
2861 cpl_propertylist_delete(comments);
2862 return cpl_error_get_code();
2863}
2864
2865
2872static char* _giraffe_calc_format_string(cpl_size maxfiles)
2873{
2874 double ndigits = 1.0;
2875 if (maxfiles > 1) {
2876 /* Figure out how many digits must be printed to support a index number
2877 * as large as 'maxfiles'. */
2878 ndigits = ceil(log10((double)maxfiles + 1.0));
2879 }
2880 return cpl_sprintf("science_spectrum_%%0%.0f"CPL_SIZE_FORMAT".fits",
2881 ndigits);
2882}
2883
2884
2901static cxint _giraffe_make_sdp_spectra(const cxchar* flux_filename,
2902 const cxchar* err_filename,
2903 cxint nassoc_keys,
2904 cpl_frameset* allframes,
2905 const cpl_parameterlist* parlist,
2906 const cxchar* recipe_id)
2907{
2908 cxint result_code = 1;
2909 cxint errorcode;
2910 cpl_error_code error = CPL_ERROR_NONE;
2911 cpl_errorstate prestate;
2912 const char* ancillary_filename = "science_ancillary.fits";
2913 GiImage* fluximage = giraffe_image_new(CPL_TYPE_DOUBLE);
2914 GiImage* errimage = giraffe_image_new(CPL_TYPE_DOUBLE);
2915 GiTable* fibertable = NULL;
2916 const cxchar* fibertable_name = NULL;
2917 irplib_sdp_spectrum* spectrum = NULL;
2918 cpl_propertylist* extrakeys = cpl_propertylist_new();
2919 cpl_propertylist* tablekeys = cpl_propertylist_new();
2920 cpl_propertylist* props;
2921 char* pipe_id = cpl_sprintf("%s/%s", PACKAGE, VERSION);
2922 const char* dict_id = PRODUCT_DID;
2923 const cpl_frame* inherit = NULL;
2924 cpl_frameset* usedframes = NULL;
2925 cpl_frameset* rawframes = NULL;
2926 cpl_frameset_iterator* iterator = NULL;
2927 cpl_size nx, ny, i, filecount;
2928 double exptime = NAN;
2929 double mjdobs = NAN;
2930 double mjdend = NAN;
2931 double wavelmin = NAN;
2932 double wavelmax = NAN;
2933 double specbin = NAN;
2934 double crpix2 = NAN;
2935 double crval2 = NAN;
2936 double cdelt2 = NAN;
2937 const char* cunit2 = NULL;
2938 double specres;
2939 const char* expmode = NULL;
2940 char strbuf[64];
2941 const int* indices = NULL;
2942 const int* fps = NULL;
2943 const char** objects = NULL;
2944 const char** spectypes = NULL;
2945 const double* ras = NULL;
2946 const double* decs = NULL;
2947 const double* gcorr = NULL;
2948 const double* hcorr = NULL;
2949 const double* bcorr = NULL;
2950 char* formatstr = NULL;
2951 char* filename = NULL;
2952 cpl_type wavecoltype;
2953 cpl_array* refwavearray = NULL;
2954 cpl_array* array = NULL;
2955 float* data_float = NULL;
2956 double* data_double = NULL;
2957 cpl_vector* fluximgcol = NULL;
2958 cpl_vector* errimgcol = NULL;
2959 cpl_boolean got_ancillary_data = CPL_FALSE;
2960 cpl_size assoc_key_offset = 1;
2961 int lamnlin = -1;
2962 double lamrms = NAN;
2963 double specerr = NAN;
2964 double specsye = NAN;
2965 int obsid = -1;
2966 const char *dateobs = NULL;
2967
2968 cpl_error_ensure(flux_filename != NULL && err_filename != NULL
2969 && allframes != NULL && parlist != NULL
2970 && recipe_id != NULL, CPL_ERROR_NULL_INPUT, goto cleanup,
2971 "NULL input parameters.");
2972
2973 error |= cpl_propertylist_append_string(extrakeys, GIALIAS_PROCATG,
2974 GIALIAS_PROCATG_RBNSPEC_IDP);
2975 error |= cpl_propertylist_set_comment(extrakeys, GIALIAS_PROCATG,
2976 GIALIAS_PROCATG_COMMENT);
2977 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
2978 "Could not set keyword '%s'.", GIALIAS_PROCATG);
2979
2980 /* Load the input flux and error data, including FITS header keywords. */
2981 errorcode = giraffe_image_load(fluximage, flux_filename, 0);
2982 cpl_error_ensure(errorcode == 0, cpl_error_get_code(), goto cleanup,
2983 "Could not load image data in primary HDU from '%s'.",
2984 flux_filename);
2985 errorcode = giraffe_image_load(errimage, err_filename, 0);
2986 cpl_error_ensure(errorcode == 0, cpl_error_get_code(), goto cleanup,
2987 "Could not load image data in primary HDU from '%s'.",
2988 err_filename);
2989
2990 giraffe_error_push();
2991 fibertable = giraffe_fiberlist_load(flux_filename, 1, GIALIAS_FIBER_SETUP);
2992 if (fibertable == NULL) {
2993 fibertable = giraffe_fiberlist_load(err_filename, 1,
2994 GIALIAS_FIBER_SETUP);
2995 fibertable_name = err_filename;
2996 } else {
2997 fibertable_name = flux_filename;
2998 }
2999 cpl_error_ensure(fibertable != NULL, CPL_ERROR_DATA_NOT_FOUND, goto cleanup,
3000 "Could not load the %s table from either '%s' or '%s'.",
3001 GIALIAS_FIBER_SETUP, flux_filename, err_filename);
3002 giraffe_error_pop();
3003
3004 /* Check that the image sizes are the same. */
3005 nx = cpl_image_get_size_x(giraffe_image_get(fluximage));
3006 ny = cpl_image_get_size_y(giraffe_image_get(fluximage));
3007 cpl_error_ensure(cpl_image_get_size_x(giraffe_image_get(errimage)) == nx
3008 && cpl_image_get_size_y(giraffe_image_get(errimage)) == ny,
3009 CPL_ERROR_INCOMPATIBLE_INPUT, goto cleanup,
3010 "The images in files '%s' and '%s' are not the same size.",
3011 flux_filename, err_filename);
3012
3013 /* Construct the used frame list from the existing list of all frames.
3014 * The frames to be included as the used frames are RAW and CALIB. */
3015 usedframes = cpl_frameset_new();
3016 rawframes = cpl_frameset_new();
3017 iterator = cpl_frameset_iterator_new(allframes);
3018 do {
3019 const cpl_frame* frame = cpl_frameset_iterator_get_const(iterator);
3020 if (frame != NULL) {
3021 switch (cpl_frame_get_group(frame)) {
3022 case CPL_FRAME_GROUP_RAW:
3023 /* Mark the first RAW frame from which to inherit keywords. */
3024 if (inherit == NULL) inherit = frame;
3025 error |= cpl_frameset_insert(rawframes,
3026 cpl_frame_duplicate(frame));
3027 error |= cpl_frameset_insert(usedframes,
3028 cpl_frame_duplicate(frame));
3029 break;
3030 case CPL_FRAME_GROUP_CALIB:
3031 error |= cpl_frameset_insert(usedframes,
3032 cpl_frame_duplicate(frame));
3033 break;
3034 default: /* Ignore all other groups */
3035 break;
3036 }
3037 }
3038 prestate = cpl_errorstate_get();
3039 error |= cpl_frameset_iterator_advance(iterator, 1);
3040 if (error == CPL_ERROR_ACCESS_OUT_OF_RANGE) {
3041 cpl_errorstate_set(prestate);
3042 break;
3043 } else if (error != CPL_ERROR_NONE) {
3044 goto cleanup;
3045 }
3046 } while (1);
3047
3048 cpl_error_ensure(inherit != NULL, CPL_ERROR_DATA_NOT_FOUND, goto cleanup,
3049 "No raw input frames found.");
3050
3051 /* Fetch the EXPTIME and MJD-OBS keywords from the product file. */
3052 props = giraffe_image_get_properties(fluximage);
3053 prestate = cpl_errorstate_get();
3054 exptime = cpl_propertylist_get_double(props, GIALIAS_EXPTIME);
3055 cpl_error_ensure(cpl_errorstate_is_equal(prestate), cpl_error_get_code(),
3056 goto cleanup, "Could not find keyword '%s' in '%s'.",
3057 GIALIAS_EXPTIME, flux_filename);
3058 mjdobs = cpl_propertylist_get_double(props, GIALIAS_MJDOBS);
3059 cpl_error_ensure(cpl_errorstate_is_equal(prestate), cpl_error_get_code(),
3060 goto cleanup, "Could not find keyword '%s' in '%s'.",
3061 GIALIAS_MJDOBS, flux_filename);
3062
3063 mjdend = mjdobs + exptime / 86400.;
3064
3065 /* Fetch the DATE-OBS and ESO OBS ID keywords from the product file. */
3066
3067 dateobs = cpl_propertylist_get_string(props, GIALIAS_DATEOBS);
3068 cpl_error_ensure(cpl_errorstate_is_equal(prestate), cpl_error_get_code(),
3069 goto cleanup, "Could not find keyword '%s' in '%s'.",
3070 GIALIAS_DATEOBS, flux_filename);
3071
3072 obsid = cpl_propertylist_get_int(props, GIALIAS_OBSID);
3073 cpl_error_ensure(cpl_errorstate_is_equal(prestate), cpl_error_get_code(),
3074 goto cleanup, "Could not find keyword '%s' in '%s'.",
3075 GIALIAS_OBSID, flux_filename);
3076
3077
3078 /* Calculate the min/max wavelength values. */
3079 crpix2 = cpl_propertylist_get_double(props, GIALIAS_CRPIX2);
3080 cpl_error_ensure(cpl_errorstate_is_equal(prestate), cpl_error_get_code(),
3081 goto cleanup, "Could not find keyword '%s' in '%s'.",
3082 GIALIAS_CRPIX2, flux_filename);
3083 crval2 = cpl_propertylist_get_double(props, GIALIAS_CRVAL2);
3084 cpl_error_ensure(cpl_errorstate_is_equal(prestate), cpl_error_get_code(),
3085 goto cleanup, "Could not find keyword '%s' in '%s'.",
3086 GIALIAS_CRVAL2, flux_filename);
3087 cdelt2 = cpl_propertylist_get_double(props, GIALIAS_CDELT2);
3088 cpl_error_ensure(cpl_errorstate_is_equal(prestate), cpl_error_get_code(),
3089 goto cleanup, "Could not find keyword '%s' in '%s'.",
3090 GIALIAS_CDELT2, flux_filename);
3091 cunit2 = cpl_propertylist_get_string(props, GIALIAS_CUNIT2);
3092 cpl_error_ensure(cunit2 != NULL, cpl_error_get_code(),
3093 goto cleanup, "Could not find keyword '%s' in '%s'.",
3094 GIALIAS_CUNIT2, flux_filename);
3095
3096 if (strcmp(cunit2, "nm") == 0) {
3097 wavelmin = (1.0 - crpix2) * cdelt2 + crval2;
3098 wavelmax = ((double)ny - crpix2) * cdelt2 + crval2;
3099 if (wavelmax < wavelmin) {
3100 double tmp = wavelmin;
3101 wavelmin = wavelmax;
3102 wavelmax = tmp;
3103 }
3104 specbin = fabs(cdelt2);
3105 } else {
3106 cpl_msg_warning(cpl_func, "Do not know how to handle keyword %s = '%s'."
3107 " Will not set WAVELMIN, WAVELMAX or SPEC_BIN.",
3108 GIALIAS_CUNIT2, cunit2);
3109 }
3110
3111 if (cpl_propertylist_has(props, GIALIAS_SETUPNAME)) {
3112 expmode = cpl_propertylist_get_string(props, GIALIAS_SETUPNAME);
3113 cpl_error_ensure(expmode != NULL, cpl_error_get_code(), goto cleanup,
3114 "Could not fetch the keyword '%s' from '%s'.",
3115 GIALIAS_SETUPNAME, flux_filename);
3116 } else if (cpl_propertylist_has(props, GIALIAS_GRATNAME)) {
3117 const char* name = cpl_propertylist_get_string(props, GIALIAS_GRATNAME);
3118 cpl_error_ensure(name != NULL, cpl_error_get_code(), goto cleanup,
3119 "Could not fetch the keyword '%s' from '%s'.",
3120 GIALIAS_GRATNAME, flux_filename);
3121 double wlen = cpl_propertylist_get_double(props, GIALIAS_GRATWLEN);
3122 cpl_error_ensure(cpl_errorstate_is_equal(prestate),
3123 cpl_error_get_code(), goto cleanup,
3124 "Could not find keyword '%s' in '%s'.",
3125 GIALIAS_GRATWLEN, flux_filename);
3126 strbuf[0] = name[0];
3127 char* numstr = cpl_sprintf("%.1f", wlen);
3128 strncpy(strbuf+1, numstr, sizeof(strbuf)-1);
3129 cpl_free(numstr);
3130 strbuf[sizeof(strbuf)-1] = '\0'; /* Ensure we have a NULL terminator. */
3131 expmode = strbuf;
3132 } else {
3133 cpl_error_set_message(cpl_func, CPL_ERROR_DATA_NOT_FOUND,
3134 "Neither '%s' nor '%s' and '%s' keywords were found in the"
3135 " file '%s'.", GIALIAS_SETUPNAME, GIALIAS_GRATNAME,
3136 GIALIAS_GRATWLEN, flux_filename);
3137 goto cleanup;
3138 }
3139
3140 specres = _giraffe_lookup_specres(expmode);
3141 if (isnan(specres)) {
3142 cpl_error_set_message(cpl_func, CPL_ERROR_ILLEGAL_INPUT,
3143 "The exposure mode '%s' is invalid or an unknown value."
3144 " Could not lookup the spectral resolution for 'SPEC_RES'.",
3145 expmode);
3146 goto cleanup;
3147 }
3148
3149 /* Add the FILTER keyword as OFILTER to the extra keywords list if it
3150 * exists. The FILTER keyword itself will be deleted. */
3151 if (cpl_propertylist_has(props, "FILTER")) {
3152 prestate = cpl_errorstate_get();
3153 cpl_propertylist_copy_property(extrakeys, props, "FILTER");
3154 cpl_property* prop = cpl_propertylist_get_property(extrakeys, "FILTER");
3155 cpl_property_set_name(prop, "OFILTER");
3156 cpl_error_ensure(cpl_errorstate_is_equal(prestate),
3157 cpl_error_get_code(), goto cleanup,
3158 "Could not rename the 'FILTER' keyword.");
3159 }
3160
3161 /* Write the ancillary data file if any ancillary data is available. */
3162 prestate = cpl_errorstate_get();
3163 got_ancillary_data = _giraffe_ancillary_data_available(flux_filename,
3164 fibertable);
3165 if (! cpl_errorstate_is_equal(prestate)) goto cleanup;
3166 if (got_ancillary_data) {
3167 error = _giraffe_make_ancillary_file(allframes, ancillary_filename,
3168 flux_filename, fluximage,
3169 fibertable);
3170 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
3171 "Failed to write the ancillary file '%s'.",
3172 ancillary_filename);
3173 }
3174
3175 /* Create a new spectrum object and setup header keywords. */
3176 spectrum = irplib_sdp_spectrum_new();
3177 error = CPL_ERROR_NONE;
3178 error |= irplib_sdp_spectrum_set_origin(spectrum, GIALIAS_ORIGIN_VALUE);
3179 error |= irplib_sdp_spectrum_set_prodlvl(spectrum, GIALIAS_PRODLVL_VALUE);
3180 error |= irplib_sdp_spectrum_copy_dispelem(spectrum,
3181 props, GIALIAS_GRATNAME);
3182 error |= irplib_sdp_spectrum_set_specsys(spectrum, GIALIAS_SPECSYS_VALUE);
3183 error |= irplib_sdp_spectrum_set_extobj(spectrum, GIALIAS_EXT_OBJ_VALUE);
3184 /* The OBJECT, RA and DEC keywords we fill now with dummy values to maintain
3185 * the order of keywords. The actual values are set later. */
3186 error |= irplib_sdp_spectrum_set_object(spectrum, "");
3187 error |= irplib_sdp_spectrum_set_ra(spectrum, 0.0);
3188 error |= irplib_sdp_spectrum_set_dec(spectrum, 0.0);
3189 error |= irplib_sdp_spectrum_copy_exptime(spectrum, props, GIALIAS_EXPTIME);
3190 error |= irplib_sdp_spectrum_copy_texptime(spectrum,
3191 props, GIALIAS_EXPTIME);
3192 error |= irplib_sdp_spectrum_copy_mjdobs(spectrum, props, GIALIAS_MJDOBS);
3193 error |= irplib_sdp_spectrum_set_mjdend(spectrum, mjdend);
3194 if (cpl_propertylist_has(props, GIALIAS_TIMESYS)) {
3195 error |= irplib_sdp_spectrum_copy_timesys(spectrum,
3196 props, GIALIAS_TIMESYS);
3197 }
3198 error |= irplib_sdp_spectrum_copy_progid(spectrum, props, GIALIAS_PROGID);
3199 error |= irplib_sdp_spectrum_copy_obid(spectrum, 1, props, GIALIAS_OBSID);
3200 error |= irplib_sdp_spectrum_set_mepoch(spectrum, GIALIAS_M_EPOCH_VALUE);
3201 error |= irplib_sdp_spectrum_append_prov(spectrum, 1, rawframes);
3202 error |= irplib_sdp_spectrum_copy_procsoft(spectrum,
3203 props, GIALIAS_PROPIPEID);
3204 error |= irplib_sdp_spectrum_copy_obstech(spectrum, props, GIALIAS_PROTECH);
3205 error |= irplib_sdp_spectrum_set_prodcatg(spectrum, GIALIAS_PRODCATG_VALUE);
3206 error |= irplib_sdp_spectrum_set_fluxcal(spectrum, GIALIAS_FLUXCAL_VALUE);
3207 error |= irplib_sdp_spectrum_set_contnorm(spectrum, GIALIAS_CONTNORM_VALUE);
3208 /* Set dummy values for WAVELMIN, WAVELMAX and SPEC_BIN to keep the order
3209 * of the keywords. Will fill these in later with actual Heliocentric
3210 * corrected values. */
3211 error |= irplib_sdp_spectrum_set_wavelmin(spectrum, 0.0);
3212 error |= irplib_sdp_spectrum_set_wavelmax(spectrum, 0.0);
3213 error |= irplib_sdp_spectrum_set_specbin(spectrum, 0.0);
3214 error |= irplib_sdp_spectrum_set_totflux(spectrum, GIALIAS_TOTFLUX_VALUE);
3215 error |= irplib_sdp_spectrum_set_fluxerr(spectrum, GIALIAS_FLUXERR_VALUE);
3216 error |= irplib_sdp_spectrum_set_ncombine(spectrum,
3217 cpl_frameset_get_size(rawframes));
3218 error |= irplib_sdp_spectrum_set_referenc(spectrum, GIALIAS_REFERENC);
3219
3220 /* Set dummy value for SNR to maintain keyword ordering. Filled later. */
3221 error |= irplib_sdp_spectrum_set_snr(spectrum, 0.0);
3222
3223 /* Copy LAMNLIN if available from flux image else try look it up. */
3224 if (cpl_propertylist_has(props, GIALIAS_LAMNLIN)) {
3225 error |= irplib_sdp_spectrum_copy_lamnlin(spectrum, props,
3226 GIALIAS_LAMNLIN);
3227 lamnlin = irplib_sdp_spectrum_get_lamnlin(spectrum);
3228 } else {
3229 lamnlin = _giraffe_lookup_lamnlin(expmode);
3230 if (lamnlin != -1) {
3231 error |= irplib_sdp_spectrum_set_lamnlin(spectrum, lamnlin);
3232 }
3233 }
3234
3235 /* Copy LAMRMS if available from flux image else try look it up.
3236 * Note: we are going to have to correct this value later and fill it in
3237 * again for each spectrum. */
3238 if (cpl_propertylist_has(props, GIALIAS_LAMRMS)) {
3239 error |= irplib_sdp_spectrum_copy_lamrms(spectrum, props,
3240 GIALIAS_LAMRMS);
3241 lamrms = irplib_sdp_spectrum_get_lamrms(spectrum);
3242 } else {
3243 lamrms = _giraffe_lookup_lamrms(expmode);
3244 if (! isnan(lamrms)) {
3245 error |= irplib_sdp_spectrum_set_lamrms(spectrum, lamrms);
3246 }
3247 }
3248
3249 /* Copy SPEC_ERR if available from the flux image, else estimate it as:
3250 * if CRDER1 is available then
3251 * SPEC_ERR = CRDER1 / sqrt(LAMNLIN)
3252 * else
3253 * SPEC_ERR = LAMRMS / sqrt(LAMNLIN)
3254 * end
3255 *
3256 * Note: we are going to have to correct this value later and fill it in
3257 * again for each spectrum.
3258 */
3259 if (cpl_propertylist_has(props, GIALIAS_SPEC_ERR)) {
3260 error |= irplib_sdp_spectrum_copy_specerr(spectrum, props,
3261 GIALIAS_SPEC_ERR);
3262 specerr = irplib_sdp_spectrum_get_specerr(spectrum);
3263 } else if (lamnlin > 0) {
3264 if (cpl_propertylist_has(props, GIALIAS_CRDER1)) {
3265 prestate = cpl_errorstate_get();
3266 double crder1 = cpl_propertylist_get_double(props, GIALIAS_CRDER1);
3267 if (cpl_errorstate_is_equal(prestate) && crder1 > 0) {
3268 specerr = crder1 / sqrt(lamnlin);
3269 } else {
3270 error = cpl_error_get_code();
3271 }
3272 } else if (! isnan(lamrms)) {
3273 specerr = lamrms / sqrt(lamnlin);
3274 }
3275 if (! isnan(specerr)) {
3276 error |= irplib_sdp_spectrum_set_specerr(spectrum, specerr);
3277 }
3278 }
3279
3280 /* Copy SPEC_SYE if available from the flux image, else estimate it as:
3281 * 0.002 nm
3282 */
3283 if (cpl_propertylist_has(props, GIALIAS_SPEC_SYE)) {
3284 error |= irplib_sdp_spectrum_copy_specsye(spectrum, props,
3285 GIALIAS_SPEC_SYE);
3286 specsye = irplib_sdp_spectrum_get_specsye(spectrum);
3287 } else {
3288 /* Don't set the local specsye variable so that it does not get
3289 * corrected later. We want it to always be 0.002 nm. Just set the
3290 * keyword in the spectrum object directly. */
3291 error |= irplib_sdp_spectrum_set_specsye(spectrum, 0.002);
3292 }
3293
3294 error |= irplib_sdp_spectrum_set_specres(spectrum, specres);
3295 error |= irplib_sdp_spectrum_copy_gain(spectrum, props, GIALIAS_CONAD);
3296 error |= irplib_sdp_spectrum_copy_detron(spectrum, props, GIALIAS_RON);
3297 if (got_ancillary_data) {
3298 /* This assumes ancillary_filename points to a fits file
3299 * this is currently always the case, otherwise the keywords
3300 * assoc and assom must also be added. */
3301 error |= irplib_sdp_spectrum_set_asson(spectrum, 1, ancillary_filename);
3302 assoc_key_offset = 2;
3303 }
3304 for (i = assoc_key_offset; i < nassoc_keys + assoc_key_offset; ++i) {
3305 /* Add extra dummy association keywords if requested. */
3306 error |= irplib_sdp_spectrum_set_asson(spectrum, i, "");
3307 error |= irplib_sdp_spectrum_set_assoc(spectrum, i, "");
3308 error |= irplib_sdp_spectrum_set_assom(spectrum, i, "");
3309 }
3310
3311 error |= irplib_sdp_spectrum_set_voclass(spectrum, GIALIAS_VOCLASS_VALUE);
3312 error |= irplib_sdp_spectrum_set_vopub(spectrum, GIALIAS_VOPUB_VALUE);
3313 error |= irplib_sdp_spectrum_set_title(spectrum, ""); /* Set dummy value */
3314 error |= cpl_propertylist_append_double(tablekeys, GIALIAS_APERTURE,
3315 GIALIAS_APERTURE_VALUE);
3316 error |= cpl_propertylist_set_comment(tablekeys, GIALIAS_APERTURE,
3317 GIALIAS_APERTURE_COMMENT);
3318 /*
3319 * Normally: telapse = (mjdend - mjdobs) * 86400.
3320 * However, doing this calculation directly leads to rounding errors that
3321 * can cause the invalid condition TELAPSE < EXPTIME. Since we have:
3322 * mjdend = mjdobs + exptime / 86400.
3323 * we can simplify the above to be: telapse = exptime
3324 * which will always satisfy the condition TELAPSE >= EXPTIME.
3325 */
3326 error |= irplib_sdp_spectrum_set_telapse(spectrum, exptime);
3327 error |= irplib_sdp_spectrum_set_tmid(spectrum, (mjdobs + mjdend) * 0.5);
3328 /* Set dummy values for the SPEC_VAL, SPEC_BW, TDMIN and TDMAX values to
3329 * keep the order of the keywords on the header. Will fill these in later
3330 * with Heliocentric corrected values. */
3331 error |= irplib_sdp_spectrum_set_specval(spectrum, 0.0);
3332 error |= irplib_sdp_spectrum_set_specbw(spectrum, 0.0);
3333 error |= irplib_sdp_spectrum_set_nelem(spectrum, ny);
3334 error |= irplib_sdp_spectrum_set_tdmin(spectrum, 0.0);
3335 error |= irplib_sdp_spectrum_set_tdmax(spectrum, 0.0);
3336
3337 /* Add the keyword FPS */
3338 error |= cpl_propertylist_append_int(extrakeys, GIALIAS_FPS, -1);
3339 error |= cpl_propertylist_set_comment(extrakeys, GIALIAS_FPS,
3340 GIALIAS_FPS_COMMENT);
3341
3342 /* Add dummy [G,H,B]CORR keywords to be updated from the fiber table. */
3343 error |= cpl_propertylist_append_double(extrakeys, GIALIAS_GEOCORR, NAN);
3344 error |= cpl_propertylist_set_comment(extrakeys, GIALIAS_GEOCORR,
3345 GIALIAS_GEOCORR_COMMENT);
3346 error |= cpl_propertylist_append_double(extrakeys, GIALIAS_HELICORR, NAN);
3347 error |= cpl_propertylist_set_comment(extrakeys, GIALIAS_HELICORR,
3348 GIALIAS_HELICORR_COMMENT);
3349 error |= cpl_propertylist_append_double(extrakeys, GIALIAS_BARYCORR, NAN);
3350 error |= cpl_propertylist_set_comment(extrakeys, GIALIAS_BARYCORR,
3351 GIALIAS_BARYCORR_COMMENT);
3352
3353 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
3354 "Could not setup the common SDP spectrum keywords.");
3355
3356 for (int qci = 0; qci < (int)(sizeof(sciqcpar) / sizeof(sciqcpar[0])); qci++) {
3357 if (cpl_propertylist_has(props, sciqcpar[qci])) {
3358 cpl_propertylist_copy_property(extrakeys, props, sciqcpar[qci]);
3359 }
3360 }
3361
3362 /* Figure out the data type required for the WAVE column to preserve the
3363 * precision of the data values. */
3364 wavecoltype = _giraffe_calc_wave_type(crval2, crpix2, cdelt2, ny);
3365
3366 /* Calculate the reference wavelength values (before corrections). */
3367 refwavearray = cpl_array_new(ny, CPL_TYPE_DOUBLE);
3368 data_double = cpl_array_get_data_double(refwavearray);
3369 assert(data_double != NULL);
3370 for (i = 1; i <= ny; ++i) {
3371 data_double[i-1] = (i-crpix2) * cdelt2 + crval2;
3372 }
3373 data_double = NULL;
3374
3375 /* Try setup the SDP table columns. */
3376 error |= irplib_sdp_spectrum_add_column(
3377 spectrum, GIALIAS_COLUMN_WAVE, wavecoltype,
3378 GIALIAS_COLUMN_WAVE_UNIT, NULL, GIALIAS_COLUMN_WAVE_TUTYP,
3379 GIALIAS_COLUMN_WAVE_TUCD, NULL);
3380
3381 /* Replace default keyword comment of the just created column */
3382 error |= irplib_sdp_spectrum_replace_column_comment(
3383 spectrum, GIALIAS_COLUMN_WAVE, "TUCD", "Air wavelength");
3384
3385 error |= irplib_sdp_spectrum_set_column_tcomm(
3386 spectrum, GIALIAS_COLUMN_WAVE, GIALIAS_COLUMN_WAVE_TCOMM);
3387 error |= irplib_sdp_spectrum_add_column(
3388 spectrum, GIALIAS_COLUMN_FLUX_REDUCED, CPL_TYPE_DOUBLE,
3389 GIALIAS_COLUMN_FLUX_REDUCED_UNIT, NULL,
3390 GIALIAS_COLUMN_FLUX_REDUCED_TUTYP,
3391 GIALIAS_COLUMN_FLUX_REDUCED_TUCD, NULL);
3392 error |= irplib_sdp_spectrum_set_column_tcomm(
3393 spectrum, GIALIAS_COLUMN_FLUX_REDUCED, "");
3394 error |= irplib_sdp_spectrum_add_column(
3395 spectrum, GIALIAS_COLUMN_ERR_REDUCED, CPL_TYPE_DOUBLE,
3396 GIALIAS_COLUMN_ERR_REDUCED_UNIT, NULL,
3397 GIALIAS_COLUMN_ERR_REDUCED_TUTYP,
3398 GIALIAS_COLUMN_ERR_REDUCED_TUCD, NULL);
3399 error |= irplib_sdp_spectrum_set_column_tcomm(
3400 spectrum, GIALIAS_COLUMN_ERR_REDUCED, "");
3401 error |= irplib_sdp_spectrum_add_column(
3402 spectrum, GIALIAS_COLUMN_SNR, CPL_TYPE_DOUBLE,
3403 GIALIAS_COLUMN_SNR_UNIT, NULL, GIALIAS_COLUMN_SNR_TUTYP,
3404 GIALIAS_COLUMN_SNR_TUCD, NULL);
3405 error |= irplib_sdp_spectrum_set_column_tcomm(
3406 spectrum, GIALIAS_COLUMN_SNR, GIALIAS_COLUMN_SNR_TCOMM);
3407 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
3408 "Could not setup the SDP spectrum columns.");
3409
3410 indices = cpl_table_get_data_int_const(giraffe_table_get(fibertable),
3411 GIALIAS_COLUMN_INDEX);
3412 fps = cpl_table_get_data_int_const(giraffe_table_get(fibertable),
3413 GIALIAS_FPS);
3414 objects = cpl_table_get_data_string_const(giraffe_table_get(fibertable),
3415 GIALIAS_COLUMN_OBJECT);
3416 spectypes = cpl_table_get_data_string_const(giraffe_table_get(fibertable),
3417 GIALIAS_COLUMN_TYPE);
3418 ras = cpl_table_get_data_double_const(giraffe_table_get(fibertable),
3419 GIALIAS_COLUMN_RA);
3420 decs = cpl_table_get_data_double_const(giraffe_table_get(fibertable),
3421 GIALIAS_COLUMN_DEC);
3422 gcorr = cpl_table_get_data_double_const(giraffe_table_get(fibertable),
3423 GIALIAS_COLUMN_GCORR);
3424 hcorr = cpl_table_get_data_double_const(giraffe_table_get(fibertable),
3425 GIALIAS_COLUMN_HCORR);
3426 bcorr = cpl_table_get_data_double_const(giraffe_table_get(fibertable),
3427 GIALIAS_COLUMN_BCORR);
3428 cpl_error_ensure(indices != NULL && fps != NULL && objects != NULL
3429 && spectypes != NULL && ras != NULL && decs != NULL
3430 && gcorr != NULL && hcorr != NULL && bcorr != NULL,
3431 cpl_error_get_code(), goto cleanup,
3432 "Could not fetch data from the fiber setup table in '%s'.",
3433 fibertable_name);
3434
3435 formatstr = _giraffe_calc_format_string(
3436 cpl_table_get_nrow(giraffe_table_get(fibertable)));
3437
3438 cpl_errorstate tempes = cpl_errorstate_get();
3439
3440 cpl_size maxx = -1;
3441 cpl_size maxy = -1;
3442 cpl_image * maxim = cpl_image_collapse_create(giraffe_image_get(fluximage),0);
3443 cpl_image_get_maxpos(maxim, &maxx, &maxy);
3444
3445 cpl_image_delete(maxim);
3446
3447 cpl_errorstate_set(tempes);
3448
3449 data_float = cpl_malloc(ny * sizeof(float));
3450 data_double = cpl_malloc(ny * sizeof(double));
3451
3452 /* Write the individual spectrum files: */
3453 filecount = 0;
3454 for (i = 0; i < cpl_table_get_nrow(giraffe_table_get(fibertable)); ++i) {
3455 const double* wave_data;
3456 double* flux_data;
3457 double* err_data;
3458 cpl_size j;
3459 double snr = 0.0;
3460 /* Keywords to remove: */
3461 const char* remregexp = "^(CDELT[0-9]+|CD[0-9]+_[0-9]+|CRPIX[0-9]+"
3462 "|CRDER[0-9]+|CSYER[0-9]+|BUNIT|BSCALE|BZERO"
3463 "|BLANK|FILTER)$";
3464 cpl_size specindex = indices[i];
3465 double vela, velb, beta;
3466 double correction_factor = 1.0;
3467
3468 /* Skip non-science spectra. */
3469 if (strcmp(spectypes[i], "M") != 0) continue;
3470
3471 filename = cpl_sprintf(formatstr, ++filecount);
3472
3473 /* Calculate the Heliocentric correction factor to apply.
3474 * The gcorr and hcorr values are in km/s, so must be converted to m/s.
3475 */
3476 vela = gcorr[i] * 1e3;
3477 velb = hcorr[i] * 1e3;
3478 beta = (vela + velb) / CPL_PHYS_C;
3479 cpl_error_ensure(-1 <= beta && beta <= 1,
3480 CPL_ERROR_ILLEGAL_OUTPUT, goto cleanup,
3481 "The velocities GCORR = %g and HCORR = %g for spectrum"
3482 "%"CPL_SIZE_FORMAT" in file '%s' give invalid"
3483 " Heliocentric correction factor values.",
3484 gcorr[i], hcorr[i], specindex, flux_filename);
3485 correction_factor = sqrt((1.0 + beta) / (1.0 - beta));
3486
3487 /* Calculate and set corrected wavelength array, remembering to cast
3488 * to the appropriate type. */
3489 wave_data = cpl_array_get_data_double_const(refwavearray);
3490 if (wavecoltype == CPL_TYPE_FLOAT) {
3491 for (j = 0; j < ny; ++j) {
3492 data_float[j] = wave_data[j] * correction_factor;
3493 }
3494 array = cpl_array_wrap_float(data_float, ny);
3495 } else {
3496 for (j = 0; j < ny; ++j) {
3497 data_double[j] = wave_data[j] * correction_factor;
3498 }
3499 array = cpl_array_wrap_double(data_double, ny);
3500 }
3501 error |= irplib_sdp_spectrum_set_column_data(
3502 spectrum, GIALIAS_COLUMN_WAVE, array);
3503 cpl_array_unwrap(array);
3504 array = NULL;
3505
3506 fluximgcol = cpl_vector_new_from_image_column(
3507 giraffe_image_get(fluximage), specindex);
3508 flux_data = cpl_vector_get_data(fluximgcol);
3509 cpl_error_ensure(flux_data != NULL, cpl_error_get_code(), goto cleanup,
3510 "Unable to extract data in column %"CPL_SIZE_FORMAT
3511 " from image in file '%s'.", specindex, flux_filename);
3512 array = cpl_array_wrap_double(flux_data, ny);
3513
3514 if(i==maxx - 1){
3515 cpl_propertylist_update_char(extrakeys, GIALIAS_QCBRIGHTFLG, 'Y');
3516 }
3517 else{
3518 cpl_propertylist_update_char(extrakeys, GIALIAS_QCBRIGHTFLG, 'N');
3519 }
3520
3521 double fluxmean = cpl_array_get_mean(array);
3522 cpl_propertylist_update_double(extrakeys, GIALIAS_QCMEANRED, fluxmean);
3523 error |= irplib_sdp_spectrum_set_column_data(
3524 spectrum, GIALIAS_COLUMN_FLUX_REDUCED, array);
3525 cpl_array_unwrap(array);
3526 array = NULL;
3527
3528 errimgcol = cpl_vector_new_from_image_column(
3529 giraffe_image_get(errimage), specindex);
3530 err_data = cpl_vector_get_data(errimgcol);
3531 cpl_error_ensure(err_data != NULL, cpl_error_get_code(), goto cleanup,
3532 "Unable to extract data in column %"CPL_SIZE_FORMAT
3533 " from image in file '%s'.", specindex, err_filename);
3534 array = cpl_array_wrap_double(err_data, ny);
3535 error |= irplib_sdp_spectrum_set_column_data(
3536 spectrum, GIALIAS_COLUMN_ERR_REDUCED, array);
3537 cpl_array_unwrap(array);
3538 array = NULL;
3539
3540 for (j = 0; j < ny; ++j) {
3541 data_double[j] = (err_data[j] != 0.0) ? flux_data[j] / err_data[j]
3542 : 0.0;
3543 }
3544 array = cpl_array_wrap_double(data_double, ny);
3545 snr = cpl_array_get_median(array);
3546 double snrmean = cpl_array_get_mean(array);
3547 cpl_propertylist_update_double(extrakeys, GIALIAS_QCSNR, snrmean);
3548 error |= irplib_sdp_spectrum_set_column_data(spectrum,
3549 GIALIAS_COLUMN_SNR, array);
3550 cpl_array_unwrap(array);
3551 array = NULL;
3552
3553 cpl_vector_delete(fluximgcol);
3554 fluximgcol = NULL;
3555 cpl_vector_delete(errimgcol);
3556 errimgcol = NULL;
3557
3558 double estmag = NAN;
3559 if(cpl_table_has_column(giraffe_table_get(fibertable), "MAGNITUDE")){
3560 estmag = cpl_table_get_double(giraffe_table_get(fibertable), "MAGNITUDE", i, NULL);
3561 }
3562 if (isfinite(estmag)) {
3563 cpl_propertylist_update_double(extrakeys, GIALIAS_QCMAG, estmag);
3564 }
3565 else {
3566 cpl_propertylist_erase(extrakeys, GIALIAS_QCMAG);
3567 cpl_msg_warning(cpl_func, "Non-finite magnitude found for spectra %"CPL_SIZE_FORMAT
3568 " in file '%s'.", specindex, flux_filename);
3569 }
3570
3571 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
3572 "Could not setup the SDP spectrum columns for '%s'.",
3573 filename);
3574
3575 error |= irplib_sdp_spectrum_set_object(spectrum, objects[i]);
3576
3577 char * tmp_string;
3578
3579 tmp_string = cpl_sprintf("%s_%d_%s",
3580 objects[i],
3581 obsid,
3582 dateobs
3583 );
3584
3585 replace_spaces_with_underscores(tmp_string);
3586
3587 error |= irplib_sdp_spectrum_set_title(spectrum, tmp_string);
3588
3589 cpl_free(tmp_string);
3590
3591 error |= irplib_sdp_spectrum_set_ra(spectrum, ras[i]);
3592 error |= irplib_sdp_spectrum_set_dec(spectrum, decs[i]);
3593 error |= irplib_sdp_spectrum_set_snr(spectrum, snr);
3594 error |= irplib_sdp_spectrum_set_column_tcomm(
3595 spectrum, GIALIAS_COLUMN_FLUX_REDUCED, flux_filename);
3596 error |= irplib_sdp_spectrum_set_column_tcomm(
3597 spectrum, GIALIAS_COLUMN_ERR_REDUCED, err_filename);
3598
3599 error |= cpl_propertylist_update_int(extrakeys, GIALIAS_FPS, fps[i]);
3600 error |= cpl_propertylist_update_double(extrakeys, GIALIAS_GEOCORR,
3601 gcorr[i]);
3602 error |= cpl_propertylist_update_double(extrakeys, GIALIAS_HELICORR,
3603 hcorr[i]);
3604 error |= cpl_propertylist_update_double(extrakeys, GIALIAS_BARYCORR,
3605 bcorr[i]);
3606
3607 /* Set corrected values that depend on wavelengths. */
3608 error |= irplib_sdp_spectrum_set_wavelmin(spectrum,
3609 wavelmin * correction_factor);
3610 error |= irplib_sdp_spectrum_set_wavelmax(spectrum,
3611 wavelmax * correction_factor);
3612 error |= irplib_sdp_spectrum_set_specval(spectrum,
3613 (wavelmax + wavelmin) * 0.5 * correction_factor);
3614 error |= irplib_sdp_spectrum_set_specbw(spectrum,
3615 (wavelmax - wavelmin) * correction_factor);
3616 error |= irplib_sdp_spectrum_set_tdmin(spectrum,
3617 wavelmin * correction_factor);
3618 error |= irplib_sdp_spectrum_set_tdmax(spectrum,
3619 wavelmax * correction_factor);
3620 error |= irplib_sdp_spectrum_set_specbin(spectrum,
3621 specbin * correction_factor);
3622 if (! isnan(lamrms)) {
3623 error |= irplib_sdp_spectrum_set_lamrms(spectrum,
3624 lamrms * correction_factor);
3625 }
3626 if (! isnan(specerr)) {
3627 error |= irplib_sdp_spectrum_set_specerr(spectrum,
3628 specerr * correction_factor);
3629 }
3630 if (! isnan(specsye)) {
3631 error |= irplib_sdp_spectrum_set_specsye(spectrum,
3632 specsye * correction_factor);
3633 }
3634
3635 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
3636 "Could not setup the SDP spectrum keywords for '%s'.",
3637 filename);
3638
3639 error |= irplib_dfs_save_spectrum(allframes, NULL, parlist, usedframes,
3640 inherit, spectrum, recipe_id,
3641 extrakeys, tablekeys, remregexp,
3642 pipe_id, dict_id, filename);
3643 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
3644 "Failed to save SDP spectrum %"CPL_SIZE_FORMAT
3645 " to file '%s'.", specindex, filename);
3646
3647 error |= irplib_fits_update_checksums(filename);
3648 cpl_error_ensure(! error, cpl_error_get_code(), goto cleanup,
3649 "Failed to save update checksums for file '%s'.",
3650 filename);
3651 cpl_free(filename);
3652 filename = NULL;
3653 }
3654
3655 if (filecount == 0) {
3656 cpl_msg_warning(cpl_func, "No science spectra found in '%s'."
3657 " No SDP spectra created.", flux_filename);
3658 }
3659
3660 result_code = 0; /* Indicate success. */
3661
3662cleanup:
3663 /* Cleanup objects and memory. Note: the delete functions already check for
3664 * NULL pointers. */
3665 cpl_vector_delete(fluximgcol);
3666 cpl_vector_delete(errimgcol);
3667 cpl_array_unwrap(array);
3668 cpl_array_delete(refwavearray);
3669 cpl_free(data_float);
3670 cpl_free(data_double);
3671 cpl_free(filename);
3672 cpl_free(formatstr);
3673 cpl_frameset_delete(usedframes);
3674 cpl_frameset_delete(rawframes);
3675 cpl_frameset_iterator_delete(iterator);
3676 cpl_free(pipe_id);
3677 cpl_propertylist_delete(tablekeys);
3678 cpl_propertylist_delete(extrakeys);
3679 giraffe_image_delete(fluximage);
3680 giraffe_image_delete(errimage);
3681 giraffe_table_delete(fibertable);
3682 irplib_sdp_spectrum_delete(spectrum);
3683 return result_code;
3684}
giraffe_add_rvcorrection
cxint giraffe_add_rvcorrection(GiTable *fibers, const GiImage *spectra)
Add the barycentric and heliocentric corrections to the given fiber setup.
Definition giastrometry.c:69
giraffe_bias_config_create
GiBiasConfig * giraffe_bias_config_create(cpl_parameterlist *list)
Creates a setup structure for a bias removal task.
Definition gibias.c:3436
giraffe_bias_config_add
void giraffe_bias_config_add(cpl_parameterlist *list)
Adds parameters for the bias removal.
Definition gibias.c:3595
giraffe_bias_remove
cxint giraffe_bias_remove(GiImage *result, const GiImage *raw, const GiImage *master_bias, const GiImage *bad_pixels, const cpl_matrix *biaslimits, const GiBiasConfig *config)
Removes the bias from an image.
Definition gibias.c:3104
giraffe_bias_config_destroy
void giraffe_bias_config_destroy(GiBiasConfig *config)
Destroys a bias removal setup structure.
Definition gibias.c:3567
giraffe_subtract_dark
cxint giraffe_subtract_dark(GiImage *image, const GiImage *dark, const GiImage *bpixel, GiDarkResults *data, const GiDarkConfig *config)
Subtract the dark current from a bias corrected image.
Definition gidark.c:480
giraffe_extract_config_add
void giraffe_extract_config_add(cpl_parameterlist *list)
Adds parameters for the spectrum extraction.
Definition giextract.c:3509
giraffe_extract_spectra
cxint giraffe_extract_spectra(GiExtraction *result, GiImage *image, GiTable *fibers, GiLocalization *sloc, GiImage *bpixel, GiImage *slight, GiExtractConfig *config)
Extracts the spectra from a preprocessed frame.
Definition giextract.c:2480
giraffe_extract_config_create
GiExtractConfig * giraffe_extract_config_create(cpl_parameterlist *list)
Creates a setup structure for the spectrum extraction.
Definition giextract.c:3405
giraffe_extract_config_destroy
void giraffe_extract_config_destroy(GiExtractConfig *config)
Destroys a spectrum extraction setup structure.
Definition giextract.c:3479
giraffe_fibers_setup
GiTable * giraffe_fibers_setup(const cpl_frame *frame, const cpl_frame *reference)
Setup a fiber list.
Definition gifibers.c:218
giraffe_fiberlist_load
GiTable * giraffe_fiberlist_load(const cxchar *filename, cxint dataset, const cxchar *tag)
Load a fiber table from a file.
Definition gifiberutils.c:754
giraffe_fiberlist_compare
cxint giraffe_fiberlist_compare(const GiTable *fibers, const GiTable *reference)
Compare two fiber lists.
Definition gifiberutils.c:952
giraffe_fiberlist_attach
cxint giraffe_fiberlist_attach(cpl_frame *frame, GiTable *fibers)
Attach a fiber table to a frame.
Definition gifiberutils.c:884
giraffe_flat_config_create
GiFlatConfig * giraffe_flat_config_create(cpl_parameterlist *list)
Creates a setup structure for the flat field correction.
Definition giflat.c:302
giraffe_flat_config_destroy
void giraffe_flat_config_destroy(GiFlatConfig *config)
Destroys a flat field setup structure.
Definition giflat.c:353
giraffe_flat_config_add
void giraffe_flat_config_add(cpl_parameterlist *list)
Adds parameters for the flat field correction.
Definition giflat.c:376
giraffe_flat_apply
cxint giraffe_flat_apply(GiExtraction *extraction, const GiTable *fibers, const GiImage *flat, const GiImage *errors, GiFlatConfig *config)
Apply the flat field correction to the given extracted spectra.
Definition giflat.c:238
giraffe_fov_config_create
GiFieldOfViewConfig * giraffe_fov_config_create(cpl_parameterlist *list)
Creates a setup structure for the field of view reconstruction.
Definition gifov.c:2002
giraffe_fov_config_destroy
void giraffe_fov_config_destroy(GiFieldOfViewConfig *config)
Destroys a field of view setup structure.
Definition gifov.c:2057
giraffe_fov_new
GiFieldOfView * giraffe_fov_new(void)
Create an empty container for the results of the field of view reconstruction.
Definition gifov.c:1383
giraffe_fov_save_cubes_eso3d
cxint giraffe_fov_save_cubes_eso3d(const GiFieldOfView *self, cpl_propertylist *properties, const cxchar *filename, cxptr data)
Write the cube components of a field-of-view object to a file.
Definition gifov.c:1663
giraffe_fov_delete
void giraffe_fov_delete(GiFieldOfView *self)
Deallocate a field of view object and its contents.
Definition gifov.c:1484
giraffe_fov_build
cxint giraffe_fov_build(GiFieldOfView *result, GiRebinning *rebinning, GiTable *fibers, GiTable *wsolution, GiTable *grating, GiTable *slitgeometry, GiFieldOfViewConfig *config)
Create and image and a data cube from extracted and rebinned spectra.
Definition gifov.c:418
giraffe_fov_save_cubes
cxint giraffe_fov_save_cubes(const GiFieldOfView *self, cpl_propertylist *properties, const cxchar *filename, cxptr data)
Write the cube components of a field-of-view object to a file.
Definition gifov.c:1520
giraffe_fov_config_add
void giraffe_fov_config_add(cpl_parameterlist *list)
Adds parameters for the image and data cube construction.
Definition gifov.c:2079
giraffe_frame_create
cpl_frame * giraffe_frame_create(const cxchar *tag, cpl_frame_level level, const cpl_propertylist *properties, cxcptr object, cxcptr data, GiFrameCreator creator)
Create a product frame using a provided frame creator.
Definition giframe.c:239
giraffe_frame_create_image
cpl_frame * giraffe_frame_create_image(GiImage *image, const cxchar *tag, cpl_frame_level level, cxbool save, cxbool update)
Create an image product frame.
Definition giframe.c:395
giraffe_get_slitgeometry
cpl_frame * giraffe_get_slitgeometry(const cpl_frameset *set)
Get the slit geometry frame from a frame set.
Definition giframe.c:786
giraffe_image_get
cpl_image * giraffe_image_get(const GiImage *self)
Gets the image data.
Definition giimage.c:218
giraffe_image_get_properties
cpl_propertylist * giraffe_image_get_properties(const GiImage *self)
Get the properties of an image.
Definition giimage.c:282
giraffe_image_delete
void giraffe_image_delete(GiImage *self)
Destroys an image.
Definition giimage.c:181
giraffe_image_create
GiImage * giraffe_image_create(cpl_type type, cxint nx, cxint ny)
Creates an image container of a given type.
Definition giimage.c:95
giraffe_image_add_info
cxint giraffe_image_add_info(GiImage *image, const GiRecipeInfo *info, const cpl_frameset *set)
Add additional frame information to an image.
Definition giimage.c:773
giraffe_image_save
cxint giraffe_image_save(GiImage *self, const cxchar *filename)
Write a Giraffe image to a file.
Definition giimage.c:570
giraffe_image_new
GiImage * giraffe_image_new(cpl_type type)
Creates an empty image container.
Definition giimage.c:65
giraffe_image_set_properties
cxint giraffe_image_set_properties(GiImage *self, cpl_propertylist *properties)
Attaches a property list to an image.
Definition giimage.c:312
giraffe_image_load
cxint giraffe_image_load(GiImage *self, const cxchar *filename, cxint position)
Gets image data and properties from a file.
Definition giimage.c:536
giraffe_rebin_config_destroy
void giraffe_rebin_config_destroy(GiRebinConfig *config)
Destroys a spectrum extraction setup structure.
Definition girebinning.c:4926
giraffe_rebin_spectra
cxint giraffe_rebin_spectra(GiRebinning *rebinning, const GiExtraction *extraction, const GiTable *fibers, const GiLocalization *localization, const GiTable *grating, const GiTable *slitgeo, const GiTable *solution, const GiRebinConfig *config)
Rebin an Extracted Spectra Frame and associated Errors Frame.
Definition girebinning.c:4052
giraffe_rebin_config_create
GiRebinConfig * giraffe_rebin_config_create(cpl_parameterlist *list)
Creates a setup structure for the rebinning.
Definition girebinning.c:4826
giraffe_rebinning_new
GiRebinning * giraffe_rebinning_new(void)
Create an empty rebinning results container.
Definition girebinning.c:4694
giraffe_rebinning_destroy
void giraffe_rebinning_destroy(GiRebinning *rebinning)
Destroys a rebinning results container and its contents.
Definition girebinning.c:4788
giraffe_rebin_config_add
void giraffe_rebin_config_add(cpl_parameterlist *list)
Adds parameters for the rebinning.
Definition girebinning.c:4950
giraffe_sgcalibration_config_create
GiSGCalConfig * giraffe_sgcalibration_config_create(cpl_parameterlist *list)
Creates a setup structure for the slit geometry calibration.
Definition gisgcalibration.c:2428
giraffe_compute_offsets
cxint giraffe_compute_offsets(GiTable *fibers, const GiRebinning *rebinning, const GiTable *grating, const GiTable *mask, const GiSGCalConfig *config)
Compute wavelength offsets for a set of rebinned input spectrum.
Definition gisgcalibration.c:2185
giraffe_sgcalibration_config_destroy
void giraffe_sgcalibration_config_destroy(GiSGCalConfig *config)
Destroys a sgcalibration field setup structure.
Definition gisgcalibration.c:2638
giraffe_sgcalibration_config_add
void giraffe_sgcalibration_config_add(cpl_parameterlist *list)
Adds parameters for the sgcalibration correction computation.
Definition gisgcalibration.c:2669
giraffe_slitgeometry_load
GiTable * giraffe_slitgeometry_load(const GiTable *fibers, const cxchar *filename, cxint pos, const cxchar *tag)
Load the slit geometry information for a given fiber setup.
Definition gislitgeometry.c:602
giraffe_table_new
GiTable * giraffe_table_new(void)
Creates a new, empty Giraffe table.
Definition gitable.c:85
giraffe_table_duplicate
GiTable * giraffe_table_duplicate(const GiTable *src)
Duplicate a Giraffe table.
Definition gitable.c:176
giraffe_table_load
cxint giraffe_table_load(GiTable *self, const cxchar *filename, cxint position, const cxchar *id)
Reads a data set from a file into a Giraffe table.
Definition gitable.c:562
giraffe_table_delete
void giraffe_table_delete(GiTable *self)
Destroys a Giraffe table.
Definition gitable.c:154
giraffe_table_get
cpl_table * giraffe_table_get(const GiTable *self)
Get the table data from a Giraffe table.
Definition gitable.c:433
giraffe_add_frameset_info
cxint giraffe_add_frameset_info(cpl_propertylist *plist, const cpl_frameset *set, cxint sequence)
Add frameset specific information to a property list.
Definition giutils.c:786
giraffe_propertylist_copy
cxint giraffe_propertylist_copy(cpl_propertylist *self, const cxchar *name, const cpl_propertylist *other, const cxchar *othername)
Copy a property from one list to another.
Definition giutils.c:1108
giraffe_get_license
const cxchar * giraffe_get_license(void)
Get the pipeline copyright and license.
Definition giutils.c:422
giraffe_get_mode
GiInstrumentMode giraffe_get_mode(cpl_propertylist *properties)
Determines the instrument mode from a property list.
Definition giutils.c:444
GiSGCalConfig
Slit geometry calibration configuration data structure.
Definition gisgcalibration.h:49
This file is part of the GIRAFFE Pipeline Reference Manual 2.19.5.
Documentation copyright © 2002-2006 European Southern Observatory.
Generated on Fri Sep 18 2026 09:09:29 by doxygen 1.13.2 written by Dimitri van Heesch, © 1997-2004