33 #include "muse_pixtable.h" 34 #include "muse_instrument.h" 36 #include "muse_cplwrappers.h" 38 #include "muse_mask.h" 39 #include "muse_pfits.h" 40 #include "muse_quadrants.h" 41 #include "muse_quality.h" 42 #include "muse_tracing.h" 43 #include "muse_wavecalib.h" 45 #include "muse_utils.h" 51 #define CREATE_MINIMAL_PIXTABLE 0 52 #define PIXTABLE_CREATE_CCDSIZED 1 53 #define DEBUG_PIXTABLE_CREATION 0 54 #define DEBUG_PIXTABLE_FEW_SLICES 0 100 #define MUSE_ORIGIN_SHIFT_XSLICE 24 101 #define MUSE_ORIGIN_SHIFT_YPIX 11 102 #define MUSE_ORIGIN_SHIFT_IFU 6 106 #define MUSE_ORIGIN_SLICE_SAFETY_OFFSET -20 110 static inline uint32_t
111 muse_pixtable_origin_encode_fast(
unsigned int aX,
unsigned int aY,
113 unsigned short aSlice,
unsigned int aOffset)
115 return ((aX - aOffset) << MUSE_ORIGIN_SHIFT_XSLICE)
116 | (aY << MUSE_ORIGIN_SHIFT_YPIX)
117 | (aIFU << MUSE_ORIGIN_SHIFT_IFU)
122 static inline unsigned int 123 muse_pixtable_origin_get_x_fast(uint32_t aOrigin, uint32_t aOffset)
125 return ((aOrigin >> MUSE_ORIGIN_SHIFT_XSLICE) & 0x7f) + aOffset;
129 static inline unsigned int 130 muse_pixtable_origin_get_y_fast(uint32_t aOrigin)
132 return (aOrigin >> MUSE_ORIGIN_SHIFT_YPIX) & 0x1fff;
136 static inline unsigned short 137 muse_pixtable_origin_get_ifu_fast(uint32_t aOrigin)
139 return (aOrigin >> MUSE_ORIGIN_SHIFT_IFU) & 0x1f;
143 static inline unsigned short 144 muse_pixtable_origin_get_slice_fast(uint32_t aOrigin)
146 return aOrigin & 0x3f;
169 unsigned short aSlice,
unsigned int aOffset)
174 cpl_ensure(aX < 8192 && aX > 0 && aY < 8192 && aY > 0 &&
175 aIFU <= kMuseNumIFUs && aIFU >= 1 &&
176 aSlice <= kMuseSlicesPerCCD && aSlice >= 1 && aOffset < 8192,
177 CPL_ERROR_ILLEGAL_INPUT, 0);
180 cpl_msg_debug(__func__,
"origin (%d, %d, %d, %d, %d) = 0x%x",
181 aX, aY, aIFU, aSlice, aOffset,
182 ((aX - aOffset) << MUSE_ORIGIN_SHIFT_XSLICE)
183 | (aY << MUSE_ORIGIN_SHIFT_YPIX)
184 | (aIFU << MUSE_ORIGIN_SHIFT_IFU)
189 return muse_pixtable_origin_encode_fast(aX, aY, aIFU, aSlice, aOffset);
215 unsigned short slice = muse_pixtable_origin_get_slice_fast(aOrigin),
216 ifu = muse_pixtable_origin_get_ifu_fast(aOrigin);
217 cpl_errorstate prestate = cpl_errorstate_get();
219 if (!cpl_errorstate_is_equal(prestate)) {
220 cpl_errorstate_set(prestate);
224 x = muse_pixtable_origin_get_x_fast(aOrigin, offset);
226 if (x > 8191 || x < 1 || !cpl_errorstate_is_equal(prestate)) {
227 cpl_msg_error(__func__,
"aOrigin=%#x x=%d (%d %d %d), %s",
228 aOrigin, x, slice, ifu, offset, cpl_error_get_message());
231 cpl_ensure(x <= 8191 && x >= 1 && cpl_errorstate_is_equal(prestate),
232 CPL_ERROR_ILLEGAL_OUTPUT, 0);
248 unsigned int y = muse_pixtable_origin_get_y_fast(aOrigin);
250 if (y > 8191 || y < 1) {
251 cpl_msg_error(__func__,
"aOrigin=%#x y=%d", aOrigin, y);
254 cpl_ensure(y <= 8191 && y >= 1, CPL_ERROR_ILLEGAL_OUTPUT, 0);
270 unsigned short ifu = muse_pixtable_origin_get_ifu_fast(aOrigin);
272 if (ifu > kMuseNumIFUs || ifu < 1) {
273 cpl_msg_error(__func__,
"aOrigin=%#x ifu=%d", aOrigin, ifu);
276 cpl_ensure(ifu <= kMuseNumIFUs && ifu >= 1, CPL_ERROR_ILLEGAL_OUTPUT, 0);
292 unsigned short slice = muse_pixtable_origin_get_slice_fast(aOrigin);
294 if (slice > kMuseSlicesPerCCD || slice < 1) {
295 cpl_msg_error(__func__,
"aOrigin=%#x slice=%d", aOrigin, slice);
298 cpl_ensure(slice <= kMuseSlicesPerCCD && slice >= 1,
299 CPL_ERROR_ILLEGAL_OUTPUT, 0);
321 cpl_polynomial *aLTrace,
322 unsigned short aIFU,
unsigned short aSlice)
324 cpl_ensure(aPixtable && aPixtable->
header, CPL_ERROR_NULL_INPUT, 0);
325 cpl_errorstate prestate = cpl_errorstate_get();
326 unsigned int offset = floor(cpl_polynomial_eval_1d(aLTrace, 1, NULL))
327 + MUSE_ORIGIN_SLICE_SAFETY_OFFSET;
328 cpl_ensure(cpl_errorstate_is_equal(prestate), cpl_error_get_code(), 0);
332 cpl_propertylist_update_int(aPixtable->
header, keyword, offset);
333 cpl_propertylist_set_comment(aPixtable->
header, keyword,
334 MUSE_HDR_PT_IFU_SLICE_OFFSET_COMMENT);
356 unsigned short aIFU,
unsigned short aSlice)
358 cpl_ensure(aPixtable && aPixtable->
header, CPL_ERROR_NULL_INPUT, 0);
361 cpl_errorstate prestate = cpl_errorstate_get();
362 unsigned int offset = cpl_propertylist_get_int(aPixtable->
header, keyword);
365 if (offset > 8191 || offset < 1) {
366 cpl_msg_error(__func__,
"aIFU=%d aSlice=%d offset=%d",
367 aIFU, aSlice, offset);
370 cpl_ensure(offset <= 8191 && offset >= 1 && cpl_errorstate_is_equal(prestate),
371 CPL_ERROR_ILLEGAL_OUTPUT, 0);
400 cpl_ensure_code(aOut && aOut->
header, CPL_ERROR_NULL_INPUT);
401 cpl_propertylist *dest = aOut->
header,
403 if (aFrom && aFrom->
header) {
406 char keyword[KEYWORD_LENGTH];
407 unsigned short nifu, nslice;
408 for (nifu = 1; nifu <= kMuseNumIFUs; nifu++) {
409 for (nslice = 1; nslice <= kMuseSlicesPerCCD; nslice++) {
413 cpl_errorstate prestate = cpl_errorstate_get();
414 unsigned int offset = cpl_propertylist_get_int(from, keyword);
415 if (!cpl_errorstate_is_equal(prestate)) {
417 cpl_errorstate_set(prestate);
421 cpl_propertylist_erase(from, keyword);
425 cpl_propertylist_update_int(dest, keyword, offset);
426 cpl_propertylist_set_comment(dest, keyword,
427 MUSE_HDR_PT_IFU_SLICE_OFFSET_COMMENT);
431 return CPL_ERROR_NONE;
453 cpl_ensure(aPixtable && aPixtable->
header, CPL_ERROR_NULL_INPUT, 0);
455 CPL_ERROR_ILLEGAL_INPUT, 0);
457 char keyword[KEYWORD_LENGTH];
458 unsigned int exposure = 0;
459 cpl_size lo = 0, hi = 0;
461 cpl_errorstate prestate = cpl_errorstate_get();
463 lo = cpl_propertylist_get_long_long(aPixtable->
header, keyword);
465 hi = cpl_propertylist_get_long_long(aPixtable->
header, keyword);
466 if (!cpl_errorstate_is_equal(prestate)) {
473 cpl_errorstate_set(prestate);
477 cpl_ensure(lo <= aRow && hi >= aRow, CPL_ERROR_ILLEGAL_OUTPUT, 0);
504 "relative x-pixel position in the output datacube", CPL_TRUE},
506 "relative y-pixel position in the output datacube", CPL_TRUE},
508 "wavelength of this pixel", CPL_TRUE},
510 "data value in this pixel", CPL_TRUE},
512 "Euro3D bad pixel status of this pixel", CPL_TRUE},
514 "data variance of this pixel", CPL_TRUE},
516 "encoded value of IFU and slice number, as well as x and y " 517 "position in the raw (trimmed) data", CPL_TRUE},
518 { NULL, 0, NULL, NULL, NULL, CPL_FALSE }
570 cpl_table *aGeoTable)
572 cpl_ensure(aImage && aWave && aTrace, CPL_ERROR_NULL_INPUT, NULL);
575 cpl_ensure(ifu >= 1 && ifu <= kMuseNumIFUs, CPL_ERROR_DATA_NOT_FOUND, NULL);
578 if (!cpl_propertylist_has(aImage->
header,
"BUNIT")) {
579 char *msg = cpl_sprintf(
"Input data of IFU %hhu does not contain a data " 580 "unit (\"BUNIT\"), cannot proceed!", ifu);
581 cpl_msg_error(__func__,
"%s", msg);
582 cpl_error_set_message(__func__, CPL_ERROR_INCOMPATIBLE_INPUT,
"%s", msg);
587 if (strncmp(unit,
"count", 6)) {
588 char *msg = cpl_sprintf(
"Input data of IFU %hhu is not in \"count\" units " 589 "but in \"%s\", not suitable for a new pixel table!",
591 cpl_msg_error(__func__,
"%s", msg);
592 cpl_error_set_message(__func__, CPL_ERROR_INCOMPATIBLE_INPUT,
"%s", msg);
600 int xsize = cpl_image_get_size_x(aImage->
data),
601 ysize = cpl_image_get_size_y(aImage->
data),
602 #if PIXTABLE_CREATE_CCDSIZED 605 isize = xsize * ysize;
609 isize = ysize * kMuseSlicesPerCCD * kMuseSliceHiLikelyWidth;
615 pt->
header = cpl_propertylist_duplicate(aImage->
header);
616 cpl_propertylist_erase_regexp(pt->
header,
617 "^SIMPLE$|^BITPIX$|^NAXIS|^EXTEND$|^XTENSION$|" 618 "^DATASUM$|^DATAMIN$|^DATAMAX$|^DATAMD5$|" 619 "^PCOUNT$|^GCOUNT$|^HDUVERS$|^BLANK$|" 620 "^BZERO$|^BSCALE$|^BUNIT$|^CHECKSUM$|^INHERIT$|" 623 aGeoTable ? MUSE_PIXTABLE_STRING_FULL
624 : MUSE_PIXTABLE_STRING_SIMPLE);
626 aGeoTable ? MUSE_PIXTABLE_COMMENT_FULL
627 : MUSE_PIXTABLE_COMMENT_SIMPLE);
647 uint32_t *cdata_origin = (uint32_t *)cpl_table_get_data_int(pt->
table,
651 const float *pixdata = cpl_image_get_data_float_const(aImage->
data),
652 *pixstat = cpl_image_get_data_float_const(aImage->
stat);
653 const int *pixdq = cpl_image_get_data_int_const(aImage->
dq);
656 unsigned short wavexorder, waveyorder;
658 cpl_msg_info(__func__,
"Creating pixel table for IFU %hhu, using order %d for " 659 "trace solution and orders %hu/%hu for wavelength solution", ifu,
663 cpl_table *geopos = NULL;
664 double *slice_x = NULL, *slice_y = NULL, *slice_angle = NULL,
668 slice_x = cpl_table_get_data_double(geopos, MUSE_GEOTABLE_X);
669 slice_y = cpl_table_get_data_double(geopos, MUSE_GEOTABLE_Y);
671 if (!geopos || !slice_x || !slice_y) {
673 ? cpl_sprintf(
"Geometry table is missing data for IFU %hhu!", ifu)
674 : cpl_sprintf(
"Geometry table is missing column%s%s%s!",
675 !slice_x && !slice_y ?
"s" :
"",
676 slice_x ?
"" :
" \""MUSE_GEOTABLE_X
"\"",
677 slice_y ?
"" :
" \""MUSE_GEOTABLE_Y
"\"");
678 cpl_error_set_message(__func__, CPL_ERROR_ILLEGAL_INPUT,
"%s", msg);
679 cpl_msg_error(__func__,
"%s", msg);
681 cpl_table_delete(geopos);
687 cpl_errorstate prestate = cpl_errorstate_get();
688 slice_angle = cpl_table_get_data_double(geopos, MUSE_GEOTABLE_ANGLE);
689 slice_width = cpl_table_get_data_double(geopos, MUSE_GEOTABLE_WIDTH);
690 if (!cpl_errorstate_is_equal(prestate)) {
691 char *msg = cpl_sprintf(
"Geometry table is missing column%s%s%s!",
692 !slice_angle && !slice_width ?
"s" :
"",
693 slice_angle ?
"" :
" \""MUSE_GEOTABLE_ANGLE
"\"",
694 slice_width ?
"" :
" \""MUSE_GEOTABLE_WIDTH
"\"");
695 cpl_error_set_message(__func__, CPL_ERROR_BAD_FILE_FORMAT,
"%s", msg);
696 cpl_msg_error(__func__,
"%s", msg);
698 if (!getenv(
"MUSE_EXPERT_USER")) {
699 cpl_table_delete(geopos);
707 cpl_size itablerow = 0;
709 double cputime = cpl_test_get_cputime(),
710 walltime = cpl_test_get_walltime();
713 #if DEBUG_PIXTABLE_FEW_SLICES > 0 714 # define SLICE_LIMIT DEBUG_PIXTABLE_FEW_SLICES 716 # define SLICE_LIMIT kMuseSlicesPerCCD 718 for (islice = 0; islice < SLICE_LIMIT; islice++) {
719 #if DEBUG_PIXTABLE_CREATION 720 cpl_msg_debug(__func__,
"Starting to process slice %2d of IFU %2hhu",
728 cpl_vector *pos = cpl_vector_new(2);
731 cpl_vector_set(pos, 0, kMuseOutputYTop / 2);
732 cpl_vector_set(pos, 1, kMuseOutputYTop / 2);
733 double ltest = cpl_polynomial_eval(pwave, pos);
734 if (!pwave || !isnormal(ltest) || fabs(ltest) < DBL_EPSILON) {
735 char *msg = cpl_sprintf(
"Wavelength calibration polynomial for slice %d " 736 "of IFU %hhu is not well defined!", islice + 1,
738 cpl_error_set_message(__func__, CPL_ERROR_ILLEGAL_INPUT,
"%s", msg);
739 cpl_msg_error(__func__,
"%s", msg);
741 cpl_polynomial_delete(pwave);
742 cpl_vector_delete(pos);
743 if (getenv(
"MUSE_EXPERT_USER")) {
746 cpl_table_delete(geopos);
756 char *msg = cpl_sprintf(
"Tracing polynomials for slice %d of IFU %hhu are" 757 " not well defined!", islice + 1, ifu);
758 cpl_error_set_message(__func__, CPL_ERROR_ILLEGAL_INPUT,
"%s", msg);
759 cpl_msg_error(__func__,
"%s", msg);
761 cpl_polynomial_delete(pwave);
762 cpl_vector_delete(pos);
763 if (getenv(
"MUSE_EXPERT_USER")) {
766 cpl_table_delete(geopos);
772 ptrace[MUSE_TRACE_LEFT],
777 for (j = 1; j <= ysize; j++) {
780 double x1 = cpl_polynomial_eval_1d(ptrace[MUSE_TRACE_LEFT], j, NULL),
781 x2 = cpl_polynomial_eval_1d(ptrace[MUSE_TRACE_RIGHT], j, NULL),
782 xcenter = (x1 + x2) / 2.,
784 if (!isnormal(x1) || !isnormal(x2) || x1 < 1 || x2 > xsize || x1 > x2) {
785 cpl_msg_warning(__func__,
"slice %2d of IFU %hhu: faulty polynomial " 786 "detected at y=%d (borders: %f ... %f)", islice + 1,
791 int ileft = ceil(x1),
793 cpl_vector_set(pos, 1, j);
797 for (i = ileft; i <= iright; i++) {
798 cpl_vector_set(pos, 0, i);
802 double lambda = cpl_polynomial_eval(pwave, pos);
803 if (lambda < 3000. || lambda > 11000. || !isfinite(lambda)) {
814 float dx = (xcenter - i) / width
815 * (slice_width ? slice_width[islice] : kMuseSliceNominalWidth);
816 #if DEBUG_PIXTABLE_CREATION 817 if (fabs(dx) > 37.5) {
818 cpl_msg_debug(__func__,
"%d - %f -> %f", i, xcenter, dx);
822 #if CREATE_MINIMAL_PIXTABLE 824 if (lambda < 6495. || lambda > 6505. || fabs(dx) > 10.) {
830 cpl_size irow = itablerow++;
831 #if !PIXTABLE_CREATE_CCDSIZED 833 #if DEBUG_PIXTABLE_CREATION 834 printf(
"pixel table before enlargement:\n");
835 cpl_table_dump(pt->
table, irow-5, 10, stdout);
841 const cpl_size nr = 1000000;
842 cpl_msg_debug(__func__,
"expand table to %"CPL_SIZE_FORMAT
" rows " 843 "(add %"CPL_SIZE_FORMAT
" lines)!", irow + nr, nr);
844 cpl_table_set_size(pt->
table, irow + nr);
862 cdata_origin = (uint32_t *)cpl_table_get_data_int(pt->
table,
864 #if DEBUG_PIXTABLE_CREATION 865 printf(
"pixel table after filling new rows:\n");
866 cpl_table_dump(pt->
table, irow-5, 10, stdout);
873 double angle = slice_angle ? slice_angle[islice] * CPL_MATH_RAD_DEG
875 cdata_xpos[irow] = aGeoTable
876 ? dx * cos(angle) + slice_x[islice]
877 : dx + islice * kMuseSliceHiLikelyWidth
878 + kMuseSliceHiLikelyWidth / 2.;
879 cdata_ypos[irow] = aGeoTable
880 ? dx * sin(angle) + slice_y[islice]
882 #if DEBUG_PIXTABLE_CREATION 883 if (!isfinite(cdata_xpos[irow]) ||
884 cdata_xpos[irow] < -200 || cdata_xpos[irow] > 4000) {
885 cpl_msg_debug(__func__,
"weird data in x: %e %e",
886 cdata_xpos[irow], cdata_ypos[irow]);
888 if (!isfinite(cdata_ypos[irow]) ||
889 cdata_ypos[irow] < -200 || cdata_ypos[irow] > 200) {
890 cpl_msg_debug(__func__,
"weird data in y: %e %e",
891 cdata_xpos[irow], cdata_ypos[irow]);
894 cdata_lambda[irow] = lambda;
897 cdata_data[irow] = pixdata[(i-1) + (j-1)*xsize];
898 cdata_dq[irow] = pixdq[(i-1) + (j-1)*xsize];
899 cdata_stat[irow] = pixstat[(i-1) + (j-1)*xsize];
901 cdata_origin[irow] = muse_pixtable_origin_encode_fast(i, j, ifu,
909 cpl_polynomial_delete(pwave);
910 cpl_vector_delete(pos);
912 cpl_table_delete(geopos);
913 cpl_msg_debug(__func__,
"IFU %hhu took %gs (CPU time), %gs (wall-clock)",
914 ifu, cpl_test_get_cputime() - cputime,
915 cpl_test_get_walltime() - walltime);
919 cpl_msg_debug(__func__,
"Trimming pixel table of IFU %hhu to %"CPL_SIZE_FORMAT
920 " of %"CPL_SIZE_FORMAT
" rows", ifu, itablerow,
922 #if DEBUG_PIXTABLE_CREATION 923 printf(
"end of used part of pixel table before trimming:\n");
924 cpl_table_dump(pt->
table, itablerow-5, 10, stdout);
927 cpl_table_set_size(pt->
table, itablerow);
932 #if DEBUG_PIXTABLE_CREATION 933 printf(
"beginning of pixel table before returning:\n");
934 cpl_table_dump(pt->
table, 0, 5, stdout);
936 printf(
"end of pixel table before returning:\n");
961 if (aPixtable == NULL) {
965 pt->
table = cpl_table_duplicate(aPixtable->
table);
966 pt->
header = cpl_propertylist_duplicate(aPixtable->
header);
989 cpl_table_delete(aPixtable->
table);
990 aPixtable->
table = NULL;
993 cpl_propertylist_delete(aPixtable->
header);
997 cpl_table_delete(aPixtable->
ffspec);
1000 cpl_free(aPixtable);
1031 cpl_table *aTrace, cpl_table *aWave,
float aSampling)
1033 cpl_ensure_code(aPixtable, CPL_ERROR_NULL_INPUT);
1037 cpl_table_delete(aPixtable->
ffspec);
1043 return cpl_error_get_code();
1052 cpl_table_erase_column(aPixtable->
ffspec,
"stat");
1053 cpl_table_erase_column(aPixtable->
ffspec,
"dq");
1055 return CPL_ERROR_NONE;
1071 static cpl_error_code
1072 muse_pixtable_save_ffspec(
muse_pixtable *aPixtable,
const char *aFilename)
1074 cpl_ensure_code(aPixtable, CPL_ERROR_NULL_INPUT);
1077 if (!aPixtable->
ffspec) {
1080 return CPL_ERROR_NONE;
1084 cpl_propertylist *hext = cpl_propertylist_new();
1086 cpl_error_code rc = cpl_table_save(aPixtable->
ffspec, NULL, hext, aFilename,
1088 cpl_propertylist_delete(hext);
1115 static cpl_error_code
1116 muse_pixtable_save_image(
muse_pixtable *aPixtable,
const char *aFilename)
1118 const char *
id =
"muse_pixtable_save";
1119 cpl_ensure_code(aPixtable, CPL_ERROR_NULL_INPUT);
1121 cpl_ensure_code(nrow > 0, CPL_ERROR_ILLEGAL_INPUT);
1123 cpl_errorstate state = cpl_errorstate_get();
1124 cpl_array *columns = cpl_table_get_column_names(aPixtable->
table);
1125 int icol, ncol = cpl_array_get_size(columns);
1126 for (icol = 0; icol < ncol; icol++) {
1127 const char *colname = cpl_array_get_string(columns, icol);
1129 cpl_type type = cpl_table_get_column_type(aPixtable->
table, colname);
1130 cpl_image *image = NULL;
1132 case CPL_TYPE_FLOAT: {
1133 float *data = cpl_table_get_data_float(aPixtable->
table, colname);
1134 image = cpl_image_wrap_float(1, nrow, data);
1137 case CPL_TYPE_INT: {
1138 int *data = cpl_table_get_data_int(aPixtable->
table, colname);
1139 image = cpl_image_wrap_int(1, nrow, data);
1143 cpl_error_set_message(
id, CPL_ERROR_UNSUPPORTED_MODE,
"type \"%s\" (of " 1144 "column %s) is not supported for MUSE pixel tables",
1145 cpl_type_get_name(type), colname);
1150 cpl_propertylist *hext = cpl_propertylist_new();
1151 cpl_propertylist_append_string(hext,
"EXTNAME", colname);
1152 const char *unit = cpl_table_get_column_unit(aPixtable->
table, colname);
1154 cpl_propertylist_append_string(hext,
"BUNIT", unit);
1156 cpl_image_save(image, aFilename, CPL_TYPE_UNSPECIFIED, hext, CPL_IO_EXTEND);
1157 cpl_image_unwrap(image);
1158 cpl_propertylist_delete(hext);
1160 cpl_array_delete(columns);
1163 muse_pixtable_save_ffspec(aPixtable, aFilename);
1165 return cpl_errorstate_is_equal(state) ? CPL_ERROR_NONE : cpl_error_get_code();
1193 cpl_ensure_code(aPixtable, CPL_ERROR_NULL_INPUT);
1199 cpl_error_code rc = cpl_propertylist_save(aPixtable->
header, aFilename,
1201 if (rc != CPL_ERROR_NONE) {
1202 cpl_error_set_message(__func__, rc,
"could not save FITS header of pixel " 1203 "table \"%s\"", aFilename);
1207 cpl_boolean astable = getenv(
"MUSE_PIXTABLE_SAVE_AS_TABLE")
1208 && atoi(getenv(
"MUSE_PIXTABLE_SAVE_AS_TABLE")) > 0;
1210 cpl_msg_debug(__func__,
"Saving pixel table \"%s\" as binary table",
1212 rc = cpl_table_save(aPixtable->
table, NULL, NULL, aFilename, CPL_IO_EXTEND);
1214 muse_pixtable_save_ffspec(aPixtable, aFilename);
1216 rc = muse_pixtable_save_image(aPixtable, aFilename);
1248 muse_pixtable_load_window_image(
const char *aFilename,
1249 cpl_size aStart, cpl_size aNRows)
1251 const char *
id =
"muse_pixtable_load";
1254 cpl_size y1 = aStart + 1,
1255 y2 = aStart + aNRows;
1258 cpl_propertylist *hdata = cpl_propertylist_load(aFilename, ext);
1264 cpl_propertylist_delete(hdata);
1268 cpl_table *table = cpl_table_new(nrow);
1269 int iext, next = cpl_fits_count_extensions(aFilename);
1270 for (iext = 1; iext <= next; iext++) {
1272 cpl_propertylist *hext = cpl_propertylist_load(aFilename, iext);
1276 cpl_propertylist_delete(hext);
1281 cpl_errorstate ps = cpl_errorstate_get();
1282 cpl_image *column = cpl_image_load_window(aFilename, CPL_TYPE_UNSPECIFIED,
1283 0, iext, 1, y1, 1, y2);
1284 if (!column || !cpl_errorstate_is_equal(ps)) {
1285 cpl_image_delete(column);
1286 cpl_error_set_message(
id, cpl_error_get_code(),
"could not load extension" 1287 " %d of pixel table \"%s\"", iext, aFilename);
1288 cpl_propertylist_delete(hext);
1291 cpl_size nrows = cpl_image_get_size_x(column)
1292 * cpl_image_get_size_y(column);
1294 cpl_table_set_size(table, nrows);
1296 }
else if (nrows != nrow) {
1297 cpl_error_set_message(
id, CPL_ERROR_INCOMPATIBLE_INPUT,
"size of column " 1298 "%s does not match", colname);
1299 cpl_propertylist_delete(hext);
1300 cpl_image_delete(column);
1303 cpl_type type = cpl_image_get_type(column);
1305 case CPL_TYPE_FLOAT:
1306 cpl_table_wrap_float(table, cpl_image_unwrap(column), colname);
1309 cpl_table_wrap_int(table, cpl_image_unwrap(column), colname);
1312 cpl_error_set_message(
id, CPL_ERROR_UNSUPPORTED_MODE,
"type \"%s\" (of " 1313 "column %s) is not supported for MUSE pixel tables",
1314 cpl_type_get_name(type), colname);
1316 cpl_errorstate state = cpl_errorstate_get();
1318 if (!cpl_errorstate_is_equal(state)) {
1319 cpl_errorstate_set(state);
1322 cpl_table_set_column_unit(table, colname, unit);
1324 cpl_propertylist_delete(hext);
1345 static cpl_error_code
1346 muse_pixtable_load_ffspec(
const char *aFilename,
muse_pixtable *aPixtable)
1348 const char *
id =
"muse_pixtable_load";
1349 cpl_ensure_code(aFilename && aPixtable, CPL_ERROR_NULL_INPUT);
1355 return CPL_ERROR_NONE;
1358 cpl_errorstate es = cpl_errorstate_get();
1359 aPixtable->
ffspec = cpl_table_load(aFilename, iext, 1);
1360 if (!cpl_errorstate_is_equal(es)) {
1361 cpl_msg_warning(
id,
"Pixel table flat-field spectrum extension %s " 1362 "exists in \"%s\", but cannot be loaded!",
1366 cpl_table_delete(aPixtable->
ffspec);
1367 aPixtable->
ffspec = NULL;
1370 cpl_errorstate_set(es);
1373 return CPL_ERROR_NONE;
1404 cpl_size aStart, cpl_size aNRows)
1410 cpl_errorstate prestate = cpl_errorstate_get();
1411 pt->
header = cpl_propertylist_load(aFilename, 0);
1412 cpl_ensure(cpl_errorstate_is_equal(prestate) && pt->
header,
1413 cpl_error_get_code(), NULL);
1415 cpl_msg_error(__func__,
"unknown pixel table type found in \"%s\"", aFilename);
1422 cpl_propertylist *hext = cpl_propertylist_load(aFilename, 1);
1423 cpl_boolean asimage = !strcmp(cpl_propertylist_get_string(hext,
"XTENSION"),
1425 cpl_propertylist_delete(hext);
1429 cpl_msg_info(__func__,
"Loading pixel table \"%s\" (image format)",
1431 pt->
table = muse_pixtable_load_window_image(aFilename, aStart, aNRows);
1433 cpl_msg_info(__func__,
"Loading pixel table \"%s\" (bintable format)",
1435 pt->
table = cpl_table_load_window(aFilename, 1, 0, NULL, aStart, aNRows);
1437 if (!cpl_errorstate_is_equal(prestate) || !pt->
table) {
1438 cpl_msg_error(__func__,
"Failed to load table part of pixel table \"%s\"",
1444 if (rc != CPL_ERROR_NONE) {
1445 cpl_error_set_message(__func__, rc,
"pixel table \"%s\" does not contain " 1446 "all expected columns", aFilename);
1450 muse_pixtable_load_ffspec(aFilename, pt);
1479 cpl_errorstate prestate = cpl_errorstate_get();
1480 cpl_propertylist *theader = cpl_propertylist_load(aFilename, 1);
1481 cpl_ensure(cpl_errorstate_is_equal(prestate) && theader,
1482 cpl_error_get_code(), NULL);
1483 cpl_size nrow = cpl_propertylist_get_long_long(theader,
"NAXIS2");
1484 cpl_propertylist_delete(theader);
1511 double aLambdaMin,
double aLambdaMax)
1519 if (rc != CPL_ERROR_NONE) {
1524 cpl_msg_error(__func__,
"Pixel table contains no entries after cutting to " 1525 "%.3f..%.3f Angstrom", aLambdaMin, aLambdaMax);
1526 cpl_error_set(__func__, CPL_ERROR_DATA_NOT_FOUND);
1583 double aLambdaMin,
double aLambdaMax)
1585 cpl_ensure(aExposureList, CPL_ERROR_NULL_INPUT, NULL);
1588 if (cpl_table_has_column(aExposureList,
"00")) {
1589 const char *filename = cpl_table_get_string(aExposureList,
"00", 0);
1600 cpl_array *lambdas = cpl_array_new((kMuseLambdaMaxX - kMuseLambdaMinX)
1601 / kMuseSpectralSamplingA + 1,
1604 int j, nwl = cpl_array_get_size(lambdas);
1605 for (j = 0; j < nwl; j++) {
1606 cpl_array_set_double(lambdas, j, kMuseLambdaMinX
1607 + j * kMuseSpectralSamplingA);
1610 cpl_boolean isfirst = CPL_TRUE;
1611 double fluxlref = 0,
1613 int i, nifu = 0, nflats = 0;
1614 for (i = 1; i <= kMuseNumIFUs; i++) {
1615 char *colname = cpl_sprintf(
"%02d", i);
1616 const char *filename = cpl_table_get_string(aExposureList, colname, 0);
1619 cpl_msg_warning(__func__,
"Channel for IFU %02d is missing", i);
1626 cpl_msg_error(__func__,
"failed to load pixel table from \"%s\"",
1628 cpl_array_delete(lambdas);
1636 cpl_msg_debug(__func__,
"loaded pixel table with %"CPL_SIZE_FORMAT
" rows",
1638 isfirst = CPL_FALSE;
1639 cpl_errorstate prestate = cpl_errorstate_get();
1640 fluxsref = cpl_propertylist_get_double(pt->
header, MUSE_HDR_FLAT_FLUX_SKY);
1641 fluxlref = cpl_propertylist_get_double(pt->
header, MUSE_HDR_FLAT_FLUX_LAMP);
1642 if (fluxsref == 0. && fluxlref == 0. && !cpl_errorstate_is_equal(prestate)) {
1646 cpl_msg_debug(__func__,
"\"%s\" was previously merged (got \"%s\" when" 1647 " asking for flat-field fluxes)", filename,
1648 cpl_error_get_message());
1649 cpl_errorstate_set(prestate);
1652 if (fluxsref == 0. && fluxlref > 0. && !cpl_errorstate_is_equal(prestate)) {
1654 cpl_msg_warning(__func__,
"only found reference lamp-flat flux (%e) in " 1655 "\"%s\", flux levels may vary between IFUs!", fluxlref,
1657 cpl_errorstate_set(prestate);
1659 cpl_msg_debug(__func__,
"reference flat fluxes sky: %e lamp: %e",
1660 fluxsref, fluxlref);
1662 cpl_propertylist_erase(pt->
header, MUSE_HDR_FLAT_FLUX_SKY);
1663 cpl_propertylist_erase(pt->
header, MUSE_HDR_FLAT_FLUX_LAMP);
1671 cpl_table_delete(pt->
ffspec);
1682 cpl_errorstate state = cpl_errorstate_get();
1683 double fluxs = cpl_propertylist_get_double(onept->
header,
1684 MUSE_HDR_FLAT_FLUX_SKY),
1685 fluxl = cpl_propertylist_get_double(onept->
header,
1686 MUSE_HDR_FLAT_FLUX_LAMP),
1688 if (fluxsref > 0. && fluxs > 0.) {
1689 scale = fluxs / fluxsref;
1690 }
else if (fluxlref > 0. && fluxl > 0.) {
1691 scale = fluxl / fluxlref;
1692 if (!cpl_errorstate_is_equal(state)) {
1693 cpl_msg_warning(__func__,
"only found relative lamp-flat flux (%e) in " 1694 "\"%s\", flux levels may vary between IFUs!", fluxl,
1696 cpl_errorstate_set(state);
1708 cpl_array_add(ffspec, ffout);
1710 cpl_array_delete(ffout);
1716 cpl_msg_debug(__func__,
"big pixel table now has %"CPL_SIZE_FORMAT
" entries," 1717 " scale was %e (flat fluxes sky: %e lamp: %e)",
1727 if (nflats > 0 && nifu != nflats) {
1728 cpl_error_set_message(__func__, CPL_ERROR_ILLEGAL_INPUT,
"Only %d of %d " 1729 "pixel tables of this exposure came with a flat-field" 1730 " spectrum, cannot continue!", nflats, nifu);
1731 cpl_array_delete(lambdas);
1732 cpl_array_delete(ffspec);
1740 cpl_array_divide_scalar(ffspec, nflats);
1741 cpl_msg_debug(__func__,
"Average of flat-field spectrum: %.4f",
1742 cpl_array_get_mean(ffspec));
1746 MUSE_PIXTABLE_OPERATION_MULTIPLY);
1751 MUSE_HDR_PT_FFCORR_COMMENT);
1755 pt->
ffspec = cpl_table_new(cpl_array_get_size(lambdas));
1760 cpl_array_delete(ffspec);
1764 cpl_table_erase_invalid(pt->
ffspec);
1766 cpl_array_delete(lambdas);
1771 cpl_error_set_message(__func__, CPL_ERROR_FILE_NOT_FOUND,
1772 "None of the pixel tables could be loaded");
1780 cpl_propertylist_erase_regexp(pt->
header,
"ESO DET (CHIP|OUT) ", 0);
1781 cpl_propertylist_erase_regexp(pt->
header,
"ESO DET2 ", 0);
1785 MUSE_HDR_PT_MERGED_COMMENT);
1806 const char *type = cpl_propertylist_get_string(aPixtable->
header,
1809 cpl_msg_debug(__func__,
"pixel table type \"%s\"", type);
1814 if (!strncmp(type, MUSE_PIXTABLE_STRING_FULL,
1815 strlen(MUSE_PIXTABLE_STRING_FULL) + 1)) {
1817 }
else if (!strncmp(type, MUSE_PIXTABLE_STRING_SIMPLE,
1818 strlen(MUSE_PIXTABLE_STRING_SIMPLE) + 1)) {
1839 cpl_ensure(aPixtable, CPL_ERROR_NULL_INPUT, 0);
1840 cpl_ensure(aPixtable->
table, CPL_ERROR_NULL_INPUT, 0);
1843 return cpl_table_get_nrow(aPixtable->
table);
1867 cpl_ensure_code(aPixtable && aPixtable->
table && aPixtable->
header,
1868 CPL_ERROR_NULL_INPUT);
1870 == CPL_ERROR_NONE, CPL_ERROR_DATA_NOT_FOUND);
1873 return CPL_ERROR_NONE;
1879 uint32_t *origin = (uint32_t *)cpl_table_get_data_int(aPixtable->
table,
1882 float xlo = FLT_MAX, xhi = -FLT_MAX,
1883 ylo = FLT_MAX, yhi = -FLT_MAX,
1884 llo = FLT_MAX, lhi = -FLT_MAX;
1885 int ifulo = INT_MAX, ifuhi = 0,
1886 slicelo = INT_MAX, slicehi = 0;
1888 for (i = 0; i < nrow; i++) {
1889 if (cdata_xpos[i] > xhi) xhi = cdata_xpos[i];
1890 if (cdata_xpos[i] < xlo) xlo = cdata_xpos[i];
1891 if (cdata_ypos[i] > yhi) yhi = cdata_ypos[i];
1892 if (cdata_ypos[i] < ylo) ylo = cdata_ypos[i];
1893 if (cdata_lambda[i] > lhi) lhi = cdata_lambda[i];
1894 if (cdata_lambda[i] < llo) llo = cdata_lambda[i];
1895 int ifu = muse_pixtable_origin_get_ifu_fast(origin[i]),
1896 slice = muse_pixtable_origin_get_slice_fast(origin[i]);
1897 if (ifu > ifuhi) ifuhi = ifu;
1898 if (ifu < ifulo) ifulo = ifu;
1899 if (slice > slicehi) slicehi = slice;
1900 if (slice < slicelo) slicelo = slice;
1902 char *dodebug = getenv(
"MUSE_DEBUG_PIXTABLE_LIMITS");
1903 if (dodebug && atoi(dodebug)) {
1904 cpl_msg_debug(__func__,
"x: %f...%f, y: %f...%f, lambda: %f...%f, " 1905 "ifu: %d...%d, slice: %d...%d",
1906 xlo, xhi, ylo, yhi, llo, lhi, ifulo, ifuhi, slicelo, slicehi);
1908 cpl_propertylist_erase_regexp(aPixtable->
header, MUSE_HDR_PT_LIMITS_REGEXP, 0);
1926 return CPL_ERROR_NONE;
1948 cpl_ensure_code(aPixtable && aPixtable->
table, CPL_ERROR_NULL_INPUT);
1949 cpl_errorstate prestate = cpl_errorstate_get();
1952 cpl_error_code rc = CPL_ERROR_NONE;
1953 if (!cpl_errorstate_is_equal(prestate)) {
1954 rc = cpl_error_get_code();
1995 const cpl_array *aLambdas,
const cpl_array *aFlux,
1998 cpl_ensure_code(aPixtable && aPixtable->
table, CPL_ERROR_NULL_INPUT);
1999 cpl_ensure_code(aLambdas && aFlux, CPL_ERROR_NULL_INPUT);
2000 cpl_ensure_code(cpl_array_get_size(aLambdas) > 0 &&
2001 cpl_array_get_size(aLambdas) == cpl_array_get_size(aFlux),
2002 CPL_ERROR_ILLEGAL_INPUT);
2003 cpl_ensure_code(cpl_array_get_type(aLambdas) == CPL_TYPE_DOUBLE &&
2004 cpl_array_get_type(aFlux) == CPL_TYPE_DOUBLE,
2005 CPL_ERROR_INCOMPATIBLE_INPUT);
2015 switch (aOperation) {
2016 case MUSE_PIXTABLE_OPERATION_SUBTRACT:
2017 cpl_msg_debug(__func__,
"Subtracting spectrum from pixel table with %" 2018 CPL_SIZE_FORMAT
" slices...", nslices);
2020 case MUSE_PIXTABLE_OPERATION_MULTIPLY:
2021 cpl_msg_debug(__func__,
"Multiplying pixel table of %"CPL_SIZE_FORMAT
2022 " slices with spectrum...", nslices);
2024 case MUSE_PIXTABLE_OPERATION_DIVIDE:
2025 cpl_msg_debug(__func__,
"Dividing pixel table of %"CPL_SIZE_FORMAT
2026 " slices with spectrum...", nslices);
2030 return cpl_error_set(__func__, CPL_ERROR_UNSUPPORTED_MODE);
2034 #pragma omp parallel for default(none) \ 2035 shared(aFlux, aLambdas, aOperation, nslices, slicepts) 2036 for (islice = 0; islice < nslices; islice++) {
2038 cpl_propertylist *order = cpl_propertylist_new();
2040 cpl_table_sort(thispt->
table, order);
2041 cpl_propertylist_delete(order);
2057 if (aOperation == MUSE_PIXTABLE_OPERATION_SUBTRACT) {
2058 cpl_array_subtract(datapt, fint);
2059 }
else if (aOperation == MUSE_PIXTABLE_OPERATION_DIVIDE) {
2060 cpl_array_divide(datapt, fint);
2062 cpl_array_divide(statpt, fint);
2063 cpl_array_divide(statpt, fint);
2065 cpl_array_multiply(datapt, fint);
2067 cpl_array_multiply(statpt, fint);
2068 cpl_array_multiply(statpt, fint);
2070 cpl_array_delete(fint);
2072 cpl_array_unwrap(datapt);
2073 cpl_array_unwrap(statpt);
2074 cpl_array_unwrap(lambdapt);
2075 cpl_table_erase_column(thispt->
table,
"lbda_d");
2079 return CPL_ERROR_NONE;
2093 static cpl_error_code
2096 cpl_ensure_code(aPixtable && aPixtable->
header && aPixtable->
table,
2097 CPL_ERROR_NULL_INPUT);
2098 if (cpl_table_count_selected(aPixtable->
table) < 1) {
2099 return CPL_ERROR_NONE;
2101 cpl_array *sel = cpl_table_where_selected(aPixtable->
table);
2102 cpl_size narray = cpl_array_get_size(sel),
2105 const cpl_size *asel = cpl_array_get_data_cplsize_const(sel);
2106 unsigned int nexp = 0;
2111 if (!cpl_propertylist_has(aPixtable->
header, kwfst) ||
2112 !cpl_propertylist_has(aPixtable->
header, kwlst)) {
2117 ifst = cpl_propertylist_get_long_long(aPixtable->
header, kwfst);
2118 ilst = cpl_propertylist_get_long_long(aPixtable->
header, kwlst);
2120 cpl_size i, nsel = 0;
2121 for (i = 0; i < narray; i++) {
2122 if (asel[i] >= ifst && asel[i] <= ilst) {
2126 cpl_size ifst2 = ifst - nselprev,
2127 ilst2 = ilst - nsel - nselprev;
2128 cpl_msg_debug(__func__,
"exp %d old %"CPL_SIZE_FORMAT
"..%"CPL_SIZE_FORMAT
2129 ", %"CPL_SIZE_FORMAT
" selected (previous: %"CPL_SIZE_FORMAT
2130 "), new %"CPL_SIZE_FORMAT
"..%"CPL_SIZE_FORMAT, nexp,
2131 ifst, ilst, nsel, nselprev, ifst2, ilst2);
2138 }
while (ilst >= ifst);
2139 cpl_array_delete(sel);
2140 return CPL_ERROR_NONE;
2166 cpl_ensure_code(aPixtable && aPixtable->
table && aPixtable->
header,
2167 CPL_ERROR_NULL_INPUT);
2171 return CPL_ERROR_NONE;
2173 #pragma omp critical(cpl_table_select) 2175 cpl_table_unselect_all(aPixtable->
table);
2177 CPL_LESS_THAN, aLow);
2179 CPL_GREATER_THAN, aHigh);
2181 cpl_table_erase_selected(aPixtable->
table);
2186 #pragma omp critical(cpl_table_select) 2188 cpl_table_unselect_all(aPixtable->
ffspec);
2191 aLow - 2. * kMuseSpectralSamplingA);
2194 aHigh + 2. * kMuseSpectralSamplingA);
2195 cpl_table_erase_selected(aPixtable->
ffspec);
2216 cpl_ensure_code(aPixtable && aPixtable->
table && aPixtable->
header,
2217 CPL_ERROR_NULL_INPUT);
2221 return CPL_ERROR_NONE;
2228 #pragma omp critical(cpl_table_select) 2230 cpl_table_unselect_all(aPixtable->
table);
2232 CPL_LESS_THAN, aLo - ptxoff);
2234 CPL_GREATER_THAN, aHi - ptxoff);
2236 cpl_table_erase_selected(aPixtable->
table);
2256 cpl_ensure_code(aPixtable && aPixtable->
table && aPixtable->
header,
2257 CPL_ERROR_NULL_INPUT);
2261 return CPL_ERROR_NONE;
2268 #pragma omp critical(cpl_table_select) 2270 cpl_table_unselect_all(aPixtable->
table);
2272 CPL_LESS_THAN, aLo - ptyoff);
2274 CPL_GREATER_THAN, aHi - ptyoff);
2276 cpl_table_erase_selected(aPixtable->
table);
2298 unsigned short aSlice)
2300 cpl_ensure_code(aPixtable, CPL_ERROR_NULL_INPUT);
2302 cpl_ensure_code(nrow > 0, CPL_ERROR_DATA_NOT_FOUND);
2304 cpl_table_unselect_all(aPixtable->
table);
2305 uint32_t *origin = (uint32_t *)cpl_table_get_data_int(aPixtable->
table,
2308 for (irow = 0; irow < nrow; irow++) {
2311 if (ifu == aIFU && slice == aSlice) {
2312 cpl_table_select_row(aPixtable->
table, irow);
2315 cpl_size nsel = cpl_table_count_selected(aPixtable->
table);
2316 cpl_error_code rc = cpl_table_erase_selected(aPixtable->
table);
2317 cpl_msg_debug(__func__,
"Erased %"CPL_SIZE_FORMAT
" rows from pixel table",
2346 cpl_ensure_code(aPixtable && aPixtable->
table, CPL_ERROR_NULL_INPUT);
2347 cpl_ensure_code(aMask && aMask->
mask, CPL_ERROR_NULL_INPUT);
2348 float *cdata_xpos = cpl_table_get_data_float(aPixtable->
table,
2350 float *cdata_ypos = cpl_table_get_data_float(aPixtable->
table,
2352 cpl_size n_rows = cpl_table_get_nrow(aPixtable->
table);
2368 cpl_size nx = cpl_mask_get_size_x(aMask->
mask);
2369 cpl_size ny = cpl_mask_get_size_y(aMask->
mask);
2370 cpl_size n_enabled = cpl_mask_count(aMask->
mask);
2371 cpl_msg_debug(__func__,
"Mask contains %"CPL_SIZE_FORMAT
" (%.2f %%) enabled " 2372 "pixels of %"CPL_SIZE_FORMAT
" total", n_enabled,
2373 100.*n_enabled/nx/ny, nx*ny);
2374 cpl_size n_sel = n_rows;
2375 cpl_size n_in_table = 0;
2376 for (i_row = 0; i_row < n_rows; i_row++) {
2377 cpl_size ix = lround((cdata_xpos[i_row] - crval_x) / cdelt_x + crpix_x);
2378 cpl_size iy = lround((cdata_ypos[i_row] - crval_y) / cdelt_y + crpix_y);
2379 if ((ix < 1) || (ix > nx) || (iy < 1) || (iy > ny)) {
2383 if (cpl_mask_get(aMask->
mask, ix, iy) != CPL_BINARY_1) {
2384 cpl_table_unselect_row(aPixtable->
table, i_row);
2388 cpl_msg_debug(__func__,
"Mask selected %"CPL_SIZE_FORMAT
" (%.2f %%/%.2f %%) " 2389 "pixels of %"CPL_SIZE_FORMAT
" total/%"CPL_SIZE_FORMAT
" in mask " 2390 "area", n_sel, 100.*n_sel/n_rows, 100.*n_sel/n_in_table,
2391 n_rows, n_in_table);
2392 return CPL_ERROR_NONE;
2425 unsigned char aDisplayHeader)
2427 cpl_ensure_code(aPixtable && aPixtable->
table && aPixtable->
header,
2428 CPL_ERROR_NULL_INPUT);
2430 cpl_ensure_code(aStart >= 0 && aStart < nrows && aCount >= 0,
2431 CPL_ERROR_ILLEGAL_INPUT);
2432 cpl_size last = aStart + aCount;
2433 if (last > nrows - 1) {
2437 double ptxoff = 0., ptyoff = 0.;
2451 cpl_errorstate es = cpl_errorstate_get();
2453 cpl_errorstate_set(es);
2455 uint32_t *cdata_origin = (uint32_t *)cpl_table_get_data_int(aPixtable->
table,
2457 cpl_ensure_code(cdata_xpos && cdata_ypos && cdata_lambda &&
2458 cdata_data && cdata_dq && cdata_stat,
2459 CPL_ERROR_BAD_FILE_FORMAT);
2462 if (aDisplayHeader) {
2463 printf(
"# xpos ypos lambda data dq stat" 2464 " weight exposure IFU xCCD yCCD xRaw yRaw slice\n");
2466 if (aDisplayHeader == 1) {
2467 printf(
"#%13s %13s %9s %11s flag %11s ---------- No No pix " 2468 " pix pix pix No\n# flux in [%s]\n# flux**2 in [%s]\n",
2472 "(flux)",
"(flux**2)",
2478 for (i = aStart; i < last; i++) {
2480 y = muse_pixtable_origin_get_y_fast(cdata_origin[i]),
2484 printf(
"%14.7e %14.7e %9.3f ", cdata_xpos[i] + ptxoff, cdata_ypos[i] + ptyoff,
2487 printf(
"%14.8f %14.8f %9.3f ", cdata_xpos[i], cdata_ypos[i], cdata_lambda[i]);
2489 printf(
"%12.5e 0x%08x %11.5e %10.4e %2d %2d %4d %4d %4d %4d %2d\n",
2490 cdata_data[i], cdata_dq[i], cdata_stat[i],
2491 cdata_weight ? cdata_weight[i] : 0.,
2493 cdata_origin ? muse_pixtable_origin_get_ifu_fast(cdata_origin[i])
2496 cdata_origin ? muse_pixtable_origin_get_slice_fast(cdata_origin[i])
2500 return CPL_ERROR_NONE;
2530 const char *unitx = cpl_table_get_column_unit(aPixtable->
table,
2532 *unity = cpl_table_get_column_unit(aPixtable->
table,
2537 cpl_ensure(!strncmp(unitx, unity, 4), CPL_ERROR_INCOMPATIBLE_INPUT,
2539 if (!strncmp(unitx,
"deg", 4)) {
2542 if (!strncmp(unitx,
"pix", 4)) {
2545 if (!strncmp(unitx,
"rad", 4)) {
2548 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT);
2564 cpl_ensure(aPixtable, CPL_ERROR_NULL_INPUT, CPL_FALSE);
2565 cpl_errorstate prestate = cpl_errorstate_get();
2566 cpl_boolean flag = cpl_propertylist_get_bool(aPixtable->
header,
2568 cpl_errorstate_set(prestate);
2584 cpl_ensure(aPixtable, CPL_ERROR_NULL_INPUT, CPL_FALSE);
2585 cpl_errorstate prestate = cpl_errorstate_get();
2586 cpl_boolean flag = cpl_propertylist_get_bool(aPixtable->
header,
2588 cpl_errorstate_set(prestate);
2604 cpl_ensure(aPixtable, CPL_ERROR_NULL_INPUT, CPL_FALSE);
2621 cpl_ensure_code(aPixtable, CPL_ERROR_NULL_INPUT);
2623 unsigned int *dq = (
unsigned int *)cpl_table_get_data_int(aPixtable->
table,
2627 #pragma omp parallel for default(none) \ 2628 shared(dq, inverse, nrow) 2629 for (i = 0; i < nrow; i++) {
2632 return CPL_ERROR_NONE;
2656 cpl_ensure(aPixtable && aPixtable->
header, CPL_ERROR_NULL_INPUT, NULL);
2660 cpl_ensure(expnum == explast, CPL_ERROR_ILLEGAL_INPUT, NULL);
2670 unsigned short ifu = 0,
2673 for (ipt = 0; ipt < npt; ipt++) {
2677 uint32_t *corigin = (uint32_t *)cpl_table_get_data_int(pts[ipt]->table,
2682 if (ifu != muse_pixtable_origin_get_ifu_fast(corigin[0])) {
2684 image->
header = cpl_propertylist_duplicate(pts[ipt]->header);
2685 cpl_propertylist_erase_regexp(image->
header,
"^ESO DRS MUSE PIXTABLE", 0);
2686 image->
data = cpl_image_new(kMuseOutputXRight, kMuseOutputYTop, CPL_TYPE_FLOAT);
2687 image->
dq = cpl_image_new(kMuseOutputXRight, kMuseOutputYTop, CPL_TYPE_INT);
2689 cpl_image_fill_noise_uniform(image->
dq, EURO3D_MISSDATA, EURO3D_MISSDATA + 0.1);
2690 image->
stat = cpl_image_new(kMuseOutputXRight, kMuseOutputYTop, CPL_TYPE_FLOAT);
2691 cpl_msg_debug(__func__,
"new image (index %hu in list)", ilist);
2695 cpl_msg_error(__func__,
"ipt = %d: no image!", ipt);
2698 float *idata = cpl_image_get_data_float(image->
data),
2699 *istat = cpl_image_get_data_float(image->
stat);
2700 int *idq = cpl_image_get_data_int(image->
dq);
2702 ifu = muse_pixtable_origin_get_ifu_fast(corigin[0]);
2703 unsigned short slice = muse_pixtable_origin_get_slice_fast(corigin[0]);
2706 x1 = INT_MAX, x2 = 0,
2708 for (irow = 0; irow < nrow; irow++) {
2710 unsigned int x = muse_pixtable_origin_get_x_fast(corigin[irow], xoff) - 1,
2711 y = muse_pixtable_origin_get_y_fast(corigin[irow]) - 1;
2712 idata[x + y*kMuseOutputXRight] = cdata[irow];
2713 idq[x + y*kMuseOutputXRight] = cdq[irow];
2714 istat[x + y*kMuseOutputXRight] = cstat[irow];
2723 char *keyword = cpl_sprintf(
"ESO DRS MUSE SLICE%hu CENTER", slice);
2724 cpl_propertylist_update_float(image->
header, keyword, (x2 + x1) / 2. + 1.);
2726 cpl_msg_debug(__func__,
"IFU %hu %s = %.1f", ifu, keyword,
2727 (x2 + x1) / 2. + 1.);
2766 cpl_ensure_code(aPixtable && aPixtable->
header && aList, CPL_ERROR_NULL_INPUT);
2770 cpl_ensure_code(expnum == explast, CPL_ERROR_ILLEGAL_INPUT);
2779 return cpl_error_set(__func__, CPL_ERROR_INCOMPATIBLE_INPUT);
2783 unsigned short ifu = 0,
2786 for (ipt = 0; ipt < npt; ipt++) {
2789 uint32_t *corigin = (uint32_t *)cpl_table_get_data_int(pts[ipt]->table,
2792 if (ifu != muse_pixtable_origin_get_ifu_fast(corigin[0])) {
2796 cpl_msg_error(__func__,
"ipt = %d: no image!", ipt);
2799 float *idata = cpl_image_get_data_float(image->
data),
2800 *istat = cpl_image_get_data_float(image->
stat);
2801 ifu = muse_pixtable_origin_get_ifu_fast(corigin[0]);
2802 unsigned short slice = muse_pixtable_origin_get_slice_fast(corigin[0]);
2806 for (irow = 0; irow < nrow; irow++) {
2808 unsigned int x = muse_pixtable_origin_get_x_fast(corigin[irow], xoff) - 1,
2809 y = muse_pixtable_origin_get_y_fast(corigin[irow]) - 1;
2810 cdata[irow] = idata[x + y*kMuseOutputXRight];
2811 cstat[irow] = istat[x + y*kMuseOutputXRight];
2816 return CPL_ERROR_NONE;
2836 cpl_ensure(aPixtable, CPL_ERROR_NULL_INPUT, NULL);
2837 cpl_size n_rows = cpl_table_get_nrow(aPixtable->
table);
2838 unsigned int ifu_slice_mask = (0x1f << MUSE_ORIGIN_SHIFT_IFU) | 0x3f;
2839 cpl_table_duplicate_column(aPixtable->
table,
"ifuslice",
2841 unsigned int *slicedata = (
unsigned int *)
2842 cpl_table_get_data_int(aPixtable->
table,
"ifuslice");
2844 unsigned int last_ifu_slice = 0;
2845 int is_sorted = CPL_TRUE;
2846 for (i_row = 0; i_row < n_rows; i_row++) {
2847 slicedata[i_row] &= ifu_slice_mask;
2848 if (is_sorted && slicedata[i_row] < last_ifu_slice) {
2849 is_sorted = CPL_FALSE;
2851 last_ifu_slice = slicedata[i_row];
2855 cpl_propertylist *order = cpl_propertylist_new();
2856 cpl_propertylist_append_bool(order,
"ifuslice", CPL_FALSE);
2858 cpl_msg_debug(__func__,
"sorting pixel table: quick sort, %"CPL_SIZE_FORMAT
2859 " entries", n_rows);
2860 cpl_table_sort(aPixtable->
table, order);
2861 cpl_propertylist_delete(order);
2863 cpl_propertylist_erase_regexp(aPixtable->
header, MUSE_HDR_PT_EXP_REGEXP, 0);
2864 cpl_msg_debug(__func__,
"pixel table sorted.");
2868 cpl_size n_col = cpl_table_get_ncol(aPixtable->
table);
2869 cpl_array *colnames = cpl_table_get_column_names(aPixtable->
table);
2871 cpl_size n_slices = 0;
2872 while (i_row < n_rows) {
2873 unsigned int ifu_slice = slicedata[i_row];
2875 for (j_row = i_row+1; j_row < n_rows && slicedata[j_row] == ifu_slice;
2878 cpl_size nrows_slice = j_row - i_row;
2880 slice_pixtable->
table = cpl_table_new(nrows_slice);
2882 for (i_col = 0; i_col < n_col; i_col++) {
2883 const char *cname = cpl_array_get_string(colnames, i_col);
2884 if (strcmp(cname,
"ifuslice") == 0)
2886 cpl_type ctype = cpl_table_get_column_type(aPixtable->
table, cname);
2887 if (ctype == CPL_TYPE_INT) {
2888 int *cdata = cpl_table_get_data_int(aPixtable->
table, cname);
2889 cpl_table_wrap_int(slice_pixtable->
table, cdata + i_row, cname);
2890 }
else if (ctype == CPL_TYPE_FLOAT) {
2891 float *cdata = cpl_table_get_data_float(aPixtable->
table, cname);
2892 cpl_table_wrap_float(slice_pixtable->
table, cdata + i_row, cname);
2893 }
else if (ctype == CPL_TYPE_DOUBLE) {
2894 double *cdata = cpl_table_get_data_double(aPixtable->
table, cname);
2895 cpl_table_wrap_double(slice_pixtable->
table, cdata + i_row, cname);
2896 }
else if (ctype == CPL_TYPE_STRING) {
2897 char **cdata = cpl_table_get_data_string(aPixtable->
table, cname);
2898 cpl_table_wrap_string(slice_pixtable->
table, cdata + i_row, cname);
2900 const char *unit = cpl_table_get_column_unit(aPixtable->
table, cname);
2901 cpl_table_set_column_unit(slice_pixtable->
table, cname, unit);
2904 slice_pixtable->
header = cpl_propertylist_duplicate(aPixtable->
header);
2906 slice_tables = cpl_realloc(slice_tables,
2908 slice_tables[n_slices] = slice_pixtable;
2910 slice_tables[n_slices] = NULL;
2913 cpl_array_delete(colnames);
2914 cpl_table_erase_column(aPixtable->
table,
"ifuslice");
2916 return slice_tables;
2934 cpl_ensure(aPixtables, CPL_ERROR_NULL_INPUT, -1);
2936 while (aPixtables[n] != NULL) {
2960 for (t = aPixtables; *t != NULL; t++) {
2961 cpl_array *colnames = cpl_table_get_column_names((*t)->table);
2962 cpl_size n_col = cpl_table_get_ncol((*t)->table);
2964 for (i_col = 0; i_col < n_col; i_col++) {
2965 const char *cname = cpl_array_get_string(colnames, i_col);
2966 cpl_table_unwrap((*t)->table, cname);
2968 cpl_array_delete(colnames);
2969 cpl_table_delete((*t)->table);
2970 cpl_propertylist_delete((*t)->header);
2973 cpl_free(aPixtables);
#define MUSE_PIXTABLE_XPOS
cpl_polynomial ** muse_trace_table_get_polys_for_slice(const cpl_table *aTable, const unsigned short aSlice)
construct polynomial from the trace table entry for the given slice
#define MUSE_PIXTABLE_FFDATA
muse_pixtable_wcs
State of the astrometric calibration of a MUSE pixel table.
unsigned int muse_pixtable_get_expnum(muse_pixtable *aPixtable, cpl_size aRow)
Get the exposure number of a given row in a pixel table.
cpl_array * muse_cplarray_interpolate_table_linear(const cpl_array *aTargetAbscissa, const cpl_table *aSrcTable, const char *aSrcAbscissa, const char *aSrcOrdinate)
Linear interpolation of a 1d column.
cpl_error_code muse_wave_table_get_orders(const cpl_table *aWave, unsigned short *aXOrder, unsigned short *aYOrder)
Determine the x- and y-order of the polynomial stored in a wavelength calibration table...
Structure definition for a collection of muse_images.
int muse_trace_table_get_order(const cpl_table *aTable)
determine order of tracing polynomial from table
#define MUSE_HDR_PT_EXP_FST
FITS header keyword defining the first row index for a given exposure.
void muse_pixtable_extracted_delete(muse_pixtable **aPixtables)
Delete a pixel table array.
double muse_pfits_get_cd(const cpl_propertylist *aHeaders, unsigned int aAxisI, unsigned int aAxisJ)
find out the WCS coordinate at the reference point
cpl_polynomial * muse_wave_table_get_poly_for_slice(const cpl_table *aTable, const unsigned short aSlice)
Construct polynomial from the wavelength calibration table entry for the given slice.
const char * muse_pfits_get_extname(const cpl_propertylist *aHeaders)
find out the extension name
muse_pixtable * muse_pixtable_load(const char *aFilename)
Load the table itself and the FITS headers of a MUSE pixel table from a file.
unsigned short muse_pixtable_origin_get_slice(uint32_t aOrigin)
Get the slice number from the encoded 32bit origin number.
muse_pixtable * muse_pixtable_duplicate(muse_pixtable *aPixtable)
Make a copy of the pixtanle.
#define MUSE_HDR_PT_XLO
FITS header keyword contains the lower limit of the data in x-direction.
cpl_size muse_pixtable_extracted_get_size(muse_pixtable **aPixtables)
Get the size of an array of extracted pixel tables.
double muse_pfits_get_crval(const cpl_propertylist *aHeaders, unsigned int aAxis)
find out the WCS coordinate at the reference point
#define MUSE_PIXTABLE_FFLAMBDA
unsigned int muse_pixtable_origin_set_offset(muse_pixtable *aPixtable, cpl_polynomial *aLTrace, unsigned short aIFU, unsigned short aSlice)
Set the slice offset from the pixel table header.
unsigned char muse_utils_get_ifu(const cpl_propertylist *aHeaders)
Find out the IFU/channel from which this header originated.
cpl_size muse_pfits_get_naxis(const cpl_propertylist *aHeaders, unsigned int aAxis)
find out the size of a given axis
#define MUSE_HDR_PT_FLUXCAL
muse_pixtable ** muse_pixtable_extracted_get_slices(muse_pixtable *aPixtable)
Extract one pixel table per IFU and slice.
cpl_image * data
the data extension
cpl_size muse_pixtable_get_nrow(const muse_pixtable *aPixtable)
get the number of rows within the pixel table
#define MUSE_HDR_PT_LHI
FITS header keyword contains the upper limit of the data in spectral direction.
#define MUSE_HDR_PT_ILO
FITS header keyword contains the lowest IFU number in the data.
int muse_pixtable_get_type(muse_pixtable *aPixtable)
Determine the type of pixel table.
#define MUSE_HDR_PT_EXP_LST
FITS header keyword defining the last row index for a given exposure.
cpl_error_code muse_pixtable_and_selected_mask(muse_pixtable *aPixtable, muse_mask *aMask)
Select all pixels where the (x,y) positions are enabled in the given mask.
cpl_error_code muse_pixtable_erase_ifu_slice(muse_pixtable *aPixtable, unsigned char aIFU, unsigned short aSlice)
Erase pixel table rows related to one slice of one IFU.
cpl_error_code muse_pixtable_dump(muse_pixtable *aPixtable, cpl_size aStart, cpl_size aCount, unsigned char aDisplayHeader)
Dump a MUSE pixel table to the screen, resolving the origin column.
cpl_image * stat
the statistics extension
#define MUSE_HDR_PT_FFCORR
#define MUSE_HDR_PT_TYPE
Pixel table "type" stored in the FITS header.
const char * muse_pfits_get_bunit(const cpl_propertylist *aHeaders)
find out the unit string
#define MUSE_HDR_PT_MERGED
Structure definition of MUSE three extension FITS file.
cpl_table * table
The pixel table.
cpl_propertylist * header
the FITS header
cpl_error_code muse_pixtable_origin_copy_offsets(muse_pixtable *aOut, muse_pixtable *aFrom, unsigned int aNum)
Copy MUSE_HDR_PT_IFU_SLICE_OFFSET keywords between pixel tables.
cpl_error_code muse_pixtable_reset_dq(muse_pixtable *aPixtable, unsigned int aDQ)
Reset a given bad pixel status (DQ flag) for all pixels in the table.
cpl_error_code muse_cpltable_check(const cpl_table *aTable, const muse_cpltable_def *aDef)
Check whether the table contains the fields of the definition.
#define MUSE_PIXTABLE_DATA
unsigned int muse_imagelist_get_size(muse_imagelist *aList)
Return the number of stored images.
muse_imagelist * muse_pixtable_to_imagelist(muse_pixtable *aPixtable)
Project a pixel table with data from one IFU back onto its image.
void muse_trace_polys_delete(cpl_polynomial *aPolys[])
Delete the multi-polynomial array created in relation to tracing.
#define MUSE_HDR_PT_SHI
FITS header keyword contains the highest slice number in the data.
cpl_image * dq
the data quality extension
cpl_error_code muse_pixtable_append_ff(muse_pixtable *aPixtable, muse_image *aFF, cpl_table *aTrace, cpl_table *aWave, float aSampling)
Create flat-field spectrum and append to pixel table.
cpl_error_code muse_pixtable_restrict_wavelength(muse_pixtable *aPixtable, double aLow, double aHigh)
Restrict a pixel table to a certain wavelength range.
cpl_table * muse_cpltable_new(const muse_cpltable_def *aDef, cpl_size aLength)
Create an empty table according to the specified definition.
double muse_pfits_get_crpix(const cpl_propertylist *aHeaders, unsigned int aAxis)
find out the WCS reference point
cpl_boolean muse_pixtable_is_skysub(muse_pixtable *aPixtable)
Determine whether the pixel table is sky subtracted.
unsigned int muse_pixtable_origin_get_x(uint32_t aOrigin, muse_pixtable *aPixtable, cpl_size aRow)
Get the horizontal coordinate from the encoded 32bit origin number.
#define MUSE_PIXTABLE_WEIGHT
Structure definition of MUSE pixel table.
#define MUSE_HDR_PT_IFU_SLICE_OFFSET
FITS header keyword for the horizontal slice offset on the CCD.
#define MUSE_WCS_KEYS
regular expression for WCS properties
muse_image * muse_imagelist_get(muse_imagelist *aList, unsigned int aIdx)
Get the muse_image of given list index.
cpl_array * muse_cpltable_extract_column(cpl_table *aTable, const char *aColumn)
Create an array from a section of a column.
#define MUSE_PIXTABLE_FF_EXT
muse_pixtable_wcs muse_pixtable_wcs_check(muse_pixtable *aPixtable)
Check the state of the world coordinate system of a pixel table.
#define MUSE_HDR_PT_YLO
FITS header keyword contains the lower limit of the data in y-direction.
#define MUSE_HDR_PT_SKYSUB
cpl_boolean muse_pixtable_is_fluxcal(muse_pixtable *aPixtable)
Determine whether the pixel table is flux calibrated.
#define MUSE_HDR_PT_ILLUM_REGEXP
cpl_table * ffspec
A flat-field spectrum.
cpl_array * muse_cplarray_interpolate_linear(const cpl_array *aTargetAbscissa, const cpl_array *aSourceAbscissa, const cpl_array *aSourceOrdinate)
Linear interpolation of a 1d array.
#define MUSE_HDR_PT_IHI
FITS header keyword contains the highest IFU number in the data.
unsigned int muse_pixtable_origin_get_offset(muse_pixtable *aPixtable, unsigned int aExpNum, unsigned short aIFU, unsigned short aSlice)
Get the slice offset from the pixel table header.
uint32_t muse_pixtable_origin_encode(unsigned int aX, unsigned int aY, unsigned short aIFU, unsigned short aSlice, unsigned int aOffset)
Encode the three CCD coordinates defining the origin of one MUSE pixel into a 32bit integer...
muse_pixtable * muse_pixtable_load_window(const char *aFilename, cpl_size aStart, cpl_size aNRows)
Load a range of rows from the table and all the FITS headers of a MUSE pixel table from a file...
muse_pixtable_operation
Type of operation to apply to a MUSE pixel table.
muse_pixtable * muse_pixtable_create(muse_image *aImage, cpl_table *aTrace, cpl_table *aWave, cpl_table *aGeoTable)
Create the pixel table for one CCD.
muse_pixtable * muse_pixtable_load_restricted_wavelength(const char *aFilename, double aLambdaMin, double aLambdaMax)
Load a pixel table from file and cut down the wavelength range.
#define MUSE_PIXTABLE_ORIGIN
cpl_error_code muse_cplpropertylist_update_long_long(cpl_propertylist *aHeader, const char *aKeyword, cpl_size aValue)
Update an integer-like property irrespective of the real type.
#define MUSE_HDR_PT_LLO
FITS header keyword contains the lower limit of the data in spectral direction.
cpl_error_code muse_pixtable_save(muse_pixtable *aPixtable, const char *aFilename)
Save a MUSE pixel table to a file on disk.
cpl_error_code muse_cpltable_copy_array(cpl_table *aTable, const char *aColumn, const cpl_array *aArray)
Copy an array into a table.
unsigned short muse_pixtable_origin_get_ifu(uint32_t aOrigin)
Get the IFU number from the encoded 32bit origin number.
muse_pixtable * muse_pixtable_load_merge_channels(cpl_table *aExposureList, double aLambdaMin, double aLambdaMax)
Load and merge the pixel tables of the 24 MUSE sub-fields.
cpl_error_code muse_pixtable_restrict_xpos(muse_pixtable *aPixtable, double aLo, double aHi)
Restrict a pixel table to a certain x coordinate range.
#define MUSE_PIXTABLE_STAT
cpl_error_code muse_pixtable_from_imagelist(muse_pixtable *aPixtable, muse_imagelist *aList)
Get pixel table values back from a per-IFU imagelist.
Handling of "mask" files.
muse_imagelist * muse_imagelist_new(void)
Create a new (empty) MUSE image list.
const muse_cpltable_def muse_pixtable_def[]
MUSE pixel table definition.
#define MUSE_PIXTABLE_LAMBDA
cpl_table * muse_geo_table_extract_ifu(const cpl_table *aTable, const unsigned char aIFU)
Extract the part of a geometry table dealing with a given IFU.
cpl_error_code muse_quadrants_coords_to_raw(cpl_propertylist *aHeader, int *aX, int *aY)
Convert coordinates of a trimmed image to raw-image coordinates.
#define MUSE_HDR_PT_SLO
FITS header keyword contains the lowest slice number in the data.
cpl_propertylist * header
the FITS header
cpl_error_code muse_pixtable_spectrum_apply(muse_pixtable *aPixtable, const cpl_array *aLambdas, const cpl_array *aFlux, muse_pixtable_operation aOperation)
Apply a spectrum given by two arrays with an operation to a pixel table.
#define MUSE_HDR_PT_RVCORR
Definition of a cpl table structure.
muse_image * muse_image_new(void)
Allocate memory for a new muse_image object.
cpl_table * muse_resampling_spectrum(muse_pixtable *aPixtable, double aBinwidth)
Resample the selected pixels of a pixel table into a spectrum.
unsigned int muse_pixtable_origin_get_y(uint32_t aOrigin)
Get the vertical coordinate from the encoded 32bit origin number.
cpl_mask * mask
The mask data.
#define MUSE_PIXTABLE_YPOS
void muse_pixtable_delete(muse_pixtable *aPixtable)
Deallocate memory associated to a pixel table object.
cpl_error_code muse_pixtable_restrict_ypos(muse_pixtable *aPixtable, double aLo, double aHi)
Restrict a pixel table to a certain y coordinate range.
#define MUSE_HDR_PT_YHI
FITS header keyword contains the upper limit of the data in y-direction.
cpl_error_code muse_imagelist_set(muse_imagelist *aList, muse_image *aImage, unsigned int aIdx)
Set the muse_image of given list index.
cpl_error_code muse_pixtable_compute_limits(muse_pixtable *aPixtable)
(Re-)Compute the limits of the coordinate columns of a pixel table.
static cpl_error_code muse_pixtable_fix_exp_headers(muse_pixtable *aPixtable)
Fix the exposure ranges in the header of a pixel table.
cpl_propertylist * header
The FITS header.
cpl_error_code muse_pixtable_flux_multiply(muse_pixtable *aPixtable, double aScale)
Scale the flux of a pixel table with correct treatment of variance.
#define MUSE_HDR_PT_XHI
FITS header keyword contains the upper limit of the data in x-ion.
cpl_boolean muse_pixtable_is_rvcorr(muse_pixtable *aPixtable)
Determine whether the pixel table is radial-velocity corrected.