35 #include "muse_instrument.h" 37 #include "muse_astro.h" 38 #include "muse_combine.h" 39 #include "muse_pfits.h" 40 #include "muse_phys.h" 41 #include "muse_quality.h" 42 #include "muse_utils.h" 48 #define DEBUG_MUSE_DARCORRECT 0 49 #define DEBUG_MUSE_DARCHECK 0 51 #define MUSE_DARCHECK_GRID_FITS 0 141 cpl_ensure_code(aPixtable && aPixtable->
header, CPL_ERROR_NULL_INPUT);
142 if (aLambdaRef > 22000. || aLambdaRef < 3500.) {
143 cpl_msg_info(__func__,
"Reference lambda %.3f Angstrom: skipping DAR " 144 "correction", aLambdaRef);
147 MUSE_HDR_PT_DAR_C_SKIP);
148 return CPL_ERROR_NONE;
152 cpl_msg_info(__func__,
"pixel table already corrected: skipping DAR " 154 return CPL_ERROR_NONE;
161 cpl_ensure_code(airmass >= 1., cpl_error_get_code());
163 double z = acos(1./airmass);
166 cpl_errorstate prestate = cpl_errorstate_get();
168 cpl_ensure_code(cpl_errorstate_is_equal(prestate), cpl_error_get_code());
171 prestate = cpl_errorstate_get();
173 cpl_ensure_code(cpl_errorstate_is_equal(prestate), cpl_error_get_code());
176 < MUSE_MODE_NFM_AO_N;
178 cpl_msg_warning(__func__,
"For NFM instrument mode DAR correction should " 183 prestate = cpl_errorstate_get();
185 if (cpl_errorstate_is_equal(prestate)) {
188 cpl_msg_warning(__func__,
"FITS header does not contain temperature: %s",
189 cpl_error_get_message());
190 temp = 11.5 + 273.15;
192 prestate = cpl_errorstate_get();
194 if (!cpl_errorstate_is_equal(prestate)) {
195 cpl_msg_warning(__func__,
"FITS header does not contain relative humidity: %s",
196 cpl_error_get_message());
200 prestate = cpl_errorstate_get();
203 if (!cpl_errorstate_is_equal(prestate)) {
204 cpl_msg_warning(__func__,
"FITS header does not contain pressure values: %s",
205 cpl_error_get_message());
213 MUSE_DAR_METHOD_FILIPPENKO = 0,
214 MUSE_DAR_METHOD_OWENS,
215 MUSE_DAR_METHOD_EDLEN,
216 MUSE_DAR_METHOD_CIDDOR,
218 int method = MUSE_DAR_METHOD_FILIPPENKO;
219 if (getenv(
"MUSE_DAR_CORRECT_METHOD") &&
220 !strncmp(getenv(
"MUSE_DAR_CORRECT_METHOD"),
"Owens", 6)) {
221 method = MUSE_DAR_METHOD_OWENS;
222 }
else if (getenv(
"MUSE_DAR_CORRECT_METHOD") &&
223 !strncmp(getenv(
"MUSE_DAR_CORRECT_METHOD"),
"Edlen", 6)) {
224 method = MUSE_DAR_METHOD_EDLEN;
225 }
else if (getenv(
"MUSE_DAR_CORRECT_METHOD") &&
226 !strncmp(getenv(
"MUSE_DAR_CORRECT_METHOD"),
"Ciddor", 7)) {
227 method = MUSE_DAR_METHOD_CIDDOR;
238 if (method == MUSE_DAR_METHOD_OWENS) {
241 cpl_msg_info(__func__,
"DAR for T=%.2f K, RH=%.2f %%, pres=%.1f mbar, " 242 "at airmass=%.4f, using ref. wavelength %.6f um (Owens)",
243 temp, rhum*100., pres, 1./cos(z), aLambdaRef / 10000.);
244 }
else if (method == MUSE_DAR_METHOD_EDLEN ||
245 method == MUSE_DAR_METHOD_CIDDOR) {
250 const double T = temp,
251 K1 = 1.16705214528E+03,
252 K2 = -7.24213167032E+05,
253 K3 = -1.70738469401E+01,
254 K4 = 1.20208247025E+04,
255 K5 = -3.23255503223E+06,
256 K6 = 1.49151086135E+01,
257 K7 = -4.82326573616E+03,
258 K8 = 4.05113405421E+05,
259 K9 = -2.38555575678E-01,
260 K10 = 6.50175348448E+02,
261 Omega = T + K9 / (T - K10),
262 A = Omega*Omega + K1 * Omega + K2,
263 B = K3 * Omega*Omega + K4 * Omega + K5,
264 C = K6 * Omega*Omega + K7 * Omega + K8,
265 X = -B + sqrt(B*B - 4 * A * C);
266 double psv_w = 1e6 * pow(2 * C / X, 4);
269 const double A1 = -13.928169,
272 Y = A1 * (1 - pow(Theta, -1.5)) + A2 * (1 - pow(Theta, -1.25));
273 double psv_i = 611.657 * exp(Y);
283 const double alpha = 1.00062,
286 double f = alpha + beta * p + gamma * t*t;
289 xv = rhum * f * psv_w / p;
291 xv = rhum * f * psv_i / p;
293 if (method == MUSE_DAR_METHOD_EDLEN) {
298 cpl_msg_info(__func__,
"DAR for T=%.2f degC, RH=%.2f %%, p=%.1f Pa, " 299 "at airmass=%.4f, using ref. wavelength %.6f um (%s)",
300 t, rhum*100., p, 1./cos(z), aLambdaRef / 10000.,
301 method == MUSE_DAR_METHOD_EDLEN ?
"Edlen" :
"Ciddor");
306 fp = rhum * ps * MUSE_PHYS_hPa_TO_mmHg;
308 pres *= MUSE_PHYS_hPa_TO_mmHg;
310 cpl_msg_info(__func__,
"DAR for T=%.2f degC, fp=%.3f mmHg, P=%.1f mmHg, " 311 "at airmass=%.4f, using ref. wavelength %.6f um (Filippenko)",
312 temp, fp, pres, 1./cos(z), aLambdaRef / 10000.);
319 double xshift = -sin((parang + posang) * CPL_MATH_RAD_DEG),
320 yshift = cos((parang + posang) * CPL_MATH_RAD_DEG);
326 double fscale = 3600.;
329 double xscale, yscale;
335 xshift /= kMuseSpaxelSizeX_WFM;
336 yshift /= kMuseSpaxelSizeY_WFM;
338 xshift /= kMuseSpaxelSizeX_NFM;
339 yshift /= kMuseSpaxelSizeY_NFM;
348 double dr0 = tan(z) * fscale * CPL_MATH_DEG_RAD;
349 #if DEBUG_MUSE_DARCORRECT 351 for (wl = 4650.; wl <= 9300; wl += 50) {
353 if (method == MUSE_DAR_METHOD_OWENS) {
355 }
else if (method == MUSE_DAR_METHOD_EDLEN) {
357 }
else if (method == MUSE_DAR_METHOD_CIDDOR) {
362 double dr = dr0 * (nr0 - nr);
363 printf(
"%.1f Angstrom: %f ==> %f %f %s (%s)\n", wl, dr,
364 xshift * dr, yshift * dr,
373 #pragma omp parallel for default(none) \ 374 shared(d1, d2, dr0, fp, lbda, method, nmax, nr0, p, pres, pv, t, \ 375 temp, xpos, ypos, xshift, xv, yshift) 376 for (i = 0; i < nmax; i++) {
377 double nr, lambda = lbda[i] * 1e-4;
378 if (method == MUSE_DAR_METHOD_OWENS) {
380 }
else if (method == MUSE_DAR_METHOD_EDLEN) {
382 }
else if (method == MUSE_DAR_METHOD_CIDDOR) {
387 double dr = dr0 * (nr0 - nr);
389 xpos[i] += xshift * dr;
390 ypos[i] += yshift * dr;
396 MUSE_HDR_PT_DAR_COMMENT);
400 cpl_propertylist_get_float(aPixtable->
header,
403 cpl_propertylist_get_float(aPixtable->
header,
406 cpl_propertylist_get_float(aPixtable->
header,
409 cpl_propertylist_get_float(aPixtable->
header,
413 return CPL_ERROR_NONE;
475 cpl_ensure_code(aMaxShift, CPL_ERROR_NULL_INPUT);
477 cpl_ensure_code(aPixtable, CPL_ERROR_NULL_INPUT);
481 cpl_msg_info(__func__,
"pixel table already checked for DAR accuracy");
483 cpl_boolean corrected = CPL_FALSE;
485 cpl_msg_info(__func__,
"pixel table already corrected for DAR residuals");
486 corrected = CPL_TRUE;
494 cpl_boolean residuals = cpl_propertylist_has(aPixtable->
header,
496 int order = residuals || corrected ? 2 : 4;
499 cpl_msg_info(__func__,
"Intermediate resampling to %s atmospheric refraction" 500 "%s (using polynomial order %d):", aCorrect ?
"correct" :
"check",
501 residuals ?
" residuals" :
"", order);
502 cpl_msg_indent_more();
508 params->dlambda = 10.;
513 cpl_msg_indent_less();
517 return cpl_error_set_message(__func__, cpl_error_get_code(),
518 "Failure while creating cube!");
526 cpl_errorstate state = cpl_errorstate_get();
527 double lambdaref = cpl_propertylist_get_double(aPixtable->
header,
529 lambda1 = (2 - crpix3) * cdelt3 + crval3,
530 lambda3 = (grid->nz - 1 - crpix3) * cdelt3 + crval3,
531 lambda2 = (lambda1 + lambda3) / 2.;
532 if (!cpl_errorstate_is_equal(state) || lambdaref < 0) {
533 cpl_errorstate_set(state);
540 cpl_msg_warning(__func__,
"Data is not DAR corrected, using reference " 541 "wavelength at %.3f Angstrom and polynomial order %d!",
544 if (lambdaref > lambda3 || lambdaref < lambda1) {
546 cpl_msg_warning(__func__,
"Reference lambda outside cube, using reference " 547 "wavelength at %.3f Angstrom!", lambdaref);
551 int iplane, irefplane = lround((lambdaref - crval3) / cdelt3 + crpix3) - 1;
556 unsigned int ilist = 0;
557 for (iplane = irefplane - 1; iplane <= irefplane + 1; iplane++) {
559 image->
data = cpl_image_duplicate(cpl_imagelist_get(cube->
data, iplane));
560 image->
dq = cpl_image_duplicate(cpl_imagelist_get(cube->
dq, iplane));
561 image->
stat = cpl_image_duplicate(cpl_imagelist_get(cube->
stat, iplane));
572 return cpl_error_set_message(__func__, cpl_error_get_code(),
573 "Failure while creating detection image!");
575 #if DEBUG_MUSE_DARCHECK > 1 576 median->
header = cpl_propertylist_new();
577 cpl_propertylist_append_string(median->
header,
"BUNIT",
578 cpl_propertylist_get_string(cube->
header,
585 double dsigmas[] = { 50., 30., 10., 8., 6., 5. };
586 int ndsigmas =
sizeof(dsigmas) /
sizeof(
double);
587 cpl_vector *vsigmas = cpl_vector_wrap(ndsigmas, dsigmas);
588 cpl_size isigma = -1;
589 state = cpl_errorstate_get();
590 cpl_apertures *apertures = cpl_apertures_extract(median->
data, vsigmas, &isigma);
591 int nx = cpl_image_get_size_x(median->
data),
592 ny = cpl_image_get_size_y(median->
data);
594 int napertures = apertures ? cpl_apertures_get_size(apertures) : 0;
595 if (napertures < 1) {
596 cpl_vector_unwrap(vsigmas);
597 cpl_apertures_delete(apertures);
598 cpl_errorstate_set(state);
604 return cpl_error_set_message(__func__, CPL_ERROR_DATA_NOT_FOUND,
"No sources " 605 "found for DAR check down to %.1f sigma " 606 "limit on plane %d", dsigmas[ndsigmas - 1],
609 cpl_msg_debug(__func__,
"The %.1f sigma threshold was used to find %d source%s," 610 " centered on plane %d", cpl_vector_get(vsigmas, isigma),
611 napertures, napertures == 1 ?
"" :
"s", iplane+1);
612 cpl_vector_unwrap(vsigmas);
613 #if DEBUG_MUSE_DARCHECK 614 cpl_apertures_dump(apertures, stdout);
616 #if DEBUG_MUSE_DARCHECK > 1 622 int xhalfsize = 2*cenhalfsize > nx ? nx/2 : cenhalfsize,
623 yhalfsize = 2*cenhalfsize > ny ? ny/2 : cenhalfsize;
626 int l, nlambda = cpl_imagelist_get_size(cube->
data);
627 cpl_vector *vlambda = cpl_vector_new(nlambda);
628 cpl_matrix *moffsets = cpl_matrix_new(2, nlambda);
629 #if !DEBUG_MUSE_DARCHECK 630 #pragma omp parallel for default(none) \ 631 shared(apertures, aPixtable, cdelt3, crpix3, crval3, cube, grid, \ 632 moffsets, napertures, nlambda, nx, ny, vlambda, xhalfsize, \ 635 for (l = 0; l < nlambda; l++) {
636 double xoffset = 0., yoffset = 0.,
637 lambda = (l+1 - crpix3) * cdelt3 + crval3;
639 for (n = 1; n <= napertures; n++) {
641 double xc = cpl_apertures_get_centroid_x(apertures, n),
642 yc = cpl_apertures_get_centroid_y(apertures, n);
643 int x1 = lround(xc) - xhalfsize,
644 x2 = lround(xc) + xhalfsize,
645 y1 = lround(yc) - yhalfsize,
646 y2 = lround(yc) + yhalfsize;
650 if (x2 > nx) x2 = nx;
651 if (y2 > ny) y2 = ny;
653 #define MUSE_DARCHECK_MAX_INTERNAL_SHIFT 5 655 #if MUSE_DARCHECK_GRID_FITS 657 #if DEBUG_MUSE_DARCHECK 658 const char *method =
"grid/moffat";
664 #define MUSE_DARCHECK_MAX_NPIX 7000 665 cpl_matrix *pos = cpl_matrix_new(MUSE_DARCHECK_MAX_NPIX, 2);
666 cpl_vector *val = cpl_vector_new(MUSE_DARCHECK_MAX_NPIX),
667 *err = cpl_vector_new(MUSE_DARCHECK_MAX_NPIX);
668 int i, npix = 0, nsize = MUSE_DARCHECK_MAX_NPIX;
669 for (i = x1 - 1; i < x2; i++) {
671 for (j = y1 - 1; j < y2; j++) {
675 for (ipx = 0; ipx < npx; ipx++) {
676 if (cdq[rows[ipx]] != EURO3D_GOODPIXEL) {
680 nsize += MUSE_DARCHECK_MAX_NPIX;
681 cpl_matrix_resize(pos, 0, nsize - cpl_matrix_get_nrow(pos), 0, 0);
682 cpl_vector_set_size(val, nsize);
683 cpl_vector_set_size(err, nsize);
685 cpl_matrix_set(pos, npix, 0, i+1);
686 cpl_matrix_set(pos, npix, 1, j+1);
687 cpl_vector_set(val, npix, cdata[rows[ipx]]);
688 cpl_vector_set(err, npix, sqrt(cstat[rows[ipx]]));
694 cpl_matrix_delete(pos);
695 cpl_vector_delete(val);
696 cpl_vector_delete(err);
699 cpl_matrix_set_size(pos, npix, 2);
700 cpl_vector_set_size(val, npix);
701 cpl_vector_set_size(err, npix);
709 cpl_array *pmoffat = cpl_array_new(8, CPL_TYPE_DOUBLE);
710 cpl_array_set(pmoffat, 2, xc1);
711 cpl_array_set(pmoffat, 3, yc1);
714 double xc2 = cpl_array_get(pmoffat, 2, NULL),
715 yc2 = cpl_array_get(pmoffat, 3, NULL);
716 cpl_array_delete(pmoffat);
717 cpl_matrix_delete(pos);
718 cpl_vector_delete(val);
719 cpl_vector_delete(err);
723 if (rc == CPL_ERROR_NONE &&
724 fabs(xc1 - xc2) < MUSE_DARCHECK_MAX_INTERNAL_SHIFT &&
725 fabs(yc1 - yc2) < MUSE_DARCHECK_MAX_INTERNAL_SHIFT) {
731 #if DEBUG_MUSE_DARCHECK 732 const char *method =
"image/gauss";
735 cpl_image *plane = cpl_imagelist_get(cube->
data, l),
736 *plerr = cpl_image_power_create(cpl_imagelist_get(cube->
stat, l), 0.5);
740 cpl_image_reject_from_mask(plane, cpl_image_get_bpm(plerr));
742 cpl_image *xtr = cpl_image_extract(plane, x1, y1, x2, y2);
743 int npix = (x2-x1+1) * (y2-y1+1)
744 - cpl_image_count_rejected(xtr);
745 cpl_image_delete(xtr);
747 cpl_image_delete(plerr);
757 cpl_array *pgauss = cpl_array_new(7, CPL_TYPE_DOUBLE);
758 cpl_array_set(pgauss, 3, xc1);
759 cpl_array_set(pgauss, 4, yc1);
760 cpl_error_code rc = cpl_fit_image_gaussian(plane, plerr,
761 lround(xc), lround(yc),
762 2*xhalfsize, 2*yhalfsize,
765 NULL, NULL, NULL, NULL);
766 double xc2 = cpl_array_get(pgauss, 3, NULL),
767 yc2 = cpl_array_get(pgauss, 4, NULL);
768 cpl_array_delete(pgauss);
769 cpl_image_delete(plerr);
773 if (rc == CPL_ERROR_NONE &&
774 fabs(xc1 - xc2) < MUSE_DARCHECK_MAX_INTERNAL_SHIFT &&
775 fabs(yc1 - yc2) < MUSE_DARCHECK_MAX_INTERNAL_SHIFT) {
781 #if DEBUG_MUSE_DARCHECK 782 printf(
"%.1f Angstrom plane %d aper %d (%f,%f): %f %f (%s) %f %f (%d, %s)\n",
783 lambda, l+1, n, xc, yc, xc - xc2, yc - yc2, __func__,
784 xc2 - xc1, yc2 - yc1, npix, method);
789 cpl_vector_set(vlambda, l, lambda);
790 cpl_matrix_set(moffsets, 0, l, xoffset / noffsets);
791 cpl_matrix_set(moffsets, 1, l, yoffset / noffsets);
792 #if DEBUG_MUSE_DARCHECK 793 printf(
"%.1f Angstrom plane %d all-apers: %f %f (%s)\n",
794 lambda, l+1, xoffset / noffsets, yoffset / noffsets, __func__);
798 cpl_apertures_delete(apertures);
801 cpl_vector *vlam1 = cpl_vector_duplicate(vlambda),
802 *vlam2 = cpl_vector_duplicate(vlambda);
803 cpl_matrix *mlam1 = cpl_matrix_wrap(1, nlambda, cpl_vector_unwrap(vlam1)),
804 *mlam2 = cpl_matrix_wrap(1, nlambda, cpl_vector_unwrap(vlam2));
805 cpl_matrix *mxoff = cpl_matrix_extract_row(moffsets, 0),
806 *myoff = cpl_matrix_extract_row(moffsets, 1);
807 cpl_vector *vxoff = cpl_vector_wrap(nlambda, cpl_matrix_unwrap(mxoff)),
808 *vyoff = cpl_vector_wrap(nlambda, cpl_matrix_unwrap(myoff));
811 double msex, msey, chisqx, chisqy;
813 NULL, NULL, order, 3.,
816 NULL, NULL, order, 3.,
818 cpl_vector_delete(vxoff);
819 cpl_vector_delete(vyoff);
820 cpl_matrix_delete(mlam1);
821 cpl_matrix_delete(mlam2);
822 #if DEBUG_MUSE_DARCHECK 823 printf(
"polynomial fit in x (order %d, rms=%f, chisq=%g):\n", order,
825 cpl_polynomial_dump(px, stdout);
826 printf(
"polynomial fit in y (order %d, rms=%f, chisq=%g):\n", order,
828 cpl_polynomial_dump(py, stdout);
834 double xref = cpl_polynomial_eval_1d(px, lambdaref, NULL),
835 yref = cpl_polynomial_eval_1d(py, lambdaref, NULL);
837 cpl_polynomial_set_coeff(px, &pows, cpl_polynomial_get_coeff(px, &pows) - xref);
838 cpl_polynomial_set_coeff(py, &pows, cpl_polynomial_get_coeff(py, &pows) - yref);
839 #if DEBUG_MUSE_DARCHECK 840 printf(
"polynomial in x (xref=%f at %f --> %f):\n", xref, lambdaref,
841 cpl_polynomial_eval_1d(px, lambdaref, NULL));
842 cpl_polynomial_dump(px, stdout);
843 printf(
"polynomial in y (yref=%f at %f --> %f):\n", yref, lambdaref,
844 cpl_polynomial_eval_1d(py, lambdaref, NULL));
845 cpl_polynomial_dump(py, stdout);
850 double xmin = DBL_MAX, xmax = -DBL_MAX,
851 ymin = DBL_MAX, ymax = -DBL_MAX;
852 for (l = 0; l < nlambda; l++) {
853 double lambda = cpl_vector_get(vlambda, l),
854 xshift = cpl_polynomial_eval_1d(px, lambda, NULL),
855 yshift = cpl_polynomial_eval_1d(py, lambda, NULL);
856 if (xshift < xmin) xmin = xshift;
857 if (xshift > xmax) xmax = xshift;
858 if (yshift < ymin) ymin = yshift;
859 if (yshift > ymax) ymax = yshift;
860 #if DEBUG_MUSE_DARCHECK 861 printf(
"%.1f Angstrom plane %d all-fitted: %f %f (%s)\n",
862 lambda, l+1, xshift, yshift, __func__);
866 cpl_vector_delete(vlambda);
867 cpl_matrix_delete(moffsets);
870 double maxdiffx = xmax - xmin,
871 maxdiffy = ymax - ymin;
874 maxdiffx *= kMuseSpaxelSizeX_WFM;
875 maxdiffy *= kMuseSpaxelSizeY_WFM;
877 maxdiffx *= kMuseSpaxelSizeX_NFM;
878 maxdiffy *= kMuseSpaxelSizeY_NFM;
880 *aMaxShift = fmax(maxdiffx, maxdiffy);
884 MUSE_HDR_PT_DAR_CHECK_C);
893 cpl_msg_debug(__func__,
"Correcting pixel table for the shift (max = %f " 894 "arcsec)", *aMaxShift);
899 #pragma omp parallel for default(none) \ 900 shared(lbda, nrow, px, py, xpos, ypos) 901 for (i = 0; i < nrow; i++) {
903 xpos[i] += cpl_polynomial_eval_1d(px, lbda[i], NULL);
904 ypos[i] += cpl_polynomial_eval_1d(py, lbda[i], NULL);
911 MUSE_HDR_PT_DAR_CORR_C);
917 cpl_propertylist_get_float(aPixtable->
header,
920 cpl_propertylist_get_float(aPixtable->
header,
923 cpl_propertylist_get_float(aPixtable->
header,
926 cpl_propertylist_get_float(aPixtable->
header,
931 cpl_polynomial_delete(px);
932 cpl_polynomial_delete(py);
937 return CPL_ERROR_NONE;
#define MUSE_PIXTABLE_XPOS
muse_pixtable_wcs
State of the astrometric calibration of a MUSE pixel table.
Structure definition of a MUSE datacube.
#define MUSE_HDR_PT_PREDAR_YLO
double muse_phys_nrindex_edlen(double l, double t, double p, double pv)
Compute the refractive index for the given wavelength following Edlén.
Structure definition for a collection of muse_images.
cpl_error_code muse_wcs_get_scales(cpl_propertylist *aHeader, double *aXScale, double *aYScale)
Compute the spatial scales (in degrees) from the FITS header WCS.
void muse_image_delete(muse_image *aImage)
Deallocate memory associated to a muse_image object.
double muse_pfits_get_cd(const cpl_propertylist *aHeaders, unsigned int aAxisI, unsigned int aAxisJ)
find out the WCS coordinate at the reference point
#define MUSE_HDR_PT_XLO
FITS header keyword contains the lower limit of the data in x-direction.
double muse_pfits_get_rhum(const cpl_propertylist *aHeaders)
find out the relavtive humidity (in %)
double muse_pfits_get_crval(const cpl_propertylist *aHeaders, unsigned int aAxis)
find out the WCS coordinate at the reference point
double muse_phys_nrindex_filippenko(double l, double T, double f, double P)
Compute the refractive index for the given wavelength following Filippenko.
cpl_image * data
the data extension
cpl_size muse_pixtable_get_nrow(const muse_pixtable *aPixtable)
get the number of rows within the pixel table
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...
void muse_phys_nrindex_owens_coeffs(double temp, double rhum, double pres, double *d1, double *d2)
Compute the two coefficients for the Owens method.
muse_image * muse_combine_median_create(muse_imagelist *aImages)
Median combine a list of input images.
muse_datacube * muse_resampling_cube(muse_pixtable *aPixtable, muse_resampling_params *aParams, muse_pixgrid **aGrid)
Resample a pixel table onto a regular grid structure representing a FITS NAXIS=3 datacube.
double muse_pfits_get_pres_start(const cpl_propertylist *aHeaders)
find out the ambient pressure at start of exposure (in mbar)
cpl_image * stat
the statistics extension
void muse_datacube_delete(muse_datacube *aCube)
Deallocate memory associated to a muse_datacube object.
void muse_imagelist_delete(muse_imagelist *aList)
Free the memory of the MUSE image list.
cpl_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...
double muse_astro_parangle(const cpl_propertylist *aHeader)
Properly average parallactic angle values of one exposure.
Structure definition of MUSE three extension FITS file.
double muse_phys_nrindex_ciddor(double l, double T, double p, double xv, double xCO2)
Compute the refractive index for the given wavelength following Ciddor.
cpl_table * table
The pixel table.
double muse_phys_nrindex_owens(double l, double d1, double d2)
Compute the refractive index for the given wavelength following Owens.
cpl_propertylist * header
the FITS header
cpl_error_code muse_pixtable_reset_dq(muse_pixtable *aPixtable, unsigned int aDQ)
Reset a given bad pixel status (DQ flag) for all pixels in the table.
cpl_error_code muse_datacube_save(muse_datacube *aCube, const char *aFilename)
Save the three cube extensions and the FITS headers of a MUSE datacube to a file. ...
#define MUSE_PIXTABLE_DATA
cpl_image * dq
the data quality extension
double muse_pfits_get_crpix(const cpl_propertylist *aHeaders, unsigned int aAxis)
find out the WCS reference point
void muse_pixgrid_delete(muse_pixgrid *aGrid)
Delete a pixel grid and remove its memory.
#define MUSE_HDR_PT_PREDAR_YHI
Structure definition of MUSE pixel table.
static const cpl_size * muse_pixgrid_get_rows(muse_pixgrid *aGrid, cpl_size aIndex)
Return a pointer to the rows stored in one pixel.
muse_pixtable_wcs muse_pixtable_wcs_check(muse_pixtable *aPixtable)
Check the state of the world coordinate system of a pixel table.
#define MUSE_HDR_PT_YLO
FITS header keyword contains the lower limit of the data in y-direction.
cpl_error_code muse_utils_fit_moffat_2d(const cpl_matrix *aPositions, const cpl_vector *aValues, const cpl_vector *aErrors, cpl_array *aParams, cpl_array *aPErrors, const cpl_array *aPFlags, double *aRMS, double *aRedChisq)
Fit a 2D Moffat function to a given set of data.
muse_resampling_params * muse_resampling_params_new(muse_resampling_type aMethod)
Create the resampling parameters structure.
cpl_polynomial * muse_utils_iterate_fit_polynomial(cpl_matrix *aPos, cpl_vector *aVal, cpl_vector *aErr, cpl_table *aExtra, const unsigned int aOrder, const double aRSigma, double *aMSE, double *aChiSq)
Iterate a polynomial fit.
#define MUSE_HDR_PT_PREDAR_XLO
cpl_imagelist * data
the cube containing the actual data values
#define MUSE_HDR_PT_PREDAR_XHI
cpl_imagelist * dq
the optional cube containing the bad pixel status
cpl_propertylist * header
the FITS header
double muse_astro_posangle(const cpl_propertylist *aHeader)
Derive the position angle of an observation from information in a FITS header.
double muse_pfits_get_pres_end(const cpl_propertylist *aHeaders)
find out the ambient pressure at end of exposure (in mbar)
#define MUSE_PIXTABLE_STAT
cpl_error_code muse_image_save(muse_image *aImage, const char *aFilename)
Save the three image extensions and the FITS headers of a MUSE image to a file.
cpl_error_code muse_image_reject_from_dq(muse_image *aImage)
Reject pixels of a muse_image depending on its DQ data.
muse_imagelist * muse_imagelist_new(void)
Create a new (empty) MUSE image list.
#define MUSE_PIXTABLE_LAMBDA
static cpl_size muse_pixgrid_get_index(muse_pixgrid *aGrid, cpl_size aX, cpl_size aY, cpl_size aZ, cpl_boolean aAllowOutside)
Get the grid index determined from all three coordinates.
#define MUSE_HDR_PT_DAR_NAME
double muse_pfits_get_temp(const cpl_propertylist *aHeaders)
find out the ambient temperature (in degrees Celsius)
cpl_error_code muse_utils_image_get_centroid_window(cpl_image *aImage, int aX1, int aY1, int aX2, int aY2, double *aX, double *aY, muse_utils_centroid_type aBgType)
Compute centroid over an image window, optionally marginalizing over the background.
double muse_phys_nrindex_owens_saturation_pressure(double temp)
Compute the saturation pressure using the Owens calibration.
void muse_resampling_params_delete(muse_resampling_params *aParams)
Delete a resampling parameters structure.
#define MUSE_HDR_PT_DAR_CORR
muse_image * muse_image_new(void)
Allocate memory for a new muse_image object.
#define MUSE_HDR_PT_DAR_CHECK
#define MUSE_PIXTABLE_YPOS
#define MUSE_HDR_PT_YHI
FITS header keyword contains the upper limit of the data in y-direction.
cpl_error_code muse_imagelist_set(muse_imagelist *aList, muse_image *aImage, unsigned int aIdx)
Set the muse_image of given list index.
cpl_error_code muse_pixtable_compute_limits(muse_pixtable *aPixtable)
(Re-)Compute the limits of the coordinate columns of a pixel table.
static cpl_size muse_pixgrid_get_count(muse_pixgrid *aGrid, cpl_size aIndex)
Return the number of rows stored in one pixel.
int muse_quality_image_reject_using_dq(cpl_image *aData, cpl_image *aDQ, cpl_image *aStat)
Reject pixels of one or two images on a DQ image.
muse_ins_mode muse_pfits_get_mode(const cpl_propertylist *aHeaders)
find out the observation mode
cpl_error_code muse_utils_get_centroid(const cpl_matrix *aPositions, const cpl_vector *aValues, const cpl_vector *aErrors, double *aX, double *aY, muse_utils_centroid_type aBgType)
Compute centroid of a two-dimensional dataset.
double muse_astro_airmass(cpl_propertylist *aHeader)
Derive the effective airmass of an observation from information in a FITS header. ...
cpl_imagelist * stat
the cube containing the data variance
cpl_propertylist * header
The FITS header.
#define MUSE_HDR_PT_XHI
FITS header keyword contains the upper limit of the data in x-ion.