33 #include "muse_instrument.h"
35 #include "muse_astro.h"
36 #include "muse_combine.h"
38 #include "muse_pfits.h"
39 #include "muse_quality.h"
40 #include "muse_resampling.h"
41 #include "muse_utils.h"
47 #define FAKE_POS_ROT 0
96 cpl_table_delete(aWCSObj->detected);
97 aWCSObj->detected = NULL;
98 cpl_propertylist_delete(aWCSObj->wcs);
117 {
"XPOS", CPL_TYPE_DOUBLE,
"pix",
"%f",
"horizontal position", CPL_TRUE },
118 {
"YPOS", CPL_TYPE_DOUBLE,
"pix",
"%f",
"vertical position", CPL_TRUE },
119 {
"XERR", CPL_TYPE_DOUBLE,
"pix",
"%f",
120 "error estimate of the horizontal position", CPL_TRUE },
121 {
"YERR", CPL_TYPE_DOUBLE,
"pix",
"%f",
122 "error estimate of the vertical position", CPL_TRUE },
123 {
"FLUX", CPL_TYPE_DOUBLE,
"count",
"%e",
"(relative) flux of the source", CPL_TRUE },
124 {
"XFWHM", CPL_TYPE_DOUBLE,
"pix",
"%f",
"horizontal FWHM", CPL_TRUE },
125 {
"YFWHM", CPL_TYPE_DOUBLE,
"pix",
"%f",
"vertical FWHM", CPL_TRUE },
126 { NULL, 0, NULL, NULL, NULL, CPL_FALSE }
142 {
"SourceID", CPL_TYPE_STRING,
"",
"%s",
"the source identification", CPL_TRUE },
143 {
"RA", CPL_TYPE_DOUBLE,
"deg",
"%f",
"right ascension (decimal degrees)", CPL_TRUE },
144 {
"DEC", CPL_TYPE_DOUBLE,
"deg",
"%f",
"declination (decimal degrees)", CPL_TRUE },
145 {
"filter", CPL_TYPE_STRING,
"",
"%s",
"the filter used for the magnitude", CPL_TRUE },
146 {
"mag", CPL_TYPE_DOUBLE,
"mag",
"%f",
"the object (Vega) magnitude", CPL_TRUE },
147 { NULL, 0, NULL, NULL, NULL, CPL_FALSE }
202 cpl_ensure_code(aPixtable && aPixtable->
header && aWCSObj,
203 CPL_ERROR_NULL_INPUT);
204 cpl_ensure_code(aSigma >0., CPL_ERROR_ILLEGAL_INPUT);
207 cpl_msg_info(__func__,
"Gaussian profile fits for object centroiding");
210 cpl_msg_info(__func__,
"Moffat profile fits for object centroiding");
213 cpl_msg_info(__func__,
"Simple square box object centroiding");
216 return cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT);
218 if (getenv(
"MUSE_DEBUG_WCS") && atoi(getenv(
"MUSE_DEBUG_WCS")) > 2) {
219 const char *fn =
"wcs__pixtable.fits";
220 cpl_msg_info(__func__,
"Saving pixel table as \"%s\"", fn);
225 cpl_boolean darcorrected = CPL_FALSE;
228 darcorrected = CPL_TRUE;
230 darcorrected = CPL_TRUE;
233 const char *message =
"Correction for differential atmospheric refraction "
234 "was not applied! Deriving astrometric correction "
235 "from such data is unlikely to give good results!";
236 cpl_msg_warning(__func__,
"%s", message);
237 cpl_error_set_message(__func__, CPL_ERROR_UNSUPPORTED_MODE,
"%s", message);
245 params->dlambda = 50.;
254 return cpl_error_set_message(__func__, cpl_error_get_code(),
255 "Failure while creating cube!");
257 aWCSObj->cube = cube;
258 if (getenv(
"MUSE_DEBUG_WCS") && atoi(getenv(
"MUSE_DEBUG_WCS")) >= 2) {
259 const char *fn =
"wcs__cube.fits";
260 cpl_msg_info(__func__,
"Saving cube as \"%s\"", fn);
265 int iplane, irefplane = cpl_imagelist_get_size(cube->
data) / 2;
270 unsigned int ilist = 0;
271 for (iplane = irefplane - 1; iplane <= irefplane + 1; iplane++) {
273 image->
data = cpl_image_duplicate(cpl_imagelist_get(cube->
data, iplane));
274 image->
dq = cpl_image_duplicate(cpl_imagelist_get(cube->
dq, iplane));
275 image->
stat = cpl_image_duplicate(cpl_imagelist_get(cube->
stat, iplane));
281 return cpl_error_set_message(__func__, cpl_error_get_code(),
282 "Failure while creating detection image!");
284 if (getenv(
"MUSE_DEBUG_WCS") && atoi(getenv(
"MUSE_DEBUG_WCS")) >= 2) {
285 const char *fn =
"wcs__image_median.fits";
286 cpl_msg_info(__func__,
"Saving median detection image as \"%s\"", fn);
287 median->
header = cpl_propertylist_new();
288 cpl_propertylist_append_string(median->
header,
"BUNIT",
289 cpl_propertylist_get_string(cube->
header,
296 cpl_apertures *apertures = cpl_apertures_extract_sigma(median->
data, aSigma);
297 int ndet = apertures ? cpl_apertures_get_size(apertures) : 0;
300 return cpl_error_set_message(__func__, CPL_ERROR_DATA_NOT_FOUND,
301 "No sources found for astrometric calibration "
302 "down to %.1f sigma limit on plane %d", aSigma,
305 cpl_msg_debug(__func__,
"The %f sigma threshold was used to find %d sources "
306 "in cube plane %d", aSigma, ndet, irefplane + 1);
308 aWCSObj->xcenter = cpl_image_get_size_x(median->
data) / 2.,
309 aWCSObj->ycenter = cpl_image_get_size_y(median->
data) / 2.;
313 cpl_msg_debug(__func__,
"image size: %d x %d --> center %f, %f",
314 (
int)cpl_image_get_size_x(median->
data),
315 (
int)cpl_image_get_size_y(median->
data),
316 aWCSObj->xcenter, aWCSObj->ycenter);
320 cpl_image_power(median->
stat, 0.5);
321 int nx = cpl_image_get_size_x(median->
data),
322 ny = cpl_image_get_size_y(median->
data);
324 const double bgguess = cpl_image_get_median(median->
data),
326 moffat_alpha_to_fwhm = 1 / (2 * sqrt(pow(2, 1/beta) - 1));
330 if (fwhmguess < FLT_EPSILON) {
334 fwhmguess /= (kMuseSpaxelSizeX_WFM + kMuseSpaxelSizeX_WFM) / 2;
336 fwhmguess /= (kMuseSpaxelSizeX_NFM + kMuseSpaxelSizeX_NFM) / 2;
342 cpl_table_unselect_all(detections);
344 for (idet = 0; idet < ndet; idet++) {
345 double xc = cpl_apertures_get_centroid_x(apertures, idet+1),
346 yc = cpl_apertures_get_centroid_y(apertures, idet+1),
347 fluxaper = cpl_apertures_get_flux(apertures, idet+1);
348 double xwaper, ywaper;
349 cpl_image_get_fwhm(median->
data, lround(xc), lround(yc), &xwaper, &ywaper);
353 if (xwaper < 0 || ywaper < 0) {
355 cpl_msg_debug(__func__,
"FWHM computation failed at %f,%f, result was "
356 "%.3f,%.3f", xc, yc, xwaper, ywaper);
358 cpl_table_select_row(detections, idet);
364 int x1 = floor(xc) - 5,
371 if (x2 > nx) x2 = nx;
372 if (y2 > ny) y2 = ny;
374 double xcen, ycen, xerr = 0., yerr = 0., xw, yw, flux;
375 cpl_errorstate state = cpl_errorstate_get();
378 cpl_array *pgauss = cpl_array_new(7, CPL_TYPE_DOUBLE),
379 *pgerr = cpl_array_new(7, CPL_TYPE_DOUBLE),
380 *pgfit = cpl_array_new(7, CPL_TYPE_INT);
382 cpl_array_set(pgfit, 3, 1);
383 cpl_array_set(pgfit, 4, 1);
384 cpl_array_set(pgfit, 0, 0);
385 cpl_array_set(pgfit, 1, 0);
386 cpl_array_set(pgfit, 2, 0);
387 cpl_array_set(pgfit, 5, 0);
388 cpl_array_set(pgfit, 6, 0);
389 cpl_array_set(pgauss, 3, xc);
390 cpl_array_set(pgauss, 4, yc);
391 cpl_array_set(pgauss, 0, bgguess);
392 cpl_array_set(pgauss, 1, fluxaper);
393 cpl_array_set(pgauss, 2, 0.);
394 cpl_array_set(pgauss, 5, xwaper * CPL_MATH_SIG_FWHM);
395 cpl_array_set(pgauss, 6, ywaper * CPL_MATH_SIG_FWHM);
396 cpl_fit_image_gaussian(median->
data, median->
stat, lround(xc), lround(yc),
397 x2-x1+1, y2-y1+1, pgauss, pgerr, pgfit, NULL, NULL,
398 NULL, NULL, NULL, NULL, NULL);
399 xcen = cpl_array_get(pgauss, 3, NULL);
400 ycen = cpl_array_get(pgauss, 4, NULL);
401 xerr = cpl_array_get(pgerr, 3, NULL);
402 yerr = cpl_array_get(pgerr, 4, NULL);
404 cpl_array_set(pgfit, 3, 0);
405 cpl_array_set(pgfit, 4, 0);
406 cpl_array_set(pgfit, 1, 1);
407 cpl_array_set(pgfit, 5, 1);
408 cpl_array_set(pgfit, 6, 1);
409 cpl_fit_image_gaussian(median->
data, median->
stat, lround(xc), lround(yc),
410 x2-x1+1, y2-y1+1, pgauss, pgerr, pgfit, NULL, NULL,
411 NULL, NULL, NULL, NULL, NULL);
412 xw = cpl_array_get(pgauss, 5, NULL) * CPL_MATH_FWHM_SIG;
413 yw = cpl_array_get(pgauss, 6, NULL) * CPL_MATH_FWHM_SIG;
414 flux = cpl_array_get(pgauss, 1, NULL);
415 cpl_array_delete(pgauss);
416 cpl_array_delete(pgerr);
417 cpl_array_delete(pgfit);
421 cpl_size nmax = (x2-x1+1) * (y2-y1+1);
422 cpl_matrix *pos = cpl_matrix_new(nmax, 2);
423 cpl_vector *val = cpl_vector_new(nmax),
424 *err = cpl_vector_new(nmax);
426 for (i = x1; i <= x2; i++) {
428 for (j = y1; j <= y2; j++) {
430 double value = cpl_image_get(median->
data, i, j, &error);
434 cpl_matrix_set(pos, npix, 0, i);
435 cpl_matrix_set(pos, npix, 1, j);
436 cpl_vector_set(val, npix, value);
438 cpl_vector_set(err, npix, cpl_image_get(median->
stat, i, j, &error));
442 cpl_matrix_set_size(pos, npix, 2);
443 cpl_vector_set_size(val, npix);
444 cpl_vector_set_size(err, npix);
445 cpl_array *pmoffat = cpl_array_new(8, CPL_TYPE_DOUBLE),
446 *pmerror = cpl_array_new(8, CPL_TYPE_DOUBLE),
447 *pmfit = cpl_array_new(8, CPL_TYPE_INT);
449 cpl_array_set(pmfit, 2, 1);
450 cpl_array_set(pmfit, 3, 1);
451 cpl_array_set(pmfit, 0, 0);
452 cpl_array_set(pmfit, 1, 0);
453 cpl_array_set(pmfit, 4, 0);
454 cpl_array_set(pmfit, 5, 0);
455 cpl_array_set(pmfit, 6, 0);
456 cpl_array_set(pmfit, 7, 0);
457 cpl_array_set(pmoffat, 2, xc);
458 cpl_array_set(pmoffat, 3, yc);
459 cpl_array_set(pmoffat, 0, bgguess);
460 cpl_array_set(pmoffat, 1, fluxaper);
461 cpl_array_set(pmoffat, 4, xwaper * moffat_alpha_to_fwhm);
462 cpl_array_set(pmoffat, 5, ywaper * moffat_alpha_to_fwhm);
463 cpl_array_set(pmoffat, 6, beta);
464 cpl_array_set(pmoffat, 7, 0.);
466 xcen = cpl_array_get(pmoffat, 2, NULL);
467 ycen = cpl_array_get(pmoffat, 3, NULL);
468 xerr = cpl_array_get(pmerror, 2, NULL);
469 yerr = cpl_array_get(pmerror, 3, NULL);
471 cpl_array_set(pmfit, 2, 0);
472 cpl_array_set(pmfit, 3, 0);
473 cpl_array_set(pmfit, 1, 1);
474 cpl_array_set(pmfit, 4, 1);
475 cpl_array_set(pmfit, 5, 1);
477 xw = cpl_array_get(pmoffat, 4, NULL) / moffat_alpha_to_fwhm;
478 yw = cpl_array_get(pmoffat, 5, NULL) / moffat_alpha_to_fwhm;
479 flux = cpl_array_get(pmoffat, 1, NULL);
480 cpl_array_delete(pmoffat);
481 cpl_array_delete(pmerror);
482 cpl_array_delete(pmfit);
483 cpl_matrix_delete(pos);
484 cpl_vector_delete(val);
485 cpl_vector_delete(err);
493 printf(
"%d apertures %f %f boxes %f %f deltas %f %f\n", idet+1, xc, yc,
494 xcen, ycen, xcen - xc, ycen - yc);
502 cpl_image_get_fwhm(median->
data, lround(xcen), lround(ycen), &xw, &yw);
507 cpl_table_set(detections,
"XPOS", idet, xcen);
508 cpl_table_set(detections,
"YPOS", idet, ycen);
510 cpl_table_set(detections,
"XERR", idet, xerr);
511 cpl_table_set(detections,
"YERR", idet, yerr);
513 cpl_table_set(detections,
"XFWHM", idet, xw);
514 cpl_table_set(detections,
"YFWHM", idet, yw);
516 cpl_table_set(detections,
"FLUX", idet, flux);
518 if (cpl_errorstate_is_equal(state) && xw > 0 && yw > 0 &&
519 xcen >= 1 && xcen <= nx && ycen >= 1 && ycen <= ny) {
523 cpl_table_select_row(detections, idet);
525 cpl_table_erase_selected(detections);
526 cpl_apertures_delete(apertures);
528 ndet = cpl_table_get_nrow(detections);
529 cpl_msg_debug(__func__,
"%d valid sources were recorded in the detections "
532 if (getenv(
"MUSE_DEBUG_WCS") && atoi(getenv(
"MUSE_DEBUG_WCS")) >= 2) {
533 const char *fn =
"wcs__detections.fits";
534 cpl_msg_info(__func__,
"Saving table of detections as \"%s\"", fn);
535 cpl_table_save(detections, NULL, NULL, fn, CPL_IO_CREATE);
538 aWCSObj->detected = detections;
539 return CPL_ERROR_NONE;
583 float aRadius,
float aFAccuracy,
int aIter,
float aThresh)
585 cpl_ensure_code(aWCSObj, CPL_ERROR_NULL_INPUT);
586 cpl_table *detected = aWCSObj->detected;
587 int ndet = cpl_table_get_nrow(detected),
588 nref = cpl_table_get_nrow(aReference);
589 cpl_ensure_code(ndet > 0 && nref > 0, CPL_ERROR_ILLEGAL_INPUT);
594 CPL_ERROR_BAD_FILE_FORMAT);
595 cpl_ensure_code(aRadius > 0.&& aFAccuracy > 0., CPL_ERROR_ILLEGAL_INPUT);
598 cpl_propertylist *order = cpl_propertylist_new();
599 cpl_propertylist_append_bool(order,
"FLUX", TRUE);
600 cpl_table_sort(detected, order);
601 cpl_propertylist_erase(order,
"FLUX");
602 cpl_propertylist_append_bool(order,
"mag", FALSE);
603 cpl_table_sort(aReference, order);
604 cpl_propertylist_delete(order);
605 if (getenv(
"MUSE_DEBUG_WCS") && atoi(getenv(
"MUSE_DEBUG_WCS")) >= 2) {
606 FILE *fp = fopen(
"wcs__detections.ascii",
"w");
607 fprintf(fp,
"%s: table of detected sources (sorted by flux):\n", __func__);
608 cpl_table_dump(detected, 0, 1000, fp);
610 fp = fopen(
"wcs__references.ascii",
"w");
611 fprintf(fp,
"%s: table of reference objects (sorted by flux):\n", __func__);
612 cpl_table_dump(aReference, 0, 1000, fp);
618 cpl_propertylist_append_double(wcsin,
"CRVAL1", aWCSObj->ra);
619 cpl_propertylist_append_double(wcsin,
"CRVAL2", aWCSObj->dec);
620 cpl_propertylist_update_double(wcsin,
"CRPIX1", aWCSObj->crpix1);
621 cpl_propertylist_update_double(wcsin,
"CRPIX2", aWCSObj->crpix2);
623 cpl_propertylist_append_int(wcsin,
"NAXIS", 2);
624 cpl_propertylist_append_int(wcsin,
"NAXIS1", aWCSObj->xcenter * 2.);
625 cpl_propertylist_append_int(wcsin,
"NAXIS2", aWCSObj->ycenter * 2.);
628 cpl_matrix *data = cpl_matrix_new(2, ndet),
629 *patt = cpl_matrix_new(2, nref);
631 for (i = 0; i < ndet; i++) {
632 cpl_matrix_set(data, 0, i, cpl_table_get(detected,
"XPOS", i, NULL));
633 cpl_matrix_set(data, 1, i, cpl_table_get(detected,
"YPOS", i, NULL));
639 for (i = 0; i < nref; i++) {
640 double ra = cpl_table_get(aReference,
"RA", i, NULL),
641 dec = cpl_table_get(aReference,
"DEC", i, NULL),
644 cpl_matrix_set(patt, 0, i, x);
645 cpl_matrix_set(patt, 1, i, y);
647 printf(
"conversion: %2d\t%.7f %.7f\t%6.2f %6.2f\n", i + 1, ra, dec, x, y);
652 if (cpl_msg_get_level() == CPL_MSG_DEBUG) {
653 printf(
"%s: data matrix:\n", __func__);
654 cpl_matrix_dump(data, stdout);
655 printf(
"%s: patt matrix:\n", __func__);
656 cpl_matrix_dump(patt, stdout);
663 double accuracy = sqrt(pow(cpl_table_get_column_mean(detected,
"XERR"), 2) +
664 pow(cpl_table_get_column_mean(detected,
"YERR"), 2));
667 accuracy *= aFAccuracy;
668 double radius = aRadius;
670 cpl_array *aid = NULL;
671 cpl_boolean dupes = CPL_FALSE;
672 double xscale, yscale;
680 int ndata = ndet < USE_DATA ? ndet : USE_DATA,
681 npatt = nref < USE_PATT ? nref : USE_PATT;
682 cpl_array_delete(aid);
684 cpl_msg_debug(__func__,
"trying pattern matching with accuracy %g and radius %g",
686 aid = cpl_ppm_match_points(data, ndata, accuracy,
688 0.1 , radius, NULL, NULL,
691 accuracy /= aid ? 1. : 1.5;
692 if (accuracy < FLT_EPSILON) {
696 cpl_errorstate state = cpl_errorstate_get();
697 nid = cpl_array_get_size(aid);
699 nid -= cpl_array_count_invalid(aid);
702 printf(
"aid (valid=%d):\n", nid);
703 cpl_array_dump(aid, 0, cpl_array_get_size(aid), stdout);
707 cpl_array_delete(aid);
708 cpl_matrix_delete(data);
709 cpl_matrix_delete(patt);
710 cpl_errorstate_set(state);
711 cpl_propertylist_delete(wcsin);
712 return cpl_error_set_message(__func__, CPL_ERROR_DATA_NOT_FOUND,
"None of "
713 "the %d detections could be identified with "
714 "the %d reference positions with radius %.3f "
715 "pix and data accuracy %.3e pix",
716 ndet, nref, radius, accuracy);
722 cpl_msg_debug(__func__,
"%d %sidentified points (scale = %g, angle = %g; "
723 "used radius = %.3f pix, data accuracy = %.3e pix)", nid,
724 dupes ?
"(non-unique!) " :
"unique ",
725 scale*1800.*(xscale+yscale), angle, radius, accuracy);
728 cpl_matrix_delete(data);
729 cpl_matrix_delete(patt);
734 cpl_matrix *mpx = cpl_matrix_new(nid, 2),
735 *msky = cpl_matrix_new(nid, 2);
737 for (i = 0; i < cpl_array_get_size(aid); i++) {
738 if (!cpl_array_is_valid(aid, i)) {
741 int idata = cpl_array_get_int(aid, i, NULL);
742 cpl_matrix_set(mpx, iid, 0, cpl_table_get(detected,
"XPOS", idata, NULL));
743 cpl_matrix_set(mpx, iid, 1, cpl_table_get(detected,
"YPOS", idata, NULL));
744 cpl_matrix_set(msky, iid, 0, cpl_table_get(aReference,
"RA", i, NULL));
745 cpl_matrix_set(msky, iid, 1, cpl_table_get(aReference,
"DEC", i, NULL));
747 printf(
"matrices: %2d\t%.7f %.7f\t%6.2f %6.2f\n", iid + 1,
748 cpl_table_get(aReference,
"RA", i, NULL),
749 cpl_table_get(aReference,
"DEC", i, NULL),
750 cpl_table_get(detected,
"XPOS", idata, NULL),
751 cpl_table_get(detected,
"YPOS", idata, NULL));
757 cpl_matrix_dump(mpx, stdout);
759 cpl_matrix_dump(msky, stdout);
762 cpl_array_delete(aid);
765 cpl_propertylist *wcsout = NULL;
766 cpl_error_code rc = cpl_wcs_platesol(wcsin, msky, mpx, aIter, aThresh,
767 CPL_WCS_PLATESOL_6, CPL_WCS_MV_CRVAL,
770 cpl_propertylist_copy_property_regexp(wcsout, aWCSObj->cube->
header,
773 cpl_matrix_delete(mpx);
774 cpl_matrix_delete(msky);
775 cpl_propertylist_delete(wcsin);
776 if (rc != CPL_ERROR_NONE) {
777 cpl_msg_warning(__func__,
"Computing the WCS solution returned an error "
778 "(%d): %s", rc, cpl_error_get_message());
782 double cd11 = cpl_propertylist_get_double(wcsout,
"CD1_1"),
783 cd22 = cpl_propertylist_get_double(wcsout,
"CD2_2"),
784 cd12 = cpl_propertylist_get_double(wcsout,
"CD1_2"),
785 cd21 = cpl_propertylist_get_double(wcsout,
"CD2_1"),
786 det = cd11 * cd22 - cd12 * cd21;
797 cpl_msg_info(__func__,
"astrometric calibration results: scales %f/%f "
798 "arcsec/spaxel, rotation %g/%g deg", xscale, yscale, xang, yang);
801 printf(
"%s: output propertylist (rc = %d):\n", __func__, rc);
803 cpl_propertylist_save(wcsout,
"astrometry_wcsout.fits", CPL_IO_CREATE);
804 system(
"fold -w 80 astrometry_wcsout.fits | grep -v \"^ \"");
805 remove(
"astrometry_wcsout.fits");
809 cpl_propertylist_update_int(wcsout, QC_ASTROMETRY_NSTARS, nid);
811 cpl_propertylist_update_float(wcsout, QC_ASTROMETRY_SCX, xscale);
812 cpl_propertylist_update_float(wcsout, QC_ASTROMETRY_SCY, yscale);
814 cpl_propertylist_update_float(wcsout, QC_ASTROMETRY_ANGX, xang);
815 cpl_propertylist_update_float(wcsout, QC_ASTROMETRY_ANGY, yang);
817 double medresx = cpl_propertylist_get_double(wcsout,
"CSYER1") * 3600.,
818 medresy = cpl_propertylist_get_double(wcsout,
"CSYER2") * 3600.;
819 cpl_propertylist_update_float(wcsout, QC_ASTROMETRY_RESX, medresx);
820 cpl_propertylist_update_float(wcsout, QC_ASTROMETRY_RESY, medresy);
822 aWCSObj->wcs = wcsout;
840 cpl_propertylist *wcs = cpl_propertylist_new();
844 cpl_propertylist_append_int(wcs,
"WCSAXES", 2);
849 double smallvalue = FLT_MIN;
850 const char *cpldesc = cpl_get_description(CPL_DESCRIPTION_DEFAULT);
852 char *pcpldesc = strstr(cpldesc,
"WCSLIB = ");
855 double wcslibversion = atof(pcpldesc);
856 if (wcslibversion >= 4.23) {
859 cpl_msg_debug(__func__,
"found wcslib %.2f, using smallvalue = %e",
860 wcslibversion, smallvalue);
865 cpl_propertylist_append_double(wcs,
"CRPIX1", smallvalue);
867 cpl_propertylist_append_double(wcs,
"CD1_1", -kMuseSpaxelSizeX_WFM / 3600);
868 cpl_propertylist_append_string(wcs,
"CTYPE1",
"RA---TAN");
869 cpl_propertylist_append_string(wcs,
"CUNIT1",
"deg");
872 cpl_propertylist_append_double(wcs,
"CRPIX2", smallvalue);
873 cpl_propertylist_append_double(wcs,
"CD2_2", kMuseSpaxelSizeY_WFM / 3600);
874 cpl_propertylist_append_string(wcs,
"CTYPE2",
"DEC--TAN");
875 cpl_propertylist_append_string(wcs,
"CUNIT2",
"deg");
878 cpl_propertylist_append_double(wcs,
"CD1_2", 0.);
879 cpl_propertylist_append_double(wcs,
"CD2_1", 0.);
904 const cpl_propertylist *aCalHeader)
906 cpl_ensure(aCalHeader && aExpHeader, CPL_ERROR_NULL_INPUT, NULL);
909 cpl_propertylist *header = cpl_propertylist_duplicate(aCalHeader);
914 xang, yang, xsc, ysc;
917 cpl_msg_debug(__func__,
"WCS solution: scales %f / %f arcsec, angles %f / %f "
918 "deg", xsc * 3600., ysc * 3600., xang, yang);
920 cpl_matrix *pc = cpl_matrix_new(2, 2);
921 cpl_matrix_set(pc, 0, 0, cpl_propertylist_get_double(header,
"CD1_1") / xsc);
922 cpl_matrix_set(pc, 0, 1, cpl_propertylist_get_double(header,
"CD1_2") / xsc);
923 cpl_matrix_set(pc, 1, 0, cpl_propertylist_get_double(header,
"CD2_1") / ysc);
924 cpl_matrix_set(pc, 1, 1, cpl_propertylist_get_double(header,
"CD2_2") / ysc);
925 cpl_matrix *pcinv = cpl_matrix_invert_create(pc);
926 cpl_matrix_delete(pc);
928 double cd11cor, cd12cor, cd21cor, cd22cor;
930 cd11cor = xsc * cpl_matrix_get(pcinv, 0, 0);
931 cd12cor = xsc * cpl_matrix_get(pcinv, 0, 1);
932 cd21cor = ysc * cpl_matrix_get(pcinv, 1, 0);
933 cd22cor = ysc * cpl_matrix_get(pcinv, 1, 1);
934 cpl_matrix_delete(pcinv);
936 cpl_msg_warning(__func__,
"CD matrix of astrometric solution could not "
937 "be inverted! Assuming negligible zeropoint rotation.");
944 double cd11 = cos(pa) * cd11cor - sin(pa) * cd21cor,
945 cd12 = cos(pa) * cd12cor - sin(pa) * cd22cor,
946 cd21 = sin(pa) * cd11cor + cos(pa) * cd21cor,
947 cd22 = sin(pa) * cd12cor + cos(pa) * cd22cor;
948 cpl_propertylist_update_double(header,
"CD1_1", cd11),
949 cpl_propertylist_update_double(header,
"CD1_2", cd12),
950 cpl_propertylist_update_double(header,
"CD2_1", cd21);
951 cpl_propertylist_update_double(header,
"CD2_2", cd22),
954 cpl_msg_debug(__func__,
"Updated CD matrix (%f deg field rotation): "
955 "%e %e %e %e (scales %f / %f arcsec, angles %f / %f deg)",
956 pa * CPL_MATH_DEG_RAD, cd11, cd12, cd21, cd22,
957 xsc * 3600., ysc * 3600., xang, yang);
996 cpl_ensure_code(nrow > 0 && aPixtable->
header && aWCS, CPL_ERROR_NULL_INPUT);
998 CPL_ERROR_INCOMPATIBLE_INPUT);
1000 const char *type1 = cpl_propertylist_get_string(aWCS,
"CTYPE1"),
1001 *type2 = cpl_propertylist_get_string(aWCS,
"CTYPE2");
1002 cpl_ensure_code(type1 && type2 && !strncmp(type1,
"RA---TAN", 9) &&
1003 !strncmp(type2,
"DEC--TAN", 9),
1004 CPL_ERROR_UNSUPPORTED_MODE);
1008 cpl_propertylist_erase_regexp(aPixtable->
header, MUSE_WCS_KEYS, 0);
1015 cpl_propertylist_erase_regexp(header,
"^CRVAL[0-9]+$|^L[OA][NT]POLE$", 0);
1016 double cd11 = cpl_propertylist_get_double(header,
"CD1_1"),
1017 cd12 = cpl_propertylist_get_double(header,
"CD1_2"),
1018 cd21 = cpl_propertylist_get_double(header,
"CD2_1"),
1019 cd22 = cpl_propertylist_get_double(header,
"CD2_2");
1024 cpl_errorstate prestate = cpl_errorstate_get();
1029 if (!cpl_errorstate_is_equal(prestate)) {
1030 cpl_errorstate_set(prestate);
1037 double wcspix1 = cpl_propertylist_get_double(header,
"CRPIX1"),
1038 wcspix2 = cpl_propertylist_get_double(header,
"CRPIX2"),
1039 crpix1 = (xhi + xlo) / 2. + wcspix1,
1040 crpix2 = (yhi + ylo) / 2. + wcspix2;
1041 cpl_propertylist_update_double(header,
"CRPIX1", crpix1),
1042 cpl_propertylist_update_double(header,
"CRPIX2", crpix2),
1043 cpl_msg_debug(__func__,
"Using reference pixel %f/%f (limits in pixel table "
1044 "%f..%f/%f..%f, WCS correction %f,%f)", crpix1, crpix2,
1045 xlo, xhi, ylo, yhi, wcspix1, wcspix2);
1048 cpl_table_set_column_unit(aPixtable->
table, MUSE_PIXTABLE_XPOS,
"");
1049 cpl_table_set_column_unit(aPixtable->
table, MUSE_PIXTABLE_YPOS,
"");
1053 float *xpos = cpl_table_get_data_float(aPixtable->
table, MUSE_PIXTABLE_XPOS),
1054 *ypos = cpl_table_get_data_float(aPixtable->
table, MUSE_PIXTABLE_YPOS);
1056 #pragma omp parallel for default(none) \
1057 shared(cd11, cd12, cd21, cd22, crpix1, crpix2, nrow, xpos, ypos)
1058 for (n = 0; n < nrow; n++) {
1061 double x = cd11 * (xpos[n] - crpix1) + cd12 * (ypos[n] - crpix2),
1062 y = cd22 * (ypos[n] - crpix2) + cd21 * (xpos[n] - crpix1);
1070 double phi = atan2(x, -y),
1071 theta = atan(CPL_MATH_DEG_RAD / sqrt(x*x + y*y));
1073 phi += CPL_MATH_2PI;
1079 ypos[n] = theta - CPL_MATH_PI_2;
1085 cpl_table_set_column_unit(aPixtable->
table, MUSE_PIXTABLE_XPOS,
"rad");
1086 cpl_table_set_column_unit(aPixtable->
table, MUSE_PIXTABLE_YPOS,
"rad");
1089 cpl_propertylist_copy_property_regexp(aPixtable->
header, header,
1091 cpl_propertylist_delete(header);
1097 MUSE_HDR_PT_WCS_COMMENT_PROJ);
1098 return CPL_ERROR_NONE;
1136 cpl_ensure_code(nrow > 0 && aPixtable->
header, CPL_ERROR_NULL_INPUT);
1138 CPL_ERROR_INCOMPATIBLE_INPUT);
1140 const char *type1 = cpl_propertylist_get_string(aPixtable->
header,
"CTYPE1"),
1141 *type2 = cpl_propertylist_get_string(aPixtable->
header,
"CTYPE2");
1142 cpl_ensure_code(type1 && type2 && !strncmp(type1,
"RA---TAN", 9) &&
1143 !strncmp(type2,
"DEC--TAN", 9),
1144 CPL_ERROR_UNSUPPORTED_MODE);
1146 cpl_msg_info(__func__,
"Adapting WCS to RA/DEC=%f/%f deg", aRA, aDEC);
1149 cpl_table_set_column_unit(aPixtable->
table, MUSE_PIXTABLE_XPOS,
"");
1150 cpl_table_set_column_unit(aPixtable->
table, MUSE_PIXTABLE_YPOS,
"");
1154 float *xpos = cpl_table_get_data_float(aPixtable->
table, MUSE_PIXTABLE_XPOS),
1155 *ypos = cpl_table_get_data_float(aPixtable->
table, MUSE_PIXTABLE_YPOS);
1156 double dp = aDEC / CPL_MATH_DEG_RAD;
1158 #pragma omp parallel for default(none) \
1159 shared(aDEC, dp, nrow, xpos, ypos)
1160 for (n = 0; n < nrow; n++) {
1166 double phi = xpos[n],
1167 theta = ypos[n] + CPL_MATH_PI_2,
1168 ra = atan2(cos(theta) * sin(phi),
1169 sin(theta) * cos(dp) + cos(theta) * sin(dp) * cos(phi))
1171 dec = asin(sin(theta) * sin(dp) - cos(theta) * cos(dp) * cos(phi))
1179 ypos[n] = dec - aDEC;
1184 cpl_table_set_column_unit(aPixtable->
table, MUSE_PIXTABLE_XPOS,
"deg");
1185 cpl_table_set_column_unit(aPixtable->
table, MUSE_PIXTABLE_YPOS,
"deg");
1187 cpl_propertylist_update_double(aPixtable->
header,
"CRVAL1", aRA);
1188 cpl_propertylist_update_double(aPixtable->
header,
"CRVAL2", aDEC);
1196 MUSE_HDR_PT_WCS_COMMENT_POSI);
1197 return CPL_ERROR_NONE;
1226 double *aRA,
double *aDEC)
1228 cpl_ensure_code(aHeader && aRA && aDEC, CPL_ERROR_NULL_INPUT);
1229 const char *type1 = cpl_propertylist_get_string(aHeader,
"CTYPE1"),
1230 *type2 = cpl_propertylist_get_string(aHeader,
"CTYPE2");
1231 cpl_ensure_code(type1 && type2 && !strncmp(type1,
"RA---TAN", 9) &&
1232 !strncmp(type2,
"DEC--TAN", 9),
1233 CPL_ERROR_UNSUPPORTED_MODE);
1239 return CPL_ERROR_NONE;
1270 double *aX,
double *aY)
1272 cpl_ensure_code(aHeader && aX && aY, CPL_ERROR_NULL_INPUT);
1274 const char *type1 = cpl_propertylist_get_string(aHeader,
"CTYPE1"),
1275 *type2 = cpl_propertylist_get_string(aHeader,
"CTYPE2");
1276 cpl_ensure_code(type1 && type2 && !strncmp(type1,
"RA---TAN", 9) &&
1277 !strncmp(type2,
"DEC--TAN", 9),
1278 CPL_ERROR_UNSUPPORTED_MODE);
1281 if (wcs->
cddet == 0.) {
1283 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT);
1285 return CPL_ERROR_ILLEGAL_INPUT;
1287 wcs->crval1 /= CPL_MATH_DEG_RAD;
1288 wcs->
crval2 /= CPL_MATH_DEG_RAD;
1290 aDEC / CPL_MATH_DEG_RAD, aX, aY);
1293 return CPL_ERROR_NONE;
1325 double aDEC,
double *aX,
double *aY)
1327 cpl_ensure_code(aHeader && aX && aY, CPL_ERROR_NULL_INPUT);
1329 const char *type1 = cpl_propertylist_get_string(aHeader,
"CTYPE1"),
1330 *type2 = cpl_propertylist_get_string(aHeader,
"CTYPE2");
1331 cpl_ensure_code(type1 && type2 && !strncmp(type1,
"RA---TAN", 9) &&
1332 !strncmp(type2,
"DEC--TAN", 9),
1333 CPL_ERROR_UNSUPPORTED_MODE);
1336 double a = aRA / CPL_MATH_DEG_RAD,
1337 d = aDEC / CPL_MATH_DEG_RAD,
1339 ap = cpl_propertylist_get_double(aHeader,
"CRVAL1") / CPL_MATH_DEG_RAD,
1340 dp = cpl_propertylist_get_double(aHeader,
"CRVAL2") / CPL_MATH_DEG_RAD,
1341 phi = atan2(-cos(d) * sin(a - ap),
1342 sin(d) * cos(dp) - cos(d) * sin(dp) * cos(a-ap))
1343 + 180 / CPL_MATH_DEG_RAD,
1344 theta = asin(sin(d) * sin(dp) + cos(d) * cos(dp) * cos(a-ap)),
1345 R_theta = CPL_MATH_DEG_RAD / tan(theta);
1347 *aX = R_theta * sin(phi),
1348 *aY = -R_theta * cos(phi);
1350 return CPL_ERROR_NONE;
1373 double *aXOut,
double *aYOut)
1375 cpl_ensure_code(aHeader && aXOut && aYOut, CPL_ERROR_NULL_INPUT);
1381 return CPL_ERROR_NONE;
1406 double *aXOut,
double *aYOut)
1408 cpl_ensure_code(aHeader && aXOut && aYOut, CPL_ERROR_NULL_INPUT);
1411 if (wcs->
cddet == 0.) {
1412 *aXOut = *aYOut = NAN;
1413 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT);
1415 return CPL_ERROR_ILLEGAL_INPUT;
1420 return CPL_ERROR_NONE;
1446 cpl_ensure_code(aHeader && aXAngle && aYAngle, CPL_ERROR_NULL_INPUT);
1448 cpl_errorstate prestate = cpl_errorstate_get();
1449 double cd11 = cpl_propertylist_get_double(aHeader,
"CD1_1"),
1450 cd22 = cpl_propertylist_get_double(aHeader,
"CD2_2"),
1451 cd12 = cpl_propertylist_get_double(aHeader,
"CD1_2"),
1452 cd21 = cpl_propertylist_get_double(aHeader,
"CD2_1"),
1453 det = cd11 * cd22 - cd12 * cd21;
1454 cpl_ensure_code(cpl_errorstate_is_equal(prestate), cpl_error_get_code());
1459 if (cd12 == 0. && cd21 == 0.) {
1462 return CPL_ERROR_NONE;
1465 *aXAngle = atan2(cd12, cd11) * CPL_MATH_DEG_RAD;
1466 *aYAngle = atan2(-cd21, cd22) * CPL_MATH_DEG_RAD;
1467 return CPL_ERROR_NONE;
1493 cpl_ensure_code(aHeader && aXScale && aYScale, CPL_ERROR_NULL_INPUT);
1495 cpl_errorstate prestate = cpl_errorstate_get();
1496 double cd11 = cpl_propertylist_get_double(aHeader,
"CD1_1"),
1497 cd22 = cpl_propertylist_get_double(aHeader,
"CD2_2"),
1498 cd12 = cpl_propertylist_get_double(aHeader,
"CD1_2"),
1499 cd21 = cpl_propertylist_get_double(aHeader,
"CD2_1"),
1500 det = cd11 * cd22 - cd12 * cd21;
1501 cpl_ensure_code(cpl_errorstate_is_equal(prestate), cpl_error_get_code());
1507 if (cd12 == 0. && cd21 == 0.) {
1510 return CPL_ERROR_NONE;
1512 *aXScale = sqrt(cd11*cd11 + cd12*cd12);
1513 *aYScale = sqrt(cd22*cd22 + cd21*cd21);
1514 return CPL_ERROR_NONE;
1540 wcs->cd11 = wcs->
cd22 = wcs->
cddet = 1.;
1544 cpl_errorstate prestate = cpl_errorstate_get();
1545 wcs->crpix1 = cpl_propertylist_get_double(aHeader,
"CRPIX1");
1546 wcs->
crpix2 = cpl_propertylist_get_double(aHeader,
"CRPIX2");
1547 wcs->crval1 = cpl_propertylist_get_double(aHeader,
"CRVAL1");
1548 wcs->
crval2 = cpl_propertylist_get_double(aHeader,
"CRVAL2");
1549 if (!cpl_errorstate_is_equal(prestate)) {
1552 cpl_errorstate_set(prestate);
1555 prestate = cpl_errorstate_get();
1556 wcs->cd11 = cpl_propertylist_get_double(aHeader,
"CD1_1");
1557 wcs->
cd22 = cpl_propertylist_get_double(aHeader,
"CD2_2");
1558 wcs->cd12 = cpl_propertylist_get_double(aHeader,
"CD1_2");
1559 wcs->cd21 = cpl_propertylist_get_double(aHeader,
"CD2_1");
1560 if (!cpl_errorstate_is_equal(prestate) &&
1561 wcs->cd11 == 0. && wcs->cd12 == 0. && wcs->cd21 == 0. && wcs->
cd22 == 0.) {
1564 wcs->cd11 = wcs->
cd22 = wcs->
cddet = 1.;
1565 cpl_errorstate_set(prestate);
1567 wcs->
cddet = wcs->cd11 * wcs->
cd22 - wcs->cd12 * wcs->cd21;
1568 if (wcs->
cddet == 0.) {
1569 cpl_error_set(__func__, CPL_ERROR_ILLEGAL_INPUT);
1573 if (getenv(
"MUSE_DEBUG_WCS") && atoi(getenv(
"MUSE_DEBUG_WCS")) > 0) {
1574 cpl_msg_debug(__func__,
"wcs: axis1 = { %f %f %e }, axis2 = { %f %f %e }, "
1575 "crossterms = { %e %e }, det = %e",
1576 wcs->crpix1, wcs->crval1, wcs->cd11,
1578 wcs->cd12, wcs->cd21, wcs->
cddet);
Structure definition of a MUSE datacube.
#define MUSE_HDR_PT_PREDAR_YLO
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.
cpl_error_code muse_wcs_pixel_from_projplane(cpl_propertylist *aHeader, double aX, double aY, double *aXOut, double *aYOut)
Convert from projection plane coordinates to pixel coordinates.
void muse_image_delete(muse_image *aImage)
Deallocate memory associated to a muse_image object.
#define MUSE_HDR_PT_XLO
FITS header keyword contains the lower limit of the data in x-direction.
cpl_size muse_pixtable_get_nrow(muse_pixtable *aPixtable)
get the number of rows within the pixel table
A structure containing a spatial two-axis WCS.
double muse_pfits_get_ra(const cpl_propertylist *aHeaders)
find out the right ascension
const muse_cpltable_def muse_wcs_detections_def[]
Definition of the table structure for the astrometric field detections.
cpl_image * data
the data extension
muse_image * muse_combine_median_create(muse_imagelist *aImages)
Median combine a list of input images.
static void muse_wcs_pixel_from_projplane_fast(muse_wcs *aWCS, double aX, double aY, double *aXOut, double *aYOut)
Convert from projection plane coordinates to pixel coordinates.
cpl_image * stat
the statistics extension
cpl_error_code muse_wcs_projplane_from_celestial(cpl_propertylist *aHeader, double aRA, double aDEC, double *aX, double *aY)
Convert from celestial spherical coordinates to projection plane coordinates.
#define MUSE_HDR_PT_WCS_PROJ
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.
WCS object to store data needed while computing the astrometric solution.
Structure definition of MUSE three extension FITS file.
cpl_table * table
The pixel table.
cpl_propertylist * header
the FITS header
cpl_boolean muse_cplarray_has_duplicate(const cpl_array *aArray)
Check, if an array contains duplicate values.
cpl_error_code muse_pixtable_reset_dq(muse_pixtable *aPixtable, unsigned int aDQ)
Reset a given bad pixel status (DQ flag) for all pixels in the table.
cpl_error_code muse_cpltable_check(const cpl_table *aTable, const muse_cpltable_def *aDef)
Check whether the table contains the fields of the definition.
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. ...
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.
muse_resampling_crstats_type crtype
cpl_image * dq
the data quality extension
cpl_table * muse_cpltable_new(const muse_cpltable_def *aDef, cpl_size aLength)
Create an empty table according to the specified definition.
const muse_cpltable_def muse_wcs_reference_def[]
Definition of the table structure for the astrometric reference sources.
muse_wcs * muse_wcs_new(cpl_propertylist *aHeader)
Create a new WCS structure from a given FITS header.
double muse_pfits_get_fwhm_end(const cpl_propertylist *aHeaders)
find out the ambient seeing at end of exposure (in arcsec)
muse_wcs_centroid_type
Type of centroiding algorithm to use.
#define MUSE_HDR_PT_PREDAR_YHI
Structure definition of MUSE pixel table.
muse_pixtable_wcs muse_pixtable_wcs_check(muse_pixtable *aPixtable)
Check the state of the world coordinate system of a pixel table.
cpl_error_code muse_wcs_celestial_from_pixel(cpl_propertylist *aHeader, double aX, double aY, double *aRA, double *aDEC)
Convert from pixel coordinates to celestial spherical coordinates.
#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.
static void muse_wcs_pixel_from_celestial_fast(muse_wcs *aWCS, double aRA, double aDEC, double *aX, double *aY)
Convert from celestial spherical coordinates to pixel coordinates.
static void muse_wcs_celestial_from_pixel_fast(muse_wcs *aWCS, double aX, double aY, double *aRA, double *aDEC)
Convert from pixel coordinates to celestial spherical coordinates.
muse_resampling_params * muse_resampling_params_new(muse_resampling_type aMethod)
Create the resampling parameters structure.
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...
muse_wcs_object * muse_wcs_object_new(void)
Allocate memory for a new muse_wcs_object object.
#define MUSE_HDR_PT_PREDAR_XLO
cpl_error_code muse_wcs_projplane_from_pixel(cpl_propertylist *aHeader, double aX, double aY, double *aXOut, double *aYOut)
Convert from pixel coordinates to projection plane coordinates.
double muse_pfits_get_fwhm_start(const cpl_propertylist *aHeaders)
find out the ambient seeing at start of exposure (in arcsec)
cpl_imagelist * data
the cube containing the actual data values
cpl_error_code muse_wcs_pixel_from_celestial(cpl_propertylist *aHeader, double aRA, double aDEC, double *aX, double *aY)
Convert from celestial spherical coordinates to pixel coordinates.
double muse_pfits_get_dec(const cpl_propertylist *aHeaders)
find out the declination
#define MUSE_HDR_PT_PREDAR_XHI
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.
cpl_imagelist * dq
the optional cube containing the bad pixel status
cpl_error_code muse_pixtable_save(muse_pixtable *aPixtable, const char *aFilename)
Save a MUSE pixel table to a file on disk.
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.
#define MUSE_HDR_PT_WCS_POSI
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.
static void muse_wcs_projplane_from_pixel_fast(muse_wcs *aWCS, double aX, double aY, double *aXOut, double *aYOut)
Convert from pixel coordinates to projection plane coordinates.
#define MUSE_HDR_PT_DAR_NAME
cpl_propertylist * muse_wcs_apply_cd(const cpl_propertylist *aExpHeader, const cpl_propertylist *aCalHeader)
Apply the CDi_j matrix of an astrometric solution to an observation.
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.
Definition of a cpl table structure.
void muse_resampling_params_delete(muse_resampling_params *aParams)
Delete a resampling parameters structure.
cpl_error_code muse_wcs_get_angles(cpl_propertylist *aHeader, double *aXAngle, double *aYAngle)
Compute the rotation angles (in degrees) from the FITS header WCS.
void muse_wcs_object_delete(muse_wcs_object *aWCSObj)
Deallocate memory associated to a muse_wcs_object.
#define MUSE_HDR_PT_DAR_CORR
muse_image * muse_image_new(void)
Allocate memory for a new muse_image object.
cpl_error_code muse_wcs_position_celestial(muse_pixtable *aPixtable, double aRA, double aDEC)
Convert native to celestial spherical coordinates in a pixel table.
#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.
muse_ins_mode muse_pfits_get_mode(const cpl_propertylist *aHeaders)
find out the observation mode
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.
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.
cpl_propertylist * muse_wcs_create_default(void)
Create FITS headers containing a default (relative) WCS.