MUSE Pipeline Reference Manual  0.18.5
muse_postproc.c
1 /* -*- Mode: C; tab-width: 2; indent-tabs-mode: nil; c-basic-offset: 2 -*- */
2 /* vim:set sw=2 sts=2 et cin: */
3 /*
4  * This file is part of the MUSE Instrument Pipeline
5  * Copyright (C) 2005-2014 European Southern Observatory
6  *
7  * This program is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU General Public License as published by
9  * the Free Software Foundation; either version 2 of the License, or
10  * (at your option) any later version.
11  *
12  * This program is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  * GNU General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with this program; if not, write to the Free Software
19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
20  */
21 
22 #ifdef HAVE_CONFIG_H
23 #include <config.h>
24 #endif
25 
26 /*----------------------------------------------------------------------------*
27  * Includes *
28  *----------------------------------------------------------------------------*/
29 #include <string.h>
30 #include <strings.h>
31 #include "muse_postproc.h"
32 #include "muse_instrument.h"
33 
34 #include "muse_astro.h"
35 #include "muse_dar.h"
36 #include "muse_flux.h"
37 #include "muse_imagelist.h"
38 #include "muse_pfits.h"
39 #include "muse_utils.h"
40 #include "muse_wcs.h"
41 
42 /*----------------------------------------------------------------------------*/
49 /*----------------------------------------------------------------------------*/
50 
53 /*----------------------------------------------------------------------------*/
64 /*----------------------------------------------------------------------------*/
67 {
68  muse_postproc_properties *prop = cpl_calloc(1, sizeof(muse_postproc_properties));
69  switch (aType) {
70  case MUSE_POSTPROC_SCIPOST:
71  case MUSE_POSTPROC_STANDARD:
72  case MUSE_POSTPROC_ASTROMETRY:
73  prop->type = aType;
74  break;
75  default:
76  cpl_msg_error(__func__, "No such setup known: %d", aType);
77  cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT);
78  cpl_free(prop);
79  return NULL;
80  }
81  return prop;
82 } /* muse_postproc_properties_new() */
83 
84 /*----------------------------------------------------------------------------*/
90 /*----------------------------------------------------------------------------*/
91 void
93 {
94  if (!aProp) {
95  return;
96  }
97  cpl_table_delete(aProp->exposures);
98  cpl_table_delete(aProp->response);
99  cpl_table_delete(aProp->extinction);
100  cpl_table_delete(aProp->telluric);
101  cpl_propertylist_delete(aProp->wcs);
102  muse_sky_master_delete(aProp->sky);
103  muse_mask_delete(aProp->sky_mask);
104  cpl_frame_delete(aProp->refframe);
105  cpl_table_delete(aProp->tellregions);
106  cpl_free(aProp);
107 } /* muse_postproc_properties_delete() */
108 
109 /*----------------------------------------------------------------------------*/
124 /*----------------------------------------------------------------------------*/
126 muse_postproc_get_resampling_type(const char *aResampleString)
127 {
128  cpl_ensure(aResampleString, CPL_ERROR_NULL_INPUT, MUSE_RESAMPLE_WEIGHTED_DRIZZLE);
129  if (!strncmp(aResampleString, "nearest", 8)) {
130  return MUSE_RESAMPLE_NEAREST;
131  }
132  if (!strncmp(aResampleString, "linear", 7)) {
134  }
135  if (!strncmp(aResampleString, "quadratic", 10)) {
137  }
138  if (!strncmp(aResampleString, "renka", 6)) {
140  }
141  if (!strncmp(aResampleString, "drizzle", 8)) {
143  }
144  if (!strncmp(aResampleString, "lanczos", 8)) {
146  }
147  return MUSE_RESAMPLE_WEIGHTED_DRIZZLE; /* catch all fallback */
148 } /* muse_postproc_get_resampling_type() */
149 
150 /*----------------------------------------------------------------------------*/
165 /*----------------------------------------------------------------------------*/
167 muse_postproc_get_cr_type(const char *aCRTypeString)
168 {
169  cpl_ensure(aCRTypeString, CPL_ERROR_NULL_INPUT, MUSE_RESAMPLING_CRSTATS_IRAF);
170  if (!strncmp(aCRTypeString, "iraf", 5)) {
172  }
173  if (!strncmp(aCRTypeString, "mean", 5)) {
175  }
176  if (!strncmp(aCRTypeString, "median", 7)) {
178  }
179  return MUSE_RESAMPLING_CRSTATS_IRAF; /* catch all fallback */
180 } /* muse_postproc_get_cr_type() */
181 
182 /*----------------------------------------------------------------------------*/
196 /*----------------------------------------------------------------------------*/
198 muse_postproc_get_cube_format(const char *aFormatString)
199 {
200  cpl_ensure(aFormatString, CPL_ERROR_NULL_INPUT, MUSE_CUBE_TYPE_FITS);
201  if (!strncmp(aFormatString, "Cube", 5)) {
202  return MUSE_CUBE_TYPE_FITS;
203  }
204  if (!strncmp(aFormatString, "Euro3D", 7)) {
205  return MUSE_CUBE_TYPE_EURO3D;
206  }
207  if (!strncmp(aFormatString, "xCube", 6)) {
208  return MUSE_CUBE_TYPE_FITS_X;
209  }
210  if (!strncmp(aFormatString, "xEuro3D", 8)) {
212  }
213  return MUSE_CUBE_TYPE_FITS; /* catch all fallback */
214 } /* muse_postproc_get_cube_format() */
215 
216 /*----------------------------------------------------------------------------*/
230 /*----------------------------------------------------------------------------*/
232 muse_postproc_get_weight_type(const char *aWeightString)
233 {
234  cpl_ensure(aWeightString, CPL_ERROR_NULL_INPUT, MUSE_XCOMBINE_EXPTIME);
235  if (!strncmp(aWeightString, "exptime", 8)) {
236  return MUSE_XCOMBINE_EXPTIME;
237  }
238  if (!strncmp(aWeightString, "fwhm", 5)) {
239  return MUSE_XCOMBINE_FWHM;
240  }
241  return MUSE_XCOMBINE_EXPTIME; /* catch all fallback */
242 } /* muse_postproc_get_weight_type() */
243 
244 /*----------------------------------------------------------------------------*/
267 /*----------------------------------------------------------------------------*/
268 static cpl_table *
269 muse_postproc_load_nearest(const cpl_propertylist *aHeader,
270  const cpl_frame *aFrame, float aWarnLimit,
271  float aErrLimit)
272 {
273  cpl_ensure(aHeader && aFrame, CPL_ERROR_NULL_INPUT, NULL);
274  cpl_errorstate state = cpl_errorstate_get();
275  double raref = muse_pfits_get_ra(aHeader),
276  decref = muse_pfits_get_dec(aHeader);
277  cpl_ensure(cpl_errorstate_is_equal(state), CPL_ERROR_DATA_NOT_FOUND, NULL);
278  cpl_msg_debug(__func__, "reference coordinates: RA = %e deg, DEC =%e deg",
279  raref, decref);
280 
281  const char *fn = cpl_frame_get_filename(aFrame);
282  cpl_propertylist *header;
283  double dmin = FLT_MAX; /* minimum distance yet found in deg */
284  int iext, inearest = -1, next = cpl_fits_count_extensions(fn);
285  for (iext = 1; iext <= next; iext++) {
286  header = cpl_propertylist_load(fn, iext);
287  const char *extname = cpl_propertylist_get_string(header, "EXTNAME");
288  double ra = muse_pfits_get_ra(header),
289  dec = muse_pfits_get_dec(header),
290  d = muse_astro_angular_distance(ra, dec, raref, decref);
291  cpl_msg_debug(__func__, "extension %d [%s]: RA = %e deg, DEC = %e deg --> "
292  "d = %e deg", iext, extname, ra, dec, d);
293  if (d < dmin) {
294  dmin = d;
295  inearest = iext;
296  }
297  cpl_propertylist_delete(header);
298  } /* for iext */
299  /* warn for distances larger than x arcsec */
300  if (dmin * 3600. > aErrLimit) {
301  cpl_error_set_message(__func__, CPL_ERROR_INCOMPATIBLE_INPUT, "Distance of "
302  "nearest reference table to observed position is %.2f "
303  "arcmin, certainly a wrong reference object!", dmin * 60.);
304  return NULL;
305  }
306  if (dmin * 3600. > aWarnLimit) {
307  cpl_msg_warning(__func__, "Distance of nearest reference table to observed "
308  "position is %.2f arcmin! (Wrong reference object?)",
309  dmin * 60.);
310  }
311  header = cpl_propertylist_load(fn, inearest);
312  const char *extname = cpl_propertylist_get_string(header, "EXTNAME");
313  cpl_msg_info(__func__, "Loading \"%s[%s]\" (distance %.1f arcsec)", fn,
314  extname, dmin * 3600.);
315  cpl_propertylist_delete(header);
316  cpl_table *table = cpl_table_load(fn, inearest, 1);
317  return table;
318 } /* muse_postproc_load_nearest() */
319 
320 /*----------------------------------------------------------------------------*/
358 /*----------------------------------------------------------------------------*/
359 void *
361  unsigned int aIndex)
362 {
363  cpl_ensure(aProp && aProp->exposures, CPL_ERROR_NULL_INPUT, NULL);
364  cpl_ensure((int)aIndex < cpl_table_get_nrow(aProp->exposures),
365  CPL_ERROR_ILLEGAL_INPUT, NULL);
366  cpl_msg_info(__func__, "Starting to process exposure %d (DATE-OBS=%s)",
367  aIndex + 1, cpl_table_get_string(aProp->exposures, "DATE-OBS",
368  aIndex));
369  cpl_ensure(aProp->exposures, CPL_ERROR_NULL_INPUT, NULL);
370 
371  cpl_table *exposure = cpl_table_extract(aProp->exposures, aIndex, 1);
373  aProp->lambdamin,
374  aProp->lambdamax);
375  cpl_table_delete(exposure);
376  if (!pt) {
377  return NULL;
378  }
379  /* erase pre-existing QC parameters */
380  cpl_propertylist_erase_regexp(pt->header, "ESO QC ", 0);
381 
382  /* do flux calibration (before sky subtraction!), if possible */
383  cpl_error_code rc = CPL_ERROR_NONE;
384  if (aProp->type != MUSE_POSTPROC_STANDARD) {
385  if (aProp->response) {
386  rc = muse_flux_calibrate(pt, aProp->response, aProp->extinction,
387  aProp->telluric);
388  cpl_msg_debug(__func__, "Flux calibration returned rc=%d: %s", rc,
389  rc != CPL_ERROR_NONE ? cpl_error_get_message_default(rc) : "");
390  } else {
391  cpl_msg_info(__func__, "Skipping flux calibration (no %s curve)",
392  MUSE_TAG_STD_RESPONSE);
393  } /* else */
394  } /* if not POSTPROC_STANDARD */
395 
396  if (aProp->skymethod == MUSE_POSTPROC_SKYMETHOD_MODEL && aProp->sky) {
397  /* Process optional sky mask */
399  if (aProp->sky_mask != NULL) {
400  cpl_table_select_all(sky_pt->table);
402  }
403  /* Create white light image and sky mask */
404  cpl_msg_info(__func__, "Intermediate resampling to produce white-light image"
405  " for sky-spectrum creation:");
406  cpl_msg_indent_more();
409  params->crsigma = 5.;
410  muse_datacube *cube = muse_resampling_cube(sky_pt, params, NULL);
411  cpl_table *fwhite = muse_table_load_filter(NULL, "white");
412  muse_image *whitelight = muse_datacube_collapse(cube, fwhite);
413  cpl_msg_indent_less();
414  muse_mask *sky_mask = muse_sky_create_skymask(whitelight, aProp->skymodel_params.fraction);
415  /* save mask and/or image ?? */
416 
417  /* Apply mask and create spectrum */
418  cpl_table_select_all(sky_pt->table);
419  muse_pixtable_and_selected_mask(sky_pt, sky_mask);
420  cpl_table *spectrum = muse_resampling_spectrum(sky_pt, aProp->skymodel_params.sampling);
421 
422  /* Fit sky lines spectrum */
423  cpl_array *lambda = muse_cpltable_extract_column(spectrum, "lambda");
424  double lambda_low = cpl_array_get_min(lambda);
425  double lambda_high = cpl_array_get_max(lambda);
426  muse_sky_lines_set_range(aProp->sky->lines, lambda_low-10, lambda_high+10);
427 
428  cpl_array *data = muse_cpltable_extract_column(spectrum, "data");
429  cpl_array *stat = muse_cpltable_extract_column(spectrum, "stat");
430  // do the fit, ignoring possible errors
431  cpl_errorstate prestate = cpl_errorstate_get();
432  muse_sky_master *master = muse_sky_master_fit(lambda, data, stat, aProp->sky->lines);
433  if (!cpl_errorstate_is_equal(prestate)) {
434  cpl_errorstate_dump(prestate, CPL_FALSE, NULL);
435  cpl_errorstate_set(prestate);
436  }
437 
438  cpl_table_delete(aProp->sky->lines);
439  aProp->sky->lines = master->lines;
440  master->lines = NULL;
441 
442  /* Create continuum */
443  if (aProp->sky->continuum == NULL) {
444  muse_sky_subtract_pixtable(sky_pt, aProp->sky, aProp->sky->lsf);
445  aProp->sky->continuum = muse_resampling_spectrum(sky_pt, aProp->skymodel_params.csampling);
446  cpl_table_erase_column(aProp->sky->continuum, "stat");
447  cpl_table_erase_column(aProp->sky->continuum, "dq");
448  cpl_table_name_column(aProp->sky->continuum, "data", "flux");
449  }
450 
451  /* clean up */
452  cpl_array_unwrap(lambda);
453  cpl_array_unwrap(data);
454  cpl_array_unwrap(stat);
455  muse_sky_master_delete(master);
456  cpl_table_delete(spectrum);
457  muse_mask_delete(sky_mask);
458  muse_image_delete(whitelight);
459  cpl_table_delete(fwhite);
461  muse_datacube_delete(cube);
462  muse_pixtable_delete(sky_pt);
463  } /* if MUSE_POSTPROC_SKYMETHOD_MODEL */
464  if (((aProp->skymethod == MUSE_POSTPROC_SKYMETHOD_NONE) ||
465  (aProp->skymethod == MUSE_POSTPROC_SKYMETHOD_MODEL)) &&
466  (aProp->sky && aProp->sky->lines && aProp->sky->continuum)) {
467  cpl_msg_debug(__func__, "doing sky subtraction by spectrum modeling");
468  rc = muse_sky_subtract_pixtable(pt, aProp->sky, aProp->sky->lsf);
469  } /* if sky modeling */
470  if (aProp->skymethod == MUSE_POSTPROC_SKYMETHOD_ROWBYROW) {
471  /* apply the row-by-row sky subtraction; needed for row-by-row */
472  cpl_msg_debug(__func__, "doing sky subtraction using row-by-row fitting");
473  /* convert pixel table into images, one per IFU */
475  unsigned int k, n = muse_imagelist_get_size(list);
476  #pragma omp parallel for default(none) \
477  shared(aProp, list, n)
478  for (k = 0; k < n; k++) {
479  muse_image *image = muse_imagelist_get(list, k);
480  muse_sky_subtract_rowbyrow_mask(image, NULL);
481  muse_sky_subtract_rowbyrow(image, NULL, 3., 4);
482  } /* for all IFUs / images */
484  muse_imagelist_delete(list);
485  } /* if row-by-row sky subtraction */
486 
487  /* XXX try to recognize a lab exposure by missing RA/DEC or weird *
488  * airmass; in the future this should be a template-based check */
489  cpl_boolean labdata = CPL_FALSE;
490  cpl_errorstate prestate = cpl_errorstate_get();
491  double airmass = muse_astro_airmass(pt->header);
492  double ra = muse_pfits_get_ra(pt->header),
493  dec = muse_pfits_get_dec(pt->header);
494  if (!cpl_errorstate_is_equal(prestate) || airmass < 1.) {
495  cpl_msg_warning(__func__, "This seems to be lab data (RA=%f DEC=%f, "
496  "airmass=%f), not doing any DAR correction!", ra, dec,
497  airmass);
498  cpl_errorstate_dump(prestate, CPL_FALSE, NULL);
499  cpl_errorstate_set(prestate);
500  labdata = CPL_TRUE;
501  }
502  if (!labdata && muse_pfits_get_mode(pt->header) <= MUSE_MODE_WFM_AO_N) {
503  cpl_msg_debug(__func__, "WFM detected: starting DAR correction");
504  rc = muse_dar_correct(pt, aProp->lambdaref);
505  cpl_msg_debug(__func__, "DAR correction returned rc=%d: %s", rc,
506  rc != CPL_ERROR_NONE ? cpl_error_get_message() : "");
507  /* check and possibly correct the DAR quality, if requested to do so */
508  if (aProp->darcheck != MUSE_POSTPROC_DARCHECK_NONE) {
509  cpl_boolean dorefine = aProp->darcheck == MUSE_POSTPROC_DARCHECK_CORRECT;
510  cpl_msg_info(__func__, "Carrying out DAR %s", dorefine ? "correction" : "check");
511  double maxshift = 0;
512  rc = muse_dar_check(pt, &maxshift, dorefine, NULL);
513  if (rc != CPL_ERROR_NONE) {
514  cpl_msg_warning(__func__, "Maximum detected shift %f\'\' (rc = %d: %s)",
515  maxshift, rc, cpl_error_get_message());
516  } else {
517  cpl_msg_info(__func__, "Maximum detected shift %f\'\'", maxshift);
518  }
519  } /* if DAR checking */
520  } /* if not labdata but WFM */
521 
522  /* do flux response computation and return */
523  if (aProp->type == MUSE_POSTPROC_STANDARD) {
524  cpl_table *reftable = muse_postproc_load_nearest(pt->header,
525  aProp->refframe, 20., 60.);
526  if (!reftable) {
527  cpl_msg_error(__func__, "Required input %s could not be loaded!",
528  MUSE_TAG_STD_FLUX_TABLE);
530  return NULL;
531  }
532  if (muse_flux_reference_table_check(reftable) != CPL_ERROR_NONE) {
533  cpl_msg_error(__func__, "%s does not have a recognized format!",
534  MUSE_TAG_STD_FLUX_TABLE);
535  cpl_table_delete(reftable);
537  return NULL;
538  }
539 
540  prestate = cpl_errorstate_get();
542  rc = muse_flux_integrate_std(pt, aProp->profile, flux);
543  if (flux->intimage) {
544  cpl_msg_debug(__func__, "Flux integration found %"CPL_SIZE_FORMAT" stars",
545  cpl_image_get_size_y(flux->intimage->data));
546  }
548  if (airmass < 0) {
549  cpl_msg_warning(__func__, "Invalid airmass derived (%.3e), resetting to "
550  " zero", airmass);
551  airmass = 0.; /* set to zero to leave out uncertain extinction term */
552  }
553  rc = muse_flux_response_compute(flux, airmass, reftable,
554  aProp->tellregions, aProp->extinction);
555  cpl_table_delete(reftable);
556  rc = muse_flux_get_response_table(flux);
557  rc = muse_flux_get_telluric_table(flux);
558 
559  if (!cpl_errorstate_is_equal(prestate)) {
560  cpl_msg_warning(__func__, "The following errors happened while computing "
561  "the flux calibration (rc = %d/%s):", rc,
562  cpl_error_get_message_default(rc));
563  cpl_errorstate_dump(prestate, CPL_FALSE, muse_cplerrorstate_dump_some);
564  }
565  return flux;
566  } /* if MUSE_POSTPROC_STANDARD */
567 
568  /* do astrometry computation and return */
569  if (aProp->type == MUSE_POSTPROC_ASTROMETRY) {
570  cpl_table *reftable = muse_postproc_load_nearest(pt->header,
571  aProp->refframe, 20., 60.);
572  if (!reftable) {
573  cpl_msg_error(__func__, "Required input %s could not be loaded!",
574  MUSE_TAG_ASTROMETRY_REFERENCE);
576  return NULL;
577  }
578  prestate = cpl_errorstate_get();
579 
581  wcs->crpix1 = aProp->crpix1; /* set center of rotation */
582  wcs->crpix2 = aProp->crpix2;
583  rc = muse_wcs_locate_sources(pt, aProp->detsigma, aProp->centroid, wcs);
584  rc = muse_wcs_solve(wcs, reftable, aProp->radius, aProp->faccuracy,
585  aProp->niter, aProp->rejsigma);
586  cpl_table_delete(reftable);
587  cpl_msg_debug(__func__, "Solving astrometry (initial radius %.1f pix, "
588  "initial relative accuracy %.1f, %d iterations with rejection"
589  " sigma %.2f) returned rc=%d: %s", aProp->radius,
590  aProp->faccuracy, aProp->niter, aProp->rejsigma,
591  rc, rc != CPL_ERROR_NONE ? cpl_error_get_message() : "");
593 
594  if (!cpl_errorstate_is_equal(prestate)) {
595  cpl_msg_warning(__func__, "The following errors happened while computing "
596  "the astrometric solution (rc = %d/%s):", rc,
597  cpl_error_get_message());
598  cpl_errorstate_dump(prestate, CPL_FALSE, muse_cplerrorstate_dump_some);
599  }
600  return wcs;
601  } /* if MUSE_POSTPROC_ASTROMETRY */
602 
603  rc = muse_wcs_project_tan(pt, aProp->wcs);
604  cpl_msg_debug(__func__, "Astrometry correction returned rc=%d: %s", rc,
605  rc != CPL_ERROR_NONE ? cpl_error_get_message() : "");
606 
607  return pt;
608 } /* muse_postproc_process_exposure() */
609 
610 /*---------------------------------------------------------------------------*/
628 /*---------------------------------------------------------------------------*/
629 cpl_propertylist *
631 {
632  cpl_ensure(aProcessing, CPL_ERROR_NULL_INPUT, NULL);
633 
634  cpl_frameset *fwcs = muse_frameset_find(aProcessing->inputFrames,
635  MUSE_TAG_OUTPUT_WCS, 0, CPL_FALSE);
636  if (!fwcs || cpl_frameset_get_size(fwcs) <= 0) {
637  /* not an error, quietly return NULL */
638  cpl_frameset_delete(fwcs);
639  return NULL;
640  }
641  cpl_frame *frame = cpl_frameset_get_position(fwcs, 0);
642  const char *fn = cpl_frame_get_filename(frame);
643  cpl_propertylist *header = NULL;
644  int iext, next = cpl_fits_count_extensions(fn);
645  for (iext = 0; iext <= next; iext++) {
646  header = cpl_propertylist_load(fn, iext);
647  cpl_errorstate es = cpl_errorstate_get();
648  cpl_wcs *wcs = cpl_wcs_new_from_propertylist(header);
649  if (!cpl_errorstate_is_equal(es)) {
650  cpl_errorstate_set(es);
651  }
652  if (!wcs) {
653  cpl_propertylist_delete(header);
654  header = NULL;
655  continue; /* try the next extension */
656  }
657  /* need 2D or 3D WCS to have something useful */
658  int naxis = cpl_wcs_get_image_naxis(wcs);
659  cpl_boolean valid = naxis == 2 || naxis == 3;
660  /* need RA---TAN and DEC--TAN for axes 1 and 2 */
661  const cpl_array *ctypes = cpl_wcs_get_ctype(wcs);
662  if (valid && cpl_array_get_string(ctypes, 0)) {
663  valid = valid && !strncmp(cpl_array_get_string(ctypes, 0), "RA---TAN", 9);
664  }
665  if (valid && cpl_array_get_string(ctypes, 1)) {
666  valid = valid && !strncmp(cpl_array_get_string(ctypes, 1), "DEC--TAN", 9);
667  }
668  /* need AWAV or AWAV-LOG for axis 3 */
669  if (valid && cpl_array_get_string(ctypes, 2)) {
670  const char *ctype3 = cpl_array_get_string(ctypes, 2);
671  valid = valid
672  && (!strncmp(ctype3, "AWAV", 5) || !strncmp(ctype3, "AWAV-LOG", 9));
673  }
674  if (valid) { /* check for unsupported PCi_j keywords */
675  /* still needs to be done on the input header not on the derives WCS! */
676  cpl_propertylist *pc = cpl_propertylist_new();
677  cpl_propertylist_copy_property_regexp(pc, header, "^PC[0-9]+_[0-9]+", 0);
678  valid = cpl_propertylist_get_size(pc) == 0;
679  cpl_propertylist_delete(pc);
680  }
681  cpl_wcs_delete(wcs);
682  if (valid) {
683  cpl_msg_debug(__func__, "Apparently valid header %s found in extension %d"
684  " of \"%s\"", MUSE_TAG_OUTPUT_WCS, iext, fn);
685  muse_processing_append_used(aProcessing, frame, CPL_FRAME_GROUP_CALIB, 1);
686  break;
687  }
688  cpl_propertylist_delete(header);
689  header = NULL;
690  } /* for iext (all extensions) */
691  if (!header) {
692  cpl_msg_warning(__func__, "No valid headers for %s found in \"%s\"",
693  MUSE_TAG_OUTPUT_WCS, fn);
694  }
695  cpl_frameset_delete(fwcs);
696  return header;
697 } /* muse_postproc_cube_load_output_wcs() */
698 
699 /*---------------------------------------------------------------------------*/
717 /*---------------------------------------------------------------------------*/
718 cpl_error_code
720  muse_pixtable *aPixtable,
721  muse_cube_type aFormat,
722  muse_resampling_params *aParams,
723  const char *aFilter)
724 {
725  cpl_ensure_code(aProcessing && aPixtable && aParams, CPL_ERROR_NULL_INPUT);
726  cpl_ensure_code(aFormat == MUSE_CUBE_TYPE_EURO3D ||
727  aFormat == MUSE_CUBE_TYPE_FITS ||
728  aFormat == MUSE_CUBE_TYPE_EURO3D_X ||
729  aFormat == MUSE_CUBE_TYPE_FITS_X, CPL_ERROR_ILLEGAL_INPUT);
730 
731  /* create cube and cast to generic pointer to save code duplication */
732  void *cube = NULL;
733  muse_pixgrid *grid = NULL;
734  if (aFormat == MUSE_CUBE_TYPE_EURO3D || aFormat == MUSE_CUBE_TYPE_EURO3D_X) {
735  cpl_msg_info(__func__, "Resampling to final cube follows, output format "
736  "\"Euro3D\"");
737  cpl_msg_indent_more();
738  muse_euro3dcube *e3d = muse_resampling_euro3d(aPixtable, aParams);
739  cpl_msg_indent_less();
740  cpl_ensure_code(e3d, cpl_error_get_code());
741  cube = e3d;
742  } else if (aFormat == MUSE_CUBE_TYPE_FITS || aFormat == MUSE_CUBE_TYPE_FITS_X) {
743  cpl_msg_info(__func__, "Resampling to final cube follows, output format "
744  "\"FITS cube\"");
745  cpl_msg_indent_more();
746  muse_datacube *fits = muse_resampling_cube(aPixtable, aParams, &grid);
747  cpl_msg_indent_less();
748  cpl_ensure_code(fits, cpl_error_get_code());
749  muse_postproc_qc_fwhm(aProcessing, fits);
750  cube = fits;
751  }
752 
753  cpl_array *filters = muse_cplarray_new_from_delimited_string(aFilter, ",");
754  int i, ipos = 0, nfilters = cpl_array_get_size(filters);
755  for (i = 0; i < nfilters; i++) {
756  /* try to find and load the filter from a table */
757  cpl_table *filter = muse_table_load_filter(aProcessing,
758  cpl_array_get_string(filters, i));
759  if (!filter) {
760  continue;
761  }
762 
763  /* if we got a filter table, do the collapsing */
764  muse_image *fov = NULL;
765  if (aFormat == MUSE_CUBE_TYPE_EURO3D || aFormat == MUSE_CUBE_TYPE_EURO3D_X) {
766  fov = muse_euro3dcube_collapse(cube, filter);
767  } else if (aFormat == MUSE_CUBE_TYPE_FITS || aFormat == MUSE_CUBE_TYPE_FITS_X) {
768  if (getenv("MUSE_COLLAPSE_PIXTABLE") &&
769  atoi(getenv("MUSE_COLLAPSE_PIXTABLE")) > 0) {
770  /* create collapsed image directly from pixel table/grid */
771  fov = muse_resampling_collapse_pixgrid(aPixtable, grid,
772  (muse_datacube *)cube,
773  filter, aParams);
774  } else {
775  fov = muse_datacube_collapse(cube, filter);
776  }
777  }
778  cpl_table_delete(filter);
779 
780  if (aFormat == MUSE_CUBE_TYPE_EURO3D_X || aFormat == MUSE_CUBE_TYPE_FITS_X) {
781  if (!((muse_datacube *)cube)->recimages) {
782  /* create empty structures before first entry */
783  ((muse_datacube *)cube)->recimages = muse_imagelist_new();
784  ((muse_datacube *)cube)->recnames = cpl_array_new(0, CPL_TYPE_STRING);
785  }
786  /* append image to cube/Euro3D structure */
787  muse_imagelist_set(((muse_datacube *)cube)->recimages, fov,
788  muse_imagelist_get_size(((muse_datacube *)cube)->recimages));
789  cpl_array_set_size(((muse_datacube *)cube)->recnames, ipos+1);
790  cpl_array_set_string(((muse_datacube *)cube)->recnames, ipos,
791  cpl_array_get_string(filters, i));
792  } else {
793  /* duplicate cube header to be able to use its WCS info as the base */
794  muse_utils_copy_modified_header(fov->header, fov->header, "OBJECT",
795  cpl_array_get_string(filters, i));
796  cpl_propertylist_update_string(fov->header, MUSE_HDR_FILTER,
797  cpl_array_get_string(filters, i));
798  muse_processing_save_image(aProcessing, -1, fov, MUSE_TAG_IMAGE_FOV);
799  muse_image_delete(fov);
800  }
801  ipos++;
802  } /* for i (all filters) */
803  cpl_array_delete(filters);
804  muse_pixgrid_delete(grid);
805 
806  cpl_error_code rc = CPL_ERROR_NONE;
807  if (aFormat == MUSE_CUBE_TYPE_EURO3D || aFormat == MUSE_CUBE_TYPE_EURO3D_X) {
808  /* save the data and clean up */
809  rc = muse_processing_save_cube(aProcessing, -1, cube, MUSE_TAG_CUBE_FINAL,
812  } else if (aFormat == MUSE_CUBE_TYPE_FITS || aFormat == MUSE_CUBE_TYPE_FITS_X) {
813  /* convert/remove the DQ cube */
815  /* save the cube and clean up */
816  rc = muse_processing_save_cube(aProcessing, -1, cube, MUSE_TAG_CUBE_FINAL,
818  muse_datacube_delete(cube);
819  }
820  return rc;
821 } /* muse_postproc_cube_resample_and_collapse() */
822 
823 /*---------------------------------------------------------------------------*/
846 /*---------------------------------------------------------------------------*/
847 cpl_error_code
849 {
850  cpl_ensure_code(aProcessing && aCube, CPL_ERROR_NULL_INPUT);
851 
852  /* determine which FITS keyword prefix to take */
853  char prefix[KEYWORD_LENGTH] = "";
854  if (!strncmp(aProcessing->recipeName, "muse_scipost", 13)) {
855  strncpy(prefix, QC_POST_PREFIX_SCIPOST, KEYWORD_LENGTH);
856  } else if (!strncmp(aProcessing->recipeName, "muse_exp_combine", 17)) {
857  strncpy(prefix, QC_POST_PREFIX_EXPCOMBINE, KEYWORD_LENGTH);
858  } else if (!strncmp(aProcessing->recipeName, "muse_standard", 14)) {
859  strncpy(prefix, QC_POST_PREFIX_STANDARD, KEYWORD_LENGTH);
860  } else if (!strncmp(aProcessing->recipeName, "muse_astrometry", 16)) {
861  strncpy(prefix, QC_POST_PREFIX_ASTROMETRY, KEYWORD_LENGTH);
862  } else {
863  cpl_msg_info(__func__, "Recipe \"%s\" found, not generating QC1 keywords",
864  aProcessing->recipeName);
865  return CPL_ERROR_NONE;
866  }
867 
868  /* get image from central plane of cube and use it to do object detection *
869  * and FWHM estimation */
870  int plane = cpl_imagelist_get_size(aCube->data) / 2;
871  cpl_image *cim = cpl_imagelist_get(aCube->data, plane);
872  double dsigmas[] = {/*50., 30., 10., 8., 6.,*/ 5., 4., 3. };
873  int ndsigmas = sizeof(dsigmas) / sizeof(double);
874  cpl_vector *vsigmas = cpl_vector_wrap(ndsigmas, dsigmas);
875  cpl_size isigma = -1;
876  cpl_errorstate prestate = cpl_errorstate_get();
877  cpl_apertures *apertures = cpl_apertures_extract(cim, vsigmas, &isigma);
878  cpl_vector_unwrap(vsigmas);
879  if (!apertures || !cpl_errorstate_is_equal(prestate)) {
880  /* just set up two fake entries, with number zero and negative FWHMs: */
881  char keyword[KEYWORD_LENGTH];
882  snprintf(keyword, KEYWORD_LENGTH, QC_POST_POSX, prefix, 0);
883  cpl_propertylist_update_float(aCube->header, keyword, -1.);
884  snprintf(keyword, KEYWORD_LENGTH, QC_POST_POSY, prefix, 0);
885  cpl_propertylist_update_float(aCube->header, keyword, -1.);
886  snprintf(keyword, KEYWORD_LENGTH, QC_POST_FWHMX, prefix, 0);
887  cpl_propertylist_update_float(aCube->header, keyword, -1.);
888  snprintf(keyword, KEYWORD_LENGTH, QC_POST_FWHMY, prefix, 0);
889  cpl_propertylist_update_float(aCube->header, keyword, -1.);
890  /* swallow the cpl_aperture error, which is too cryptic for users */
891  cpl_errorstate_set(prestate);
892  cpl_msg_warning(__func__, "No sources found for FWHM measurement down to "
893  "%.1f sigma limit on plane %d, QC parameters will not "
894  "contain useful information", dsigmas[ndsigmas-1], plane+1);
895  return CPL_ERROR_DATA_NOT_FOUND;
896  }
897  int ndet = cpl_apertures_get_size(apertures);
898  cpl_msg_info(__func__, "Computing FWHM QC parameters for %d source%s found "
899  "down to the %.1f sigma threshold on plane %d", ndet,
900  ndet == 1 ? "" : "s", dsigmas[isigma], plane + 1);
901 
902  /* get some kind of WCS for conversion of FWHM from pixels to arcsec, *
903  * by default assume WFM and x=RA and y=DEC */
904  double cd11 = kMuseSpaxelSizeX_WFM, cd12 = 0., cd21 = 0.,
905  cd22 = kMuseSpaxelSizeY_WFM;
906  cpl_errorstate es = cpl_errorstate_get();
907  cpl_wcs *wcs = cpl_wcs_new_from_propertylist(aCube->header);
908  if (!cpl_errorstate_is_equal(es)) {
909  cpl_errorstate_set(es); /* swallow possible errors from creating the WCS */
910  }
911  if (!wcs ||
912  !strncasecmp(cpl_propertylist_get_string(aCube->header, "CTYPE1"),
913  "PIXEL", 5)) {
914  cpl_msg_debug(__func__, "QC1 FWHM parameter estimation (%d sources): "
915  "simple conversion to arcsec (CTYPE=%s/%s)!", ndet,
916  cpl_propertylist_get_string(aCube->header, "CTYPE1"),
917  cpl_propertylist_get_string(aCube->header, "CTYPE2"));
918  if (muse_pfits_get_mode(aCube->header) > MUSE_MODE_WFM_AO_N) { /* NFM scaling */
919  cd11 = kMuseSpaxelSizeX_WFM;
920  cd22 = kMuseSpaxelSizeY_WFM;
921  }
922  } else {
923  const cpl_matrix *cd = cpl_wcs_get_cd(wcs);
924  /* take the absolute and scale by 3600 to get positive arcseconds */
925  cd11 = fabs(cpl_matrix_get(cd, 0, 0)) * 3600.,
926  cd12 = fabs(cpl_matrix_get(cd, 0, 1)) * 3600.,
927  cd21 = fabs(cpl_matrix_get(cd, 1, 0)) * 3600.,
928  cd22 = fabs(cpl_matrix_get(cd, 1, 1)) * 3600.;
929  cpl_msg_debug(__func__, "QC1 FWHM parameter estimation (%d sources): "
930  "full conversion to arcsec (CD=%f,%f,%f,%f)", ndet,
931  cd11, cd12, cd21, cd22);
932  }
933  cpl_wcs_delete(wcs); /* rely on internal NULL check */
934 
935  /* Compute FWHM of all apertures */
936 #if 0
937  printf("index RA_FWHM DEC_FWHM\n");
938 #endif
939  int n;
940  for (n = 1; n <= ndet; n++) {
941  cpl_size xcen = lround(cpl_apertures_get_centroid_x(apertures, n)),
942  ycen = lround(cpl_apertures_get_centroid_y(apertures, n));
943  double x, y;
944  cpl_image_get_fwhm(cim, xcen, ycen, &x, &y);
945 
946  /* linear WCS scaling */
947  x = cd11 * x + cd12 * y;
948  y = cd22 * y + cd21 * x;
949 #if 0
950  printf("%4d %f %f\n", n, x, y);
951 #endif
952 
953  /* generate and write FITS keywords */
954  char keyword[KEYWORD_LENGTH];
955  snprintf(keyword, KEYWORD_LENGTH, QC_POST_POSX, prefix, n);
956  cpl_propertylist_update_float(aCube->header, keyword, xcen);
957  snprintf(keyword, KEYWORD_LENGTH, QC_POST_POSY, prefix, n);
958  cpl_propertylist_update_float(aCube->header, keyword, ycen);
959  snprintf(keyword, KEYWORD_LENGTH, QC_POST_FWHMX, prefix, n);
960  cpl_propertylist_update_float(aCube->header, keyword, x);
961  snprintf(keyword, KEYWORD_LENGTH, QC_POST_FWHMY, prefix, n);
962  cpl_propertylist_update_float(aCube->header, keyword, y);
963  } /* for n (aperture number) */
964 #if 0
965  fflush(stdout);
966 #endif
967  cpl_apertures_delete(apertures);
968 
969  return CPL_ERROR_NONE;
970 } /* muse_postproc_qc_fwhm() */
971 
Structure definition of a MUSE datacube.
Definition: muse_datacube.h:48
cpl_propertylist * muse_postproc_cube_load_output_wcs(muse_processing *aProcessing)
Find a file with a usable output WCS in the input frameset.
Structure definition for a collection of muse_images.
muse_postproc_properties * muse_postproc_properties_new(muse_postproc_type aType)
Create a post-processing properties object.
Definition: muse_postproc.c:66
void muse_image_delete(muse_image *aImage)
Deallocate memory associated to a muse_image object.
Definition: muse_image.c:84
muse_pixtable * muse_pixtable_duplicate(muse_pixtable *aPixtable)
Make a copy of the pixtanle.
muse_wcs_centroid_type centroid
double muse_pfits_get_ra(const cpl_propertylist *aHeaders)
find out the right ascension
Definition: muse_pfits.c:206
cpl_image * data
the data extension
Definition: muse_image.h:47
muse_xcombine_types muse_postproc_get_weight_type(const char *aWeightString)
Select correct weighting type for weight string.
cpl_error_code muse_dar_correct(muse_pixtable *aPixtable, double aLambdaRef)
Correct the pixel coordinates of all pixels of a given pixel table for differential atmospheric refra...
Definition: muse_dar.c:306
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.
Structure to hold the MASTER SKY result.
Definition: muse_sky.h:59
muse_resampling_crstats_type muse_postproc_get_cr_type(const char *aCRTypeString)
Select correct cosmic ray rejection type for crtype string.
cpl_error_code muse_sky_subtract_rowbyrow(muse_image *, cpl_table *, float, unsigned int)
Subtract the sky row-by-row from a CCD-based image.
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.
const char * recipeName
muse_image * muse_datacube_collapse(muse_datacube *aCube, cpl_table *aFilter)
Integrate a FITS NAXIS=3 datacube along the wavelength direction.
cpl_error_code muse_dar_check(muse_pixtable *aPixtable, double *aMaxShift, cpl_boolean aCorrect, muse_datacube **aCube)
check that the continuum of objects in the frame is well-aligned perpendicular to the spatial axis...
Definition: muse_dar.c:641
cpl_error_code muse_flux_reference_table_check(cpl_table *aTable)
Check and/or adapt the standard flux reference table format.
Definition: muse_flux.c:162
cpl_error_code muse_postproc_qc_fwhm(muse_processing *aProcessing, muse_datacube *aCube)
Compute QC1 parameters for datacubes and save them in the FITS header.
muse_mask * muse_sky_create_skymask(muse_image *, double)
Select spaxels to be considered as sky.
WCS object to store data needed while computing the astrometric solution.
Definition: muse_wcs.h:58
The pixel grid.
Definition: muse_pixgrid.h:56
Structure definition of MUSE three extension FITS file.
Definition: muse_image.h:41
muse_flux_object * muse_flux_object_new(void)
Allocate memory for a new muse_flux_object object.
Definition: muse_flux.c:66
cpl_error_code muse_postproc_cube_resample_and_collapse(muse_processing *aProcessing, muse_pixtable *aPixtable, muse_cube_type aFormat, muse_resampling_params *aParams, const char *aFilter)
High level function to resample to a datacube and collapse that to an image of the field of view and ...
cpl_table * table
The pixel table.
muse_flux_profile_type profile
muse_postproc_skymethod skymethod
Definition: muse_postproc.h:98
cpl_propertylist * header
the FITS header
Definition: muse_image.h:73
cpl_table * muse_table_load_filter(muse_processing *aProcessing, const char *aFilterName)
Load a table for a given filter name.
Definition: muse_utils.c:585
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.
cpl_error_code muse_wcs_project_tan(muse_pixtable *aPixtable, const cpl_propertylist *aWCS)
Carry out a gnomonic projection of a pixel table into native spherical coordinates.
Definition: muse_wcs.c:992
muse_resampling_crstats_type crtype
cpl_error_code muse_flux_get_telluric_table(muse_flux_object *aFluxObj)
Get the table of the telluric correction.
Definition: muse_flux.c:2051
cpl_array * muse_cplarray_new_from_delimited_string(const char *aString, const char *aDelim)
Convert a delimited string into an array of strings.
cpl_error_code muse_flux_get_response_table(muse_flux_object *aFluxObj)
Get the table of the standard star response function.
Definition: muse_flux.c:2009
Structure definition of MUSE pixel table.
Flux object to store data needed while computing the flux calibration.
Definition: muse_flux.h:57
cpl_error_code muse_datacube_convert_dq(muse_datacube *aCube)
Convert the DQ extension of a datacube to NANs in DATA and STAT.
muse_cube_type muse_postproc_get_cube_format(const char *aFormatString)
Select correct cube format for format string.
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.
cpl_table * lines
Table of Atmospheric emission lines and their intensities.
Definition: muse_sky.h:63
static cpl_table * muse_postproc_load_nearest(const cpl_propertylist *aHeader, const cpl_frame *aFrame, float aWarnLimit, float aErrLimit)
Load the calibration from a multi-table FITS file that is nearest on the sky.
muse_resampling_crstats_type
Cosmic ray rejection statistics type.
muse_resampling_params * muse_resampling_params_new(muse_resampling_type aMethod)
Create the resampling parameters structure.
cpl_error_code muse_processing_save_cube(muse_processing *aProcessing, int aIFU, void *aCube, const char *aTag, muse_cube_type aType)
Save a MUSE datacube to disk.
cpl_error_code muse_wcs_solve(muse_wcs_object *aWCSObj, cpl_table *aReference, float aRadius, float aFAccuracy, int aIter, float aThresh)
Find the world coordinate solution for a given field using coordinates of detected sources and coordi...
Definition: muse_wcs.c:582
muse_wcs_object * muse_wcs_object_new(void)
Allocate memory for a new muse_wcs_object object.
Definition: muse_wcs.c:71
muse_cube_type
cpl_propertylist * wcs
Definition: muse_postproc.h:97
Structure definition of the post-processing properties.
Definition: muse_postproc.h:87
void muse_cplerrorstate_dump_some(unsigned aCurrent, unsigned aFirst, unsigned aLast)
Dump some CPL errors.
void muse_euro3dcube_delete(muse_euro3dcube *aEuro3D)
Deallocate memory associated to a muse_euro3dcube object.
void muse_postproc_properties_delete(muse_postproc_properties *aProp)
Free memory taken by a post-processing properties object and all its components.
Definition: muse_postproc.c:92
cpl_error_code muse_flux_calibrate(muse_pixtable *aPixtable, const cpl_table *aResponse, const cpl_table *aExtinction, const cpl_table *aTelluric)
Convert the input pixel table from counts to fluxes.
Definition: muse_flux.c:2161
cpl_imagelist * data
the cube containing the actual data values
Definition: muse_datacube.h:76
cpl_error_code muse_sky_lines_set_range(cpl_table *, double, double)
Limit the lines in the table to a wavelength range.
Structure definition of a Euro3D datacube.
Definition: muse_datacube.h:97
double muse_pfits_get_dec(const cpl_propertylist *aHeaders)
find out the declination
Definition: muse_pfits.c:224
void muse_processing_append_used(muse_processing *aProcessing, cpl_frame *aFrame, cpl_frame_group aGroup, int aDuplicate)
Add a frame to the set of used frames.
void muse_sky_master_delete(muse_sky_master *)
Delete a MASTER SKY structure.
muse_datacube * muse_resampling_cube(muse_pixtable *aPixtable, muse_resampling_params *aParams, muse_pixgrid **aPixgrid)
Resample a pixel table onto a regular grid structure representing a FITS NAXIS=3 datacube.
muse_sky_master * muse_sky_master_fit(const cpl_array *, const cpl_array *, const cpl_array *, const cpl_table *)
Fit all entries of the pixel table to the master sky.
muse_image * muse_resampling_collapse_pixgrid(muse_pixtable *aPixtable, muse_pixgrid *aPixgrid, muse_datacube *aCube, cpl_table *aFilter, muse_resampling_params *aParams)
Integrate a pixel table / pixel grid along the wavelength direction.
void muse_pixgrid_delete(muse_pixgrid *aPixels)
Delete a pixgrid and remove its memory.
Definition: muse_pixgrid.c:398
muse_postproc_type
Type of per-exposure processing to run.
Definition: muse_postproc.h:54
cpl_error_code muse_flux_integrate_std(muse_pixtable *aPixtable, muse_flux_profile_type aProfile, muse_flux_object *aFluxObj)
Integrate the flux of the standard star(s) in the field over all wavelengths.
Definition: muse_flux.c:1006
cpl_propertylist * header
the FITS header
Definition: muse_datacube.h:57
int muse_processing_save_image(muse_processing *aProcessing, int aIFU, muse_image *aImage, const char *aTag)
Save a computed MUSE image to disk.
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_flux_response_compute(muse_flux_object *aFluxObj, double aAirmass, const cpl_table *aReference, const cpl_table *aTellBands, const cpl_table *aExtinct)
Compare measured flux distribution over wavelength with calibrated stellar fluxes and derive instrume...
Definition: muse_flux.c:1833
muse_resampling_type muse_postproc_get_resampling_type(const char *aResampleString)
Select correct resampling type for resample string.
muse_image * muse_euro3dcube_collapse(muse_euro3dcube *aEuro3D, cpl_table *aFilter)
Integrate a Euro3D datacube along the wavelength direction.
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.
Definition: muse_mask.h:43
muse_imagelist * muse_imagelist_new(void)
Create a new (empty) MUSE image list.
muse_sky_master * sky
Definition: muse_postproc.h:99
muse_postproc_type type
Definition: muse_postproc.h:88
muse_resampling_type
Resampling types.
cpl_error_code muse_sky_subtract_rowbyrow_mask(muse_image *, cpl_table *)
Prepare an (object) mask for the sky row-by-row fitting.
muse_euro3dcube * muse_resampling_euro3d(muse_pixtable *aPixtable, muse_resampling_params *aParams)
Resample a pixel table onto a regular grid structure representing a Euro3D format file...
Resampling parameters.
void muse_resampling_params_delete(muse_resampling_params *aParams)
Delete a resampling parameters structure.
void muse_mask_delete(muse_mask *aMask)
Deallocate memory associated to a muse_mask object.
Definition: muse_mask.c:70
void * muse_postproc_process_exposure(muse_postproc_properties *aProp, unsigned int aIndex)
Merge and process pixel tables from one exposure.
cpl_table * continuum
Continuum flux table
Definition: muse_sky.h:65
cpl_table * muse_resampling_spectrum(muse_pixtable *aPixtable, double aBinwidth)
Resample the selected pixels of a pixel table into a spectrum.
cpl_frameset * muse_frameset_find(const cpl_frameset *aFrames, const char *aTag, unsigned char aIFU, cpl_boolean aInvert)
return frameset containing data from an IFU/channel with a certain tag
Definition: muse_utils.c:154
muse_sky_params skymodel_params
cpl_frameset * inputFrames
void muse_pixtable_delete(muse_pixtable *aPixtable)
Deallocate memory associated to a pixel table object.
muse_postproc_darcheck darcheck
Definition: muse_postproc.h:93
muse_xcombine_types
Xposure combination types.
Definition: muse_xcombine.h:45
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_sky_subtract_pixtable(muse_pixtable *a, muse_sky_master *, muse_lsf_params **)
Subtract the sky spectrum from the "data" column of a pixel table.
double muse_astro_angular_distance(double aRA1, double aDEC1, double aRA2, double aDEC2)
Compute angular distance in the sky between two positions.
Definition: muse_astro.c:487
muse_lsf_params ** lsf
LSF parameter for the resampled spectrum.
Definition: muse_sky.h:67
muse_ins_mode muse_pfits_get_mode(const cpl_propertylist *aHeaders)
find out the observation mode
Definition: muse_pfits.c:1003
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.
Definition: muse_utils.c:2140
cpl_error_code muse_wcs_locate_sources(muse_pixtable *aPixtable, float aSigma, muse_wcs_centroid_type aCentroid, muse_wcs_object *aWCSObj)
Determine the centers of stars on an astrometric exposure.
Definition: muse_wcs.c:198
double muse_astro_airmass(cpl_propertylist *aHeader)
Derive the effective airmass of an observation from information in a FITS header. ...
Definition: muse_astro.c:322
cpl_propertylist * header
The FITS header.