CR2RE Pipeline Reference Manual 1.6.10
hdrl_spectrum.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_spectrum.h"
28
29#include "hdrl_DER_SNR.h"
30#include "hdrl_utils.h"
31
32#include <math.h>
33
42/*----------------------------------------------------------------------------*/
43
45/*-----------------------------------------------------------------------------
46 Private Functions
47 -----------------------------------------------------------------------------*/
48typedef cpl_error_code (* operate_image_mutable) (hdrl_image * self,
49 const hdrl_image * other);
50
51typedef cpl_error_code (* operate_image_scalar_mutable) (hdrl_image * self,
52 hdrl_value scalar);
53
54static inline cpl_error_code
55operate_spectra_flux_mutate(hdrl_spectrum1D * self,
56 const hdrl_spectrum1D * other,
57 operate_image_mutable func);
58
59static inline hdrl_spectrum1D *
60operate_spectra_flux_create(const hdrl_spectrum1D * self,
61 const hdrl_spectrum1D * other,
62 operate_image_mutable func);
63
64static inline cpl_error_code
65operate_spectra_scalar_flux_mutate(hdrl_spectrum1D * self,
66 hdrl_value scalar,
67 operate_image_scalar_mutable func);
68
69static inline hdrl_spectrum1D *
70operate_spectra_scalar_flux_create(const hdrl_spectrum1D * self,
71 hdrl_value scalar,
72 operate_image_scalar_mutable func);
73
74static inline
75hdrl_spectrum1D * hdrl_spectrum1D_wrap(hdrl_image * arg_flux,
76 cpl_array * wavelength,
77 hdrl_spectrum1D_wave_scale scale);
78
79static inline int compare_double(const void * a, const void * b);
80
81static inline cpl_boolean is_uniformly_sampled(const double * v, cpl_size sz,
82 double * bin);
83
84
85static inline cpl_boolean
86is_wlen_selected(const cpl_bivector * windows, const cpl_boolean is_internal,
87 const hdrl_data_t w);
88
89/*-----------------------------------------------------------------------------
90 Functions
91 -----------------------------------------------------------------------------*/
92
93/* ---------------------------------------------------------------------------*/
116/* ---------------------------------------------------------------------------*/
117hdrl_spectrum1D * hdrl_spectrum1D_create(const cpl_image * arg_flux,
118 const cpl_image * arg_flux_e,
119 const cpl_array * wavelength,
120 hdrl_spectrum1D_wave_scale wave_scale){
121
122 cpl_ensure(arg_flux != NULL && wavelength != NULL && arg_flux_e != NULL,
123 CPL_ERROR_NULL_INPUT, NULL);
124
125 cpl_ensure(cpl_image_get_size_y(arg_flux) == 1
126 && cpl_image_get_size_y(arg_flux_e) == 1,
127 CPL_ERROR_INCOMPATIBLE_INPUT, NULL);
128
129 cpl_ensure(cpl_image_get_size_x(arg_flux) == cpl_array_get_size(wavelength)
130 && cpl_image_get_size_x(arg_flux_e) == cpl_array_get_size(wavelength),
131 CPL_ERROR_INCOMPATIBLE_INPUT, NULL);
132
133 cpl_image * flux_error = cpl_image_cast(arg_flux_e, HDRL_TYPE_ERROR);
134 cpl_image * flux = cpl_image_cast(arg_flux, HDRL_TYPE_DATA);
135
136 hdrl_image * flux_img = hdrl_image_wrap(flux, flux_error, NULL, CPL_TRUE);
137 cpl_array * lambda = cpl_array_cast(wavelength, HDRL_TYPE_DATA);
138
139 return hdrl_spectrum1D_wrap(flux_img, lambda, wave_scale);
140}
141
142/* ---------------------------------------------------------------------------*/
159/* ---------------------------------------------------------------------------*/
161 calculate_analytic_spectrum_point func,
162 const cpl_array * wavelength,
163 hdrl_spectrum1D_wave_scale scale){
164
165 cpl_ensure(wavelength != NULL, CPL_ERROR_NULL_INPUT, NULL);
166 cpl_ensure(func != NULL, CPL_ERROR_NULL_INPUT, NULL);
167
168 cpl_size sx = cpl_array_get_size(wavelength);
169 cpl_image * flux = cpl_image_new(sx, 1, HDRL_TYPE_DATA);
170 cpl_image * flux_e = cpl_image_new(sx, 1, HDRL_TYPE_ERROR);
171
172 for(cpl_size i = 0; i < sx; ++i){
173 hdrl_data_t lambda = cpl_array_get(wavelength, i, NULL);
174 hdrl_value v = func(lambda);
175 cpl_image_set(flux, i + 1, 1, v.data);
176 cpl_image_set(flux_e, i + 1, 1, v.error);
177 }
178
179 hdrl_spectrum1D * to_ret =
180 hdrl_spectrum1D_create(flux, flux_e, wavelength, scale);
181
182 cpl_image_delete(flux);
183 cpl_image_delete(flux_e);
184
185 return to_ret;
186}
187
188/* ---------------------------------------------------------------------------*/
203/* ---------------------------------------------------------------------------*/
205 (const cpl_image * arg_flux,
206 const cpl_array * wavelength,
207 hdrl_spectrum1D_wave_scale scale){
208
209 cpl_ensure(arg_flux != NULL, CPL_ERROR_NULL_INPUT, NULL);
210
211 cpl_size sx = cpl_image_get_size_x(arg_flux);
212 cpl_size sy = cpl_image_get_size_y(arg_flux);
213
214 cpl_ensure(sy == 1 && sx > 0, CPL_ERROR_INCOMPATIBLE_INPUT, NULL);
215
216 cpl_image * zero_errors = cpl_image_new(sx, sy, HDRL_TYPE_ERROR);
217 cpl_image_fill_window(zero_errors, 1, 1, sx, sy, 0.0);
218
219 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_create
220 (arg_flux, zero_errors, wavelength, scale);
221
222 cpl_image_delete(zero_errors);
223
224 return to_ret;
225
226}
227
228/* ---------------------------------------------------------------------------*/
249/* ---------------------------------------------------------------------------*/
251 (const cpl_image * arg_flux,
252 cpl_size half_window,
253 const cpl_array * wavelength,
254 hdrl_spectrum1D_wave_scale scale){
255 cpl_ensure(arg_flux != NULL, CPL_ERROR_NULL_INPUT, NULL);
256 cpl_ensure(wavelength != NULL, CPL_ERROR_NULL_INPUT, NULL);
257
258 cpl_size sx = cpl_image_get_size_x(arg_flux);
259 cpl_size sy = cpl_image_get_size_y(arg_flux);
260
261 cpl_ensure(sy == 1 && sx > 0, CPL_ERROR_INCOMPATIBLE_INPUT, NULL);
262
263 cpl_image * flux = cpl_image_cast(arg_flux, HDRL_TYPE_DATA);
264
265 const hdrl_data_t * flux_data =
266 (const hdrl_data_t*)cpl_image_get_data_const(flux);
267 const cpl_mask * mask_in = cpl_image_get_bpm_const(flux);
268 const cpl_binary * msk_in =
269 mask_in == NULL ? NULL : cpl_mask_get_data_const(mask_in);
270
271 cpl_image * DER_SNR_errors = estimate_noise_DER_SNR
272 (flux_data, msk_in, wavelength, sx, half_window);
273
274 if(!DER_SNR_errors || cpl_error_get_code() != CPL_ERROR_NONE){
275 cpl_image_delete(flux);
276 cpl_image_delete(DER_SNR_errors);
277 return NULL;
278 }
279
280 /* DER_SNR_errors might contain more bad pixels than flux, e.g. 1 good pixel
281 * in position i that is surrounded by bad pixels
282 */
283 cpl_mask * msk = cpl_image_unset_bpm(DER_SNR_errors);
284 cpl_mask_delete(cpl_image_set_bpm(flux, msk));
285
286 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_create
287 (flux , DER_SNR_errors, wavelength, scale);
288
289 cpl_image_delete(DER_SNR_errors);
290 cpl_image_delete(flux);
291
292 return to_ret;
293}
294
295/* ---------------------------------------------------------------------------*/
301/* ---------------------------------------------------------------------------*/
302hdrl_spectrum1D * hdrl_spectrum1D_duplicate(const hdrl_spectrum1D* self)
303{
304 if(!self) return NULL;
305
306 hdrl_image * flux = hdrl_image_duplicate(self->flux);
307 cpl_array * lambdas = cpl_array_duplicate(self->wavelength);
308
309
310 return hdrl_spectrum1D_wrap(flux, lambdas, self->wave_scale);
311}
312
313
314/* ---------------------------------------------------------------------------*/
322/* ---------------------------------------------------------------------------*/
323void hdrl_spectrum1D_delete(hdrl_spectrum1D ** p_self){
324
325 if(!p_self) return;
326
327 hdrl_spectrum1D * spectrum = *p_self;
328
329 if(!spectrum) return;
330
331 cpl_array_delete(spectrum->wavelength);
332 hdrl_image_delete(spectrum->flux);
333 cpl_free(spectrum);
334
335 *p_self = NULL;
336}
337
343/* ---------------------------------------------------------------------------*/
344cpl_size hdrl_spectrum1D_get_size(const hdrl_spectrum1D * self){
345
346 if(!self) return 0;
347
348 return cpl_array_get_size(self->wavelength);
349}
350
356/* ---------------------------------------------------------------------------*/
357const hdrl_image*
358hdrl_spectrum1D_get_flux(const hdrl_spectrum1D * self){
359
360 if(!self) return NULL;
361
362 return self->flux;
363}
364
371/* ---------------------------------------------------------------------------*/
372hdrl_spectrum1D_wavelength
373hdrl_spectrum1D_get_wavelength(const hdrl_spectrum1D * self){
374
375 hdrl_spectrum1D_wavelength to_ret = {NULL, NULL,
376 hdrl_spectrum1D_wave_scale_linear};
377
378 cpl_ensure(self != NULL,
379 CPL_ERROR_NULL_INPUT, to_ret);
380
381 const hdrl_image * flux = hdrl_spectrum1D_get_flux(self);
382
383 const cpl_array * lambdas = self->wavelength;
384
385 const cpl_mask *
386 mask = cpl_image_get_bpm_const(hdrl_image_get_image_const(flux));
387
388 to_ret.bpm = mask;
389 to_ret.wavelength = lambdas;
390 to_ret.scale = hdrl_spectrum1D_get_scale(self);
391 return to_ret;
392}
393
399/* ---------------------------------------------------------------------------*/
400hdrl_spectrum1D_wave_scale
401hdrl_spectrum1D_get_scale(const hdrl_spectrum1D * self){
402 cpl_ensure(self != NULL,
403 CPL_ERROR_NULL_INPUT, hdrl_spectrum1D_wave_scale_linear);
404
405 return self->wave_scale;
406}
407
415/* ---------------------------------------------------------------------------*/
416hdrl_value
417hdrl_spectrum1D_get_flux_value(const hdrl_spectrum1D * self, int idx, int * rej){
418
419 cpl_ensure(self != NULL, CPL_ERROR_NULL_INPUT, ((hdrl_value){0.0, 0.0}));
420
421 const hdrl_image * flx = hdrl_spectrum1D_get_flux(self);
422
423 return hdrl_image_get_pixel(flx, idx + 1, 1, rej);
424}
425
433/* ---------------------------------------------------------------------------*/
434hdrl_data_t
435hdrl_spectrum1D_get_wavelength_value(const hdrl_spectrum1D * self, int idx, int * rej){
436
437 cpl_ensure(self != NULL, CPL_ERROR_NULL_INPUT, 0.0);
438
439 const hdrl_spectrum1D_wavelength lambdas =
441
442 const cpl_array * wlengs = lambdas.wavelength;
443 const hdrl_data_t to_ret = (hdrl_data_t) cpl_array_get(wlengs, idx, NULL);
444
445 if(rej){
446 const cpl_mask * msk = lambdas.bpm;
447 if(msk)
448 *rej = cpl_mask_get(msk, idx + 1, 1);
449 else
450 *rej = CPL_BINARY_0;
451 }
452
453 return to_ret;
454}
455
465/* ---------------------------------------------------------------------------*/
466hdrl_spectrum1D *
467hdrl_spectrum1D_div_spectrum_create(const hdrl_spectrum1D * num,
468 const hdrl_spectrum1D * den){
469 return operate_spectra_flux_create(num, den, hdrl_image_div_image);
470}
471
481/* ---------------------------------------------------------------------------*/
482hdrl_spectrum1D *
483hdrl_spectrum1D_mul_spectrum_create(const hdrl_spectrum1D * f1,
484 const hdrl_spectrum1D * f2){
485 return operate_spectra_flux_create(f1, f2, hdrl_image_mul_image);
486}
487
497/* ---------------------------------------------------------------------------*/
498hdrl_spectrum1D *
499hdrl_spectrum1D_add_spectrum_create(const hdrl_spectrum1D * f1,
500 const hdrl_spectrum1D * f2){
501 return operate_spectra_flux_create(f1, f2, hdrl_image_add_image);
502}
503
504
514/* ---------------------------------------------------------------------------*/
515hdrl_spectrum1D *
516hdrl_spectrum1D_sub_spectrum_create(const hdrl_spectrum1D * f1,
517 const hdrl_spectrum1D * f2){
518 return operate_spectra_flux_create(f1, f2, hdrl_image_sub_image);
519}
520
521
528/* ---------------------------------------------------------------------------*/
529cpl_error_code hdrl_spectrum1D_div_spectrum(hdrl_spectrum1D * self,
530 const hdrl_spectrum1D * other){
531 return operate_spectra_flux_mutate(self, other, hdrl_image_div_image);
532}
533
540/* ---------------------------------------------------------------------------*/
541cpl_error_code hdrl_spectrum1D_mul_spectrum(hdrl_spectrum1D * self,
542 const hdrl_spectrum1D * other){
543 return operate_spectra_flux_mutate(self, other, hdrl_image_mul_image);
544}
545
552/* ---------------------------------------------------------------------------*/
553cpl_error_code hdrl_spectrum1D_add_spectrum(hdrl_spectrum1D * self,
554 const hdrl_spectrum1D * other){
555 return operate_spectra_flux_mutate(self, other, hdrl_image_add_image);
556}
557
558
565/* ---------------------------------------------------------------------------*/
566cpl_error_code hdrl_spectrum1D_sub_spectrum(hdrl_spectrum1D * self,
567 const hdrl_spectrum1D * other){
568 return operate_spectra_flux_mutate(self, other, hdrl_image_sub_image);
569}
570
580/* ---------------------------------------------------------------------------*/
581hdrl_spectrum1D *
582hdrl_spectrum1D_div_scalar_create(const hdrl_spectrum1D * self,
583 hdrl_value scalar_operator){
584 return
585 operate_spectra_scalar_flux_create(self, scalar_operator, hdrl_image_div_scalar);
586}
587
597/* ---------------------------------------------------------------------------*/
598hdrl_spectrum1D *
599hdrl_spectrum1D_mul_scalar_create(const hdrl_spectrum1D * self,
600 hdrl_value scalar_operator){
601 return
602 operate_spectra_scalar_flux_create(self, scalar_operator, hdrl_image_mul_scalar);
603}
604
614/* ---------------------------------------------------------------------------*/
615hdrl_spectrum1D *
616hdrl_spectrum1D_add_scalar_create(const hdrl_spectrum1D * self,
617 hdrl_value scalar_operator){
618 return
619 operate_spectra_scalar_flux_create(self, scalar_operator, hdrl_image_add_scalar);
620}
621
631/* ---------------------------------------------------------------------------*/
632hdrl_spectrum1D *
633hdrl_spectrum1D_sub_scalar_create(const hdrl_spectrum1D * self,
634 hdrl_value scalar_operator){
635 return
636 operate_spectra_scalar_flux_create(self, scalar_operator, hdrl_image_sub_scalar);
637}
638
648/* ---------------------------------------------------------------------------*/
649hdrl_spectrum1D *
650hdrl_spectrum1D_pow_scalar_create(const hdrl_spectrum1D * self,
651 hdrl_value scalar_operator){
652 return
653 operate_spectra_scalar_flux_create(self, scalar_operator, hdrl_image_pow_scalar);
654}
655
665/* ---------------------------------------------------------------------------*/
666hdrl_spectrum1D *
667hdrl_spectrum1D_exp_scalar_create(const hdrl_spectrum1D * self,
668 hdrl_value scalar_operator){
669 return
670 operate_spectra_scalar_flux_create(self, scalar_operator, hdrl_image_exp_scalar);
671}
672
680/* ---------------------------------------------------------------------------*/
681cpl_error_code hdrl_spectrum1D_div_scalar(hdrl_spectrum1D * self,
682 hdrl_value scalar_operator){
683 return
684 operate_spectra_scalar_flux_mutate(self, scalar_operator, hdrl_image_div_scalar);
685}
686
687
695/* ---------------------------------------------------------------------------*/
696cpl_error_code hdrl_spectrum1D_mul_scalar(hdrl_spectrum1D * self,
697 hdrl_value scalar_operator){
698 return
699 operate_spectra_scalar_flux_mutate(self, scalar_operator, hdrl_image_mul_scalar);
700}
701
709/* ---------------------------------------------------------------------------*/
710cpl_error_code hdrl_spectrum1D_add_scalar(hdrl_spectrum1D * self,
711 hdrl_value scalar_operator){
712 return
713 operate_spectra_scalar_flux_mutate(self, scalar_operator, hdrl_image_add_scalar);
714}
715
723/* ---------------------------------------------------------------------------*/
724cpl_error_code hdrl_spectrum1D_sub_scalar(hdrl_spectrum1D * self,
725 hdrl_value scalar_operator){
726 return
727 operate_spectra_scalar_flux_mutate(self, scalar_operator, hdrl_image_sub_scalar);
728}
729
737/* ---------------------------------------------------------------------------*/
738cpl_error_code hdrl_spectrum1D_pow_scalar(hdrl_spectrum1D * self,
739 hdrl_value scalar_operator){
740 return
741 operate_spectra_scalar_flux_mutate(self, scalar_operator, hdrl_image_pow_scalar);
742}
743
751/* ---------------------------------------------------------------------------*/
752cpl_error_code hdrl_spectrum1D_exp_scalar(hdrl_spectrum1D * self,
753 hdrl_value scalar_operator){
754 return
755 operate_spectra_scalar_flux_mutate(self, scalar_operator, hdrl_image_exp_scalar);
756}
757
767/* ---------------------------------------------------------------------------*/
768cpl_error_code hdrl_spectrum1D_wavelength_mult_scalar_linear(hdrl_spectrum1D * self,
769 hdrl_data_t scale_linear)
770{
771 cpl_ensure_code(scale_linear > 0, CPL_ERROR_INCOMPATIBLE_INPUT);
772
773 if(self == NULL) return CPL_ERROR_NONE;
774
775 if(self->wave_scale == hdrl_spectrum1D_wave_scale_linear)
776 return cpl_array_multiply_scalar(self->wavelength, scale_linear);
777
778 return cpl_array_add_scalar(self->wavelength, log(scale_linear));
779}
780
792/* ---------------------------------------------------------------------------*/
793hdrl_spectrum1D *
795 hdrl_data_t scale_linear)
796{
797 if(!self) return NULL;
798
799 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_duplicate(self);
800
801 cpl_error_code fail =
803 if(fail){
804 hdrl_spectrum1D_delete(&to_ret);
805 }
806 return to_ret;
807}
808
816/* ---------------------------------------------------------------------------*/
817cpl_error_code hdrl_spectrum1D_wavelength_shift(hdrl_spectrum1D * self,
818 hdrl_data_t shift)
819{
820 if(!self) return CPL_ERROR_NONE;
821
822 return cpl_array_add_scalar(self->wavelength, shift);
823}
824
834/* ---------------------------------------------------------------------------*/
835hdrl_spectrum1D *
836hdrl_spectrum1D_wavelength_shift_create(const hdrl_spectrum1D * self,
837 hdrl_data_t shift){
838
839 if(!self) return NULL;
840
841 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_duplicate(self);
842
843 cpl_error_code fail = hdrl_spectrum1D_wavelength_shift(to_ret, shift);
844
845 if(fail){
846 hdrl_spectrum1D_delete(&to_ret);
847 }
848
849 return to_ret;
850}
851
859/* ---------------------------------------------------------------------------*/
860cpl_error_code hdrl_spectrum1D_wavelength_convert_to_linear(hdrl_spectrum1D * self){
861
862 if(self == NULL || self->wave_scale == hdrl_spectrum1D_wave_scale_linear)
863 return CPL_ERROR_NONE;
864
865 cpl_error_code fail = cpl_array_exponential(self->wavelength, CPL_MATH_E);
866 self->wave_scale = hdrl_spectrum1D_wave_scale_linear;
867 return fail;
868}
869
879/* ---------------------------------------------------------------------------*/
880hdrl_spectrum1D *
882
883 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_duplicate(self);
884 cpl_error_code fail = hdrl_spectrum1D_wavelength_convert_to_linear(to_ret);
885
886 if(fail){
887 hdrl_spectrum1D_delete(&to_ret);
888 }
889 return to_ret;
890}
891
898/* ---------------------------------------------------------------------------*/
899cpl_error_code hdrl_spectrum1D_wavelength_convert_to_log(hdrl_spectrum1D * self){
900
901 if(self == NULL || self->wave_scale == hdrl_spectrum1D_wave_scale_log)
902 return CPL_ERROR_NONE;
903
904 cpl_error_code fail = cpl_array_logarithm(self->wavelength, CPL_MATH_E);
905 self->wave_scale = hdrl_spectrum1D_wave_scale_log;
906 return fail;
907}
908
920/* ---------------------------------------------------------------------------*/
921hdrl_spectrum1D *
923
924 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_duplicate(self);
925 cpl_error_code fail = hdrl_spectrum1D_wavelength_convert_to_log(to_ret);
926
927 if(fail){
928 hdrl_spectrum1D_delete(&to_ret);
929 }
930 return to_ret;
931}
932
953/* ---------------------------------------------------------------------------*/
954hdrl_spectrum1D *
955hdrl_spectrum1D_select_wavelengths(const hdrl_spectrum1D * self,
956 const cpl_bivector * windows, const cpl_boolean is_internal){
957
958 cpl_ensure(self != NULL, CPL_ERROR_NULL_INPUT, NULL);
959 cpl_ensure(windows != NULL, CPL_ERROR_NULL_INPUT, NULL);
960
961 const cpl_size sz = hdrl_spectrum1D_get_size(self);
962 cpl_size num_selected = 0;
963
964 for(cpl_size i = 0; i < sz; ++i){
965 const hdrl_data_t w = hdrl_spectrum1D_get_wavelength_value(self, i, NULL);
966 if(is_wlen_selected(windows, is_internal, w))
967 num_selected++;
968 }
969
970 if(num_selected == sz)
971 return hdrl_spectrum1D_duplicate(self);
972
973 cpl_ensure(num_selected > 0, CPL_ERROR_ILLEGAL_OUTPUT, NULL);
974
975 cpl_image * flux = cpl_image_new(num_selected, 1, HDRL_TYPE_DATA);
976 cpl_image * flux_e = cpl_image_new(num_selected, 1, HDRL_TYPE_ERROR);
977 cpl_array * wavs = cpl_array_new(num_selected, HDRL_TYPE_DATA);
978
979 cpl_size idx_this = 0;
980 for(cpl_size i = 0; i < sz; ++i){
981
982 int rej = 0;
983 const hdrl_data_t w = hdrl_spectrum1D_get_wavelength_value(self, i, NULL);
984 if(!is_wlen_selected(windows, is_internal, w)) continue;
985
986 hdrl_value v = hdrl_spectrum1D_get_flux_value(self, i, &rej);
987
988 if(!rej){
989 cpl_image_set(flux, idx_this + 1, 1, v.data);
990 cpl_image_set(flux_e, idx_this + 1, 1, v.error);
991 }
992 else{
993 cpl_image_reject(flux, idx_this + 1, 1);
994 cpl_image_reject(flux_e, idx_this + 1, 1);
995 }
996
997 cpl_array_set(wavs, idx_this, w);
998
999 idx_this++;
1000 }
1001
1002 const hdrl_spectrum1D_wave_scale scale = hdrl_spectrum1D_get_scale(self);
1003 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_create(flux,flux_e, wavs, scale);
1004
1005 cpl_image_delete(flux);
1006 cpl_image_delete(flux_e);
1007 cpl_array_delete(wavs);
1008
1009 return to_ret;
1010}
1011
1025/* ---------------------------------------------------------------------------*/
1026hdrl_spectrum1D *
1027hdrl_spectrum1D_reject_pixels(const hdrl_spectrum1D * self,
1028 const cpl_array * bad_samples){
1029
1030 const cpl_size sz = cpl_array_get_size(bad_samples);
1031
1032 cpl_ensure(self != NULL, CPL_ERROR_NULL_INPUT, NULL);
1033 cpl_ensure(bad_samples != NULL, CPL_ERROR_NULL_INPUT, NULL);
1034 cpl_ensure(sz == hdrl_spectrum1D_get_size(self), CPL_ERROR_ILLEGAL_INPUT, NULL);
1035
1036 hdrl_image * flux = hdrl_image_duplicate(hdrl_spectrum1D_get_flux(self));
1037
1038 for(cpl_size i = 0; i < sz; ++i){
1039
1040 if(cpl_array_get_int(bad_samples, i, NULL))
1041 hdrl_image_reject(flux, i + 1, 1);
1042 }
1043
1044 const cpl_image * flx_cpl = hdrl_image_get_image(flux);
1045 const cpl_image * flx_e_cpl = hdrl_image_get_error(flux);
1046 hdrl_spectrum1D_wavelength wlen = hdrl_spectrum1D_get_wavelength(self);
1047
1048 hdrl_spectrum1D * to_ret =
1049 hdrl_spectrum1D_create(flx_cpl, flx_e_cpl,
1050 wlen.wavelength, wlen.scale);
1051
1052 hdrl_image_delete(flux);
1053 return to_ret;
1054}
1055
1073/* ---------------------------------------------------------------------------*/
1075 const hdrl_spectrum1D * self, const char * flux_col_name,
1076 const char* wavelength_col_name, const char * flux_e_col_name,
1077 const char * flux_bpm_col_name){
1078
1079 cpl_ensure(self != NULL, CPL_ERROR_NULL_INPUT, NULL);
1080 cpl_ensure(flux_col_name != NULL || wavelength_col_name != NULL,
1081 CPL_ERROR_NULL_INPUT, NULL);
1082
1083 cpl_size sz = hdrl_spectrum1D_get_size(self);
1084
1085 cpl_table * tb = cpl_table_new(sz);
1086 cpl_ensure(tb != NULL, CPL_ERROR_NULL_INPUT, NULL);
1087
1088 cpl_error_code fail = hdrl_spectrum1D_append_to_table
1089 (self, tb, flux_col_name, wavelength_col_name, flux_e_col_name,
1090 flux_bpm_col_name);
1091
1092 if(fail){
1093 cpl_table_delete(tb);
1094 tb = NULL;
1095 }
1096 return tb;
1097}
1098
1114/* ---------------------------------------------------------------------------*/
1116(const hdrl_spectrum1D * self, cpl_table * dest,
1117const char * flux_col_name, const char* wavelength_col_name,
1118const char * flux_e_col_name, const char * flux_bpm_col_name){
1119
1120 cpl_ensure_code(self != NULL && dest != NULL,CPL_ERROR_NULL_INPUT);
1121 cpl_ensure_code(flux_col_name != NULL || wavelength_col_name != NULL,
1122 CPL_ERROR_NULL_INPUT);
1123
1124 cpl_size cl_sz = cpl_table_get_nrow(dest);
1125 cpl_size sz = hdrl_spectrum1D_get_size(self);
1126
1127 cpl_ensure_code(sz == cl_sz, CPL_ERROR_INCOMPATIBLE_INPUT);
1128
1129 if(wavelength_col_name){
1130 double * lambdas = cpl_calloc(sz ,sizeof(double));
1131 for(cpl_size i = 0; i < sz; i++){
1132 lambdas[i] = hdrl_spectrum1D_get_wavelength_value(self, i, NULL);
1133 }
1134 cpl_error_code fail =
1135 cpl_table_wrap_double(dest, lambdas, wavelength_col_name);
1136
1137 if(fail){
1138 cpl_free(lambdas);
1139 return fail;
1140 }
1141 }
1142
1143 if(flux_col_name){
1144 double * flux = cpl_calloc(sz ,sizeof(double));
1145 for(cpl_size i = 0; i < sz; i++){
1146 flux[i] = hdrl_spectrum1D_get_flux_value(self, i, NULL).data;
1147 }
1148
1149 cpl_error_code fail = cpl_table_wrap_double(dest, flux, flux_col_name);
1150 if(fail){
1151 cpl_free(flux);
1152 return fail;
1153 }
1154 }
1155
1156 if(flux_e_col_name){
1157 double * e_flux = cpl_calloc(sz ,sizeof(double));
1158 for(cpl_size i = 0; i < sz; i++){
1159 e_flux[i] = hdrl_spectrum1D_get_flux_value(self, i, NULL).error;
1160 }
1161 cpl_error_code fail =
1162 cpl_table_wrap_double(dest, e_flux, flux_e_col_name);
1163
1164 if(fail){
1165 cpl_free(e_flux);
1166 return fail;
1167 }
1168 }
1169
1170 if(flux_bpm_col_name){
1171 int * bpm_flux = cpl_calloc(sz ,sizeof(int));
1172 for(cpl_size i = 0; i < sz; i++){
1173 hdrl_spectrum1D_get_flux_value(self, i, &bpm_flux[i]);
1174 }
1175 cpl_error_code fail =
1176 cpl_table_wrap_int(dest, bpm_flux, flux_bpm_col_name);
1177 if(fail){
1178 cpl_free(bpm_flux);
1179 return fail;
1180 }
1181 }
1182
1183 return CPL_ERROR_NONE;
1184}
1185
1199/* ---------------------------------------------------------------------------*/
1201(const cpl_table * self, const char * flux_col_name,
1202const char* wavelength_col_name, const char * flux_e_col_name,
1203const char * flux_bpm_col_name, hdrl_spectrum1D_wave_scale scale){
1204
1205 cpl_ensure(self != NULL, CPL_ERROR_NULL_INPUT, NULL);
1206 cpl_ensure(flux_col_name != NULL, CPL_ERROR_NULL_INPUT, NULL);
1207 cpl_ensure(wavelength_col_name != NULL, CPL_ERROR_NULL_INPUT, NULL);
1208
1209 cpl_size sz = cpl_table_get_nrow(self);
1210
1211 cpl_ensure(sz > 1, CPL_ERROR_ILLEGAL_INPUT, NULL);
1212
1213 cpl_image * flux = cpl_image_new(sz, 1, HDRL_TYPE_DATA);
1214 cpl_image * flux_e = cpl_image_new(sz, 1, HDRL_TYPE_ERROR);
1215 cpl_array * lambdas = cpl_array_new(sz, HDRL_TYPE_DATA);
1216
1217
1218 for(cpl_size i = 0; i < sz; ++i){
1219 const double fx = cpl_table_get(self, flux_col_name, i, NULL);
1220 const double l = cpl_table_get(self, wavelength_col_name, i, NULL);
1221
1222 double fx_e = 0;
1223 if(flux_e_col_name)
1224 fx_e = cpl_table_get(self, flux_e_col_name, i, NULL);
1225
1226 int rej = 0;
1227 if(flux_bpm_col_name)
1228 rej = cpl_table_get_int(self, flux_bpm_col_name, i, NULL);
1229
1230 cpl_image_set(flux, i + 1, 1, fx);
1231
1232 if(rej)
1233 cpl_image_reject(flux, i + 1, 1);
1234
1235 cpl_image_set(flux_e, i + 1, 1, fx_e);
1236 cpl_array_set(lambdas, i, l);
1237 }
1238
1239 hdrl_spectrum1D * sp = hdrl_spectrum1D_create(flux, flux_e, lambdas, scale);
1240
1241 cpl_image_delete(flux);
1242 cpl_image_delete(flux_e);
1243 cpl_array_delete(lambdas);
1244 return sp;
1245}
1246
1247void hdrl_spectrum1D_save(const hdrl_spectrum1D * s, const char * fname){
1248 if(s == NULL) return;
1249 cpl_table * tb = hdrl_spectrum1D_convert_to_table(s, "FLX", "WLN", "FLX_E",
1250 "FLX_BPM");
1251 cpl_table_save(tb, NULL, NULL, fname, CPL_IO_CREATE);
1252 cpl_table_delete(tb);
1253}
1254
1261/* ---------------------------------------------------------------------------*/
1262cpl_boolean
1264 const cpl_array * w2){
1265 if(w1 == NULL && w2 == NULL) return CPL_TRUE;
1266
1267 if(w1 == NULL) return CPL_FALSE;
1268
1269 if(w2 == NULL) return CPL_FALSE;
1270
1271 cpl_size sz = cpl_array_get_size(w1);
1272 if(sz != cpl_array_get_size(w2)) return CPL_FALSE;
1273 for(cpl_size i = 0; i < sz; i++){
1274 const double wa = cpl_array_get(w1, i, NULL);
1275 const double wb = cpl_array_get(w2, i, NULL);
1276
1277 const double d = wa - wb;
1278 if(fabs(d) > 1e-10 * CPL_MIN(wa, wb)){
1279 return CPL_FALSE;
1280 }
1281 }
1282
1283 return CPL_TRUE;
1284}
1285
1292/* ---------------------------------------------------------------------------*/
1293cpl_boolean
1294hdrl_spectrum1D_are_spectra_compatible(const hdrl_spectrum1D_wavelength* s1,
1295 const hdrl_spectrum1D_wavelength* s2){
1296
1297 if(s1 == NULL && s2 == NULL) return CPL_TRUE;
1298
1299 if(s1 == NULL) return CPL_FALSE;
1300
1301 if(s2 == NULL) return CPL_FALSE;
1302
1303 if(s1->scale != s2->scale) return CPL_FALSE;
1304
1305 return hdrl_spectrum1D_are_wavelengths_compatible(s1->wavelength,
1306 s2->wavelength);
1307}
1308
1317/* ---------------------------------------------------------------------------*/
1318cpl_boolean hdrl_spectrum1D_is_uniformly_sampled(const hdrl_spectrum1D * self,
1319 double * bin){
1320
1321 *bin = 0.0;
1322
1323 if(self == NULL) return CPL_FALSE;
1324
1325 const cpl_size sz = hdrl_spectrum1D_get_size(self);
1326
1327 if(sz <= 2) return CPL_TRUE;
1328
1329 double * vd = cpl_calloc(sz, sizeof(double));
1330
1331 for(cpl_size i = 0; i < sz; ++i){
1332 vd[i] = hdrl_spectrum1D_get_wavelength_value(self, i, NULL);
1333 }
1334
1335 qsort(vd, sz, sizeof(double), compare_double);
1336
1337 cpl_boolean to_ret = is_uniformly_sampled(vd, sz, bin);
1338
1339 cpl_free(vd);
1340
1341 return to_ret;
1342
1343}
1344
1345/*-----------------------------------------------------------------------------
1346 Private Functions Implementation
1347 -----------------------------------------------------------------------------*/
1357/* ---------------------------------------------------------------------------*/
1358static inline cpl_error_code operate_spectra_flux_mutate(hdrl_spectrum1D * self,
1359 const hdrl_spectrum1D * other,
1360 operate_image_mutable func){
1361 cpl_ensure_code(self != NULL && other != NULL,
1362 CPL_ERROR_NULL_INPUT);
1363
1364 hdrl_spectrum1D_wavelength w_self = hdrl_spectrum1D_get_wavelength(self);
1365 hdrl_spectrum1D_wavelength w_other = hdrl_spectrum1D_get_wavelength(other);
1366 cpl_ensure_code(hdrl_spectrum1D_are_spectra_compatible(&w_self, &w_other),
1367 CPL_ERROR_INCOMPATIBLE_INPUT);
1368
1369 hdrl_image * f1 = self->flux;
1370 const hdrl_image * f2 = other->flux;
1371
1372 cpl_ensure_code(f1 != NULL && f2 != NULL,
1373 CPL_ERROR_NULL_INPUT);
1374
1375 func(f1, f2);
1376
1377 return CPL_ERROR_NONE;
1378}
1379
1389/* ---------------------------------------------------------------------------*/
1390static inline hdrl_spectrum1D *
1391operate_spectra_flux_create(const hdrl_spectrum1D * self,
1392 const hdrl_spectrum1D * other,
1393 operate_image_mutable func){
1394
1395 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_duplicate(self);
1396
1397 cpl_error_code fail = operate_spectra_flux_mutate(to_ret, other, func);
1398
1399 if(fail){
1400 hdrl_spectrum1D_delete(&to_ret);
1401 }
1402 return to_ret;
1403}
1404
1414/* ---------------------------------------------------------------------------*/
1415static inline cpl_error_code
1416operate_spectra_scalar_flux_mutate(hdrl_spectrum1D * self,
1417 hdrl_value scalar,
1418 operate_image_scalar_mutable func){
1419
1420 if(!self) return CPL_ERROR_NONE;
1421
1422 hdrl_image * f1 = self->flux;
1423
1424 cpl_ensure_code(f1 != NULL, CPL_ERROR_NULL_INPUT);
1425
1426 func(f1, scalar);
1427
1428 return CPL_ERROR_NONE;
1429}
1430
1440/* ---------------------------------------------------------------------------*/
1441static inline hdrl_spectrum1D *
1442operate_spectra_scalar_flux_create(const hdrl_spectrum1D * self,
1443 hdrl_value scalar,
1444 operate_image_scalar_mutable func){
1445
1446 hdrl_spectrum1D * to_ret = hdrl_spectrum1D_duplicate(self);
1447
1448 cpl_error_code fail =
1449 operate_spectra_scalar_flux_mutate(to_ret, scalar, func);
1450
1451 if(fail){
1452 hdrl_spectrum1D_delete(&to_ret);
1453 }
1454
1455 return to_ret;
1456}
1457
1458static inline
1459hdrl_spectrum1D * hdrl_spectrum1D_wrap(hdrl_image * arg_flux,
1460 cpl_array * wavelength,
1461 hdrl_spectrum1D_wave_scale scale){
1462
1463 hdrl_spectrum1D * to_ret = cpl_calloc(1, sizeof(*to_ret));
1464
1465 to_ret->flux = arg_flux;
1466 to_ret->wavelength = wavelength;
1467 to_ret->wave_scale = scale;
1468 return to_ret;
1469}
1470
1471static inline
1472int compare_double(const void * a, const void * b){
1473 const double ad = *((const double *)a);
1474 const double bd = *((const double *)b);
1475 const double delta = ad - bd;
1476 if(delta > 0.0) return 1;
1477 if(delta < 0.0) return -1;
1478 return 0;
1479}
1480
1481#if HDRL_SIZEOF_DATA == 4
1482 const double wave_delta = 1.e-5;
1483#else
1484 const double wave_delta = 1.e-6;
1485#endif
1486
1487static inline
1488cpl_boolean is_uniformly_sampled(const double * v, cpl_size sz, double * bin){
1489
1490 const double d = v[1] - v[0];
1491 *bin = d;
1492 for(cpl_size i = 1; i < sz - 1; ++i){
1493 const double d2 = (v[i + 1] - v[i]);
1494 const double eps = (fabs(d2 - d) / d);
1495 if(eps > wave_delta){
1496 return CPL_FALSE;
1497 }
1498 }
1499
1500 return CPL_TRUE;
1501}
1502
1503
1504static inline cpl_boolean
1505is_contained_in_at_least_one_window(const cpl_bivector * windows, const hdrl_data_t w){
1506 const cpl_size sz = cpl_bivector_get_size(windows);
1507
1508 for(cpl_size i = 0; i < sz; ++i){
1509 const double wmin = cpl_vector_get(cpl_bivector_get_x_const(windows), i);
1510 const double wmax = cpl_vector_get(cpl_bivector_get_y_const(windows), i);
1511
1512 if(w >= wmin && w <= wmax) return CPL_TRUE;
1513 }
1514
1515 return CPL_FALSE;
1516}
1517
1518static inline cpl_boolean
1519is_wlen_selected(const cpl_bivector * windows, const cpl_boolean is_internal,
1520 const hdrl_data_t w){
1521
1522 if(is_internal) return is_contained_in_at_least_one_window(windows, w);
1523
1524 return !is_contained_in_at_least_one_window(windows, w);
1525
1526}
1527
hdrl_value hdrl_image_get_pixel(const hdrl_image *self, cpl_size xpos, cpl_size ypos, int *pis_rejected)
get pixel values of hdrl_image
Definition: hdrl_image.c:559
cpl_error_code hdrl_image_sub_image(hdrl_image *self, const hdrl_image *other)
Subtract two images, store the result in the first image.
cpl_error_code hdrl_image_div_scalar(hdrl_image *self, hdrl_value value)
Elementwise division of an image with a scalar.
cpl_error_code hdrl_image_pow_scalar(hdrl_image *self, const hdrl_value exponent)
computes the power of an image by a scalar
cpl_error_code hdrl_image_add_image(hdrl_image *self, const hdrl_image *other)
Add two images, store the result in the first image.
cpl_error_code hdrl_image_div_image(hdrl_image *self, const hdrl_image *other)
Divide two images, store the result in the first image.
cpl_error_code hdrl_image_mul_scalar(hdrl_image *self, hdrl_value value)
Elementwise multiplication of an image with a scalar.
hdrl_image * hdrl_image_duplicate(const hdrl_image *himg)
copy hdrl_image
Definition: hdrl_image.c:391
cpl_error_code hdrl_image_add_scalar(hdrl_image *self, hdrl_value value)
Elementwise addition of a scalar to an image.
cpl_error_code hdrl_image_mul_image(hdrl_image *self, const hdrl_image *other)
Multiply two images, store the result in the first image.
cpl_error_code hdrl_image_exp_scalar(hdrl_image *self, const hdrl_value base)
computes the exponential of an image by a scalar
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
const cpl_image * hdrl_image_get_image_const(const hdrl_image *himg)
get data as cpl image
Definition: hdrl_image.c:118
void hdrl_image_delete(hdrl_image *himg)
delete hdrl_image
Definition: hdrl_image.c:379
cpl_error_code hdrl_image_sub_scalar(hdrl_image *self, hdrl_value value)
Elementwise subtraction of a scalar from an image.
hdrl_spectrum1D * hdrl_spectrum1D_pow_scalar_create(const hdrl_spectrum1D *self, hdrl_value scalar_operator)
subtract a scalar from a spectrum
const hdrl_image * hdrl_spectrum1D_get_flux(const hdrl_spectrum1D *self)
hdrl_spectrum1D getter flux
cpl_error_code hdrl_spectrum1D_mul_scalar(hdrl_spectrum1D *self, hdrl_value scalar_operator)
computes the elementwise multiplication of a spectrum by a scalar, the self parameter is modified
cpl_error_code hdrl_spectrum1D_wavelength_mult_scalar_linear(hdrl_spectrum1D *self, hdrl_data_t scale_linear)
computes the elementwise multiplication of the scalar for the wavelength. The scalar is assumed to be...
cpl_error_code hdrl_spectrum1D_wavelength_convert_to_linear(hdrl_spectrum1D *self)
converts the wavelength scale to linear.
hdrl_spectrum1D * hdrl_spectrum1D_sub_spectrum_create(const hdrl_spectrum1D *f1, const hdrl_spectrum1D *f2)
subtract two spectra
cpl_error_code hdrl_spectrum1D_wavelength_shift(hdrl_spectrum1D *self, hdrl_data_t shift)
computes the elementwise shift of the wavelength by the shift parameter. The self parameter is modifi...
cpl_error_code hdrl_spectrum1D_div_spectrum(hdrl_spectrum1D *self, const hdrl_spectrum1D *other)
divide one spectrum by another spectrum
cpl_error_code hdrl_spectrum1D_pow_scalar(hdrl_spectrum1D *self, hdrl_value scalar_operator)
computes the elementwise power of the flux to the scalar, the self parameter is modified
hdrl_spectrum1D * hdrl_spectrum1D_mul_spectrum_create(const hdrl_spectrum1D *f1, const hdrl_spectrum1D *f2)
multiply one spectrum by another spectrum
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
void hdrl_spectrum1D_delete(hdrl_spectrum1D **p_self)
hdrl_spectrum1D destructor
cpl_error_code hdrl_spectrum1D_mul_spectrum(hdrl_spectrum1D *self, const hdrl_spectrum1D *other)
multiply one spectrum by another spectrum
cpl_error_code hdrl_spectrum1D_wavelength_convert_to_log(hdrl_spectrum1D *self)
converts the wavelength scale to log. If the spectrum is already in log scale nothing is done.
hdrl_spectrum1D * hdrl_spectrum1D_div_spectrum_create(const hdrl_spectrum1D *num, const hdrl_spectrum1D *den)
divide one spectrum by another spectrum
cpl_error_code hdrl_spectrum1D_exp_scalar(hdrl_spectrum1D *self, hdrl_value scalar_operator)
computes the elementwise power of the scalar to the flux, the self parameter is modified
hdrl_spectrum1D * hdrl_spectrum1D_mul_scalar_create(const hdrl_spectrum1D *self, hdrl_value scalar_operator)
multiply a spectrum by a scalar
cpl_error_code hdrl_spectrum1D_sub_spectrum(hdrl_spectrum1D *self, const hdrl_spectrum1D *other)
subtract two spectra
hdrl_spectrum1D * hdrl_spectrum1D_wavelength_shift_create(const hdrl_spectrum1D *self, hdrl_data_t shift)
computes the elementwise shift of the wavelength by the shift parameter.
hdrl_spectrum1D * hdrl_spectrum1D_div_scalar_create(const hdrl_spectrum1D *self, hdrl_value scalar_operator)
divide a spectrum by a scalar
hdrl_spectrum1D * hdrl_spectrum1D_create_error_DER_SNR(const cpl_image *arg_flux, cpl_size half_window, const cpl_array *wavelength, hdrl_spectrum1D_wave_scale scale)
hdrl_spectrum1D constructor when no error information is available, in this case we use DER_SNR to es...
hdrl_spectrum1D * hdrl_spectrum1D_convert_from_table(const cpl_table *self, const char *flux_col_name, const char *wavelength_col_name, const char *flux_e_col_name, const char *flux_bpm_col_name, hdrl_spectrum1D_wave_scale scale)
convert a table to a spectrum
hdrl_spectrum1D * hdrl_spectrum1D_wavelength_convert_to_log_create(const hdrl_spectrum1D *self)
converts the wavelength scale to log. It returns a modified version of self. self is not modified....
cpl_boolean hdrl_spectrum1D_are_wavelengths_compatible(const cpl_array *w1, const cpl_array *w2)
checks if two wavelengths array are defined on the same wavelengths.
cpl_error_code hdrl_spectrum1D_sub_scalar(hdrl_spectrum1D *self, hdrl_value scalar_operator)
computes the elementwise subtraction of a spectrum by a scalar, the self parameter is modified
cpl_boolean hdrl_spectrum1D_is_uniformly_sampled(const hdrl_spectrum1D *self, double *bin)
checks if the spectrum is defined on uniformly sampled wavelengths.
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
cpl_error_code hdrl_spectrum1D_add_scalar(hdrl_spectrum1D *self, hdrl_value scalar_operator)
computes the elementwise addition of a spectrum by a scalar, the self parameter is modified
cpl_error_code hdrl_spectrum1D_add_spectrum(hdrl_spectrum1D *self, const hdrl_spectrum1D *other)
sum two spectra
hdrl_spectrum1D * hdrl_spectrum1D_exp_scalar_create(const hdrl_spectrum1D *self, hdrl_value scalar_operator)
subtract a scalar from a spectrum
hdrl_spectrum1D * hdrl_spectrum1D_sub_scalar_create(const hdrl_spectrum1D *self, hdrl_value scalar_operator)
subtract a scalar from a spectrum
hdrl_spectrum1D * hdrl_spectrum1D_wavelength_mult_scalar_linear_create(const hdrl_spectrum1D *self, hdrl_data_t scale_linear)
computes the elementwise multiplication of the scalar for the wavelength. The scalar is assumed to be...
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_spectrum1D * hdrl_spectrum1D_add_scalar_create(const hdrl_spectrum1D *self, hdrl_value scalar_operator)
add a scalar to a spectrum
hdrl_data_t hdrl_spectrum1D_get_wavelength_value(const hdrl_spectrum1D *self, int idx, int *rej)
hdrl_spectrum1D getter for a wavelength value
hdrl_spectrum1D * hdrl_spectrum1D_create_analytic(calculate_analytic_spectrum_point func, const cpl_array *wavelength, hdrl_spectrum1D_wave_scale scale)
hdrl_spectrum1D constructor in the case of a spectrum defined by an analytical function
cpl_error_code hdrl_spectrum1D_append_to_table(const hdrl_spectrum1D *self, cpl_table *dest, const char *flux_col_name, const char *wavelength_col_name, const char *flux_e_col_name, const char *flux_bpm_col_name)
append a spectrum to a table.
hdrl_spectrum1D * hdrl_spectrum1D_wavelength_convert_to_linear_create(const hdrl_spectrum1D *self)
converts the wavelength scale to linear.
cpl_table * hdrl_spectrum1D_convert_to_table(const hdrl_spectrum1D *self, const char *flux_col_name, const char *wavelength_col_name, const char *flux_e_col_name, const char *flux_bpm_col_name)
converts a spectrum in a table.
hdrl_spectrum1D * hdrl_spectrum1D_select_wavelengths(const hdrl_spectrum1D *self, const cpl_bivector *windows, const cpl_boolean is_internal)
the function selects or discards flux values according to whether the value of the corresponding wave...
hdrl_spectrum1D * hdrl_spectrum1D_add_spectrum_create(const hdrl_spectrum1D *f1, const hdrl_spectrum1D *f2)
sum one spectrum to another spectrum
cpl_error_code hdrl_spectrum1D_div_scalar(hdrl_spectrum1D *self, hdrl_value scalar_operator)
computes the elementwise division of a spectrum by a scalar, the self parameter is modified
cpl_image * estimate_noise_DER_SNR(const hdrl_data_t *flux_in, const cpl_binary *msk_in, const cpl_array *wavelengths, const cpl_size length, const cpl_size half_window)
For every pixel in position i in img_arg, the function estimates the noise using the pixels in the wi...
Definition: hdrl_DER_SNR.c:160
cpl_boolean hdrl_spectrum1D_are_spectra_compatible(const hdrl_spectrum1D_wavelength *s1, const hdrl_spectrum1D_wavelength *s2)
checks if two spectrum wavelengths are equal.
hdrl_spectrum1D * hdrl_spectrum1D_reject_pixels(const hdrl_spectrum1D *self, const cpl_array *bad_samples)
For every i-th element in bad_samples having value CPL_TRUE, the i-th pixel in the 1D spectrum is mar...
hdrl_value hdrl_spectrum1D_get_flux_value(const hdrl_spectrum1D *self, int idx, int *rej)
hdrl_spectrum1D getter for a flux value