CR2RE Pipeline Reference Manual 1.6.10
hdrl_spectrumlist.c
1/*
2 * This file is part of the HDRL
3 * Copyright (C) 2017 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/*-----------------------------------------------------------------------------
25 Includes
26 -----------------------------------------------------------------------------*/
27#include "hdrl_spectrumlist.h"
28#include "hdrl_spectrum_resample.h"
29#include "hdrl_imagelist.h"
30
31#include <math.h>
32
38/*----------------------------------------------------------------------------*/
39
42/*-----------------------------------------------------------------------------
43 Private Functions
44 -----------------------------------------------------------------------------*/
45static inline
46cpl_error_code check_getter(const hdrl_spectrum1Dlist * s, const cpl_size idx);
47
48static inline
49void push_back(hdrl_spectrum1Dlist *self, hdrl_spectrum1D * s);
50
51static inline hdrl_spectrum1D **
52safe_realloc(hdrl_spectrum1D ** current, const cpl_size old_capacity,
53 const cpl_size new_capacity);
54
55static inline hdrl_imagelist *
56create_list(hdrl_spectrum1D ** resampled_spectra,
57 const hdrl_spectrum1Dlist * ori_spectra,
58 const cpl_boolean mark_bpm_in_interpolation);
59
60static inline void
61delete_spectrum1D_array(hdrl_spectrum1D ** resampled_spectra,
62 const cpl_size num_spectra);
63
64static inline hdrl_image *
65get_padded_flux(const hdrl_spectrum1D * resampled_spectrum,
66 const hdrl_spectrum1D * ori_spectrum,
67 const cpl_boolean mark_bpm_in_interpolation);
68
69static inline cpl_boolean
70are_spectra_valid(const hdrl_spectrum1Dlist * list);
71
72static inline cpl_error_code
73get_first_error_code(const cpl_error_code * cds, const cpl_size length);
74
75static inline cpl_boolean
76check_scales_are_same(const hdrl_spectrum1Dlist * list);
77
78static inline double
79get_wmin_valid(const hdrl_spectrum1D * s);
80
81static inline double
82get_wmax_valid(const hdrl_spectrum1D * s);
83
84static inline void
85remove_if_neighbors_are_rejected(hdrl_image * flx,
86 const hdrl_spectrum1D * ori_spectrum, const cpl_array * wlens);
87
88static inline hdrl_spectrum1D *
89get_interp_bpm(const hdrl_spectrum1D * s, const cpl_array * wlens);
90
91static inline cpl_boolean
92contains(const hdrl_spectrum1Dlist * list, const hdrl_spectrum1D * s);
93/* ---------------------------------------------------------------------------*/
101/* ---------------------------------------------------------------------------*/
102hdrl_spectrum1Dlist * hdrl_spectrum1Dlist_new(void){
103 hdrl_spectrum1Dlist * to_ret = cpl_calloc(1, sizeof(*to_ret));
104
105 to_ret->length = 0;
106 to_ret->capacity= 0;
107 to_ret->spectra = NULL;
108 return to_ret;
109}
110
111/* ---------------------------------------------------------------------------*/
121/* ---------------------------------------------------------------------------*/
122hdrl_spectrum1Dlist * hdrl_spectrum1Dlist_duplicate(const hdrl_spectrum1Dlist * l){
123
124 if(l == NULL) return NULL;
125
126 hdrl_spectrum1Dlist * to_ret = hdrl_spectrum1Dlist_new();
127
128 for(cpl_size i = 0; i < hdrl_spectrum1Dlist_get_size(l); ++i){
129 const hdrl_spectrum1D * s_old = hdrl_spectrum1Dlist_get_const(l, i);
130 hdrl_spectrum1D * s = hdrl_spectrum1D_duplicate(s_old);
131 hdrl_spectrum1Dlist_set(to_ret, s, i);
132 }
133 return to_ret;
134}
135
147/* ---------------------------------------------------------------------------*/
148hdrl_spectrum1Dlist *
149hdrl_spectrum1Dlist_wrap(hdrl_spectrum1D ** self, const cpl_size sz){
150 hdrl_spectrum1Dlist * to_ret = hdrl_spectrum1Dlist_new();
151 to_ret->spectra = self;
152 to_ret->capacity = to_ret->length = sz;
153 return to_ret;
154}
155
166/* ---------------------------------------------------------------------------*/
167hdrl_spectrum1D *
168hdrl_spectrum1Dlist_get(hdrl_spectrum1Dlist * self, const cpl_size idx){
169
170 cpl_error_code fail = check_getter(self, idx);
171 cpl_ensure(fail == CPL_ERROR_NONE, fail, NULL);
172
173 return self->spectra[idx];
174
175}
176
187/* ---------------------------------------------------------------------------*/
188const hdrl_spectrum1D *
189hdrl_spectrum1Dlist_get_const(const hdrl_spectrum1Dlist *self , const cpl_size idx){
190
191 cpl_error_code fail = check_getter(self, idx);
192 cpl_ensure(fail == CPL_ERROR_NONE, fail, NULL);
193
194 return self->spectra[idx];
195}
196
214/* ---------------------------------------------------------------------------*/
215cpl_error_code
216hdrl_spectrum1Dlist_set(hdrl_spectrum1Dlist * self,
217 hdrl_spectrum1D * s, const cpl_size idx){
218
219 cpl_ensure_code(self != NULL, CPL_ERROR_NULL_INPUT);
220 cpl_ensure_code(self->length >= 0, CPL_ERROR_ACCESS_OUT_OF_RANGE);
221 cpl_ensure_code(idx <= self->length, CPL_ERROR_ACCESS_OUT_OF_RANGE);
222 cpl_ensure_code(!contains(self, s), CPL_ERROR_ILLEGAL_INPUT);
223
224 if(idx == self->length){
225 push_back(self, s);
226 return CPL_ERROR_NONE;
227 }
228
229 hdrl_spectrum1D * old = self->spectra[idx];
231
232 self->spectra[idx] = s;
233
234 return CPL_ERROR_NONE;
235}
236
252/* ---------------------------------------------------------------------------*/
253hdrl_spectrum1D *
254hdrl_spectrum1Dlist_unset(hdrl_spectrum1Dlist * self, const cpl_size idx){
255
256 cpl_error_code fail = check_getter(self, idx);
257 cpl_ensure(fail == CPL_ERROR_NONE, fail, NULL);
258
259 hdrl_spectrum1D * to_ret = self->spectra[idx];
260
261 const cpl_size sz = hdrl_spectrum1Dlist_get_size(self);
262
263 for(cpl_size i = idx; i < sz - 1; ++i){
264 self->spectra[i] = self->spectra[i + 1];
265 }
266 self->length--;
267
268 const cpl_size new_capacity = self->capacity / 2;
269
270 if(hdrl_spectrum1Dlist_get_size(self) <= new_capacity){
271 self->spectra = safe_realloc(self->spectra, self->capacity, new_capacity);
272 self->capacity = new_capacity;
273 }
274
275 return to_ret;
276}
277
285/* ---------------------------------------------------------------------------*/
286void hdrl_spectrum1Dlist_delete(hdrl_spectrum1Dlist * l){
287
288 if(l == NULL) return;
289
290 for(cpl_size i = 0; i < hdrl_spectrum1Dlist_get_size(l); ++i){
291 hdrl_spectrum1D_delete(&l->spectra[i]);
292 }
293 cpl_free(l->spectra);
294 cpl_free(l);
295}
296
307/* ---------------------------------------------------------------------------*/
308cpl_size
309hdrl_spectrum1Dlist_get_size(const hdrl_spectrum1Dlist * l){
310
311 cpl_ensure(l != NULL, CPL_ERROR_NULL_INPUT, 0);
312
313 return l->length;
314}
315
316/* ---------------------------------------------------------------------------*/
348/* ---------------------------------------------------------------------------*/
349cpl_error_code
350hdrl_spectrum1Dlist_collapse (const hdrl_spectrum1Dlist *list,
351 const hdrl_parameter * stacking_par,
352 const cpl_array * wlengths, const hdrl_parameter * resample_par,
353 const cpl_boolean mark_bpm_in_interpolation,
354 hdrl_spectrum1D ** result, cpl_image ** contrib,
355 hdrl_imagelist ** resampled_and_aligned_fluxes){
356
357 cpl_ensure_code(are_spectra_valid(list), CPL_ERROR_NULL_INPUT);
358 cpl_ensure_code(wlengths != NULL, CPL_ERROR_NULL_INPUT);
359 cpl_ensure_code(check_scales_are_same(list), CPL_ERROR_ILLEGAL_INPUT);
360 cpl_ensure_code(result != NULL, CPL_ERROR_NULL_INPUT);
361 cpl_ensure_code(resampled_and_aligned_fluxes != NULL, CPL_ERROR_NULL_INPUT);
362 *result = NULL;
363 *contrib = NULL;
364
365 const cpl_size num_spectra = hdrl_spectrum1Dlist_get_size(list);
366 hdrl_spectrum1D ** resampled_spectra =
367 cpl_calloc(num_spectra, sizeof(hdrl_spectrum1D*));
368
369 cpl_ensure_code(num_spectra > 0, CPL_ERROR_ILLEGAL_INPUT);
370 cpl_error_code * errors = cpl_calloc(num_spectra, sizeof(cpl_error_code));
371
372 HDRL_OMP(omp parallel for)
373 for(cpl_size i = 0; i < num_spectra; ++i){
374 const hdrl_spectrum1D * this_s = hdrl_spectrum1Dlist_get_const(list, i);
375 resampled_spectra[i] =
376 hdrl_spectrum1D_resample_on_array(this_s, wlengths, resample_par);
377 errors[i] = cpl_error_get_code();
378 }
379
380
381 cpl_error_code err = get_first_error_code(errors, num_spectra);
382 cpl_free(errors);
383
384 if(!err){
385 hdrl_imagelist * stack_list =
386 create_list(resampled_spectra, list, mark_bpm_in_interpolation);
387
388 hdrl_image * stacked_img = NULL;
389
390 err = hdrl_imagelist_collapse(stack_list, stacking_par,
391 &stacked_img, contrib);
392
393 *resampled_and_aligned_fluxes = stack_list;
394
395 if(!err){
396 const hdrl_spectrum1D_wave_scale scale =
398
399 *result = hdrl_spectrum1D_create(hdrl_image_get_image(stacked_img),
400 hdrl_image_get_error(stacked_img),
401 wlengths,
402 scale);
403 }
404
405 hdrl_image_delete(stacked_img);
406 }
407
408 delete_spectrum1D_array(resampled_spectra, num_spectra);
409 return err;
410}
411
412/*-----------------------------------------------------------------------------
413 Private Functions Implementation
414 -----------------------------------------------------------------------------*/
415/*Gets the minimum wavelengths ignoring bad pixels*/
416static inline double
417get_wmin_valid(const hdrl_spectrum1D * s){
418 const cpl_size sz = hdrl_spectrum1D_get_size(s);
419
420 double wmin = INFINITY;
421
422 for(cpl_size i = 0; i < sz; ++i){
423 int rej = 0;
424 double w = hdrl_spectrum1D_get_wavelength_value(s, i, &rej);
425 if(rej) continue;
426 if(w < wmin) wmin = w;
427 }
428 return wmin;
429}
430/*Gets the maximum wavelengths ignoring bad pixels*/
431static inline double
432get_wmax_valid(const hdrl_spectrum1D * s){
433 const cpl_size sz = hdrl_spectrum1D_get_size(s);
434
435 double wmax = -INFINITY;
436
437 for(cpl_size i = 0; i < sz; ++i){
438 int rej = 0;
439 double w = hdrl_spectrum1D_get_wavelength_value(s, i, &rej);
440 if(rej) continue;
441 if(w > wmax) wmax = w;
442 }
443 return wmax;
444}
445
446static inline cpl_image *
447get_img_from_bpm(hdrl_spectrum1D_wavelength * w){
448 if(w->bpm)
449 return cpl_image_new_from_mask(w->bpm);
450
451 const cpl_size sz = cpl_array_get_size(w->wavelength);
452 return cpl_image_new(sz, 1, CPL_TYPE_INT);
453}
454
455static inline hdrl_spectrum1D *
456get_interp_bpm(const hdrl_spectrum1D * s, const cpl_array * wlens){
457
458 hdrl_spectrum1D_wavelength wlen_ori = hdrl_spectrum1D_get_wavelength(s);
459 cpl_image * flx = get_img_from_bpm(&wlen_ori);
460 hdrl_spectrum1D * bpm = hdrl_spectrum1D_create_error_free(flx,
461 wlen_ori.wavelength, wlen_ori.scale);
462 cpl_image_delete(flx);
463
464 hdrl_parameter * par =
466 hdrl_spectrum1D_interp_linear);
467 hdrl_spectrum1D * interp_bpm = hdrl_spectrum1D_resample_on_array(bpm, wlens,
468 par);
471 return interp_bpm;
472}
473
474static inline void
475remove_if_neighbors_are_rejected(hdrl_image * flx,
476 const hdrl_spectrum1D * ori_spectrum, const cpl_array * wlens){
477 /*We interpolate the bpm to see if the elements in wlens have a bad pixel
478 * close by. If that is the case we reject them.*/
479 hdrl_spectrum1D * interp_bpm = get_interp_bpm(ori_spectrum, wlens);
480
481 for(cpl_size i = 0; i < hdrl_spectrum1D_get_size(interp_bpm); ++i){
482 const double v = hdrl_spectrum1D_get_flux_value(interp_bpm, i, NULL).data;
483 if(v > HDRL_EPS_DATA){
484 hdrl_image_reject(flx, i + 1, 1);
485 }
486 }
487 hdrl_spectrum1D_delete(&interp_bpm);
488}
489
490static inline hdrl_image *
491get_padded_flux(const hdrl_spectrum1D * resampled_spectrum,
492 const hdrl_spectrum1D * ori_spectrum,
493 const cpl_boolean mark_bpm_in_interpolation){
494
495 if(resampled_spectrum == NULL) return NULL;
496
497 const double wmin = get_wmin_valid(ori_spectrum);
498
499 const double wmax = get_wmax_valid(ori_spectrum);
500
501 if(isinf(wmin) || isinf(wmax)) return NULL;
502
503 hdrl_image * flx =
505
506 const cpl_array * wlens =
507 hdrl_spectrum1D_get_wavelength(resampled_spectrum).wavelength;
508
509 for(cpl_size i = 0; i < hdrl_spectrum1D_get_size(resampled_spectrum); ++i){
510 const double wlen = cpl_array_get(wlens, i, NULL);
511 if(wlen < wmin || wlen > wmax) {
512 hdrl_image_reject(flx, i + 1, 1);
513 }
514 }
515
516 if(mark_bpm_in_interpolation)
517 remove_if_neighbors_are_rejected(flx, ori_spectrum, wlens);
518
519 return flx;
520}
521
522static inline cpl_boolean
523check_scales_are_same(const hdrl_spectrum1Dlist * list){
524 if(list == NULL) return CPL_TRUE;
525 const cpl_size sz = hdrl_spectrum1Dlist_get_size(list);
526
527 if(sz <= 1) return CPL_TRUE;
528
529 const hdrl_spectrum1D_wave_scale scale =
531
532 for(cpl_size i = 1; i < sz; ++i){
533 const hdrl_spectrum1D_wave_scale scale_this =
535
536 if(scale_this != scale) return CPL_FALSE;
537 }
538
539 return CPL_TRUE;
540}
541/*Given a array of hdrl_spectrum1D this function creates an hdrl_imagelist containing
542 * the fluxes of the spectra. The current implementation of resample extrapolates
543 * points outside the interval where values are available. We reject those points.*/
544static inline hdrl_imagelist *
545create_list(hdrl_spectrum1D ** resampled_spectra,
546 const hdrl_spectrum1Dlist * ori_spectra,
547 const cpl_boolean mark_bpm_in_interpolation){
548
549 const cpl_size num_spectra = hdrl_spectrum1Dlist_get_size(ori_spectra);
550 hdrl_image ** images = cpl_calloc(num_spectra, sizeof(hdrl_image*));
551 cpl_error_code * errs = cpl_calloc(num_spectra, sizeof(cpl_error_code));
552 HDRL_OMP(omp parallel for)
553 for(cpl_size i = 0; i < num_spectra; i++){
554 hdrl_image * img = get_padded_flux(resampled_spectra[i],
555 hdrl_spectrum1Dlist_get_const(ori_spectra, i),
556 mark_bpm_in_interpolation);
557 images[i] = img;
558 errs[i] = cpl_error_get_code();
559 }
560
561 cpl_error_code fail = get_first_error_code(errs, num_spectra);
562 cpl_free(errs);
563
564 hdrl_imagelist * list = NULL;
565 if(!fail){
566 list = hdrl_imagelist_new();
567
568 /* hdrl_imagelist_set has a race condition when growing the internal buffer
569 * it cannot be parallelized */
570 for(cpl_size i = 0; i < num_spectra; i++){
571 if(images[i] != NULL)
572 hdrl_imagelist_set(list, images[i], i);
573 }
574 }
575 cpl_free(images);
576 return list;
577}
578/*Function used in the getters: check that idx is in the correct ranges.*/
579static inline
580cpl_error_code check_getter(const hdrl_spectrum1Dlist *s, const cpl_size idx){
581
582 if(s == NULL) return CPL_ERROR_NULL_INPUT;
583
584 if(idx < 0) return CPL_ERROR_ACCESS_OUT_OF_RANGE;
585
586 if(idx >= s->length) return CPL_ERROR_ACCESS_OUT_OF_RANGE;
587
588 return CPL_ERROR_NONE;
589}
590/*wrapper around realloc, works are realloc, except if old_capacity == 0 (in this
591 * case it works like calloc) and if new_capacity == 0 (in this case it works like free)*/
592static inline hdrl_spectrum1D **
593safe_realloc(hdrl_spectrum1D ** current, const cpl_size old_capacity,
594 const cpl_size new_capacity){
595
596 /* Going from empty list to list with elements */
597 if(old_capacity == 0 && new_capacity > 0){
598 return cpl_calloc(new_capacity, sizeof(hdrl_spectrum1D * ));
599 }
600
601 /* Clean capacity */
602 if(new_capacity == 0){
603 cpl_free(current);
604 return NULL;
605 }
606
607 hdrl_spectrum1D ** to_ret =
608 cpl_realloc(current, sizeof(hdrl_spectrum1D * ) * new_capacity);
609
610 for(cpl_size i = old_capacity + 1; i < new_capacity; ++i){
611 to_ret[i] = NULL;
612 }
613
614 return to_ret;
615}
616
617/*increase the capacity of the list if length >= capacity*/
618static inline
619void resize_if_needed(hdrl_spectrum1Dlist *self){
620
621 if(self->length < self->capacity) return;
622
623 const cpl_size new_capacity = self->capacity == 0 ? 1 : 2 * self->capacity;
624
625 self->spectra = safe_realloc(self->spectra, self->capacity,
626 new_capacity);
627 self->capacity = new_capacity;
628}
629
630/*add an element in the back of a list, it can increase the capacity if needed*/
631static inline
632void push_back(hdrl_spectrum1Dlist *self, hdrl_spectrum1D * s){
633
634 resize_if_needed(self);
635 self->spectra[self->length] = s;
636 self->length++;
637}
638
639/*frees an array of spectrum1D and frees every spectrum1D*/
640static inline void
641delete_spectrum1D_array(hdrl_spectrum1D ** resampled_spectra,
642 const cpl_size num_spectra){
643
644 hdrl_spectrum1Dlist * l =
645 hdrl_spectrum1Dlist_wrap(resampled_spectra, num_spectra);
647}
648
649/*check that all the spectra in the list are != NULL*/
650static inline cpl_boolean
651are_spectra_valid(const hdrl_spectrum1Dlist * list){
652 if(list == NULL) return CPL_FALSE;
653
654 const cpl_size length = hdrl_spectrum1Dlist_get_size(list);
655 for(cpl_size i = 0; i < length; ++i){
656 const hdrl_spectrum1D * s = hdrl_spectrum1Dlist_get_const(list, i);
657 if(s == NULL) return CPL_FALSE;
658 }
659 return CPL_TRUE;
660}
661/*get the first element in the array different from CPL_ERROR_NONE*/
662static inline cpl_error_code
663get_first_error_code(const cpl_error_code * cds, const cpl_size length){
664 for(cpl_size i = 0; i < length; ++i){
665 if(cds[i]) return cds[i];
666 }
667 return CPL_ERROR_NONE;
668}
669
670static inline cpl_boolean
671contains(const hdrl_spectrum1Dlist * list, const hdrl_spectrum1D * s){
672 const cpl_size sz = hdrl_spectrum1Dlist_get_size(list);
673 for(cpl_size i = 0; i < sz; ++i){
674 const hdrl_spectrum1D * s_in = hdrl_spectrum1Dlist_get_const(list, i);
675 if(s_in == s) return CPL_TRUE;
676 }
677 return CPL_FALSE;
678}
679
hdrl_image * hdrl_image_duplicate(const hdrl_image *himg)
copy hdrl_image
Definition: hdrl_image.c:391
cpl_image * hdrl_image_get_error(hdrl_image *himg)
get error as cpl image
Definition: hdrl_image.c:131
cpl_image * hdrl_image_get_image(hdrl_image *himg)
get data as cpl image
Definition: hdrl_image.c:105
cpl_error_code hdrl_image_reject(hdrl_image *self, cpl_size xpos, cpl_size ypos)
mark pixel as bad
Definition: hdrl_image.c:427
void hdrl_image_delete(hdrl_image *himg)
delete hdrl_image
Definition: hdrl_image.c:379
cpl_error_code hdrl_imagelist_set(hdrl_imagelist *himlist, hdrl_image *himg, cpl_size pos)
Insert an image into an imagelist.
cpl_error_code hdrl_imagelist_collapse(const hdrl_imagelist *himlist, const hdrl_parameter *param, hdrl_image **out, cpl_image **contrib)
collapsing of image list
hdrl_imagelist * hdrl_imagelist_new(void)
Create an empty imagelist.
void hdrl_parameter_delete(hdrl_parameter *obj)
shallow delete of a parameter
const hdrl_image * hdrl_spectrum1D_get_flux(const hdrl_spectrum1D *self)
hdrl_spectrum1D getter flux
cpl_size hdrl_spectrum1Dlist_get_size(const hdrl_spectrum1Dlist *l)
hdrl_spectrum1Dlist getter for size
hdrl_spectrum1Dlist * hdrl_spectrum1Dlist_new(void)
hdrl_spectrum1Dlist default constructor
hdrl_parameter * hdrl_spectrum1D_resample_interpolate_parameter_create(const hdrl_spectrum1D_interpolation_method method)
constructor for the hdrl_parameter in the case of interpolation
hdrl_spectrum1D * hdrl_spectrum1D_duplicate(const hdrl_spectrum1D *self)
hdrl_spectrum1D copy constructor
cpl_size hdrl_spectrum1D_get_size(const hdrl_spectrum1D *self)
hdrl_spectrum1D getter for size
hdrl_spectrum1D * hdrl_spectrum1Dlist_unset(hdrl_spectrum1Dlist *self, const cpl_size idx)
hdrl_spectrum1Dlist remove of the i-th element
void hdrl_spectrum1D_delete(hdrl_spectrum1D **p_self)
hdrl_spectrum1D destructor
cpl_error_code hdrl_spectrum1Dlist_set(hdrl_spectrum1Dlist *self, hdrl_spectrum1D *s, const cpl_size idx)
hdrl_spectrum1Dlist setter of the i-th element
void hdrl_spectrum1Dlist_delete(hdrl_spectrum1Dlist *l)
hdrl_spectrum1Dlist destructor
hdrl_spectrum1D * hdrl_spectrum1Dlist_get(hdrl_spectrum1Dlist *self, const cpl_size idx)
hdrl_spectrum1Dlist getter of the i-th element
const hdrl_spectrum1D * hdrl_spectrum1Dlist_get_const(const hdrl_spectrum1Dlist *self, const cpl_size idx)
hdrl_spectrum1Dlist getter of the i-th element
hdrl_spectrum1D * hdrl_spectrum1D_create(const cpl_image *arg_flux, const cpl_image *arg_flux_e, const cpl_array *wavelength, hdrl_spectrum1D_wave_scale wave_scale)
hdrl_spectrum1D default constructor
hdrl_spectrum1D * hdrl_spectrum1D_resample_on_array(const hdrl_spectrum1D *self, const cpl_array *waves, const hdrl_parameter *par)
resample a hdrl_spectrum1D on the wavelengths contained in waves
hdrl_spectrum1D_wave_scale hdrl_spectrum1D_get_scale(const hdrl_spectrum1D *self)
hdrl_spectrum1D getter for scale
hdrl_spectrum1D_wavelength hdrl_spectrum1D_get_wavelength(const hdrl_spectrum1D *self)
hdrl_spectrum1D getter for wavelengths
hdrl_spectrum1D * hdrl_spectrum1D_create_error_free(const cpl_image *arg_flux, const cpl_array *wavelength, hdrl_spectrum1D_wave_scale scale)
hdrl_spectrum1D constructor in the case of error-free spectrum (i.e. the error on the flux is zero fo...
hdrl_data_t hdrl_spectrum1D_get_wavelength_value(const hdrl_spectrum1D *self, int idx, int *rej)
hdrl_spectrum1D getter for a wavelength value
hdrl_spectrum1Dlist * hdrl_spectrum1Dlist_duplicate(const hdrl_spectrum1Dlist *l)
hdrl_spectrum1Dlist copy-constructor
hdrl_spectrum1Dlist * hdrl_spectrum1Dlist_wrap(hdrl_spectrum1D **self, const cpl_size sz)
hdrl_spectrum1Dlist wrapper
cpl_error_code hdrl_spectrum1Dlist_collapse(const hdrl_spectrum1Dlist *list, const hdrl_parameter *stacking_par, const cpl_array *wlengths, const hdrl_parameter *resample_par, const cpl_boolean mark_bpm_in_interpolation, hdrl_spectrum1D **result, cpl_image **contrib, hdrl_imagelist **resampled_and_aligned_fluxes)
collapsing a hdrl_spectrum1Dlist. The spectra in list are first resampled on the wavelengths wlengths...
hdrl_value hdrl_spectrum1D_get_flux_value(const hdrl_spectrum1D *self, int idx, int *rej)
hdrl_spectrum1D getter for a flux value