Index: trunk/psModules/src/astrom/Makefile.am =================================================================== --- trunk/psModules/src/astrom/Makefile.am (revision 7151) +++ trunk/psModules/src/astrom/Makefile.am (revision 7152) @@ -4,11 +4,9 @@ libpsmoduleastrom_la_LDFLAGS = -release $(PACKAGE_VERSION) libpsmoduleastrom_la_SOURCES = \ - pmAstrometryObjects.c - -# pmFPAAstrometry.c + pmAstrometryObjects.c \ + pmAstrometryDistortion.c psmoduleincludedir = $(includedir) psmoduleinclude_HEADERS = \ - pmAstrometryObjects.h - -# pmFPAAstrometry.h + pmAstrometryObjects.h \ + pmAstrometryDistortion.h Index: trunk/psModules/src/astrom/pmAstrometryDistortion.c =================================================================== --- trunk/psModules/src/astrom/pmAstrometryDistortion.c (revision 7152) +++ trunk/psModules/src/astrom/pmAstrometryDistortion.c (revision 7152) @@ -0,0 +1,189 @@ +/** @file pmAstrometryDistortion.c +* +* @brief This file defines the basic types for measuring the focal-plane distortion. +* +* @ingroup AstroImage +* +* @author EAM, IfA +* +* @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ +* @date $Date: 2006-05-19 15:39:53 $ +* +* Copyright 2006 Institute for Astronomy, University of Hawaii +*/ + +/******************************************************************************/ +/* INCLUDE FILES */ +/******************************************************************************/ + +//#include +//#include +//#include +//#include // for unlink +#include "pslib.h" +#include "pmFPA.h" +#include "pmAstrometryObjects.h" +#include "pmAstrometryDistortion.h" + +static void pmAstromGradientFree (pmAstromGradient *grad) +{ + + if (grad == NULL) + return; + + return; +} + +pmAstromGradient *pmAstromGradientAlloc () +{ + + pmAstromGradient *gradient = psAlloc (sizeof(pmAstromGradient)); + psMemSetDeallocator(gradient, (psFreeFunc) pmAstromGradientFree); + + // XXX init any of the data values? + + return (gradient); +} + +psArray *pmAstromMeasureGradients(psArray *grads, psArray *rawstars, psArray *refstars, psArray *matches, psMetadata *config) +{ + + if (grads == NULL) { + grads = psArrayAlloc (100); + } + + // determine range + // XXX for the moment, I'm hard-wiring this for megacam + // XXX we need to get the chip dimensions from the metadata somewhere + int Nx = 2; + int Ny = 5; + int DX = 2080 / Nx; + int DY = 4600 / Ny; + + psPolynomial2D *local = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 1, 1); + local->mask[1][1] = 1; + + // measure gradient for fractional chip regions + for (int nx = 0; nx < Nx; nx++) { + for (int ny = 0; ny < Ny; ny++) { + int Xmin = nx*DX; + int Xmax = Xmin + DX; + int Ymin = ny*DY; + int Ymax = Ymin + DY; + + psVector *L = psVectorAlloc (100, PS_TYPE_F32); + psVector *M = psVectorAlloc (100, PS_TYPE_F32); + psVector *dP = psVectorAlloc (100, PS_TYPE_F32); + psVector *dQ = psVectorAlloc (100, PS_TYPE_F32); + int Npts = 0; + + // XXX this is a bit inefficient: first sorting by X or Y could speed this up. + // XXX or assigning to a segment in a single pass first + // select the stars within this chip region + for (int i = 0; i < matches->n; i++) { + + pmAstromMatch *match = matches->data[i]; + + pmAstromObj *raw = rawstars->data[match->i1]; + + if (raw->chip->x < Xmin) + continue; + if (raw->chip->x > Xmax) + continue; + if (raw->chip->y < Ymin) + continue; + if (raw->chip->y > Ymax) + continue; + + pmAstromObj *ref = refstars->data[match->i2]; + + L->data.F32[Npts] = ref->FP->x; + M->data.F32[Npts] = ref->FP->y; + + dP->data.F32[Npts] = ref->TP->x - raw->TP->x; + dQ->data.F32[Npts] = ref->TP->y - raw->TP->y; + + psVectorExtend (L, 100, 1); + psVectorExtend (M, 100, 1); + psVectorExtend (dP, 100, 1); + psVectorExtend (dQ, 100, 1); + } + + // fit the collection of positions and offsets with a local 1st order gradient + // 3hi/3lo sigma clipping on the rflux vs metric fit + psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN | PS_STAT_SAMPLE_STDEV); + psVector *mask = psVectorAlloc (Npts, PS_TYPE_MASK); + + pmAstromGradient *grad = pmAstromGradientAlloc (); + + psVectorStats (stats, L, NULL, NULL, 0); + grad->FP.x = stats->sampleMedian; + + psVectorStats (stats, M, NULL, NULL, 0); + grad->FP.y = stats->sampleMedian; + + psVectorClipFitPolynomial2D (local, stats, mask, 0xff, dP, NULL, L, M); + + grad->dTPdL.x = local->coeff[1][0]; + grad->dTPdM.x = local->coeff[0][1]; + + psVectorClipFitPolynomial2D (local, stats, mask, 0xff, dQ, NULL, L, M); + + grad->dTPdL.y = local->coeff[1][0]; + grad->dTPdM.y = local->coeff[0][1]; + + psArrayAdd (grads, 100, grad); + } + } + return grads; +} + +bool pmAstromFitDistortion(pmFPA *fpa, psArray *grads, psMetadata *config) +{ + + // assign the gradient elements to psVectors for fitting + psVector *dPdL = psVectorAlloc (grads->n, PS_TYPE_F32); + psVector *dQdL = psVectorAlloc (grads->n, PS_TYPE_F32); + psVector *dPdM = psVectorAlloc (grads->n, PS_TYPE_F32); + psVector *dQdM = psVectorAlloc (grads->n, PS_TYPE_F32); + psVector *L = psVectorAlloc (grads->n, PS_TYPE_F32); + psVector *M = psVectorAlloc (grads->n, PS_TYPE_F32); + + for (int i = 0; i < grads->n; i++) { + + pmAstromGradient *grad = grads->data[i]; + + dPdL->data.F32[i] = grad->dTPdL.x; + dQdL->data.F32[i] = grad->dTPdL.y; + + dPdM->data.F32[i] = grad->dTPdM.x; + dQdM->data.F32[i] = grad->dTPdM.y; + + L->data.F32[i] = grad->FP.x; + M->data.F32[i] = grad->FP.y; + } + dPdL->n = dQdL->n = dPdM->n = dQdM->n = grads->n; + L->n = M->n = grads->n; + + psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN | PS_STAT_SAMPLE_STDEV); + psVector *mask = psVectorAlloc (grads->n, PS_TYPE_MASK); + + // XXX the order of the gradient fits need to be 1 less than the distortion term + // XXX determine the distortion order(s) from the fpa->toTP structure + + psPolynomial2D *localX = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 1, 1); + psPolynomial2D *localY = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, 1, 1); + // XXX set masks based on fpa->toTP + + psVectorClipFitPolynomial2D (localX, stats, mask, 0xff, dPdL, NULL, L, M); + psVectorClipFitPolynomial2D (localY, stats, mask, 0xff, dPdM, NULL, L, M); + + // XXX update fpa->toTP distortion terms + + psVectorClipFitPolynomial2D (localX, stats, mask, 0xff, dQdL, NULL, L, M); + psVectorClipFitPolynomial2D (localY, stats, mask, 0xff, dQdM, NULL, L, M); + + // XXX free unneeded structures + + return true; +} Index: trunk/psModules/src/astrom/pmAstrometryDistortion.h =================================================================== --- trunk/psModules/src/astrom/pmAstrometryDistortion.h (revision 7152) +++ trunk/psModules/src/astrom/pmAstrometryDistortion.h (revision 7152) @@ -0,0 +1,58 @@ +/** @file pmAstrometryDistortion.h +* +* @brief This file defines the basic types for measuring and fitting the focal-plane distortion. +* +* @ingroup AstroImage +* +* @author EAM, IfA +* +* @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ +* @date $Date: 2006-05-19 15:39:53 $ +* +* Copyright 2006 Institute for Astronomy, University of Hawaii +*/ + +#ifndef PM_ASTROMETRY_DISTORTION_H +#define PM_ASTROMETRY_DISTORTION_H + +/* The following data structure carries the information about the residual + * gradient of source positions in the tangent plane (pmAstromObj.TP) as a + * function of position in the focal plane (pmAstromObj.FP). + */ +typedef struct +{ + psPlane FP; + psPlane dTPdL; + psPlane dTPdM; +} +pmAstromGradient; + +pmAstromGradient *pmAstromGradientAlloc (); + +/* The following function determines the position residual, in the tangent + * plane, as a function of position in the focal plane, for a collection of raw + * measurements and matched reference stars. The configuration data must include + * the bin size over which the gradient is measured (keyword: ASTROM.GRAD.BOX). + * The function returns an array of pmAstromGradient structures, defined below. + */ +psArray *pmAstromMeasureGradients( + psArray *gradients, + psArray *rawstars, + psArray *refstars, + psArray *matches, + psMetadata *config +); + +/* The gradient set measured above can be fitted with a pair of 2D + * polynomials. The resulting fits can then be related back to the implied + * polynomials which represent the distortion. The following function performs + * the fit and applies the result to the distortion transformation of the + * supplied pmFPA structure. The configuration variable supplies the polynomial + * order (keyword: ASTROM.DISTORT.ORDER). + */ +bool pmAstromFitDistortion( + pmFPA *fpa, + psArray *gradients, + psMetadata *config); + +#endif // PM_ASTROMETRY_DISTORTION_H Index: trunk/psModules/src/astrom/pmAstrometryObjects.c =================================================================== --- trunk/psModules/src/astrom/pmAstrometryObjects.c (revision 7151) +++ trunk/psModules/src/astrom/pmAstrometryObjects.c (revision 7152) @@ -8,6 +8,6 @@ * @author EAM, IfA * -* @version $Revision: 1.4 $ $Name: not supported by cvs2svn $ -* @date $Date: 2006-04-30 22:03:58 $ +* @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ +* @date $Date: 2006-05-19 15:39:53 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -17,6 +17,8 @@ /* INCLUDE FILES */ /******************************************************************************/ +#include #include #include +#include // for unlink #include "pslib.h" #include "psVectorSmooth.h" Index: trunk/psModules/src/astrom/pmAstrometryObjects.h =================================================================== --- trunk/psModules/src/astrom/pmAstrometryObjects.h (revision 7151) +++ trunk/psModules/src/astrom/pmAstrometryObjects.h (revision 7152) @@ -8,18 +8,12 @@ * @author EAM, IfA * -* @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ -* @date $Date: 2006-04-30 22:03:58 $ +* @version $Revision: 1.6 $ $Name: not supported by cvs2svn $ +* @date $Date: 2006-05-19 15:39:53 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii */ -#if !defined(PM_ASTROMETRYOBJECTS_H) -#define PM_PM_ASTROMETRYOBJECTS_H_H - -# include -# include // for strcasecmp -# include // for unlink -# include -# include +#ifndef PM_ASTROMETRY_OBJECTS_H +#define PM_ASTROMETRY_OBJECTS_H /* @@ -46,5 +40,5 @@ psPlane *FP; ///< the position in the pmFPA frame psPlane *TP; ///< the position in the tangent plane - psSphere *sky; ///< the position on the Celestial Sphere. + psSphere *sky; ///< the position on the Celestial Sphere. double Mag; ///< object magnitude XXX what filter? double dMag; ///< error on object magnitude @@ -255,61 +249,4 @@ -/* - * - * The following function determines the position residual, in the tangent - * plane, as a function of position in the focal plane, for a collection of raw - * measurements and matched reference stars. The configuration data must include - * the bin size over which the gradient is measured (keyword: ASTROM.GRAD.BOX). - * The function returns an array of pmAstromGradient structures, defined below. - * - * XXX: No prototype code. - * - */ -psArray pmAstromMeasureGradients( - psArray *starlist1, - psArray *starlist2, - psArray *match, - psMetadata *config -); - - - -/* - * - * The following data structure carries the information about the residual - * gradient of source positions in the tangent plane (pmAstromObj.TP) as a - * function of position in the focal plane (pmAstromObj.FP). - * - */ -typedef struct -{ - psPlane FP; - psPlane dTPdL; - psPlane dTPdM; -} -pmAstromGradient; - - -/* - * - * The gradient set measured above can be fitted with a pair of 2D - * polynomials. The resulting fits can then be related back to the implied - * polynomials which represent the distortion. The following function performs - * the fit and applies the result to the distortion transformation of the - * supplied pmFPA structure. The configuration variable supplies the polynomial - * order (keyword: ASTROM.DISTORT.ORDER). - * - * XXX: No prototype code. - * - */ -psArray pmAstromFitDistortion( - pmFPA *fpa, - psArray *gradients, - psMetadata *config); - - - - - /******************************************************************************* The following functions and structs were in the prototype code, but not the @@ -394,4 +331,3 @@ ); -#endif - +#endif // PM_ASTROMETRY_OBJECTS_H Index: trunk/psModules/src/psmodules.h =================================================================== --- trunk/psModules/src/psmodules.h (revision 7151) +++ trunk/psModules/src/psmodules.h (revision 7152) @@ -41,8 +41,9 @@ #include #include +#include // the following headers are from psModule:astrom #include -#include +#include // the following headers are from psModule:detrend