27 #include <cxmessages.h>
28 #include <cxstrutils.h>
30 #include <cpl_error.h>
31 #include <cpl_parameterlist.h>
32 #include <cpl_propertylist.h>
33 #include <cpl_matrix.h>
34 #include <cpl_image.h>
35 #include <cpl_table.h>
40 #include "gifiberutils.h"
41 #include "gigrating.h"
43 #include "gilocalization.h"
44 #include "giextraction.h"
45 #include "girebinning.h"
46 #include "gisgcalibration.h"
57 struct GiMeasurement {
62 typedef struct GiMeasurement GiMeasurement;
80 typedef struct GiSGSetup GiSGSetup;
83 struct GiCPFitParams {
98 typedef struct GiCPFitParams GiCPFitParams;
102 GiMeasurement amplitude;
103 GiMeasurement background;
104 GiMeasurement center;
110 typedef struct GiCPeakFit GiCPeakFit;
123 cpl_matrix* wavelength;
128 typedef struct GiSGMask GiSGMask;
131 inline static GiSGMask*
132 _giraffe_sgmask_new(cxsize size)
135 GiSGMask*
self = cx_calloc(1,
sizeof *
self);
137 self->wavelength = cpl_matrix_new(1, size);
138 self->flux = cpl_matrix_new(1, size);
143 self->scale = GIREBIN_SCALE_LINEAR;
154 _giraffe_sgmask_delete(GiSGMask*
self)
159 if (self->wavelength != NULL) {
160 cpl_matrix_delete(self->wavelength);
161 self->wavelength = NULL;
164 if (self->flux != NULL) {
165 cpl_matrix_delete(self->flux);
178 inline static GiSGMask*
179 _giraffe_sgmask_create(cxsize size, cxdouble start, cxdouble step,
180 GiRebinScale scale,
const GiTable* mask)
187 cxdouble wlstep = 0.;
189 cpl_table* _mask = NULL;
191 GiSGMask*
self = NULL;
194 cx_assert(mask != NULL);
197 cx_assert(_mask != NULL);
199 self = _giraffe_sgmask_new(size);
210 for (i = 0; i <
self->size; i++) {
211 cpl_matrix_set(self->wavelength, 0, i, self->start + i * self->step);
215 wlmin = cpl_matrix_get(self->wavelength, 0, 0);
216 wlmax = cpl_matrix_get(self->wavelength, 0, self->size - 1);
219 if (self->scale == GIREBIN_SCALE_LOG) {
223 wlstep = exp(wlstep);
233 cpl_table_select_all(_mask);
235 cpl_table_and_selected_double(_mask,
"WLEN1", CPL_GREATER_THAN, wlmin);
236 cpl_table_and_selected_double(_mask,
"WLEN2", CPL_LESS_THAN, wlmax);
238 _mask = cpl_table_extract_selected(_mask);
240 if (_mask == NULL || cpl_table_get_nrow(_mask) <= 0) {
241 _giraffe_sgmask_delete(
self);
248 self->nholes = cpl_table_get_nrow(_mask);
250 for (i = 0; i <
self->nholes; i++) {
254 cxdouble hstart = cpl_table_get(_mask,
"WLEN1", i, NULL) - wlmin;
255 cxdouble hend = cpl_table_get(_mask,
"WLEN2", i, NULL) - wlmin;
261 for (j = (cxsize)(hstart + 0.5); j <= (cxsize)(hend + 0.5); j++) {
263 cpl_matrix_set(self->flux, 0, j, 1.);
269 cpl_table_delete(_mask);
277 _giraffe_sgmask_size(
const GiSGMask*
self)
280 cx_assert(
self != NULL);
288 _giraffe_sgmask_holes(
const GiSGMask*
self)
291 cx_assert(
self != NULL);
299 _giraffe_sgmask_set_flux(GiSGMask*
self, cxsize position, cxdouble value)
302 cx_assert(
self != NULL);
304 if (position >= (cxsize)cpl_matrix_get_ncol(self->flux)) {
308 cpl_matrix_set(self->flux, 0, position, value);
315 inline static cxdouble
316 _giraffe_sgmask_get_flux(GiSGMask*
self, cxsize position)
319 const cxchar*
const fctid =
"_giraffe_sgmask_get_flux";
322 cx_assert(
self != NULL);
324 if (position >= (cxsize)cpl_matrix_get_ncol(self->flux)) {
325 cpl_error_set(fctid, CPL_ERROR_ILLEGAL_INPUT);
329 return cpl_matrix_get(self->flux, 0, position);
334 inline static const cpl_matrix*
335 _giraffe_sgmask_get(GiSGMask*
self)
338 cx_assert(
self != NULL);
346 _giraffe_sgmask_crop(GiSGMask*
self, cxsize begin, cxsize end)
351 cpl_matrix* buffer = NULL;
354 cx_assert(
self != NULL);
355 cx_assert(end > begin);
356 cx_assert(cpl_matrix_get_nrow(self->wavelength) == 1);
357 cx_assert(cpl_matrix_get_nrow(self->flux) == 1);
359 if (begin >= (cxsize)cpl_matrix_get_ncol(self->flux)) {
363 if (end > (cxsize)cpl_matrix_get_ncol(self->flux)) {
364 end = cpl_matrix_get_ncol(self->flux);
367 if (begin == 0 && end == (cxsize)cpl_matrix_get_ncol(self->flux)) {
373 buffer = cpl_matrix_extract(self->wavelength, 0, begin, 1, 1, 1, size);
374 cpl_matrix_delete(self->wavelength);
375 self->wavelength = buffer;
377 buffer = cpl_matrix_extract(self->flux, 0, begin, 1, 1, 1, size);
378 cpl_matrix_delete(self->flux);
381 cx_assert(cpl_matrix_get_nrow(self->flux) == 1);
382 cx_assert((cxsize)cpl_matrix_get_ncol(self->flux) == size);
391 inline static cxdouble
392 _giraffe_clip_value(cxdouble value, cxdouble low, cxdouble high,
396 cxbool status = FALSE;
423 inline static cpl_image*
424 _giraffe_resample_image(cpl_image* signal, cxdouble step1, cxdouble step2)
431 cxint step = CX_MAX(1, (cxint)(step1/step2));
433 cpl_image* _signal = NULL;
436 cx_assert(signal != NULL);
438 ny = cpl_image_get_size_x(signal);
439 nx1 = cpl_image_get_size_y(signal);
441 nx2 = (nx1 - 1) * step + 1;
443 _signal = cpl_image_new(ny, nx2, CPL_TYPE_DOUBLE);
445 for (i = 0; i < ny; i++) {
449 register cxdouble* data = cpl_image_get_data(signal);
450 register cxdouble* _data = cpl_image_get_data(_signal);
453 for (j = 0; j < nx1 - 1; j++) {
456 register cxint l = j * ny + i;
457 register cxint m = j * ny * step + i;
459 for (k = 0; k < step; k++) {
461 cxdouble f = (cxdouble)k / (cxdouble)step;
463 _data[m + k * ny] = (1. - f) * data[l] + f * data[l + ny];
469 _data[nx2 - 1] = data[nx1 - 1];
482 inline static cpl_matrix*
483 _giraffe_compute_cross_correlation(
const cpl_matrix* signal,
484 const cpl_matrix*
template,
485 cxint start, cxint end)
488 const cxchar*
const fctid =
"_giraffe_compute_cross_correlation";
499 cpl_matrix* _signal = (cpl_matrix*)signal;
500 cpl_matrix* _template = (cpl_matrix*)
template;
501 cpl_matrix* ccf = NULL;
502 cpl_matrix* _ccf = NULL;
505 cx_assert(_signal != NULL);
506 cx_assert(cpl_matrix_get_nrow(_signal) == 1);
508 cx_assert(_template != NULL);
509 cx_assert(cpl_matrix_get_nrow(_template) == 1);
511 ns = cpl_matrix_get_ncol(_signal);
512 cx_assert(ns == cpl_matrix_get_ncol(_template));
514 cx_assert(start <= end);
521 nmax = cpl_matrix_get_ncol(_signal) / 2;
523 start = CX_MAX(CX_MIN(start, nmax), -nmax);
524 end = CX_MAX(CX_MIN(end, nmax), -nmax);
528 cpl_msg_debug(fctid,
"Cross-correlation function: signal size = %"
529 CPL_SIZE_FORMAT
", template size = %" CPL_SIZE_FORMAT
530 ", window start = %d, window end = %d",
531 cpl_matrix_get_ncol(_signal), cpl_matrix_get_ncol(_template),
535 ccf = cpl_matrix_new(1, nccf);
537 for (i = start; i < end; i++) {
550 for (j = 0; j < ns + i; j++) {
552 cxdouble s = cpl_matrix_get(_signal, 0, j);
553 cxdouble t = cpl_matrix_get(_template, 0, j - i);
559 sum /= (cxdouble)(ns + i);
561 cpl_matrix_set(ccf, 0, i - start, sum);
575 for (j = i; j < ns; j++) {
577 cxdouble s = cpl_matrix_get(_signal, 0, j);
578 cxdouble t = cpl_matrix_get(_template, 0, j - i);
584 sum /= (cxdouble)(ns - i);
586 cpl_matrix_set(ccf, 0, i - start, sum);
600 for (j = 0; j < ns; j++) {
602 cxdouble t = cpl_matrix_get(_template, 0, j);
603 cxdouble s = cpl_matrix_get(_signal, 0, j);
611 cpl_matrix_set(ccf, 0, -start, sum);
625 n = CX_MAX(1, nccf / 10);
627 _ccf = cpl_matrix_duplicate(ccf);
632 for (i = nccf - n; i < nccf; i++) {
633 sum += cpl_matrix_get(_ccf, 0, i);
638 cpl_matrix_delete(_ccf);
643 for (i = 0; i < nccf; i++) {
644 cpl_matrix_set(ccf, 0, i, cpl_matrix_get(ccf, 0, i) / sum);
655 _giraffe_create_setup(GiSGSetup* setup,
const GiImage* spectra)
658 cpl_propertylist* properties = NULL;
660 cpl_image* _spectra = NULL;
663 cx_assert(setup != NULL);
664 cx_assert(spectra != NULL);
667 cx_assert(properties != NULL);
670 cx_assert(_spectra != NULL);
677 setup->nx = cpl_image_get_size_y(_spectra);
680 if (!cpl_propertylist_has(properties, GIALIAS_EXT_NX)) {
685 setup->nex = cpl_propertylist_get_int(properties, GIALIAS_EXT_NX);
689 if (!cpl_propertylist_has(properties, GIALIAS_BINSCALE)) {
694 const cxchar* s = cpl_propertylist_get_string(properties,
698 if (cx_strncasecmp(s,
"log", 3) == 0) {
699 setup->scale = GIREBIN_SCALE_LOG;
702 setup->scale = GIREBIN_SCALE_LINEAR;
707 if (!cpl_propertylist_has(properties, GIALIAS_BINWLMIN)) {
711 setup->wlmin = cpl_propertylist_get_double(properties,
715 if (!cpl_propertylist_has(properties, GIALIAS_BINSTEP)) {
719 setup->wlstep = cpl_propertylist_get_double(properties,
723 setup->wlmax = setup->wlmin + (setup->nx - 1) * setup->wlstep;
726 if (!cpl_propertylist_has(properties, GIALIAS_PIXSIZY)) {
730 setup->pixelsize = cpl_propertylist_get_double(properties,
740 _giraffe_peak_fit(GiCPeakFit* peak,
const cpl_matrix* lambda,
741 const cpl_matrix* ccf,
const GiGrating* grating,
742 const GiCPFitParams* setup)
745 const cxchar*
const fctid =
"_giraffe_peak_fit";
754 cxdouble amplitude = 0.;
755 cxdouble background = 0.;
756 cxdouble center = 0.;
766 } initial = {0., 0., 0., 0.};
771 GiModel* model = giraffe_model_new(
"gaussian");
775 cx_assert(model != NULL);
776 cx_assert(strcmp(giraffe_model_get_name(model),
"gaussian") == 0);
777 cx_assert(lambda != NULL);
778 cx_assert(ccf != NULL);
779 cx_assert(grating != NULL);
780 cx_assert(setup != NULL);
792 amplitude = cpl_matrix_get_max((cpl_matrix*)ccf) - background;
794 cpl_matrix_get_maxpos((cpl_matrix*)ccf, &nr, &nc);
797 center = cpl_matrix_get((cpl_matrix*)lambda, 0, nc);
800 if (setup->scale == GIREBIN_SCALE_LOG) {
801 width = 0.5 / grating->
resol;
804 width = 0.5 / grating->
resol * grating->
wlen0;
807 giraffe_model_set_parameter(model,
"Background", background);
808 giraffe_model_set_parameter(model,
"Amplitude", amplitude);
809 giraffe_model_set_parameter(model,
"Center", center);
810 giraffe_model_set_parameter(model,
"Width1", width);
812 giraffe_model_thaw(model);
814 giraffe_model_set_iterations(model, setup->fit.iterations);
815 giraffe_model_set_tests(model, setup->fit.tests);
816 giraffe_model_set_delta(model, setup->fit.delta);
823 initial.amplitude = amplitude;
824 initial.background = background;
825 initial.center = center;
826 initial.width = width;
830 while (i < setup->iterations && !stop) {
839 cpl_matrix* tlambda = (cpl_matrix*)lambda;
840 cpl_matrix* tccf = (cpl_matrix*)ccf;
851 const cxdouble da = 0.2;
852 const cxdouble dc = 1.;
853 const cxdouble db = 1.;
854 const cxdouble dw = 0.2;
858 value = giraffe_model_get_parameter(model,
"Amplitude") * da;
859 value += (1. - da) * initial.amplitude;
861 giraffe_model_set_parameter(model,
"Amplitude", value);
864 value = giraffe_model_get_parameter(model,
"Center") * dc;
865 value += (1. - dc) * initial.center;
867 giraffe_model_set_parameter(model,
"Center", value);
870 value = giraffe_model_get_parameter(model,
"Background") * db;
871 value += (1. - db) * initial.background;
873 giraffe_model_set_parameter(model,
"Background", value);
876 value = giraffe_model_get_parameter(model,
"Width1") * dw;
877 value += (1. - dw) * initial.width;
879 giraffe_model_set_parameter(model,
"Width1", value);
889 dwd = 2. * giraffe_model_get_parameter(model,
"Width1") *
891 dwc = giraffe_model_get_parameter(model,
"Center");
893 dwd = CX_MAX(setup->dnmin, 2. * dwd / setup->step) * setup->step / 2.;
895 lower = dwc + 0.5 * setup->step - dwd;
896 upper = dwc + 0.5 * setup->step + dwd;
905 for (j = 0; j < cpl_matrix_get_ncol(tlambda); j++) {
907 if (cpl_matrix_get(tlambda, 0, j) > lower) {
914 for (j = cpl_matrix_get_ncol(tlambda) - 1; j >= 0; j--) {
916 if (cpl_matrix_get(tlambda, 0, j) < upper) {
924 if (i > 0 && dn1 == _dn1 && dn2 == _dn2) {
926 cxdouble _width = giraffe_model_get_parameter(model,
"Width1");
933 dwd = CX_MAX(setup->dnmin, 4. * _width * setup->wfactor /
934 setup->step) * setup->step / 2.;
936 lower = dwc + 0.5 * setup->step - dwd;
937 upper = dwc + 0.5 * setup->step + dwd;
939 for (j = 0; j < cpl_matrix_get_ncol(tlambda); j++) {
941 if (cpl_matrix_get(tlambda, 0, j) > lower) {
948 for (j = cpl_matrix_get_ncol(tlambda) - 1; j <= 0; j--) {
950 if (cpl_matrix_get(tlambda, 0, j) < upper) {
967 if (i <= 1 || dn1 != _dn1 || dn2 != _dn2) {
971 const cxint pflag = 1;
974 cxdouble damplitude = 0.;
975 cxdouble dbackground = 0.;
976 cxdouble dcenter = 0.;
977 cxdouble dwidth = 0.;
979 cpl_matrix* x = NULL;
980 cpl_matrix* y = NULL;
981 cpl_matrix* sigma = NULL;
987 x = cpl_matrix_new(dn2 - dn1, 1);
988 y = cpl_matrix_new(dn2 - dn1, 1);
989 sigma = cpl_matrix_new(dn2 - dn1, 1);
991 for (j = 0; j < cpl_matrix_get_nrow(y); j++) {
993 cpl_matrix_set(x, j, 0, cpl_matrix_get(tlambda, 0, dn1 + j));
994 cpl_matrix_set(y, j, 0, cpl_matrix_get(tccf, 0, dn1 + j));
995 cpl_matrix_set(sigma, j, 0, setup->sigma);
1004 status = giraffe_model_fit(model, x, y, sigma);
1008 peak->amplitude.value = initial.amplitude;
1009 peak->background.value = initial.background;
1010 peak->center.value = initial.center;
1011 peak->width.value = initial.width;
1013 peak->amplitude.sigma = 1.;
1014 peak->background.sigma = 1.;
1015 peak->center.sigma = 1.;
1016 peak->width.sigma = 1.;
1020 cpl_matrix_delete(x);
1021 cpl_matrix_delete(y);
1022 cpl_matrix_delete(sigma);
1024 giraffe_model_delete(model);
1036 amplitude = giraffe_model_get_parameter(model,
"Amplitude");
1037 damplitude = giraffe_model_get_sigma(model,
"Amplitude");
1039 center = giraffe_model_get_parameter(model,
"Center");
1040 dcenter = giraffe_model_get_sigma(model,
"Center");
1042 background = giraffe_model_get_parameter(model,
"Background");
1043 dbackground = giraffe_model_get_sigma(model,
"Background");
1045 width = giraffe_model_get_parameter(model,
"Width1");
1046 dwidth = giraffe_model_get_sigma(model,
"Width1");
1054 lower = -9. * (1 - pflag) + 1.e-5 * pflag;
1055 upper = 9. * pflag - 1.e-5 * (1 - pflag);
1057 peak->amplitude.value = _giraffe_clip_value(amplitude, lower,
1059 peak->amplitude.sigma = _giraffe_clip_value(damplitude, 0.,
1062 stop = stop == FALSE ? flag == TRUE ? TRUE : FALSE : stop;
1066 lower = cpl_matrix_get(x, 1, 0);
1067 upper = cpl_matrix_get(x, cpl_matrix_get_nrow(x) - 2, 0);
1069 peak->center.value = _giraffe_clip_value(center, lower,
1071 peak->center.sigma = _giraffe_clip_value(dcenter, 0.,
1072 initial.width, NULL);
1074 stop = stop == FALSE ? flag == TRUE ? TRUE : FALSE : stop;
1081 peak->background.value = _giraffe_clip_value(background, lower,
1083 peak->background.sigma = _giraffe_clip_value(dbackground, 0.,
1086 stop = stop == FALSE ? flag == TRUE ? TRUE : FALSE : stop;
1090 lower = 0.5 * initial.width;
1091 upper = 2. * (cpl_matrix_get(x, cpl_matrix_get_nrow(x) - 2, 0) -
1092 cpl_matrix_get(x, 0, 0));
1094 peak->width.value = _giraffe_clip_value(width, lower,
1096 peak->width.sigma = _giraffe_clip_value(dwidth, 0.,
1099 stop = stop == FALSE ? flag == TRUE ? TRUE : FALSE : stop;
1101 cpl_matrix_delete(x);
1102 cpl_matrix_delete(y);
1103 cpl_matrix_delete(sigma);
1106 cpl_msg_debug(fctid,
"Cross-correlation peak fit "
1107 "parameter out of bounds!");
1126 giraffe_model_delete(model);
1134 _giraffe_compute_fiber_offsets(cpl_table* offsets,
1136 const GiSGSetup* setup)
1141 const cxdouble ccdfactor = 1.1;
1143 cxdouble gcamera = 1.;
1144 cxdouble cfactor = 1.;
1145 cxdouble lincorr = 1.;
1146 cxdouble wlen0 = 0.;
1149 cx_assert(offsets != NULL);
1151 if (!cpl_table_has_column(offsets,
"WAVELENGTH")) {
1155 if (!cpl_table_has_column(offsets,
"DWF")) {
1156 cpl_table_new_column(offsets,
"DWF", CPL_TYPE_DOUBLE);
1159 if (!cpl_table_has_column(offsets,
"DXF")) {
1160 cpl_table_new_column(offsets,
"DXF", CPL_TYPE_DOUBLE);
1169 if (setup->scale == GIREBIN_SCALE_LOG) {
1170 wlen0 = 0.5 * (exp(setup->wlmin) + exp(setup->wlmax));
1173 wlen0 = 0.5 * (setup->wlmin + setup->wlmax);
1185 gcamera = 0.3894 - 5. * (1. / wlen0 - 1. / 550.) -
1186 0.00025 * pow(1. / wlen0 - 1. / 550., 2.);
1196 cfactor = (setup->nex * setup->pixelsize / 1000. * ccdfactor) /
1204 if (setup->scale == GIREBIN_SCALE_LOG) {
1208 lincorr = 0.5 * (setup->wlmin + setup->wlmax) /
1209 exp(0.5 * (log(setup->wlmin) + log(setup->wlmax)));
1217 for (i = 0; i < cpl_table_get_nrow(offsets); i++) {
1220 cxdouble dwf = cpl_table_get_double(offsets,
"WAVELENGTH", i, NULL);
1225 dxf = dwf * cfactor;
1227 cpl_table_set_double(offsets,
"DWF", i, dwf);
1228 cpl_table_set_double(offsets,
"DXF", i, dxf);
1237 inline static cpl_table*
1238 _giraffe_compute_offsets(
const GiImage* spectra,
const GiTable* mask,
1239 const cpl_table* fibers,
const GiGrating* grating,
1243 const cxchar*
const fctid =
"_giraffe_compute_offsets";
1245 const cxint dnmin = 7;
1259 const cxdouble clight = 299702.547;
1261 cxdouble cstep = 0.;
1262 cxdouble wlen0 = 0.;
1263 cxdouble nm2km = clight;
1264 cxdouble hpixels = 0.;
1272 cpl_matrix* spectrum = NULL;
1274 cpl_image* _spectra = NULL;
1275 cpl_image* tspectra = NULL;
1277 cpl_table* peakdata = NULL;
1279 GiSGMask* _mask = NULL;
1282 cx_assert(spectra != NULL);
1283 cx_assert(mask != NULL);
1284 cx_assert(fibers != NULL);
1285 cx_assert(grating != NULL);
1286 cx_assert(setup != NULL);
1287 cx_assert(config != NULL);
1299 if (setup->scale == GIREBIN_SCALE_LOG) {
1301 cxdouble wlstep = (exp(setup->wlmax) - exp(setup->wlmin)) /
1304 sampling = (cxint)(0.5 + wlstep / config->
cc_step);
1309 sampling = (cxint)(0.5 + setup->wlstep / config->
cc_step);
1315 cstep = setup->wlstep / sampling;
1322 _mask = _giraffe_sgmask_create((setup->nx - 1) * sampling + 1,
1323 setup->wlmin, cstep, setup->scale,
1326 if (_mask == NULL) {
1335 pixel0 = setup->nx / 2;
1337 if (setup->scale == GIREBIN_SCALE_LOG) {
1343 wlen0 = 0.5 * (exp(setup->wlmin) + exp(setup->wlmax));
1353 wlen0 = 0.5 * (setup->wlmin + setup->wlmax);
1354 nm2km = clight / wlen0;
1370 cpl_msg_debug(fctid,
"Cross-correlation limits: RVlow = %.4f km/s "
1371 "(%.4f nm), RVhigh = %.4f km/s (%.4f nm)", dv1, dw1,
1374 dwd = (dw2 - dw1) / 2.;
1375 dwc = (dw2 + dw1) / 2.;
1377 dnd = CX_MIN(pixel0, CX_MAX(dnmin, (cxint)(dwd / cstep + 0.5)));
1378 dnc = CX_MIN(pixel0, CX_MAX(-pixel0, (cxint)(dwc / cstep + 0.5)));
1380 dn1 = CX_MIN(pixel0 + 1, CX_MAX(-pixel0, dnc - dnd));
1381 dn2 = CX_MIN(pixel0 + 1, CX_MAX(-pixel0, dnc + dnd + 1));
1383 cpl_msg_debug(fctid,
"Cross-correlation window: center = %.4f nm "
1384 "(%d pxl) half-width = %.4f nm (%d pxl)", dwc, dnc,
1396 if (xc1 > 0 || xc2 > 0) {
1397 _giraffe_sgmask_crop(_mask, xc1, xc2);
1400 for (i = 0; (cxsize)i < _giraffe_sgmask_size(_mask); i++) {
1402 cxdouble value = _giraffe_sgmask_get_flux(_mask, i);
1410 hpixels /= _giraffe_sgmask_holes(_mask);
1420 k = CX_MAX(0, -dn1);
1422 while (i < k || _giraffe_sgmask_get_flux(_mask, i) > 0.) {
1424 _giraffe_sgmask_set_flux(_mask, i, 0.);
1429 cpl_msg_debug(fctid,
"Mask cleared from 0 to %d", i - 1);
1431 i = _giraffe_sgmask_size(_mask);
1432 k = _giraffe_sgmask_size(_mask) - CX_MAX(0, dn2);
1434 while (i > k || _giraffe_sgmask_get_flux(_mask, i) > 0.) {
1436 _giraffe_sgmask_set_flux(_mask, i, 0.);
1441 cpl_msg_debug(fctid,
"Mask cleared from %d to %ld", k,
1442 (
long)(_giraffe_sgmask_size(_mask) - 1));
1452 if (_spectra == NULL) {
1454 _giraffe_sgmask_delete(_mask);
1461 if (config->
zmax > 0.) {
1463 cpl_image_threshold(_spectra, CX_MINDOUBLE, config->
zmax,
1469 tspectra = _giraffe_resample_image(_spectra, setup->wlstep, cstep);
1471 if (tspectra == NULL) {
1473 cpl_image_delete(_spectra);
1475 _giraffe_sgmask_delete(_mask);
1481 cpl_image_delete(_spectra);
1484 if (xc1 > 0 || xc2 > 0) {
1486 _spectra = cpl_image_extract(tspectra, 1, xc1 + 1,
1487 cpl_image_get_size_x(tspectra), xc2 + 1);
1489 if (_spectra == NULL) {
1491 cpl_image_delete(tspectra);
1493 _giraffe_sgmask_delete(_mask);
1499 cpl_image_delete(tspectra);
1505 _spectra = tspectra;
1516 peakdata = cpl_table_new(cpl_table_get_nrow(fibers));
1518 cpl_table_duplicate_column(peakdata,
"INDEX", (cpl_table*)fibers,
1520 cpl_table_duplicate_column(peakdata,
"FPS", (cpl_table*)fibers,
1523 cpl_table_new_column(peakdata,
"WAVELENGTH", CPL_TYPE_DOUBLE);
1524 cpl_table_new_column(peakdata,
"FWHM", CPL_TYPE_DOUBLE);
1525 cpl_table_new_column(peakdata,
"AMPLITUDE", CPL_TYPE_DOUBLE);
1526 cpl_table_new_column(peakdata,
"BACKGROUND", CPL_TYPE_DOUBLE);
1527 cpl_table_new_column(peakdata,
"RV", CPL_TYPE_DOUBLE);
1528 cpl_table_new_column(peakdata,
"RVERR", CPL_TYPE_DOUBLE);
1529 cpl_table_new_column(peakdata,
"RESOLUTION", CPL_TYPE_DOUBLE);
1530 cpl_table_new_column(peakdata,
"STATUS", CPL_TYPE_INT);
1538 cpl_msg_debug(fctid,
"Computing cross-correlation: central wavelength = "
1539 "%.4f, window = [%.4f, %.4f] [km/s]", wlen0, dv1, dv2);
1541 spectrum = cpl_matrix_new(1, cpl_image_get_size_y(_spectra));
1543 for (i = 0; i < cpl_table_get_nrow(fibers); i++) {
1546 cxint ns = cpl_image_get_size_x(_spectra);
1547 cxint fiber = cpl_table_get_int(fibers,
"FPS", i, NULL);
1548 cxint idx = cpl_table_get_int(fibers,
"INDEX", i, NULL) - 1;
1550 const cxdouble fwhm_ratio = 2. * sqrt(2. * log(2.));
1552 cxdouble avsigma = 0.;
1554 cxdouble fxtotal = 0.;
1556 cxdouble fxmask = 0.;
1558 cxdouble position = 0.;
1560 cxdouble width = 0.;
1561 cxdouble resolution = 0.;
1563 cxdouble rverr = 0.;
1564 cxdouble* data = cpl_image_get_data(_spectra);
1566 const cpl_matrix*
template = NULL;
1567 cpl_matrix* ccf = NULL;
1568 cpl_matrix* lambda = NULL;
1570 GiCPFitParams peak_setup;
1580 for (j = 0; j < cpl_matrix_get_ncol(spectrum); j++) {
1582 cxdouble flux = data[j * ns + idx];
1585 cpl_matrix_set(spectrum, 0, j, flux);
1588 fxmask += _giraffe_sgmask_get_flux(_mask, j);
1589 fx += flux * _giraffe_sgmask_get_flux(_mask, j);
1597 avsigma = 1. / sqrt(fx);
1600 cpl_msg_debug(fctid,
"Cross-correlation of spectrum %d in window "
1601 "from %d pxl to %d pxl (%.4f nm to %.4f nm)", fiber,
1602 dn1, dn2, dw1, dw2);
1609 lambda = cpl_matrix_new(1, dn2 - dn1);
1611 for (j = dn1; j < dn2; j++) {
1612 cpl_matrix_set(lambda, 0, j - dn1, j * cstep);
1620 template = _giraffe_sgmask_get(_mask);
1622 ccf = _giraffe_compute_cross_correlation(spectrum,
template, dn1, dn2);
1626 cpl_matrix_delete(lambda);
1627 cpl_matrix_delete(spectrum);
1629 cpl_image_delete(_spectra);
1631 cpl_table_delete(peakdata);
1633 _giraffe_sgmask_delete(_mask);
1641 for (j = 0; j < cpl_matrix_get_ncol(ccf); j++) {
1642 sum += cpl_matrix_get(ccf, 0, j);
1646 cpl_msg_debug(fctid,
"Cross-correlation failed: Skipping "
1647 "spectrum %d.", fiber);
1649 cpl_matrix_delete(lambda);
1660 peak_setup.dnmin = dnmin;
1661 peak_setup.iterations = config->
rv_niter;
1662 peak_setup.step = cstep;
1664 peak_setup.sigma = avsigma;
1665 peak_setup.scale = setup->scale;
1667 peak_setup.fit.iterations = config->
pf_niter;
1668 peak_setup.fit.tests = config->
pf_ntest;
1669 peak_setup.fit.delta = config->
pf_dchisq;
1671 status = _giraffe_peak_fit(&peak, lambda, ccf, grating, &peak_setup);
1675 cpl_matrix_delete(ccf);
1676 cpl_matrix_delete(lambda);
1678 cpl_matrix_delete(spectrum);
1679 cpl_image_delete(_spectra);
1681 cpl_table_delete(peakdata);
1683 _giraffe_sgmask_delete(_mask);
1694 if (setup->scale == GIREBIN_SCALE_LOG) {
1695 position = peak.center.value * wlen0;
1696 fwhm = (exp(peak.width.value) - 1.) * wlen0;
1699 position = peak.center.value;
1700 fwhm = peak.width.value;
1703 width = pow(fwhm_ratio * fwhm, 2.) - pow(0.6 * hpixels * cstep, 2.);
1704 resolution = width > 0. ? wlen0 / sqrt(width) : 0.;
1708 rv = CX_MAX(dv1, CX_MIN(dv2, peak.center.value * nm2km));
1709 rverr = CX_MIN(dv2 - dv1, peak.center.sigma * nm2km);
1711 cpl_table_set_double(peakdata,
"WAVELENGTH", i, position);
1712 cpl_table_set_double(peakdata,
"FWHM", i, fwhm);
1713 cpl_table_set_double(peakdata,
"AMPLITUDE", i, peak.amplitude.value);
1714 cpl_table_set_double(peakdata,
"BACKGROUND", i,
1715 peak.background.value);
1716 cpl_table_set_double(peakdata,
"RESOLUTION", i, resolution);
1717 cpl_table_set_double(peakdata,
"RV", i, rv);
1718 cpl_table_set_double(peakdata,
"RVERR", i, rverr);
1719 cpl_table_set_int(peakdata,
"STATUS", i, peak.status);
1721 cpl_matrix_delete(lambda);
1722 cpl_matrix_delete(ccf);
1726 cpl_matrix_delete(spectrum);
1727 cpl_image_delete(_spectra);
1729 _giraffe_sgmask_delete(_mask);
1736 inline static cpl_table*
1737 _giraffe_compute_slitgeometry(
const GiImage* spectra,
const GiTable* mask,
1738 const GiTable* slitgeometry,
1746 cpl_table* peakdata = NULL;
1755 status = _giraffe_create_setup(&setup, spectra);
1773 peakdata = _giraffe_compute_offsets(spectra, mask, _slitgeometry,
1774 grating, &setup, config);
1776 if (peakdata == NULL) {
1786 status = _giraffe_compute_fiber_offsets(peakdata, grating, &setup);
1789 cpl_table_delete(peakdata);
1798 cpl_table_duplicate_column(peakdata,
"XF", _slitgeometry,
"XF");
1799 cpl_table_add_columns(peakdata,
"XF",
"DXF");
1806 inline static GiTable*
1807 _giraffe_slitgeometry_table(
const cpl_table* offsets,
1808 const GiImage* spectra,
1809 const GiTable* fibers,
1810 const GiTable* slitgeometry,
1814 const cxchar* idx = NULL;
1818 cpl_propertylist* properties = NULL;
1819 cpl_propertylist* _properties = NULL;
1821 cpl_table* _slit = NULL;
1822 cpl_table* _fibers = NULL;
1823 cpl_table* _slitgeometry = NULL;
1825 GiTable* slit = NULL;
1828 cx_assert(spectra != NULL);
1829 cx_assert(fibers != NULL);
1832 cx_assert(_fibers != NULL);
1835 cx_assert(_slitgeometry != NULL);
1837 if (offsets == NULL) {
1845 cx_assert(properties != NULL);
1847 giraffe_error_push();
1849 _properties = cpl_propertylist_new();
1852 GIALIAS_INSTRUMENT);
1884 cpl_propertylist_update_double(_properties, GIALIAS_SCAL_CUTOFF,
1886 cpl_propertylist_set_comment(_properties, GIALIAS_SCAL_CUTOFF,
1887 "Cutoff pixel value.");
1889 cpl_propertylist_update_string(_properties, GIALIAS_GIRFTYPE,
"SLITGEOTAB");
1890 cpl_propertylist_set_comment(_properties, GIALIAS_GIRFTYPE,
1891 "Giraffe frame type.");
1894 _slit = cpl_table_new(cpl_table_get_nrow(_fibers));
1896 cpl_table_new_column(_slit,
"INDEX", CPL_TYPE_INT);
1898 for (i = 0; i < cpl_table_get_nrow(_slit); i++) {
1899 cpl_table_set_int(_slit,
"INDEX", i, i + 1);
1902 cpl_table_duplicate_column(_slit,
"FPS", (cpl_table*)_fibers,
"FPS");
1903 cpl_table_duplicate_column(_slit,
"SSN", (cpl_table*)_fibers,
"SSN");
1904 cpl_table_duplicate_column(_slit,
"XF", (cpl_table*)offsets,
"XF");
1905 cpl_table_duplicate_column(_slit,
"YF", (cpl_table*)_slitgeometry,
"YF");
1907 if (cpl_table_has_column(_slitgeometry,
"ZF")) {
1908 cpl_table_duplicate_column(_slit,
"ZF",
1909 (cpl_table*)_slitgeometry,
"ZF");
1912 if (cpl_table_has_column(_slitgeometry,
"ZDEFOCUS")) {
1913 cpl_table_duplicate_column(_slit,
"ZDEFOCUS",
1914 (cpl_table*)_slitgeometry,
"ZDEFOCUS");
1917 cpl_table_duplicate_column(_slit,
"RV", (cpl_table*)offsets,
"RV");
1918 cpl_table_duplicate_column(_slit,
"RVERR", (cpl_table*)offsets,
"RVERR");
1919 cpl_table_duplicate_column(_slit,
"RESOLUTION", (cpl_table*)offsets,
1923 cpl_table_duplicate_column(_slit,
"RINDEX", _fibers, idx);
1925 if (cpl_error_get_code() != CPL_ERROR_NONE) {
1926 cpl_propertylist_delete(_properties);
1927 cpl_table_delete(_slit);
1934 giraffe_error_pop();
1937 cpl_propertylist_delete(_properties);
1940 cpl_table_delete(_slit);
1954 const GiLocalization* localization,
1955 const GiTable* fibers,
const GiTable* wlsolution,
1956 const GiTable* slitgeometry,
const GiTable* grating,
1960 const cxchar*
const fctid =
"giraffe_calibrate_slit";
1965 cpl_table* _fibers = NULL;
1966 cpl_table* _slitgeometry = NULL;
1967 cpl_table* offsets = NULL;
1969 GiTable* slit = NULL;
1973 GiExtraction* _extraction = NULL;
1976 if (result == NULL) {
1980 if (extraction == NULL) {
1984 if (extraction->spectra == NULL) {
1988 if (localization == NULL) {
1992 if (fibers == NULL) {
1996 if (wlsolution == NULL) {
2000 if (slitgeometry == NULL) {
2004 if (grating == NULL) {
2012 if (config == NULL) {
2018 cx_assert(_fibers != NULL);
2021 cx_assert(_slitgeometry != NULL);
2023 if (cpl_table_get_nrow(_fibers) != cpl_table_get_nrow(_slitgeometry)) {
2034 if (setup == NULL) {
2044 _extraction = giraffe_extraction_new();
2046 _extraction->spectra = extraction->spectra;
2047 _extraction->error = NULL;
2052 for (i = 0; i < config->
repeat; i++) {
2054 cxint fps_rvmin = 0;
2055 cxint fps_rvmax = 0;
2057 cxdouble rvmin = 0.;
2058 cxdouble rvmax = 0.;
2059 cxdouble rvmean = 0.;
2065 GiRebinConfig rebin_config = {
2066 GIREBIN_METHOD_LINEAR,
2069 GIREBIN_SCALE_LINEAR,
2076 localization, grating, slit,
2077 wlsolution, &rebin_config);
2082 giraffe_extraction_delete(_extraction);
2090 offsets = _giraffe_compute_slitgeometry(rebinning->spectra, mask,
2091 slit, setup, config);
2093 if (offsets == NULL) {
2096 giraffe_extraction_delete(_extraction);
2109 cx_assert(cpl_table_get_nrow(offsets) == cpl_table_get_nrow(_fibers));
2112 slit = _giraffe_slitgeometry_table(offsets, rebinning->spectra,
2113 fibers, slitgeometry, config);
2117 cpl_table_delete(offsets);
2119 giraffe_extraction_delete(_extraction);
2127 cpl_table_delete(offsets);
2142 cpl_msg_info(fctid,
"Iteration %d: Fiber offsets [km/s]: minimum = "
2143 "%.6e (fps %d), maximum = %.6e (fps %d), mean = %.6e",
2144 i + 1, rvmin, fps_rvmin, rvmax, fps_rvmax, rvmean);
2151 giraffe_extraction_delete(_extraction);
2154 cx_assert(slit != NULL);
2192 const GiTable* grating,
const GiTable* mask,
2206 cpl_table* _fibers = NULL;
2207 cpl_table* peakdata = NULL;
2214 if ((rebinning == NULL) || (rebinning->spectra == NULL)) {
2218 if (fibers == NULL) {
2222 if (grating == NULL) {
2230 if (config == NULL) {
2236 cx_assert(_fibers != NULL);
2244 cpl_table_unselect_all(_fibers);
2245 cpl_table_or_selected_int(_fibers,
"RP", CPL_EQUAL_TO, -1);
2247 _fibers = cpl_table_extract_selected(_fibers);
2249 if (_fibers == NULL) {
2260 if (_grating == NULL) {
2261 cpl_table_delete(_fibers);
2269 status = _giraffe_create_setup(&setup, rebinning->spectra);
2274 cpl_table_delete(_fibers);
2293 peakdata = _giraffe_compute_offsets(rebinning->spectra, mask, _fibers,
2294 _grating, &setup, config);
2296 if (peakdata == NULL) {
2299 cpl_table_delete(_fibers);
2311 status = _giraffe_compute_fiber_offsets(peakdata, _grating, &setup);
2314 cpl_table_delete(peakdata);
2315 cpl_table_delete(_fibers);
2321 cpl_table_delete(_fibers);
2337 giraffe_error_push();
2339 cpl_table_new_column(_fibers,
"WLRES", CPL_TYPE_DOUBLE);
2340 cpl_table_set_column_unit(_fibers,
"WLRES",
"nm");
2342 if (cpl_error_get_code() != CPL_ERROR_NONE) {
2343 cpl_table_delete(peakdata);
2347 giraffe_error_pop();
2350 giraffe_error_push();
2352 fps0 = cpl_table_get_int(peakdata,
"FPS", 0, NULL);
2353 dwf0 = cpl_table_get_double(peakdata,
"DWF", 0, NULL);
2355 nfibers = cpl_table_get_nrow(_fibers);
2357 fps = cpl_table_get_int(_fibers,
"FPS", 0, NULL);
2359 while (fps != fps0) {
2361 cpl_table_set_double(_fibers,
"WLRES", fiber, dwf0);
2364 fps = cpl_table_get_int(_fibers,
"FPS", fiber, NULL);
2368 for (peak = 1; peak < cpl_table_get_nrow(peakdata); ++peak) {
2370 cxdouble dwf1 = cpl_table_get_double(peakdata,
"DWF", peak, NULL);
2371 cxdouble slope = 0.;
2374 fps1 = cpl_table_get_int(peakdata,
"FPS", peak, NULL);
2376 slope = (dwf1 - dwf0) / ((cxdouble)(fps1 - fps0));
2378 while (fps != fps1) {
2380 cpl_table_set_double(_fibers,
"WLRES", fiber,
2381 slope * (fps - fps0) + dwf0);
2384 fps = cpl_table_get_int(_fibers,
"FPS", fiber, NULL);
2392 fps1 = cpl_table_get_int(_fibers,
"FPS", nfibers - 1, NULL);
2394 while (fps != fps1) {
2396 cpl_table_set_double(_fibers,
"WLRES", fiber, dwf0);
2399 fps = cpl_table_get_int(_fibers,
"FPS", fiber, NULL);
2403 cpl_table_set_double(_fibers,
"WLRES", fiber, dwf0);
2405 cx_assert(fiber == nfibers - 1);
2407 if (cpl_error_get_code() != CPL_ERROR_NONE) {
2408 cpl_table_delete(peakdata);
2412 cpl_table_delete(peakdata);
2414 giraffe_error_pop();
2437 const cxchar* s = NULL;
2439 cpl_parameter* p = NULL;
2448 config = cx_calloc(1,
sizeof *config);
2452 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.iterations");
2453 config->
repeat = cpl_parameter_get_int(p);
2455 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.zmax");
2456 config->
zmax = cpl_parameter_get_double(p);
2458 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.cc.step");
2459 config->
cc_step = cpl_parameter_get_double(p);
2461 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.cc.domain");
2462 s = cpl_parameter_get_string(p);
2466 cxchar** values = cx_strsplit(s,
",", 3);
2468 if (values == NULL) {
2478 cxdouble lower = 0.;
2479 cxdouble upper = 0.;
2482 lower = strtod(values[0], &last);
2484 if (*last !=
'\0') {
2486 cx_strfreev(values);
2493 lower = lower >= 0. ? lower : 0.;
2496 if (values[1] != NULL) {
2498 upper = strtod(values[1], &last);
2500 if (*last !=
'\0') {
2502 cx_strfreev(values);
2509 upper = upper > lower ? upper : 0.;
2518 cx_strfreev(values);
2524 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.rv.limits");
2525 s = cpl_parameter_get_string(p);
2529 cxchar** values = cx_strsplit(s,
",", 3);
2531 if (values == NULL) {
2541 cxdouble lower = 0.;
2542 cxdouble upper = 0.;
2545 lower = strtod(values[0], &last);
2547 if (*last !=
'\0') {
2549 cx_strfreev(values);
2556 if (values[1] != NULL) {
2558 upper = strtod(values[1], &last);
2560 if (*last !=
'\0') {
2562 cx_strfreev(values);
2569 if (lower > 0 || upper < lower) {
2571 cx_strfreev(values);
2595 cx_assert(lower <= 0);
2596 cx_assert(lower < upper);
2603 cx_strfreev(values);
2609 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.rv.iterations");
2610 config->
rv_niter = cpl_parameter_get_int(p);
2612 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.rv.wfactor");
2613 config->
rv_wfactor = cpl_parameter_get_double(p);
2615 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.peak.iterations");
2616 config->
pf_niter = cpl_parameter_get_int(p);
2618 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.peak.tests");
2619 config->
pf_ntest = cpl_parameter_get_int(p);
2621 p = cpl_parameterlist_find(list,
"giraffe.sgcalibration.peak.dchisquare");
2622 config->
pf_dchisq = cpl_parameter_get_double(p);
2685 p = cpl_parameter_new_value(
"giraffe.sgcalibration.iterations",
2687 "Slit geometry calibration maximum number "
2689 "giraffe.sgcalibration",
2691 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-cniter");
2692 cpl_parameterlist_append(list, p);
2695 p = cpl_parameter_new_value(
"giraffe.sgcalibration.zmax",
2697 "Maximum allowed pixel value. To be "
2698 "effective it must be larger than 0.",
2699 "giraffe.sgcalibration",
2701 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-zmax");
2702 cpl_parameterlist_append(list, p);
2705 p = cpl_parameter_new_value(
"giraffe.sgcalibration.cc.step",
2707 "Cross-correlation step.",
2708 "giraffe.sgcalibration",
2710 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-cstep");
2711 cpl_parameterlist_append(list, p);
2714 p = cpl_parameter_new_value(
"giraffe.sgcalibration.cc.domain",
2716 "Restricts the cross-correlation to the "
2718 "giraffe.sgcalibration",
2720 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-cdomain");
2721 cpl_parameterlist_append(list, p);
2724 p = cpl_parameter_new_value(
"giraffe.sgcalibration.rv.limits",
2726 "Delta RV limits of the cross-correlation "
2728 "giraffe.sgcalibration",
2730 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-rvlimits");
2731 cpl_parameterlist_append(list, p);
2734 p = cpl_parameter_new_value(
"giraffe.sgcalibration.rv.iterations",
2736 "Maximum number of iterations used for the "
2737 "RV determination.",
2738 "giraffe.sgcalibration",
2740 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-rvniter");
2741 cpl_parameterlist_append(list, p);
2744 p = cpl_parameter_new_value(
"giraffe.sgcalibration.rv.wfactor",
2746 "Data window width factor. The FWHM times "
2747 "this value determines the window width.",
2748 "giraffe.sgcalibration",
2750 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-rvwfactor");
2751 cpl_parameterlist_append(list, p);
2754 p = cpl_parameter_new_value(
"giraffe.sgcalibration.peak.iterations",
2756 "Peak model fit maximum number of "
2758 "giraffe.sgcalibration",
2761 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-pfniter");
2762 cpl_parameterlist_append(list, p);
2765 p = cpl_parameter_new_value(
"giraffe.sgcalibration.peak.tests",
2767 "Cross-correlation peak fit maximum number "
2769 "giraffe.sgcalibration",
2772 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-pfntest");
2773 cpl_parameterlist_append(list, p);
2776 p = cpl_parameter_new_value(
"giraffe.sgcalibration.peak.dchisquare",
2778 "Cross-correlation peak fit minimum "
2779 "chi-square difference.",
2780 "giraffe.sgcalibration",
2783 cpl_parameter_set_alias(p, CPL_PARAMETER_MODE_CLI,
"scal-pfdchisq");
2784 cpl_parameterlist_append(list, p);
const cxchar * giraffe_fiberlist_query_index(const cpl_table *fibers)
Query a fiber list for the name of the fiber reference index column.
void giraffe_grating_delete(GiGrating *self)
Destroys an GiGrating object.
GiGrating * giraffe_grating_create(const GiImage *spectra, const GiTable *grating)
Create a GiGrating from a reference image.
cpl_propertylist * giraffe_image_get_properties(const GiImage *self)
Get the properties of an image.
cpl_image * giraffe_image_get(const GiImage *self)
Gets the image data.
cxint giraffe_matrix_sort(cpl_matrix *mA)
Sort in place the matrix elements in ascending order.
GiRange * giraffe_range_create(cxdouble min, cxdouble max)
Creates a new range from the given minimum and maximum values.
void giraffe_range_delete(GiRange *self)
Destroys a range object.
cxdouble giraffe_range_get_min(const GiRange *const self)
Get the minimum of a range.
cxdouble giraffe_range_get_max(const GiRange *const self)
Get the maximum of a range.
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.
void giraffe_rebinning_destroy(GiRebinning *rebinning)
Destroys a rebinning results container and its contents.
GiRebinning * giraffe_rebinning_new(void)
Create an empty rebinning results container.
GiSGCalConfig * giraffe_sgcalibration_config_create(cpl_parameterlist *list)
Creates a setup structure for the slit geometry calibration.
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.
void giraffe_sgcalibration_config_destroy(GiSGCalConfig *config)
Destroys a sgcalibration field setup structure.
void giraffe_sgcalibration_config_add(cpl_parameterlist *list)
Adds parameters for the sgcalibration correction computation.
cxint giraffe_calibrate_slit(GiTable *result, const GiExtraction *extraction, const GiLocalization *localization, const GiTable *fibers, const GiTable *wlsolution, const GiTable *slitgeometry, const GiTable *grating, const GiTable *mask, const GiSGCalConfig *config)
Compute a slit geometry corresponding to the given rebinned spectrum.
cxint giraffe_table_set(GiTable *self, cpl_table *table)
Sets the table data.
cpl_propertylist * giraffe_table_get_properties(const GiTable *self)
Gets the table properties.
GiTable * giraffe_table_new(void)
Creates a new, empty Giraffe table.
void giraffe_table_delete(GiTable *self)
Destroys a Giraffe table.
cpl_table * giraffe_table_get(const GiTable *self)
Get the table data from a Giraffe table.
GiTable * giraffe_table_duplicate(const GiTable *src)
Duplicate a Giraffe table.
cxint giraffe_table_set_properties(GiTable *self, cpl_propertylist *properties)
Attaches a property list to an table.
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.
Structure to handle Grating Information.
Slit geometry calibration configuration data structure.