CR2RE Pipeline Reference Manual 1.6.8
hdrl_spectrum1dlist-test.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
28#include "hdrl.h"
29#include "hdrl_spectrum_resample.h"
30
31#include <math.h>
32
33
34hdrl_spectrum1D * create_spectrum(double fx){
35 cpl_image * flx = cpl_image_new(1, 1, HDRL_TYPE_DATA);
36 cpl_array * wln = cpl_array_new(1, HDRL_TYPE_DATA);
37 cpl_image_set(flx, 1, 1 , fx);
38 cpl_array_set(wln, 0, 100);
39 hdrl_spectrum1D * s = hdrl_spectrum1D_create_error_free(flx, wln,
40 hdrl_spectrum1D_wave_scale_linear);
41
42 cpl_array_delete(wln);
43 cpl_image_delete(flx);
44
45 return s;
46}
47
48void check_equal(const hdrl_spectrum1D * s1, const hdrl_spectrum1D * s2){
49
50 const cpl_size sz = hdrl_spectrum1D_get_size(s1);
51 cpl_test_eq(sz, hdrl_spectrum1D_get_size(s2));
52
53 for(cpl_size i = 0; i < sz; ++i){
54 int rej1 = 0, rej2 = 0;
55 hdrl_value v1 = hdrl_spectrum1D_get_flux_value(s1, i, &rej1);
56 hdrl_value v2 = hdrl_spectrum1D_get_flux_value(s2, i, &rej2);
57
58 cpl_test_eq(rej1, rej2);
59 cpl_test_rel(v1.data, v2.data, 1e-5);
60 cpl_test_rel(v1.error, v2.error, 1e-5);
61 }
62}
63
64void check_list_sequential(const hdrl_spectrum1Dlist * list,
65 cpl_size idx, double mlx){
66
67 for(cpl_size i = 0; i < hdrl_spectrum1Dlist_get_size(list); ++i){
68 const hdrl_spectrum1D * s = hdrl_spectrum1Dlist_get_const(list, i);
69 hdrl_value v = hdrl_spectrum1D_get_flux_value(s, 0, NULL);
70 cpl_test_rel(v.data, idx * mlx, 1e-10);
71 idx++;
72 }
73
74}
75
76void test_spectrum1dlist_wrap(void){
77 hdrl_spectrum1D ** sl = cpl_calloc(6, sizeof(hdrl_spectrum1D*));
78
79 sl[0] = create_spectrum(1);
80 sl[1] = create_spectrum(2);
81 sl[2] = create_spectrum(3);
82 sl[3] = create_spectrum(4);
83 sl[4] = create_spectrum(5);
84 sl[5] = create_spectrum(6);
85
86 hdrl_spectrum1Dlist * list = hdrl_spectrum1Dlist_wrap(sl, 6);
87
88 cpl_test_eq(hdrl_spectrum1Dlist_get_size(list), 6);
89
90 for(cpl_size i = 0; i < hdrl_spectrum1Dlist_get_size(list); ++i){
91 const hdrl_spectrum1D * s = hdrl_spectrum1Dlist_get_const(list, i);
92 cpl_test_eq_ptr(s, sl[i]);
93 }
94
96}
97
98void test_spectrum1dlist(void){
99
100 hdrl_spectrum1Dlist * list1 = hdrl_spectrum1Dlist_new();
101
102 cpl_test_eq(0, hdrl_spectrum1Dlist_get_size(list1));
103
104 hdrl_spectrum1D * s4 = create_spectrum(4);
105
106 {
107 hdrl_spectrum1D * s1 = create_spectrum(1);
108 hdrl_spectrum1D * s11 = create_spectrum(1);
109 hdrl_spectrum1D * s2 = create_spectrum(2);
110 hdrl_spectrum1D * s3 = create_spectrum(3);
111 hdrl_spectrum1D * s5 = create_spectrum(5);
112 hdrl_spectrum1D * s6 = create_spectrum(6);
113
114 hdrl_spectrum1Dlist_set(list1, s1, 0);
115 /*Removed s1 when setting*/
116 hdrl_spectrum1Dlist_set(list1, s11, 0);
117 hdrl_spectrum1Dlist_set(list1, s2, 1);
118 hdrl_spectrum1Dlist_set(list1, s3, 2);
119 hdrl_spectrum1Dlist_set(list1, s4, 3);
120 hdrl_spectrum1Dlist_set(list1, s5, 4);
121 hdrl_spectrum1Dlist_set(list1, s6, 5);
122 }
123
124 cpl_test_eq(6, hdrl_spectrum1Dlist_get_size(list1));
125
126 /* Check duplication */
127 hdrl_spectrum1Dlist * list2 = hdrl_spectrum1Dlist_duplicate(list1);
128
129 cpl_test_eq(hdrl_spectrum1Dlist_get_size(list1),
131
132 for(cpl_size i = 0; i < hdrl_spectrum1Dlist_get_size(list1); ++i){
133 const hdrl_spectrum1D * s1 = hdrl_spectrum1Dlist_get_const(list1, i);
134 const hdrl_spectrum1D * s2 = hdrl_spectrum1Dlist_get_const(list2, i);
135
136 check_equal(s1, s2);
137
138 /*check that list1 and list2 point to different memory spaces*/
139 cpl_test_noneq_ptr(s1, s2);
140 }
141
142 /*Check getters (mutability)*/
143 for(cpl_size i = 0; i < hdrl_spectrum1Dlist_get_size(list1); ++i){
144 hdrl_spectrum1D * s = hdrl_spectrum1Dlist_get(list2, i);
145 hdrl_value v = hdrl_spectrum1D_get_flux_value(s, 0, NULL);
146 cpl_test_rel(v.data, (double)i + 1.0 , 1e-5);
147 hdrl_spectrum1D_mul_scalar(s, (hdrl_value){5.0, 0.0});
148 }
149
150 for(cpl_size i = 0; i < hdrl_spectrum1Dlist_get_size(list1); ++i){
151 hdrl_spectrum1D * s = hdrl_spectrum1Dlist_get(list2, i);
152 hdrl_value v = hdrl_spectrum1D_get_flux_value(s, 0, NULL);
153 cpl_test_rel(v.data, ((double)i + 1.0) * 5.0 , 1e-5);
154 }
155
156 hdrl_spectrum1D * new_s4 = hdrl_spectrum1Dlist_unset(list1, 3);
157
158 cpl_test_eq_ptr(new_s4, s4);
159
160 cpl_test_eq(hdrl_spectrum1Dlist_get_size(list1), 5);
161
162 const double flx_value[] = {1, 2, 3, 5, 6};
163 for(cpl_size i = 0; i < hdrl_spectrum1Dlist_get_size(list1); ++i){
164
165 const hdrl_spectrum1D * s = hdrl_spectrum1Dlist_get_const(list1, i);
166 hdrl_value v = hdrl_spectrum1D_get_flux_value(s, 0, NULL);
167 cpl_test_rel(v.data, flx_value[i], 1e-5);
168 }
169
170 cpl_size i = 1;
171 while(hdrl_spectrum1Dlist_get_size(list2)){
172 hdrl_spectrum1D * s = hdrl_spectrum1Dlist_unset(list2, 0);
173 check_list_sequential(list2, i + 1, 5.0);
174 hdrl_value v = hdrl_spectrum1D_get_flux_value(s, 0, NULL);
175 cpl_test_rel(v.data, i * 5.0, 1e-5);
177 i++;
178 }
179
183}
184 /*data[i] is both wavelength and flux. if data[i] the i-th pixel is rejected*/
185hdrl_spectrum1D *
186create_spectrum_long(const double * data, const cpl_size length){
187 cpl_array * wlens = cpl_array_new(length, HDRL_TYPE_DATA);
188 cpl_image * flx = cpl_image_new(length, 1, HDRL_TYPE_DATA);
189
190 for(cpl_size i = 0; i < length; ++i){
191 cpl_array_set(wlens, i, fabs(data[i]));
192
193 if(data[i] >= 0)
194 cpl_image_set(flx, i + 1, 1, data[i]);
195 else
196 cpl_image_reject(flx, i + 1, 1);
197 }
198
199 hdrl_spectrum1D * s = hdrl_spectrum1D_create_error_free(flx, wlens,
200 hdrl_spectrum1D_wave_scale_linear);
201
202 cpl_image_delete(flx);
203 cpl_array_delete(wlens);
204 return s;
205}
206
207/*check that bad pixels at the beginning or at the end of spectrum
208 * do not contribute to the collapsed spectrum*/
209void test_spectrum1dlist_collapse_badpix(void){
210
211 hdrl_spectrum1Dlist * l = hdrl_spectrum1Dlist_new();
212
213 double aValues1[4] = {1, 2, 3, 4};
214 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues1, 4), 0);
215
216 double aValues2[3] = {-1, 2, 4};
217 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues2, 3), 1);
218
219 double aValues3[3] = {1, 3, -4};
220 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues3, 3), 2);
221
222 cpl_array * wlenghts = cpl_array_new(6, HDRL_TYPE_DATA);
223 for(cpl_size i = 0; i < cpl_array_get_size(wlenghts); ++i){
224 cpl_array_set(wlenghts, i, i);
225 }
226
227 hdrl_parameter * stacking_par = hdrl_collapse_mean_parameter_create();
228 hdrl_parameter * resampling_par =
229 hdrl_spectrum1D_resample_interpolate_parameter_create(hdrl_spectrum1D_interp_linear);
230
231 hdrl_spectrum1D * res = NULL;
232 cpl_image * contrib = NULL;
233 hdrl_imagelist * aligned_fluxes = NULL;
234
235 hdrl_spectrum1Dlist_collapse(l, stacking_par, wlenghts, resampling_par, CPL_FALSE,
236 &res, &contrib, &aligned_fluxes);
237
238 int contribs[] = {2, 3, 3, 2};
239 for(cpl_size i = 1; i < 5; ++i){
240 int rej = 0;
241 const double el = cpl_image_get(contrib, i + 1, 1, &rej);
242 cpl_test_eq(rej, 0);
243 cpl_test_eq(el, contribs[i - 1]);
244 }
245
246 int reject = 0;
247 double el = cpl_image_get(contrib, 6, 1, &reject);
248 cpl_test_eq(reject, 0);
249 cpl_test_eq(el, 0);
250
251 reject = 0;
252 el = cpl_image_get(contrib, 1, 1, &reject);
253 cpl_test_eq(reject, 0);
254 cpl_test_eq(el, 0);
255
256 const cpl_size sz = hdrl_spectrum1D_get_size(res);
257 cpl_test_eq(cpl_array_get_size(wlenghts), sz);
258
259 for(cpl_size i = 0; i < sz; ++i){
260 int rej = 0;
261 const hdrl_value v = hdrl_spectrum1D_get_flux_value(res, i, &rej);
262 if(i == 0 || i == (sz - 1)){
263 cpl_test_eq(rej, 1);
264 continue;
265 }
266 cpl_test_rel(v.data, i, 1e-5);
267 cpl_test_eq(rej, 0);
268 }
269
272 cpl_image_delete(contrib);
273 hdrl_parameter_delete(stacking_par);
274 hdrl_parameter_delete(resampling_par);
275 hdrl_imagelist_delete(aligned_fluxes);
276 cpl_array_delete(wlenghts);
277}
278
279/*Check that resamplking pixels having rejected neighbors in the original spectrum
280 * do not contribute to the stacking*/
281void test_spectrum1dlist_collapse_mark_rej_in_interpolation(void){
282
283 hdrl_spectrum1Dlist * l = hdrl_spectrum1Dlist_new();
284
285 double aValues1[4] = {1, 2, 3, 4};
286 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues1, 4), 0);
287
288 double aValues2[4] = {-1, 2, -3, 4};
289 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues2, 4), 1);
290
291 double aValues3[4] = {1,-2, 3, -4};
292 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues3, 4), 2);
293
294 cpl_array * wlenghts = cpl_array_new(6, HDRL_TYPE_DATA);
295 for(cpl_size i = 0; i < cpl_array_get_size(wlenghts); ++i){
296 cpl_array_set(wlenghts, i, i);
297 }
298
299 hdrl_parameter * stacking_par = hdrl_collapse_mean_parameter_create();
300 hdrl_parameter * resampling_par =
301 hdrl_spectrum1D_resample_interpolate_parameter_create(hdrl_spectrum1D_interp_linear);
302
303 hdrl_spectrum1D * res = NULL;
304 cpl_image * contrib = NULL;
305 hdrl_imagelist * aligned_fluxes = NULL;
306
307 hdrl_spectrum1Dlist_collapse(l, stacking_par, wlenghts, resampling_par, CPL_TRUE,
308 &res, &contrib, &aligned_fluxes);
309
310 int contribs[] = {2, 2, 2, 2};
311 for(cpl_size i = 1; i < 5; ++i){
312 int reject = 0;
313 const double el = cpl_image_get(contrib, i + 1, 1, &reject);
314 cpl_test_eq(reject, 0);
315 cpl_test_eq(el, contribs[i - 1]);
316 }
317
318 int rej_contr = 0;
319 double el = cpl_image_get(contrib, 6, 1, &rej_contr);
320 cpl_test_eq(rej_contr, 0);
321 cpl_test_eq(el, 0);
322
323 rej_contr = 0;
324 el = cpl_image_get(contrib, 1, 1, &rej_contr);
325 cpl_test_eq(rej_contr, 0);
326 cpl_test_eq(el, 0);
327
328 const cpl_size sz = hdrl_spectrum1D_get_size(res);
329 cpl_test_eq(cpl_array_get_size(wlenghts), sz);
330
331 for(cpl_size i = 0; i < sz; ++i){
332 int rej = 0;
333 const hdrl_value v = hdrl_spectrum1D_get_flux_value(res, i, &rej);
334 if(i == 0 || i == (sz - 1)){
335 cpl_test_eq(rej, 1);
336 continue;
337 }
338 cpl_test_rel(v.data, i, 1e-5);
339 cpl_test_eq(rej, 0);
340 }
341
344 cpl_image_delete(contrib);
345 hdrl_parameter_delete(stacking_par);
346 hdrl_parameter_delete(resampling_par);
347 hdrl_imagelist_delete(aligned_fluxes);
348 cpl_array_delete(wlenghts);
349}
350
351/*check that shorted spectra are handled correctly*/
352void test_spectrum1dlist_collapse_holes(void){
353
354 hdrl_spectrum1Dlist * l = hdrl_spectrum1Dlist_new();
355
356 double aValues1[4] = {1, 2, 3, 4};
357 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues1, 4), 0);
358
359 double aValues2[2] = {2, 4};
360 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues2, 2), 1);
361
362 double aValues3[2] = {1, 3};
363 hdrl_spectrum1Dlist_set(l, create_spectrum_long(aValues3, 2), 2);
364
365 cpl_array * wlenghts = cpl_array_new(6, HDRL_TYPE_DATA);
366 for(cpl_size i = 0; i < cpl_array_get_size(wlenghts); ++i){
367 cpl_array_set(wlenghts, i, i);
368 }
369
370 hdrl_parameter * stacking_par = hdrl_collapse_mean_parameter_create();
371 hdrl_parameter * resampling_par =
372 hdrl_spectrum1D_resample_interpolate_parameter_create(hdrl_spectrum1D_interp_linear);
373
374 hdrl_spectrum1D * res = NULL;
375 cpl_image * contrib = NULL;
376 hdrl_imagelist * aligned_fluxes = NULL;
377 hdrl_spectrum1Dlist_collapse(l, stacking_par, wlenghts, resampling_par,
378 CPL_FALSE, &res, &contrib, &aligned_fluxes);
379
380 int contribs[] = {2, 3, 3, 2};
381 for(cpl_size i = 1; i < 5; ++i){
382 int rej_contr = 0;
383 const double el = cpl_image_get(contrib, i + 1, 1, &rej_contr);
384 cpl_test_eq(rej_contr, 0);
385 cpl_test_eq(el, contribs[i - 1]);
386 }
387
388 int rej_contrib = 0;
389 double el = cpl_image_get(contrib, 6, 1, &rej_contrib);
390 cpl_test_eq(rej_contrib, 0);
391 cpl_test_eq(el, 0);
392
393 rej_contrib = 0;
394 el = cpl_image_get(contrib, 1, 1, &rej_contrib);
395 cpl_test_eq(rej_contrib, 0);
396 cpl_test_eq(el, 0);
397
398 const cpl_size sz = hdrl_spectrum1D_get_size(res);
399 cpl_test_eq(cpl_array_get_size(wlenghts), sz);
400
401 for(cpl_size i = 0; i < sz; ++i){
402 int rej = 0;
403 const hdrl_value v = hdrl_spectrum1D_get_flux_value(res, i, &rej);
404 if(i == 0 || i == (sz - 1)){
405 cpl_test_eq(rej, 1);
406 continue;
407 }
408 cpl_test_rel(v.data, i, 1e-5);
409 cpl_test_eq(rej, 0);
410 }
411
414 cpl_image_delete(contrib);
415 hdrl_imagelist_delete(aligned_fluxes);
416 hdrl_parameter_delete(stacking_par);
417 hdrl_parameter_delete(resampling_par);
418 cpl_array_delete(wlenghts);
419}
420
421void test_error_and_reset(cpl_error_code expected){
422 cpl_test_eq_error(expected, cpl_error_get_code());
423 cpl_error_reset();
424}
425
426void test_spectrum1dlist_insert_duplication(void){
427
428 hdrl_spectrum1Dlist * list1 = hdrl_spectrum1Dlist_new();
429
430 cpl_test_eq(0, hdrl_spectrum1Dlist_get_size(list1));
431
432 hdrl_spectrum1D * s1 = create_spectrum(1);
433 hdrl_spectrum1D * s2 = create_spectrum(2);
434 hdrl_spectrum1D * s3 = create_spectrum(3);
435 hdrl_spectrum1D * s4 = create_spectrum(4);
436 hdrl_spectrum1D * s5 = create_spectrum(5);
437 hdrl_spectrum1D * s6 = create_spectrum(6);
438
439 hdrl_spectrum1Dlist_set(list1, s1, 0);
440 hdrl_spectrum1Dlist_set(list1, s2, 1);
441 hdrl_spectrum1Dlist_set(list1, s3, 2);
442 hdrl_spectrum1Dlist_set(list1, s4, 3);
443 hdrl_spectrum1Dlist_set(list1, s5, 4);
444 hdrl_spectrum1Dlist_set(list1, s6, 5);
445
446 hdrl_spectrum1Dlist_set(list1, s1, 4);
447 test_error_and_reset(CPL_ERROR_ILLEGAL_INPUT);
448
449 hdrl_spectrum1Dlist_set(list1, s2, 3);
450 test_error_and_reset(CPL_ERROR_ILLEGAL_INPUT);
451
452 hdrl_spectrum1Dlist_set(list1, s3, 4);
453 test_error_and_reset(CPL_ERROR_ILLEGAL_INPUT);
454
455 hdrl_spectrum1Dlist_set(list1, s4, 5);
456 test_error_and_reset(CPL_ERROR_ILLEGAL_INPUT);
457
458 hdrl_spectrum1Dlist_set(list1, s5, 0);
459 test_error_and_reset(CPL_ERROR_ILLEGAL_INPUT);
460
461 hdrl_spectrum1Dlist_set(list1, s6, 2);
462 test_error_and_reset(CPL_ERROR_ILLEGAL_INPUT);
463
465}
466
467/*----------------------------------------------------------------------------*/
471/*----------------------------------------------------------------------------*/
472int main(void)
473{
474 cpl_test_init(PACKAGE_BUGREPORT, CPL_MSG_WARNING);
475
476 test_spectrum1dlist();
477 test_spectrum1dlist_wrap();
478 test_spectrum1dlist_collapse_holes();
479 test_spectrum1dlist_collapse_badpix();
480 test_spectrum1dlist_collapse_mark_rej_in_interpolation();
481 test_spectrum1dlist_insert_duplication();
482
483 cpl_test_error(CPL_ERROR_NONE);
484
485 return cpl_test_end(0);
486}
487
hdrl_parameter * hdrl_collapse_mean_parameter_create(void)
create a parameter object for mean
void hdrl_imagelist_delete(hdrl_imagelist *himlist)
Free all memory used by a hdrl_imagelist object including the images.
void hdrl_parameter_delete(hdrl_parameter *obj)
shallow delete of a parameter
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_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
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_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_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