MUSE Pipeline Reference Manual  1.0.2
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_rvcorrect.h"
40 #include "muse_utils.h"
41 #include "muse_wcs.h"
42 
43 /*----------------------------------------------------------------------------*/
50 /*----------------------------------------------------------------------------*/
51 
54 /*----------------------------------------------------------------------------*/
65 /*----------------------------------------------------------------------------*/
68 {
69  muse_postproc_properties *prop = cpl_calloc(1, sizeof(muse_postproc_properties));
70  switch (aType) {
71  case MUSE_POSTPROC_SCIPOST:
72  case MUSE_POSTPROC_STANDARD:
73  case MUSE_POSTPROC_ASTROMETRY:
74  prop->type = aType;
75  break;
76  default:
77  cpl_msg_error(__func__, "No such setup known: %d", aType);
78  cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT);
79  cpl_free(prop);
80  return NULL;
81  }
82  /* rvtype of 0 is already MUSE_RVCORRECT_NONE */
83  return prop;
84 } /* muse_postproc_properties_new() */
85 
86 /*----------------------------------------------------------------------------*/
92 /*----------------------------------------------------------------------------*/
93 void
95 {
96  if (!aProp) {
97  return;
98  }
99  cpl_table_delete(aProp->exposures);
100  cpl_table_delete(aProp->response);
101  cpl_table_delete(aProp->extinction);
102  cpl_table_delete(aProp->telluric);
103  cpl_propertylist_delete(aProp->wcs);
104  muse_sky_master_delete(aProp->sky);
105  muse_mask_delete(aProp->sky_mask);
106  cpl_frame_delete(aProp->refframe);
107  cpl_table_delete(aProp->tellregions);
108  cpl_free(aProp);
109 } /* muse_postproc_properties_delete() */
110 
111 /*----------------------------------------------------------------------------*/
118 /*----------------------------------------------------------------------------*/
119 cpl_boolean
120 muse_postproc_check_save_param(const char *aSave, const char *aValid)
121 {
122  cpl_ensure(aSave, CPL_ERROR_NULL_INPUT, CPL_FALSE);
123  if (strlen(aSave) < 4) { /* less than "cube"! */
124  cpl_error_set_message(__func__, CPL_ERROR_DATA_NOT_FOUND,
125  "no (valid) save option given!");
126  return CPL_FALSE;
127  }
128  /* now compare all given output products against the allowed options */
129  cpl_boolean allvalid = CPL_TRUE;
130  cpl_array *asave = muse_cplarray_new_from_delimited_string(aSave, ","),
131  *avalid = muse_cplarray_new_from_delimited_string(aValid, ",");
132  int i, j,
133  ns = cpl_array_get_size(asave),
134  nv = cpl_array_get_size(avalid);
135  for (i = 0; i < ns; i++) {
136  cpl_boolean valid = CPL_FALSE;
137  for (j = 0; j < nv; j++) { /* go through all valid options */
138  if (!strcmp(cpl_array_get_string(asave, i),
139  cpl_array_get_string(avalid, j))) {
140  valid = CPL_TRUE;
141  }
142  } /* for j (all valid strings) */
143  if (!valid) { /* this string was apparently not valid! */
144  cpl_error_set_message(__func__, CPL_ERROR_ILLEGAL_INPUT,
145  "save option %d (%s) is not valid!", i + 1,
146  cpl_array_get_string(asave, i));
147  allvalid = CPL_FALSE;
148  }
149  } /* for i (all given strings) */
150  cpl_array_delete(asave);
151  cpl_array_delete(avalid);
152  return allvalid;
153 } /* muse_scipost_check_save_param() */
154 
155 /*----------------------------------------------------------------------------*/
170 /*----------------------------------------------------------------------------*/
172 muse_postproc_get_resampling_type(const char *aResampleString)
173 {
174  cpl_ensure(aResampleString, CPL_ERROR_NULL_INPUT, MUSE_RESAMPLE_WEIGHTED_DRIZZLE);
175  if (!strncmp(aResampleString, "nearest", 8)) {
176  return MUSE_RESAMPLE_NEAREST;
177  }
178  if (!strncmp(aResampleString, "linear", 7)) {
180  }
181  if (!strncmp(aResampleString, "quadratic", 10)) {
183  }
184  if (!strncmp(aResampleString, "renka", 6)) {
186  }
187  if (!strncmp(aResampleString, "drizzle", 8)) {
189  }
190  if (!strncmp(aResampleString, "lanczos", 8)) {
192  }
193  return MUSE_RESAMPLE_WEIGHTED_DRIZZLE; /* catch all fallback */
194 } /* muse_postproc_get_resampling_type() */
195 
196 /*----------------------------------------------------------------------------*/
211 /*----------------------------------------------------------------------------*/
213 muse_postproc_get_cr_type(const char *aCRTypeString)
214 {
215  cpl_ensure(aCRTypeString, CPL_ERROR_NULL_INPUT, MUSE_RESAMPLING_CRSTATS_IRAF);
216  if (!strncmp(aCRTypeString, "iraf", 5)) {
218  }
219  if (!strncmp(aCRTypeString, "mean", 5)) {
221  }
222  if (!strncmp(aCRTypeString, "median", 7)) {
224  }
225  return MUSE_RESAMPLING_CRSTATS_MEDIAN; /* catch all fallback */
226 } /* muse_postproc_get_cr_type() */
227 
228 /*----------------------------------------------------------------------------*/
242 /*----------------------------------------------------------------------------*/
244 muse_postproc_get_cube_format(const char *aFormatString)
245 {
246  cpl_ensure(aFormatString, CPL_ERROR_NULL_INPUT, MUSE_CUBE_TYPE_FITS);
247  if (!strncmp(aFormatString, "Cube", 5)) {
248  return MUSE_CUBE_TYPE_FITS;
249  }
250  if (!strncmp(aFormatString, "Euro3D", 7)) {
251  return MUSE_CUBE_TYPE_EURO3D;
252  }
253  if (!strncmp(aFormatString, "xCube", 6)) {
254  return MUSE_CUBE_TYPE_FITS_X;
255  }
256  if (!strncmp(aFormatString, "xEuro3D", 8)) {
258  }
259  return MUSE_CUBE_TYPE_FITS; /* catch all fallback */
260 } /* muse_postproc_get_cube_format() */
261 
262 /*----------------------------------------------------------------------------*/
276 /*----------------------------------------------------------------------------*/
278 muse_postproc_get_weight_type(const char *aWeightString)
279 {
280  cpl_ensure(aWeightString, CPL_ERROR_NULL_INPUT, MUSE_XCOMBINE_EXPTIME);
281  if (!strncmp(aWeightString, "exptime", 8)) {
282  return MUSE_XCOMBINE_EXPTIME;
283  }
284  if (!strncmp(aWeightString, "fwhm", 5)) {
285  return MUSE_XCOMBINE_FWHM;
286  }
287  if (!strncmp(aWeightString, "none", 5)) {
288  return MUSE_XCOMBINE_NONE;
289  }
290  return MUSE_XCOMBINE_EXPTIME; /* catch all fallback */
291 } /* muse_postproc_get_weight_type() */
292 
293 /*----------------------------------------------------------------------------*/
318 /*----------------------------------------------------------------------------*/
319 static cpl_table *
320 muse_postproc_load_nearest(const cpl_propertylist *aHeader,
321  const cpl_frame *aFrame, float aWarnLimit,
322  float aErrLimit, double *aRA, double *aDEC)
323 {
324  cpl_ensure(aHeader && aFrame, CPL_ERROR_NULL_INPUT, NULL);
325  cpl_errorstate state = cpl_errorstate_get();
326  double raref = muse_pfits_get_ra(aHeader),
327  decref = muse_pfits_get_dec(aHeader);
328  cpl_ensure(cpl_errorstate_is_equal(state), CPL_ERROR_DATA_NOT_FOUND, NULL);
329  cpl_msg_debug(__func__, "reference coordinates: RA = %e deg, DEC =%e deg",
330  raref, decref);
331  if (aRA) {
332  *aRA = raref;
333  }
334  if (aDEC) {
335  *aDEC = decref;
336  }
337 
338  const char *fn = cpl_frame_get_filename(aFrame);
339  cpl_propertylist *header;
340  double dmin = FLT_MAX; /* minimum distance yet found in deg */
341  int iext, inearest = -1, next = cpl_fits_count_extensions(fn);
342  for (iext = 1; iext <= next; iext++) {
343  header = cpl_propertylist_load(fn, iext);
344  const char *extname = cpl_propertylist_get_string(header, "EXTNAME");
345  double ra = muse_pfits_get_ra(header),
346  dec = muse_pfits_get_dec(header),
347  d = muse_astro_angular_distance(ra, dec, raref, decref);
348  cpl_msg_debug(__func__, "extension %d [%s]: RA = %e deg, DEC = %e deg --> "
349  "d = %e deg", iext, extname, ra, dec, d);
350  if (d < dmin) {
351  dmin = d;
352  inearest = iext;
353  }
354  cpl_propertylist_delete(header);
355  } /* for iext */
356  /* warn for distances larger than x arcsec */
357  if (dmin * 3600. > aErrLimit) {
358  cpl_error_set_message(__func__, CPL_ERROR_INCOMPATIBLE_INPUT, "Distance of "
359  "nearest reference table to observed position is %.2f "
360  "arcmin, certainly a wrong reference object!", dmin * 60.);
361  return NULL;
362  }
363  if (dmin * 3600. > aWarnLimit) {
364  cpl_msg_warning(__func__, "Distance of nearest reference table to observed "
365  "position is %.2f arcmin! (Wrong reference object?)",
366  dmin * 60.);
367  }
368  header = cpl_propertylist_load(fn, inearest);
369  const char *extname = cpl_propertylist_get_string(header, "EXTNAME");
370  cpl_msg_info(__func__, "Loading \"%s[%s]\" (distance %.1f arcsec)", fn,
371  extname, dmin * 3600.);
372  cpl_propertylist_delete(header);
373  cpl_table *table = cpl_table_load(fn, inearest, 1);
374  return table;
375 } /* muse_postproc_load_nearest() */
376 
377 /*----------------------------------------------------------------------------*/
420 /*----------------------------------------------------------------------------*/
421 void *
423  unsigned int aIndex,
424  muse_postproc_sky_outputs *aSkyOut)
425 {
426  cpl_ensure(aProp && aProp->exposures, CPL_ERROR_NULL_INPUT, NULL);
427  cpl_ensure((int)aIndex < cpl_table_get_nrow(aProp->exposures),
428  CPL_ERROR_ILLEGAL_INPUT, NULL);
429  cpl_msg_info(__func__, "Starting to process exposure %d (DATE-OBS=%s)",
430  aIndex + 1, cpl_table_get_string(aProp->exposures, "DATE-OBS",
431  aIndex));
432  cpl_ensure(aProp->exposures, CPL_ERROR_NULL_INPUT, NULL);
433 
434  cpl_table *exposure = cpl_table_extract(aProp->exposures, aIndex, 1);
436  aProp->lambdamin,
437  aProp->lambdamax);
438  cpl_table_delete(exposure);
439  if (!pt) {
440  return NULL;
441  }
442  /* erase pre-existing QC parameters */
443  cpl_propertylist_erase_regexp(pt->header, "ESO QC ", 0);
444 
445  /* XXX try to recognize a lab exposure by missing RA/DEC or weird *
446  * airmass; in the future this should be a template-based check */
447  cpl_boolean labdata = CPL_FALSE;
448  cpl_errorstate prestate = cpl_errorstate_get();
449  double airmass = muse_astro_airmass(pt->header);
450  double ra = muse_pfits_get_ra(pt->header),
451  dec = muse_pfits_get_dec(pt->header);
452  if (!cpl_errorstate_is_equal(prestate) || airmass < 1.) {
453  cpl_msg_warning(__func__, "This seems to be lab data (RA=%f DEC=%f, "
454  "airmass=%f), not doing any DAR correction!", ra, dec,
455  airmass);
456  cpl_errorstate_dump(prestate, CPL_FALSE, NULL);
457  cpl_errorstate_set(prestate);
458  labdata = CPL_TRUE;
459  }
460  cpl_error_code rc = CPL_ERROR_NONE;
461  if (!labdata && muse_pfits_get_mode(pt->header) <= MUSE_MODE_WFM_AO_N) {
462  cpl_msg_debug(__func__, "WFM detected: starting DAR correction");
463  rc = muse_dar_correct(pt, aProp->lambdaref);
464  cpl_msg_debug(__func__, "DAR correction returned rc=%d: %s", rc,
465  rc != CPL_ERROR_NONE ? cpl_error_get_message() : "");
466  /* check and possibly correct the DAR quality, if requested to do so */
467  if (aProp->darcheck != MUSE_POSTPROC_DARCHECK_NONE) {
468  cpl_boolean dorefine = aProp->darcheck == MUSE_POSTPROC_DARCHECK_CORRECT;
469  cpl_msg_info(__func__, "Carrying out DAR %s", dorefine ? "correction" : "check");
470  double maxshift = 0;
471  rc = muse_dar_check(pt, &maxshift, dorefine, NULL);
472  if (rc != CPL_ERROR_NONE) {
473  cpl_msg_warning(__func__, "Maximum detected shift %f\'\' (rc = %d: %s)",
474  maxshift, rc, cpl_error_get_message());
475  } else {
476  cpl_msg_info(__func__, "Maximum detected shift %f\'\'", maxshift);
477  }
478  } /* if DAR checking */
479  } /* if not labdata but WFM */
480 
481  /* do flux calibration (before sky subtraction!), if possible */
482  if (aProp->type != MUSE_POSTPROC_STANDARD) {
483  if (aProp->response) {
484  rc = muse_flux_calibrate(pt, aProp->response, aProp->extinction,
485  aProp->telluric);
486  cpl_msg_debug(__func__, "Flux calibration returned rc=%d: %s", rc,
487  rc != CPL_ERROR_NONE ? cpl_error_get_message_default(rc) : "");
488  } else {
489  cpl_msg_info(__func__, "Skipping flux calibration (no %s curve)",
490  MUSE_TAG_STD_RESPONSE);
491  } /* else */
492  } /* if not POSTPROC_STANDARD */
493 
494  cpl_table *sky_cont = NULL,
495  *sky_lines = NULL;
496  if (muse_pixtable_is_fluxcal(pt) &&
497  aProp->skymethod == MUSE_POSTPROC_SKYMETHOD_MODEL &&
498  aProp->sky && aProp->sky->lines) {
499  /* Process optional sky mask */
501  if (aProp->sky_mask) {
502  cpl_table_select_all(sky_pt->table);
504  }
505  /* Create white light image and sky mask */
506  cpl_msg_info(__func__, "Intermediate resampling to produce white-light image"
507  " for sky-spectrum creation:");
508  cpl_msg_indent_more();
511  params->crsigma = 15.;
512  muse_datacube *cube = muse_resampling_cube(sky_pt, params, NULL);
513  if (!cube) {
514  cpl_msg_error(__func__, "Could not create cube for whitelight image");
515  muse_pixtable_delete(sky_pt);
517  return NULL;
518  }
519  cpl_table *fwhite = muse_table_load_filter(NULL, "white");
520  muse_image *whitelight = muse_datacube_collapse(cube, fwhite);
521  cpl_msg_indent_less();
522  muse_mask *sky_mask = muse_sky_create_skymask(whitelight,
523  aProp->skymodel_params.fraction,
524  /* XXX hardcode prefix for now... */ "ESO QC SCIPOST");
525 
526  /* Apply mask and create spectrum */
527  cpl_table_select_all(sky_pt->table);
528  muse_pixtable_and_selected_mask(sky_pt, sky_mask);
529  cpl_table *spectrum = muse_resampling_spectrum(sky_pt, aProp->skymodel_params.sampling);
530 
531  /* Fit sky lines spectrum */
532  cpl_array *lambda = muse_cpltable_extract_column(spectrum, "lambda");
533  double lambda_low = cpl_array_get_min(lambda);
534  double lambda_high = cpl_array_get_max(lambda);
535  /* duplicate input sky lines, to crop and save them, for this one exposure */
536  sky_lines = cpl_table_duplicate(aProp->sky->lines);
537  muse_sky_lines_set_range(sky_lines, lambda_low-5, lambda_high+5);
538 
539  cpl_array *data = muse_cpltable_extract_column(spectrum, "data");
540  cpl_array *stat = muse_cpltable_extract_column(spectrum, "stat");
541  /* do the fit, ignoring possible errors */
542  prestate = cpl_errorstate_get();
543  muse_sky_master *master = muse_sky_master_fit(lambda, data, stat, sky_lines);
544  if (!cpl_errorstate_is_equal(prestate)) {
545  cpl_errorstate_dump(prestate, CPL_FALSE, NULL);
546  cpl_errorstate_set(prestate);
547  }
548  cpl_table_delete(sky_lines);
549  sky_lines = cpl_table_duplicate(master->lines); /* keep the modified lines */
550 
551  /* Create continuum */
552  if (!aProp->sky->continuum) {
553  cpl_msg_info(__func__, "No sky continuum given, create it from the data");
554  /* remove the continuum from the master sky object, *
555  * otherwise it messes up the continuum creation here! */
556  cpl_table_delete(master->continuum);
557  master->continuum = NULL;
558  muse_sky_subtract_pixtable(sky_pt, master, aProp->sky->lsf);
559  sky_cont = muse_resampling_spectrum(sky_pt, aProp->skymodel_params.csampling);
560  cpl_table_erase_column(sky_cont, "stat");
561  cpl_table_erase_column(sky_cont, "dq");
562  cpl_table_name_column(sky_cont, "data", "flux");
563  } else {
564  cpl_msg_info(__func__, "Using sky continuum given as input");
565  sky_cont = cpl_table_duplicate(aProp->sky->continuum);
566  }
567 
568  /* clean up */
569  cpl_array_unwrap(lambda);
570  cpl_array_unwrap(data);
571  cpl_array_unwrap(stat);
572  muse_sky_master_delete(master);
573  if (aSkyOut) {
574  aSkyOut->mask = sky_mask;
575  aSkyOut->spectrum = spectrum;
576  } else {
577  muse_mask_delete(sky_mask);
578  cpl_table_delete(spectrum);
579  }
580  muse_image_delete(whitelight);
581  cpl_table_delete(fwhite);
583  muse_datacube_delete(cube);
584  muse_pixtable_delete(sky_pt);
585  } /* if MUSE_POSTPROC_SKYMETHOD_MODEL */
586  /* In the following we do the sky subtraction, if we have the products *
587  * for it. The SKYMETHOD_NONE here reflects the subtract-model option *
588  * of the muse_scipost recipe. */
589  if (muse_pixtable_is_fluxcal(pt) && aProp->sky &&
590  ((aProp->skymethod == MUSE_POSTPROC_SKYMETHOD_NONE) ||
591  (aProp->skymethod == MUSE_POSTPROC_SKYMETHOD_MODEL))) {
592  /* create a new (partially duplicate) master sky that we are going to use */
594  /* leave sky->lsf NULL, since that is not used for subtraction */
595  if (aProp->sky->continuum) {
596  cpl_msg_debug(__func__, "doing sky continuum subtraction with input "
597  "spectrum");
598  sky->continuum = cpl_table_duplicate(aProp->sky->continuum);
599  } else if (sky_cont) {
600  cpl_msg_debug(__func__, "doing sky continuum subtraction with created "
601  "spectrum");
602  sky->continuum = cpl_table_duplicate(sky_cont);
603  }
604  if (sky_lines) {
605  cpl_msg_debug(__func__, "doing sky line subtraction with created list");
606  sky->lines = cpl_table_duplicate(sky_lines);
607  } else if (aProp->sky->lines) {
608  cpl_msg_debug(__func__, "doing sky line subtraction with input list");
609  sky->lines = cpl_table_duplicate(aProp->sky->lines);
610  }
611  cpl_msg_debug(__func__, "doing sky subtraction by spectrum modeling");
612  rc = muse_sky_subtract_pixtable(pt, sky, aProp->sky->lsf);
614  } /* if sky modeling */
615  /* either transfer the sky lines and continuum out or clean *
616  * them up; both should still work, if one or both are NULL */
617  if (aSkyOut) {
618  aSkyOut->lines = sky_lines;
619  aSkyOut->continuum = sky_cont;
620  } else {
621  cpl_table_delete(sky_lines);
622  cpl_table_delete(sky_cont);
623  }
624  if (aProp->skymethod == MUSE_POSTPROC_SKYMETHOD_ROWBYROW) {
625  /* apply the row-by-row sky subtraction; needed for row-by-row */
626  cpl_msg_debug(__func__, "doing sky subtraction using row-by-row fitting");
627  /* convert pixel table into images, one per IFU */
629  unsigned int k, n = muse_imagelist_get_size(list);
630  #pragma omp parallel for default(none) \
631  shared(aProp, list, n)
632  for (k = 0; k < n; k++) {
633  muse_image *image = muse_imagelist_get(list, k);
634  muse_sky_subtract_rowbyrow_mask(image, NULL);
635  muse_sky_subtract_rowbyrow(image, NULL, 3., 4);
636  } /* for all IFUs / images */
638  muse_imagelist_delete(list);
639  } /* if row-by-row sky subtraction */
640 
641  /* carry out radial velocity correction, if specified */
642  muse_rvcorrect(pt, aProp->rvtype);
643 
644  /* do flux response computation and return */
645  if (aProp->type == MUSE_POSTPROC_STANDARD) {
646  double raref, decref;
647  cpl_table *reftable = muse_postproc_load_nearest(pt->header,
648  aProp->refframe, 20., 60.,
649  &raref, &decref);
650  if (!reftable) {
651  cpl_msg_error(__func__, "Required input %s could not be loaded!",
652  MUSE_TAG_STD_FLUX_TABLE);
654  return NULL;
655  }
656  if (muse_flux_reference_table_check(reftable) != CPL_ERROR_NONE) {
657  cpl_msg_error(__func__, "%s does not have a recognized format!",
658  MUSE_TAG_STD_FLUX_TABLE);
659  cpl_table_delete(reftable);
661  return NULL;
662  }
663 
664  prestate = cpl_errorstate_get();
666  rc = muse_flux_integrate_std(pt, aProp->profile, flux);
667  if (flux->intimage) {
668  cpl_msg_debug(__func__, "Flux integration found %"CPL_SIZE_FORMAT" stars",
669  cpl_image_get_size_y(flux->intimage->data));
670  }
672  if (airmass < 0) {
673  cpl_msg_warning(__func__, "Invalid airmass derived (%.3e), resetting to "
674  " zero", airmass);
675  airmass = 0.; /* set to zero to leave out uncertain extinction term */
676  }
677  flux->raref = raref;
678  flux->decref = decref;
679  rc = muse_flux_response_compute(flux, aProp->select, airmass, reftable,
680  aProp->tellregions, aProp->extinction);
681  cpl_table_delete(reftable);
682  rc = muse_flux_get_response_table(flux, aProp->smooth);
683  rc = muse_flux_get_telluric_table(flux);
684 
685  if (!cpl_errorstate_is_equal(prestate)) {
686  cpl_msg_warning(__func__, "The following errors happened while computing "
687  "the flux calibration (rc = %d/%s):", rc,
688  cpl_error_get_message_default(rc));
689  cpl_errorstate_dump(prestate, CPL_FALSE, muse_cplerrorstate_dump_some);
690  }
691  return flux;
692  } /* if MUSE_POSTPROC_STANDARD */
693 
694  /* do astrometry computation and return */
695  if (aProp->type == MUSE_POSTPROC_ASTROMETRY) {
696  cpl_table *reftable = muse_postproc_load_nearest(pt->header,
697  aProp->refframe, 20., 60.,
698  NULL, NULL);
699  if (!reftable) {
700  cpl_msg_error(__func__, "Required input %s could not be loaded!",
701  MUSE_TAG_ASTROMETRY_REFERENCE);
703  return NULL;
704  }
705  prestate = cpl_errorstate_get();
706 
708  wcs->crpix1 = aProp->crpix1; /* set center of rotation */
709  wcs->crpix2 = aProp->crpix2;
710  rc = muse_wcs_locate_sources(pt, aProp->detsigma, aProp->centroid, wcs);
711  rc = muse_wcs_solve(wcs, reftable, aProp->radius, aProp->faccuracy,
712  aProp->niter, aProp->rejsigma);
713  cpl_table_delete(reftable);
714  cpl_msg_debug(__func__, "Solving astrometry (initial radius %.1f pix, "
715  "initial relative accuracy %.1f, %d iterations with rejection"
716  " sigma %.2f) returned rc=%d: %s", aProp->radius,
717  aProp->faccuracy, aProp->niter, aProp->rejsigma,
718  rc, rc != CPL_ERROR_NONE ? cpl_error_get_message() : "");
720 
721  if (!cpl_errorstate_is_equal(prestate)) {
722  cpl_msg_warning(__func__, "The following errors happened while computing "
723  "the astrometric solution (rc = %d/%s):", rc,
724  cpl_error_get_message());
725  cpl_errorstate_dump(prestate, CPL_FALSE, muse_cplerrorstate_dump_some);
726  }
727  return wcs;
728  } /* if MUSE_POSTPROC_ASTROMETRY */
729 
730  rc = muse_wcs_project_tan(pt, aProp->wcs);
731  cpl_msg_debug(__func__, "Astrometry correction returned rc=%d: %s", rc,
732  rc != CPL_ERROR_NONE ? cpl_error_get_message() : "");
733 
734  return pt;
735 } /* muse_postproc_process_exposure() */
736 
737 /*---------------------------------------------------------------------------*/
755 /*---------------------------------------------------------------------------*/
756 cpl_propertylist *
758 {
759  cpl_ensure(aProcessing, CPL_ERROR_NULL_INPUT, NULL);
760 
761  cpl_frameset *fwcs = muse_frameset_find(aProcessing->inframes,
762  MUSE_TAG_OUTPUT_WCS, 0, CPL_FALSE);
763  if (!fwcs || cpl_frameset_get_size(fwcs) <= 0) {
764  /* not an error, quietly return NULL */
765  cpl_frameset_delete(fwcs);
766  return NULL;
767  }
768  cpl_frame *frame = cpl_frameset_get_position(fwcs, 0);
769  const char *fn = cpl_frame_get_filename(frame);
770  cpl_propertylist *header = NULL;
771  int iext, next = cpl_fits_count_extensions(fn);
772  for (iext = 0; iext <= next; iext++) {
773  header = cpl_propertylist_load(fn, iext);
774  cpl_errorstate es = cpl_errorstate_get();
775  cpl_wcs *wcs = cpl_wcs_new_from_propertylist(header);
776  if (!cpl_errorstate_is_equal(es)) {
777  cpl_errorstate_set(es);
778  }
779  if (!wcs) {
780  cpl_propertylist_delete(header);
781  header = NULL;
782  continue; /* try the next extension */
783  }
784  /* need 2D or 3D WCS to have something useful */
785  int naxis = cpl_wcs_get_image_naxis(wcs);
786  cpl_boolean valid = naxis == 2 || naxis == 3;
787  /* need RA---TAN and DEC--TAN for axes 1 and 2 */
788  const cpl_array *ctypes = cpl_wcs_get_ctype(wcs);
789  if (valid && cpl_array_get_string(ctypes, 0)) {
790  valid = valid && !strncmp(cpl_array_get_string(ctypes, 0), "RA---TAN", 9);
791  }
792  if (valid && cpl_array_get_string(ctypes, 1)) {
793  valid = valid && !strncmp(cpl_array_get_string(ctypes, 1), "DEC--TAN", 9);
794  }
795  /* need AWAV or AWAV-LOG for axis 3 */
796  if (valid && cpl_array_get_string(ctypes, 2)) {
797  const char *ctype3 = cpl_array_get_string(ctypes, 2);
798  valid = valid
799  && (!strncmp(ctype3, "AWAV", 5) || !strncmp(ctype3, "AWAV-LOG", 9));
800  }
801  if (valid) { /* check for unsupported PCi_j keywords */
802  /* still needs to be done on the input header not on the derives WCS! */
803  cpl_propertylist *pc = cpl_propertylist_new();
804  cpl_propertylist_copy_property_regexp(pc, header, "^PC[0-9]+_[0-9]+", 0);
805  valid = cpl_propertylist_get_size(pc) == 0;
806  cpl_propertylist_delete(pc);
807  }
808  cpl_wcs_delete(wcs);
809  if (valid) {
810  cpl_msg_debug(__func__, "Apparently valid header %s found in extension %d"
811  " of \"%s\"", MUSE_TAG_OUTPUT_WCS, iext, fn);
812  muse_processing_append_used(aProcessing, frame, CPL_FRAME_GROUP_CALIB, 1);
813  break;
814  }
815  cpl_propertylist_delete(header);
816  header = NULL;
817  } /* for iext (all extensions) */
818  if (!header) {
819  cpl_msg_warning(__func__, "No valid headers for %s found in \"%s\"",
820  MUSE_TAG_OUTPUT_WCS, fn);
821  }
822  cpl_frameset_delete(fwcs);
823  return header;
824 } /* muse_postproc_cube_load_output_wcs() */
825 
826 /*---------------------------------------------------------------------------*/
844 /*---------------------------------------------------------------------------*/
845 cpl_error_code
847  muse_pixtable *aPixtable,
848  muse_cube_type aFormat,
849  muse_resampling_params *aParams,
850  const char *aFilter)
851 {
852  cpl_ensure_code(aProcessing && aPixtable && aParams, CPL_ERROR_NULL_INPUT);
853  cpl_ensure_code(aFormat == MUSE_CUBE_TYPE_EURO3D ||
854  aFormat == MUSE_CUBE_TYPE_FITS ||
855  aFormat == MUSE_CUBE_TYPE_EURO3D_X ||
856  aFormat == MUSE_CUBE_TYPE_FITS_X, CPL_ERROR_ILLEGAL_INPUT);
857 
858  /* create cube and cast to generic pointer to save code duplication */
859  void *cube = NULL;
860  muse_pixgrid *grid = NULL;
861  if (aFormat == MUSE_CUBE_TYPE_EURO3D || aFormat == MUSE_CUBE_TYPE_EURO3D_X) {
862  cpl_msg_info(__func__, "Resampling to final cube follows, output format "
863  "\"Euro3D\"");
864  cpl_msg_indent_more();
865  muse_euro3dcube *e3d = muse_resampling_euro3d(aPixtable, aParams);
866  cpl_msg_indent_less();
867  cpl_ensure_code(e3d, cpl_error_get_code());
868  cube = e3d;
869  } else if (aFormat == MUSE_CUBE_TYPE_FITS || aFormat == MUSE_CUBE_TYPE_FITS_X) {
870  cpl_msg_info(__func__, "Resampling to final cube follows, output format "
871  "\"FITS cube\"");
872  cpl_msg_indent_more();
873  muse_datacube *fits = muse_resampling_cube(aPixtable, aParams, &grid);
874  cpl_msg_indent_less();
875  cpl_ensure_code(fits, cpl_error_get_code());
876  muse_postproc_qc_fwhm(aProcessing, fits);
877  cube = fits;
878  }
879 
880  cpl_array *filters = muse_cplarray_new_from_delimited_string(aFilter, ",");
881  int i, ipos = 0, nfilters = cpl_array_get_size(filters);
882  for (i = 0; i < nfilters; i++) {
883  /* try to find and load the filter from a table */
884  cpl_table *filter = muse_table_load_filter(aProcessing,
885  cpl_array_get_string(filters, i));
886  if (!filter) {
887  continue;
888  }
889 
890  /* if we got a filter table, do the collapsing */
891  muse_image *fov = NULL;
892  if (aFormat == MUSE_CUBE_TYPE_EURO3D || aFormat == MUSE_CUBE_TYPE_EURO3D_X) {
893  fov = muse_euro3dcube_collapse(cube, filter);
894  } else if (aFormat == MUSE_CUBE_TYPE_FITS || aFormat == MUSE_CUBE_TYPE_FITS_X) {
895  if (getenv("MUSE_COLLAPSE_PIXTABLE") &&
896  atoi(getenv("MUSE_COLLAPSE_PIXTABLE")) > 0) {
897  /* create collapsed image directly from pixel table/grid */
898  fov = muse_resampling_collapse_pixgrid(aPixtable, grid,
899  (muse_datacube *)cube,
900  filter, aParams);
901  } else {
902  fov = muse_datacube_collapse(cube, filter);
903  }
904  }
905  cpl_table_delete(filter);
906 
907  if (aFormat == MUSE_CUBE_TYPE_EURO3D_X || aFormat == MUSE_CUBE_TYPE_FITS_X) {
908  if (!((muse_datacube *)cube)->recimages) {
909  /* create empty structures before first entry */
910  ((muse_datacube *)cube)->recimages = muse_imagelist_new();
911  ((muse_datacube *)cube)->recnames = cpl_array_new(0, CPL_TYPE_STRING);
912  }
913  /* append image to cube/Euro3D structure */
914  muse_imagelist_set(((muse_datacube *)cube)->recimages, fov,
915  muse_imagelist_get_size(((muse_datacube *)cube)->recimages));
916  cpl_array_set_size(((muse_datacube *)cube)->recnames, ipos+1);
917  cpl_array_set_string(((muse_datacube *)cube)->recnames, ipos,
918  cpl_array_get_string(filters, i));
919  } else {
920  /* we want NANs instead of DQ in the output images */
922  /* duplicate cube header to be able to use its WCS info as the base */
923  muse_utils_copy_modified_header(fov->header, fov->header, "OBJECT",
924  cpl_array_get_string(filters, i));
925  cpl_propertylist_update_string(fov->header, MUSE_HDR_FILTER,
926  cpl_array_get_string(filters, i));
927  muse_processing_save_image(aProcessing, -1, fov, MUSE_TAG_IMAGE_FOV);
928  muse_image_delete(fov);
929  }
930  ipos++;
931  } /* for i (all filters) */
932  cpl_array_delete(filters);
933  muse_pixgrid_delete(grid);
934 
935  cpl_error_code rc = CPL_ERROR_NONE;
936  if (aFormat == MUSE_CUBE_TYPE_EURO3D || aFormat == MUSE_CUBE_TYPE_EURO3D_X) {
937  /* save the data and clean up */
938  rc = muse_processing_save_cube(aProcessing, -1, cube, MUSE_TAG_CUBE_FINAL,
941  } else if (aFormat == MUSE_CUBE_TYPE_FITS || aFormat == MUSE_CUBE_TYPE_FITS_X) {
942  /* convert/remove the DQ cube */
944  /* save the cube and clean up */
945  rc = muse_processing_save_cube(aProcessing, -1, cube, MUSE_TAG_CUBE_FINAL,
947  muse_datacube_delete(cube);
948  }
949  return rc;
950 } /* muse_postproc_cube_resample_and_collapse() */
951 
952 /*---------------------------------------------------------------------------*/
975 /*---------------------------------------------------------------------------*/
976 cpl_error_code
978 {
979  cpl_ensure_code(aProcessing && aCube, CPL_ERROR_NULL_INPUT);
980 
981  /* determine which FITS keyword prefix to take */
982  char prefix[KEYWORD_LENGTH] = "";
983  if (!strncmp(aProcessing->name, "muse_scipost", 13)) {
984  strncpy(prefix, QC_POST_PREFIX_SCIPOST, KEYWORD_LENGTH);
985  } else if (!strncmp(aProcessing->name, "muse_exp_combine", 17)) {
986  strncpy(prefix, QC_POST_PREFIX_EXPCOMBINE, KEYWORD_LENGTH);
987  } else if (!strncmp(aProcessing->name, "muse_standard", 14)) {
988  strncpy(prefix, QC_POST_PREFIX_STANDARD, KEYWORD_LENGTH);
989  } else if (!strncmp(aProcessing->name, "muse_astrometry", 16)) {
990  strncpy(prefix, QC_POST_PREFIX_ASTROMETRY, KEYWORD_LENGTH);
991  } else {
992  cpl_msg_info(__func__, "Recipe \"%s\" found, not generating QC1 keywords",
993  aProcessing->name);
994  return CPL_ERROR_NONE;
995  }
996 
997  /* get image from central plane of cube and use it to do object detection *
998  * and FWHM estimation */
999  int plane = cpl_imagelist_get_size(aCube->data) / 2;
1000  cpl_image *cim = cpl_imagelist_get(aCube->data, plane);
1001  double dsigmas[] = {/*50., 30., 10., 8., 6.,*/ 5., 4., 3. };
1002  int ndsigmas = sizeof(dsigmas) / sizeof(double);
1003  cpl_vector *vsigmas = cpl_vector_wrap(ndsigmas, dsigmas);
1004  cpl_size isigma = -1;
1005  cpl_errorstate prestate = cpl_errorstate_get();
1006  cpl_apertures *apertures = cpl_apertures_extract(cim, vsigmas, &isigma);
1007  cpl_vector_unwrap(vsigmas);
1008  if (!apertures || !cpl_errorstate_is_equal(prestate)) {
1009  /* just set up two fake entries, with number zero and negative FWHMs: */
1010  char keyword[KEYWORD_LENGTH];
1011  snprintf(keyword, KEYWORD_LENGTH, QC_POST_POSX, prefix, 0);
1012  cpl_propertylist_update_float(aCube->header, keyword, -1.);
1013  snprintf(keyword, KEYWORD_LENGTH, QC_POST_POSY, prefix, 0);
1014  cpl_propertylist_update_float(aCube->header, keyword, -1.);
1015  snprintf(keyword, KEYWORD_LENGTH, QC_POST_FWHMX, prefix, 0);
1016  cpl_propertylist_update_float(aCube->header, keyword, -1.);
1017  snprintf(keyword, KEYWORD_LENGTH, QC_POST_FWHMY, prefix, 0);
1018  cpl_propertylist_update_float(aCube->header, keyword, -1.);
1019  /* swallow the cpl_aperture error, which is too cryptic for users */
1020  cpl_errorstate_set(prestate);
1021  cpl_msg_warning(__func__, "No sources found for FWHM measurement down to "
1022  "%.1f sigma limit on plane %d, QC parameters will not "
1023  "contain useful information", dsigmas[ndsigmas-1], plane+1);
1024  return CPL_ERROR_DATA_NOT_FOUND;
1025  }
1026  int ndet = cpl_apertures_get_size(apertures);
1027  cpl_msg_info(__func__, "Computing FWHM QC parameters for %d source%s found "
1028  "down to the %.1f sigma threshold on plane %d", ndet,
1029  ndet == 1 ? "" : "s", dsigmas[isigma], plane + 1);
1030 
1031  /* get some kind of WCS for conversion of FWHM from pixels to arcsec, *
1032  * by default assume WFM and x=RA and y=DEC */
1033  double cd11 = kMuseSpaxelSizeX_WFM, cd12 = 0., cd21 = 0.,
1034  cd22 = kMuseSpaxelSizeY_WFM;
1035  cpl_errorstate es = cpl_errorstate_get();
1036  cpl_wcs *wcs = cpl_wcs_new_from_propertylist(aCube->header);
1037  if (!cpl_errorstate_is_equal(es)) {
1038  cpl_errorstate_set(es); /* swallow possible errors from creating the WCS */
1039  }
1040  if (!wcs ||
1041  !strncasecmp(cpl_propertylist_get_string(aCube->header, "CTYPE1"),
1042  "PIXEL", 5)) {
1043  cpl_msg_debug(__func__, "QC1 FWHM parameter estimation (%d sources): "
1044  "simple conversion to arcsec (CTYPE=%s/%s)!", ndet,
1045  cpl_propertylist_get_string(aCube->header, "CTYPE1"),
1046  cpl_propertylist_get_string(aCube->header, "CTYPE2"));
1047  if (muse_pfits_get_mode(aCube->header) > MUSE_MODE_WFM_AO_N) { /* NFM scaling */
1048  cd11 = kMuseSpaxelSizeX_WFM;
1049  cd22 = kMuseSpaxelSizeY_WFM;
1050  }
1051  } else {
1052  const cpl_matrix *cd = cpl_wcs_get_cd(wcs);
1053  /* take the absolute and scale by 3600 to get positive arcseconds */
1054  cd11 = fabs(cpl_matrix_get(cd, 0, 0)) * 3600.,
1055  cd12 = fabs(cpl_matrix_get(cd, 0, 1)) * 3600.,
1056  cd21 = fabs(cpl_matrix_get(cd, 1, 0)) * 3600.,
1057  cd22 = fabs(cpl_matrix_get(cd, 1, 1)) * 3600.;
1058  cpl_msg_debug(__func__, "QC1 FWHM parameter estimation (%d sources): "
1059  "full conversion to arcsec (CD=%f,%f,%f,%f)", ndet,
1060  cd11, cd12, cd21, cd22);
1061  }
1062  cpl_wcs_delete(wcs); /* rely on internal NULL check */
1063 
1064  /* Compute FWHM of all apertures */
1065 #if 0
1066  printf("index RA_FWHM DEC_FWHM\n");
1067 #endif
1068  int n;
1069  for (n = 1; n <= ndet; n++) {
1070  cpl_size xcen = lround(cpl_apertures_get_centroid_x(apertures, n)),
1071  ycen = lround(cpl_apertures_get_centroid_y(apertures, n));
1072  double x, y;
1073  cpl_image_get_fwhm(cim, xcen, ycen, &x, &y);
1074 
1075  /* linear WCS scaling */
1076  x = cd11 * x + cd12 * y;
1077  y = cd22 * y + cd21 * x;
1078 #if 0
1079  printf("%4d %f %f\n", n, x, y);
1080 #endif
1081 
1082  /* generate and write FITS keywords */
1083  char keyword[KEYWORD_LENGTH];
1084  snprintf(keyword, KEYWORD_LENGTH, QC_POST_POSX, prefix, n);
1085  cpl_propertylist_update_float(aCube->header, keyword, xcen);
1086  snprintf(keyword, KEYWORD_LENGTH, QC_POST_POSY, prefix, n);
1087  cpl_propertylist_update_float(aCube->header, keyword, ycen);
1088  snprintf(keyword, KEYWORD_LENGTH, QC_POST_FWHMX, prefix, n);
1089  cpl_propertylist_update_float(aCube->header, keyword, x);
1090  snprintf(keyword, KEYWORD_LENGTH, QC_POST_FWHMY, prefix, n);
1091  cpl_propertylist_update_float(aCube->header, keyword, y);
1092  } /* for n (aperture number) */
1093 #if 0
1094  fflush(stdout);
1095 #endif
1096  cpl_apertures_delete(apertures);
1097 
1098  return CPL_ERROR_NONE;
1099 } /* muse_postproc_qc_fwhm() */
1100 
Structure definition of a MUSE datacube.
Definition: muse_datacube.h:47
muse_rvcorrect_type rvtype
Definition: muse_postproc.h:96
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:67
muse_flux_selection_type select
void muse_image_delete(muse_image *aImage)
Deallocate memory associated to a muse_image object.
Definition: muse_image.c:85
muse_image * intimage
Definition: muse_flux.h:70
muse_pixtable * muse_pixtable_duplicate(muse_pixtable *aPixtable)
Make a copy of the pixtanle.
muse_wcs_centroid_type centroid
const char * name
double muse_pfits_get_ra(const cpl_propertylist *aHeaders)
find out the right ascension
Definition: muse_pfits.c:232
cpl_image * data
the data extension
Definition: muse_image.h:46
cpl_error_code muse_rvcorrect(muse_pixtable *aPixtable, muse_rvcorrect_type aType)
Correct the wavelengths of all pixels of a given pixel table for radial velocity shift.
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:57
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.
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:168
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, const char *)
Select spaxels to be considered as sky.
void * muse_postproc_process_exposure(muse_postproc_properties *aProp, unsigned int aIndex, muse_postproc_sky_outputs *aSkyOut)
Merge and process pixel tables from one exposure.
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:40
muse_flux_object * muse_flux_object_new(void)
Allocate memory for a new muse_flux_object object.
Definition: muse_flux.c:68
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
cpl_propertylist * header
the FITS header
Definition: muse_image.h:72
cpl_table * muse_table_load_filter(muse_processing *aProcessing, const char *aFilterName)
Load a table for a given filter name.
Definition: muse_utils.c:853
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:993
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:2458
cpl_array * muse_cplarray_new_from_delimited_string(const char *aString, const char *aDelim)
Convert a delimited string into an array of strings.
Structure definition of MUSE pixel table.
Flux object to store data needed while computing the flux calibration.
Definition: muse_flux.h:65
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_boolean muse_postproc_check_save_param(const char *aSave, const char *aValid)
Check the –save parameter contents against allowed options.
cpl_table * lines
Table of Atmospheric emission lines and their intensities.
Definition: muse_sky.h:61
cpl_error_code muse_flux_response_compute(muse_flux_object *aFluxObj, muse_flux_selection_type aSelect, 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:1930
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:580
cpl_boolean muse_pixtable_is_fluxcal(muse_pixtable *aPixtable)
Determine whether the pixel table is flux calibrated.
muse_wcs_object * muse_wcs_object_new(void)
Allocate memory for a new muse_wcs_object object.
Definition: muse_wcs.c:71
muse_flux_smooth_type smooth
muse_cube_type
cpl_propertylist * wcs
Structure definition of the post-processing properties.
Definition: muse_postproc.h:89
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:94
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:2568
cpl_imagelist * data
the cube containing the actual data values
Definition: muse_datacube.h:75
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:96
double muse_pfits_get_dec(const cpl_propertylist *aHeaders)
find out the declination
Definition: muse_pfits.c:250
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.
cpl_error_code muse_flux_get_response_table(muse_flux_object *aFluxObj, muse_flux_smooth_type aSmooth)
Get the table of the standard star response function.
Definition: muse_flux.c:2397
void muse_sky_master_delete(muse_sky_master *)
Delete a MASTER SKY structure.
static cpl_table * muse_postproc_load_nearest(const cpl_propertylist *aHeader, const cpl_frame *aFrame, float aWarnLimit, float aErrLimit, double *aRA, double *aDEC)
Load the calibration from a multi-table FITS file that is nearest on the sky.
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:56
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:1020
cpl_propertylist * header
the FITS header
Definition: muse_datacube.h:56
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.
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:42
muse_imagelist * muse_imagelist_new(void)
Create a new (empty) MUSE image list.
muse_sky_master * sky
cpl_frameset * inframes
muse_postproc_type type
Definition: muse_postproc.h:90
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.
Structure definition of the post-processing output sky data.
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:68
cpl_table * continuum
Continuum flux table
Definition: muse_sky.h:63
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:158
muse_sky_params skymodel_params
void muse_pixtable_delete(muse_pixtable *aPixtable)
Deallocate memory associated to a pixel table object.
muse_postproc_darcheck darcheck
Definition: muse_postproc.h:95
muse_xcombine_types
Xposure combination types.
Definition: muse_xcombine.h:44
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:496
muse_lsf_params ** lsf
LSF parameter for the resampled spectrum.
Definition: muse_sky.h:65
muse_ins_mode muse_pfits_get_mode(const cpl_propertylist *aHeaders)
find out the observation mode
Definition: muse_pfits.c:1097
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:2489
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
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).
Definition: muse_image.c:897
double muse_astro_airmass(cpl_propertylist *aHeader)
Derive the effective airmass of an observation from information in a FITS header. ...
Definition: muse_astro.c:331
cpl_propertylist * header
The FITS header.
muse_sky_master * muse_sky_master_new(void)
Create a muse_sky_master structure.