33 #include "muse_datacube.h" 35 #include "muse_flux.h" 36 #include "muse_pfits.h" 37 #include "muse_quality.h" 38 #include "muse_utils.h" 72 muse_datacube_collapse_filter_buffer(
const muse_table *aFilter,
73 double aCRVAL,
double aCRPIX,
double aCDELT,
74 cpl_boolean aIsLog,
int *aL1,
int *aL2,
77 if (!aFilter || !aL1 || !aL2) {
80 if (!aFilter->
table) {
83 double *fdata = (
double *)cpl_calloc(*aL2,
sizeof(
double));
86 for (l = 0; l < *aL2; l++) {
87 double lambda = (l + 1. - aCRPIX) * aCDELT + aCRVAL;
89 lambda = aCRVAL * exp((l + 1. - aCRPIX) * aCDELT / aCRVAL);
97 while (l < *aL2 && fabs(fdata[l]) < DBL_EPSILON) {
101 while (l > *aL1 && fabs(fdata[l]) < DBL_EPSILON) {
104 double lbda1 = (*aL1 + 1. - aCRPIX) * aCDELT + aCRVAL,
105 lbda2 = (*aL2 + 1. - aCRPIX) * aCDELT + aCRVAL;
107 lbda1 = aCRVAL * exp((*aL1 + 1. - aCRPIX) * aCDELT / aCRVAL);
108 lbda2 = aCRVAL * exp((*aL2 + 1. - aCRPIX) * aCDELT / aCRVAL);
111 cpl_msg_debug(__func__,
"Filter wavelength range %.1f..%.1fA (cube planes " 112 "%d..%d), %.2f%% of filter area covered by data.", lbda1, lbda2,
113 *aL1, *aL2, fraction * 100.);
115 *aFraction = fraction;
147 cpl_ensure(aEuro3D && aEuro3D->
dtable && aEuro3D->
hdata, CPL_ERROR_NULL_INPUT,
154 const char *unitx = cpl_table_get_column_unit(aEuro3D->
dtable,
"XPOS"),
155 *unity = cpl_table_get_column_unit(aEuro3D->
dtable,
"YPOS");
156 cpl_ensure(unitx && unity, CPL_ERROR_DATA_NOT_FOUND, NULL);
158 if (!strncmp(unitx, unity, 4) && !strncmp(unitx,
"deg", 4)) {
162 double xmin = cpl_table_get_column_min(aEuro3D->
dtable,
"XPOS"),
163 xmax = cpl_table_get_column_max(aEuro3D->
dtable,
"XPOS"),
164 ymin = cpl_table_get_column_min(aEuro3D->
dtable,
"YPOS"),
165 ymax = cpl_table_get_column_max(aEuro3D->
dtable,
"YPOS");
166 double x1, x2, y1, y2;
168 wcs->crval1 /= CPL_MATH_DEG_RAD;
169 wcs->
crval2 /= CPL_MATH_DEG_RAD;
171 ymin / CPL_MATH_DEG_RAD, &x1, &y1);
173 ymax / CPL_MATH_DEG_RAD, &x2, &y2);
180 int zmin = cpl_table_get_column_min(aEuro3D->
dtable,
"SPEC_STA"),
181 zmax = cpl_table_get_column_max(aEuro3D->
dtable,
"SPEC_STA"),
182 nx = lround(fabs(x2 - x1)) + 1,
183 ny = lround(fabs(y2 - y1)) + 1,
184 nz = zmax - zmin + 1;
187 cpl_size zmaxpos = -1;
188 cpl_table_get_column_maxpos(aEuro3D->
dtable,
"SPEC_STA", &zmaxpos);
189 const cpl_array *amax = cpl_table_get_array(aEuro3D->
dtable,
"DATA_SPE",
191 int l, nsize = cpl_array_get_size(amax);
192 for (l = nsize - 1; l > 0; l--) {
194 if (cpl_array_is_valid(amax, l) == 1) {
199 int nspec = cpl_table_get_nrow(aEuro3D->
dtable);
200 cpl_msg_debug(__func__,
"Euro3D dimensions: %dx%dx%d (z = %d - %d, valid %d)," 201 " %d spectra", nx, ny, nz, zmax, zmin, l + 1, nspec);
204 double crvals = cpl_propertylist_get_double(aEuro3D->
hdata,
"CRVALS"),
205 cdelts = cpl_propertylist_get_double(aEuro3D->
hdata,
"CDELTS");
206 const char *ctypes = cpl_propertylist_get_string(aEuro3D->
hdata,
"CTYPES");
207 cpl_boolean loglambda = ctypes && (!strncmp(ctypes,
"AWAV-LOG", 9) ||
208 !strncmp(ctypes,
"WAVE-LOG", 9));
209 cpl_msg_debug(__func__,
"spectral WCS: %f / %f %s", crvals, cdelts,
210 loglambda ?
"log" :
"linear");
212 double *fdata = NULL;
213 double ffraction = 1.;
215 fdata = muse_datacube_collapse_filter_buffer(aFilter, crvals, zmin,
216 cdelts, loglambda, &l1, &l2,
222 image->
header = cpl_propertylist_duplicate(aEuro3D->
header);
227 image->
data = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
228 float *outdata = cpl_image_get_data_float(image->
data);
229 image->
dq = cpl_image_new(nx, ny, CPL_TYPE_INT);
231 cpl_image_add_scalar(image->
dq, EURO3D_MISSDATA);
232 int *outdq = cpl_image_get_data_int(image->
dq);
236 cpl_boolean usevariance = CPL_FALSE;
237 if (getenv(
"MUSE_COLLAPSE_USE_VARIANCE") &&
238 atoi(getenv(
"MUSE_COLLAPSE_USE_VARIANCE")) > 0) {
239 usevariance = CPL_TRUE;
244 #pragma omp parallel for default(none) private(l) \ 245 shared(aEuro3D, fdata, l1, l2, noutside, nspec, nx, ny, outdata, \ 246 outdq, usevariance, wcs) 247 for (k = 0; k < nspec; k++) {
249 double xpos = cpl_table_get(aEuro3D->
dtable,
"XPOS", k, &err),
250 ypos = cpl_table_get(aEuro3D->
dtable,
"YPOS", k, &err);
252 cpl_msg_warning(__func__,
"spectrum %d in Euro3D table does not have " 253 "position information!", k + 1);
260 ypos * CPL_MATH_RAD_DEG, &xpx, &ypx);
264 int i = lround(xpx) - 1,
266 if (i >= nx || j >= ny) {
272 int nstart = cpl_table_get_int(aEuro3D->
dtable,
"SPEC_STA", k, &err);
273 const cpl_array *adata = cpl_table_get_array(aEuro3D->
dtable,
"DATA_SPE", k),
274 *adq = cpl_table_get_array(aEuro3D->
dtable,
"QUAL_SPE", k),
277 astat = cpl_table_get_array(aEuro3D->
dtable,
"STAT_SPE", k);
280 double sumdata = 0., sumweight = 0.;
281 for (l = l1; l < l2; l++) {
283 int idx = l - nstart + 1;
286 double fweight = fdata ? fdata[l] : 1.;
287 cpl_errorstate prestate = cpl_errorstate_get();
288 int dq = cpl_array_get_int(adq, idx, &err);
290 cpl_errorstate_set(prestate);
294 double variance = 1.;
296 variance = astat ? cpl_array_get(astat, idx, &err) : 1.;
298 variance *= variance;
302 if (err || !isnormal(variance)) {
306 double data = cpl_array_get(adata, idx, &err);
308 sumdata += data * fweight / variance;
309 sumweight += fweight / variance;
311 if (isnormal(sumweight)) {
312 outdata[i + j*nx] = sumdata / sumweight;
313 outdq[i + j*nx] = EURO3D_GOODPIXEL;
315 outdq[i + j*nx] = EURO3D_MISSDATA;
321 cpl_msg_warning(__func__,
"Skipped %d spaxels, due to their location " 322 "outside the recostructed image!", noutside);
353 cpl_ensure(aCube && aCube->
data && aCube->
header,
354 CPL_ERROR_NULL_INPUT, NULL);
355 cpl_ensure(cpl_image_get_type(cpl_imagelist_get(aCube->
data, 0)) == CPL_TYPE_FLOAT,
356 CPL_ERROR_INCOMPATIBLE_INPUT, NULL);
358 cpl_ensure(cpl_image_get_type(cpl_imagelist_get(aCube->
dq, 0)) == CPL_TYPE_INT,
359 CPL_ERROR_INCOMPATIBLE_INPUT, NULL);
362 cpl_ensure(cpl_image_get_type(cpl_imagelist_get(aCube->
stat, 0)) == CPL_TYPE_FLOAT,
363 CPL_ERROR_INCOMPATIBLE_INPUT, NULL);
367 int nx = cpl_image_get_size_x(cpl_imagelist_get(aCube->
data, 0)),
368 ny = cpl_image_get_size_y(cpl_imagelist_get(aCube->
data, 0)),
369 nz = cpl_imagelist_get_size(aCube->
data);
375 cpl_boolean loglambda = ctype3 && (!strncmp(ctype3,
"AWAV-LOG", 9) ||
376 !strncmp(ctype3,
"WAVE-LOG", 9));
378 double *fdata = NULL;
379 double ffraction = 1.;
381 fdata = muse_datacube_collapse_filter_buffer(aFilter, crval3, crpix3,
382 cd33, loglambda, &l1, &l2,
388 image->
header = cpl_propertylist_duplicate(aCube->
header);
389 cpl_propertylist_erase_regexp(image->
header,
"^C...*3$|^CD3_.$", 0);
393 image->
data = cpl_image_new(nx, ny, CPL_TYPE_FLOAT);
394 float *outdata = cpl_image_get_data_float(image->
data);
395 image->
dq = cpl_image_new(nx, ny, CPL_TYPE_INT);
396 int *outdq = cpl_image_get_data_int(image->
dq);
400 cpl_boolean usevariance = CPL_FALSE;
401 if (getenv(
"MUSE_COLLAPSE_USE_VARIANCE") &&
402 atoi(getenv(
"MUSE_COLLAPSE_USE_VARIANCE")) > 0) {
403 usevariance = CPL_TRUE;
408 #pragma omp parallel for default(none) \ 409 shared(aCube, nx, ny, l1, l2, fdata, outdata, outdq, usevariance) 410 for (i = 0; i < nx; i++) {
412 for (j = 0; j < ny; j++) {
413 double sumdata = 0., sumweight = 0.;
415 for (l = l1; l < l2; l++) {
417 double fweight = fdata ? fdata[l] : 1.;
418 float *pdata = cpl_image_get_data_float(cpl_imagelist_get(aCube->
data, l)),
420 if (!isfinite(pdata[i + j*nx])) {
424 int *pdq = cpl_image_get_data_int(cpl_imagelist_get(aCube->
dq, l));
425 if (pdq && pdq[i + j*nx]) {
430 double variance = 1.;
432 pstat = cpl_image_get_data_float(cpl_imagelist_get(aCube->
stat, l));
433 variance = pstat ? pstat[i + j*nx] : 1.;
434 if (!isnormal(variance)) {
439 sumdata += pdata[i + j*nx] * fweight / variance;
440 sumweight += fweight / variance;
442 if (isnormal(sumweight)) {
443 outdata[i + j*nx] = sumdata / sumweight;
444 outdq[i + j*nx] = EURO3D_GOODPIXEL;
446 outdq[i + j*nx] = EURO3D_MISSDATA;
476 cpl_ensure_code(aFilename, CPL_ERROR_NULL_INPUT);
477 cpl_error_code rc = CPL_ERROR_NONE;
478 if (!aImages || !aNames) {
482 for (i = 0; i < nimages; i++) {
486 cpl_propertylist *header = cpl_propertylist_new();
487 cpl_errorstate es = cpl_errorstate_get();
489 *ucomment = cpl_propertylist_get_comment(image->
header,
"BUNIT");
490 if (!cpl_errorstate_is_equal(es) && !unit) {
491 cpl_errorstate_set(es);
493 char dataext[KEYWORD_LENGTH], obj[KEYWORD_LENGTH],
494 *dqext = NULL, *statext = NULL;
495 snprintf(dataext, KEYWORD_LENGTH,
"%s", cpl_array_get_string(aNames, i));
497 dqext = cpl_sprintf(
"%s_%s", cpl_array_get_string(aNames, i), EXTNAME_DQ);
500 statext = cpl_sprintf(
"%s_%s", cpl_array_get_string(aNames, i), EXTNAME_STAT);
502 snprintf(obj, KEYWORD_LENGTH,
"%s", cpl_array_get_string(aNames, i));
503 cpl_propertylist_append_string(header,
"EXTNAME", dataext);
504 cpl_propertylist_set_comment(header,
"EXTNAME",
"reconstructed image (data values)");
507 cpl_propertylist_append_string(header,
"BUNIT", unit);
508 cpl_propertylist_set_comment(header,
"BUNIT", ucomment);
512 cpl_propertylist_copy_property_regexp(header, image->
header,
516 rc = cpl_image_save(image->
data, aFilename, CPL_TYPE_UNSPECIFIED, header,
520 cpl_propertylist_update_string(header,
"EXTNAME", dqext);
521 cpl_propertylist_set_comment(header,
"EXTNAME",
"reconstructed image (bad pixel status values)");
522 cpl_propertylist_erase(header,
"BUNIT");
523 snprintf(obj, KEYWORD_LENGTH,
"%s, %s", cpl_array_get_string(aNames, i),
527 rc = cpl_image_save(image->
dq, aFilename, CPL_TYPE_INT, header,
531 cpl_propertylist_update_string(header,
"EXTNAME", statext);
532 cpl_propertylist_set_comment(header,
"EXTNAME",
"reconstructed image (variance)");
534 char *ustat = cpl_sprintf(
"(%s)**2", unit);
535 cpl_propertylist_append_string(header,
"BUNIT", ustat);
538 snprintf(obj, KEYWORD_LENGTH,
"%s, %s", cpl_array_get_string(aNames, i),
542 rc = cpl_image_save(image->
stat, aFilename, CPL_TYPE_UNSPECIFIED, header,
546 cpl_propertylist_delete(header);
573 cpl_ensure_code(aEuro3D && aFilename, CPL_ERROR_NULL_INPUT);
576 cpl_error_code rc = cpl_table_save(aEuro3D->
dtable, aEuro3D->
header,
577 aEuro3D->
hdata, aFilename, CPL_IO_CREATE);
578 if (rc != CPL_ERROR_NONE) {
579 cpl_msg_warning(__func__,
"failed to save data part of the Euro3D table: " 580 "%s", cpl_error_get_message());
583 rc = cpl_table_save(aEuro3D->
gtable, NULL, aEuro3D->
hgroup, aFilename,
585 if (rc != CPL_ERROR_NONE) {
586 cpl_msg_warning(__func__,
"failed to save group part of the Euro3D table: " 587 "%s", cpl_error_get_message());
616 cpl_table_delete(aEuro3D->
dtable);
617 cpl_table_delete(aEuro3D->
gtable);
618 cpl_propertylist_delete(aEuro3D->
header);
619 cpl_propertylist_delete(aEuro3D->
hdata);
620 cpl_propertylist_delete(aEuro3D->
hgroup);
622 cpl_array_delete(aEuro3D->
recnames);
637 static cpl_error_code
640 cpl_ensure_code(aCube, CPL_ERROR_NULL_INPUT);
642 return CPL_ERROR_NONE;
645 for (k = 0; k < nimages; k++) {
653 return CPL_ERROR_NONE;
673 cpl_ensure_code(aCube && aCube->
data && aCube->
stat && aCube->
dq,
674 CPL_ERROR_NULL_INPUT);
677 int l, nx = cpl_image_get_size_x(cpl_imagelist_get(aCube->
data, 0)),
678 ny = cpl_image_get_size_y(cpl_imagelist_get(aCube->
data, 0)),
679 nz = cpl_imagelist_get_size(aCube->
data);
680 #pragma omp parallel for default(none) \ 681 shared(aCube, nx, ny, nz) 682 for (l = 0; l < nz; l++) {
684 for (i = 0; i < nx; i++) {
686 for (j = 0; j < ny; j++) {
687 float *pdata = cpl_image_get_data_float(cpl_imagelist_get(aCube->
data, l)),
688 *pstat = cpl_image_get_data_float(cpl_imagelist_get(aCube->
stat, l));
689 int *pdq = cpl_image_get_data_int(cpl_imagelist_get(aCube->
dq, l));
690 if (pdq[i + j*nx] == EURO3D_GOODPIXEL) {
694 pdata[i + j*nx] = NAN;
695 pstat[i + j*nx] = NAN;
701 cpl_imagelist_delete(aCube->
dq);
705 muse_datacube_convert_dq_recimages(aCube);
707 return CPL_ERROR_NONE;
743 cpl_ensure_code(aCube && aCube->
header && aFilename, CPL_ERROR_NULL_INPUT);
746 cpl_propertylist *header = cpl_propertylist_new();
748 cpl_propertylist_copy_property_regexp(header, aCube->
header,
750 cpl_error_code rc = cpl_propertylist_save(header, aFilename, CPL_IO_CREATE);
751 cpl_propertylist_delete(header);
754 header = cpl_propertylist_new();
755 cpl_propertylist_append_string(header,
"EXTNAME", EXTNAME_DATA);
756 cpl_propertylist_set_comment(header,
"EXTNAME", EXTNAME_DATA_COMMENT);
759 cpl_propertylist_copy_property_regexp(header, aCube->
header,
762 aCube->
dq ?
"DQ" : NULL, aCube->
stat ?
"STAT" : NULL);
763 rc = cpl_imagelist_save(aCube->
data, aFilename, CPL_TYPE_FLOAT, header,
765 cpl_propertylist_delete(header);
767 if (rc == CPL_ERROR_NONE && aCube->
dq) {
768 header = cpl_propertylist_new();
769 cpl_propertylist_append_string(header,
"EXTNAME", EXTNAME_DQ);
770 cpl_propertylist_set_comment(header,
"EXTNAME", EXTNAME_DQ_COMMENT);
773 cpl_propertylist_copy_property_regexp(header, aCube->
header,
776 aCube->
stat ?
"STAT" : NULL);
777 rc = cpl_imagelist_save(aCube->
dq, aFilename, CPL_TYPE_INT, header,
779 cpl_propertylist_delete(header);
782 if (rc == CPL_ERROR_NONE && aCube->
stat) {
783 header = cpl_propertylist_new();
784 cpl_propertylist_append_string(header,
"EXTNAME", EXTNAME_STAT);
785 cpl_propertylist_set_comment(header,
"EXTNAME", EXTNAME_STAT_COMMENT);
788 kMuseFluxUnitString, strlen(kMuseFluxUnitString) + 1)) {
790 cpl_propertylist_append_string(header,
"BUNIT", kMuseFluxStatString);
794 cpl_propertylist_copy_property_regexp(header, aCube->
header,
798 rc = cpl_imagelist_save(aCube->
stat, aFilename, CPL_TYPE_FLOAT, header,
800 cpl_propertylist_delete(header);
816 static cpl_propertylist *
817 muse_datacube_load_header(
const char *aFilename)
819 cpl_ensure(aFilename, CPL_ERROR_NULL_INPUT, NULL);
820 int extdata = cpl_fits_find_extension(aFilename,
"DATA");
821 cpl_ensure(extdata >= 0, CPL_ERROR_ILLEGAL_INPUT, NULL);
822 cpl_ensure(extdata > 0, CPL_ERROR_DATA_NOT_FOUND, NULL);
824 cpl_propertylist *header = cpl_propertylist_load(aFilename, 0);
825 cpl_propertylist *hdata = cpl_propertylist_load(aFilename, extdata);
826 cpl_propertylist_copy_property_regexp(header, hdata,
828 cpl_propertylist_delete(hdata);
866 cpl_ensure(aFilename, CPL_ERROR_NULL_INPUT, NULL);
869 cpl_errorstate state = cpl_errorstate_get();
870 cube->
header = muse_datacube_load_header(aFilename);
871 if (!cpl_errorstate_is_equal(state) || !cube->
header) {
872 cpl_msg_error(__func__,
"Loading cube-like headers from \"%s\" failed!",
877 int ext = cpl_fits_find_extension(aFilename,
"DATA");
878 cube->
data = cpl_imagelist_load(aFilename, CPL_TYPE_FLOAT, ext);
880 ext = cpl_fits_find_extension(aFilename,
"DQ");
882 cube->
stat = cpl_imagelist_load(aFilename, CPL_TYPE_INT, ext);
884 ext = cpl_fits_find_extension(aFilename,
"STAT");
886 cube->
stat = cpl_imagelist_load(aFilename, CPL_TYPE_FLOAT, ext);
888 int next = cpl_fits_count_extensions(aFilename);
889 while (++ext <= next) {
892 image->
header = cpl_propertylist_load(aFilename, ext);
893 image->
data = cpl_image_load(aFilename, CPL_TYPE_UNSPECIFIED, 0, ext);
896 char *extname2 = cpl_sprintf(
"%s_DQ", extname);
897 int ext2 = cpl_fits_find_extension(aFilename, extname2);
900 image->
dq = cpl_image_load(aFilename, CPL_TYPE_INT, 0, ext2);
903 extname2 = cpl_sprintf(
"%s_STAT", extname);
904 ext2 = cpl_fits_find_extension(aFilename, extname2);
907 image->
stat = cpl_image_load(aFilename, CPL_TYPE_UNSPECIFIED, 0, ext2);
912 cube->
recnames = cpl_array_new(1, CPL_TYPE_STRING);
917 cpl_array_set_string(cube->
recnames, cpl_array_get_size(cube->
recnames) - 1,
953 cpl_ensure_code(aCube && aAppend, CPL_ERROR_NULL_INPUT);
955 cpl_size ndata = cpl_imagelist_get_size(aCube->
data),
956 nstat = cpl_imagelist_get_size(aCube->
stat);
957 cpl_ensure_code(ndata == nstat, CPL_ERROR_ILLEGAL_INPUT);
958 cpl_size i, n = cpl_imagelist_get_size(aAppend->
data);
959 cpl_ensure_code(n == cpl_imagelist_get_size(aAppend->
stat),
960 CPL_ERROR_ILLEGAL_INPUT);
963 cpl_size nx1 = cpl_image_get_size_x(cpl_imagelist_get(aCube->
data, ndata - 1)),
964 ny1 = cpl_image_get_size_y(cpl_imagelist_get(aCube->
data, ndata - 1)),
965 nx2 = cpl_image_get_size_x(cpl_imagelist_get(aAppend->
data, 0)),
966 ny2 = cpl_image_get_size_y(cpl_imagelist_get(aAppend->
data, 0));
967 cpl_ensure_code(nx1 == nx2 && ny1 == ny2, CPL_ERROR_ILLEGAL_INPUT);
970 cpl_ensure_code(ctype1 && ctype2, CPL_ERROR_UNSUPPORTED_MODE);
971 cpl_ensure_code((!strncmp(ctype1,
"AWAV", 5) && !strncmp(ctype2,
"AWAV", 5)) ||
972 (!strncmp(ctype1,
"WAVE", 5) && !strncmp(ctype2,
"WAVE", 5)),
973 CPL_ERROR_UNSUPPORTED_MODE);
982 lbda1 = ((ndata - 1) + 1. - pix1) * cd1 + val1,
983 lbda2 = (0 + 1. - pix2) * cd2 + val2;
985 cpl_msg_debug(__func__,
"lambdas: %f %f (%f %f)", lbda1, lbda2, cd1, cd2);
986 cpl_ensure_code(fabs(cd1 - cd2) < DBL_EPSILON &&
987 fabs(lbda2 - cd2 - lbda1) < DBL_EPSILON,
988 CPL_ERROR_ILLEGAL_INPUT);
999 cpl_boolean usedq = CPL_FALSE;
1000 if (aCube->
dq && cpl_imagelist_get_size(aCube->
dq) == ndata &&
1001 aAppend->
dq && cpl_imagelist_get_size(aAppend->
dq) == n) {
1004 cpl_imagelist_delete(aCube->
dq);
1009 for (i = 0; i < n; i++) {
1010 cpl_imagelist_set(aCube->
data,
1011 cpl_image_duplicate(cpl_imagelist_get(aAppend->
data, i)),
1012 cpl_imagelist_get_size(aCube->
data));
1014 cpl_imagelist_set(aCube->
dq,
1015 cpl_image_duplicate(cpl_imagelist_get(aAppend->
dq, i)),
1016 cpl_imagelist_get_size(aCube->
dq));
1018 cpl_imagelist_set(aCube->
stat,
1019 cpl_image_duplicate(cpl_imagelist_get(aAppend->
stat, i)),
1020 cpl_imagelist_get_size(aCube->
stat));
1022 return CPL_ERROR_NONE;
1046 cpl_imagelist_delete(aCube->
data);
1048 cpl_imagelist_delete(aCube->
dq);
1050 cpl_imagelist_delete(aCube->
stat);
1054 cpl_propertylist_delete(aCube->
header);
Structure definition of a MUSE datacube.
Structure definition for a collection of muse_images.
double muse_pfits_get_cd(const cpl_propertylist *aHeaders, unsigned int aAxisI, unsigned int aAxisJ)
find out the WCS coordinate at the reference point
const char * muse_pfits_get_extname(const cpl_propertylist *aHeaders)
find out the extension name
cpl_error_code muse_datacube_concat(muse_datacube *aCube, const muse_datacube *aAppend)
Concatenate one datacube at the end of another one.
double muse_pfits_get_crval(const cpl_propertylist *aHeaders, unsigned int aAxis)
find out the WCS coordinate at the reference point
A structure containing a spatial two-axis WCS.
cpl_error_code muse_euro3dcube_save(muse_euro3dcube *aEuro3D, const char *aFilename)
Save a Euro3D cube object to a file.
cpl_image * data
the data extension
muse_datacube * muse_datacube_load(const char *aFilename)
Load header, DATA and optionally STAT and DQ extensions as well as the reconstructed images of a MUSE...
cpl_error_code muse_utils_filter_copy_properties(cpl_propertylist *aHeader, const muse_table *aFilter, double aFraction)
Add/propagate filter properties to header of collapsed image.
static void muse_wcs_pixel_from_projplane_fast(muse_wcs *aWCS, double aX, double aY, double *aXOut, double *aYOut)
Convert from projection plane coordinates to pixel coordinates.
cpl_image * stat
the statistics extension
void muse_datacube_delete(muse_datacube *aCube)
Deallocate memory associated to a muse_datacube object.
void muse_imagelist_delete(muse_imagelist *aList)
Free the memory of the MUSE image list.
cpl_propertylist * hgroup
the group FITS header
cpl_array * recnames
the reconstructed image filter names
const char * muse_pfits_get_bunit(const cpl_propertylist *aHeaders)
find out the unit string
muse_imagelist * recimages
the reconstructed image data
cpl_propertylist * header
the primary FITS header
Structure definition of MUSE three extension FITS file.
cpl_array * recnames
the reconstructed image filter names
cpl_propertylist * header
the FITS header
cpl_error_code muse_utils_set_hduclass(cpl_propertylist *aHeader, const char *aClass2, const char *aExtData, const char *aExtDQ, const char *aExtStat)
Set HDU headers for the ESO FITS data format.
cpl_error_code muse_datacube_save(muse_datacube *aCube, const char *aFilename)
Save the three cube extensions and the FITS headers of a MUSE datacube to a file. ...
unsigned int muse_imagelist_get_size(muse_imagelist *aList)
Return the number of stored images.
cpl_image * dq
the data quality extension
double muse_pfits_get_crpix(const cpl_propertylist *aHeaders, unsigned int aAxis)
find out the WCS reference point
double muse_flux_response_interpolate(const cpl_table *aResponse, double aLambda, double *aError, muse_flux_interpolation_type aType)
Compute linearly interpolated response of some kind at given wavelength.
muse_wcs * muse_wcs_new(cpl_propertylist *aHeader)
Create a new WCS structure from a given FITS header.
cpl_error_code muse_datacube_convert_dq(muse_datacube *aCube)
Convert the DQ extension of a datacube to NANs in DATA and STAT.
#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.
static void muse_wcs_pixel_from_celestial_fast(muse_wcs *aWCS, double aRA, double aDEC, double *aX, double *aY)
Convert from celestial spherical coordinates to pixel coordinates.
cpl_table * table
The table.
cpl_error_code muse_datacube_save_recimages(const char *aFilename, muse_imagelist *aImages, cpl_array *aNames)
Save reconstructed images of a cube in extra extensions.
void muse_euro3dcube_delete(muse_euro3dcube *aEuro3D)
Deallocate memory associated to a muse_euro3dcube object.
cpl_imagelist * data
the cube containing the actual data values
Structure definition of a Euro3D datacube.
Structure to store a table together with a property list.
double muse_utils_filter_fraction(const muse_table *aFilter, double aLambda1, double aLambda2)
Compute fraction of filter area covered by the data range.
cpl_imagelist * dq
the optional cube containing the bad pixel status
cpl_table * dtable
the table containing the actual Euro3D data
const char * muse_pfits_get_ctype(const cpl_propertylist *aHeaders, unsigned int aAxis)
find out the WCS axis type string
cpl_propertylist * hdata
the data FITS header
cpl_table * gtable
the table containing the Euro3D groups
cpl_propertylist * header
the FITS header
muse_imagelist * muse_imagelist_new(void)
Create a new (empty) MUSE image list.
muse_image * muse_datacube_collapse(muse_datacube *aCube, const muse_table *aFilter)
Integrate a FITS NAXIS=3 datacube along the wavelength direction.
muse_image * muse_image_new(void)
Allocate memory for a new muse_image object.
cpl_error_code muse_imagelist_set(muse_imagelist *aList, muse_image *aImage, unsigned int aIdx)
Set the muse_image of given list index.
muse_imagelist * recimages
the reconstructed image data
cpl_error_code muse_utils_copy_modified_header(cpl_propertylist *aH1, cpl_propertylist *aH2, const char *aKeyword, const char *aString)
Copy a modified header keyword from one header to another.
muse_image * muse_euro3dcube_collapse(muse_euro3dcube *aEuro3D, const muse_table *aFilter)
Integrate a Euro3D datacube along the wavelength direction.
cpl_error_code muse_image_dq_to_nan(muse_image *aImage)
Convert pixels flagged in the DQ extension to NANs in DATA (and STAT, if present).
cpl_imagelist * stat
the cube containing the data variance