Changeset 10829 for trunk/psModules/src/astrom/pmAstrometryWCS.c

Timestamp:

Dec 24, 2006, 3:49:50 PM (20 years ago)

Author:

magnier

Message:

some API work to reduce dependencies on the recipes, support for mosaic astrometry

File:

• trunk/psModules/src/astrom/pmAstrometryWCS.c (modified) (8 diffs)

trunk/psModules/src/astrom/pmAstrometryWCS.c

@@ -7 +7 @@
  *  @author EAM, IfA
  *
- *  @version $Revision: 1.11 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2006-12-22 21:23:06 $
+ *  @version $Revision: 1.12 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-12-25 01:49:49 $
  *
  *  Copyright 2006 Institute for Astronomy, University of Hawaii
@@ -32 +32 @@
 
 // interpret header WCS (only handles traditional WCS for the moment)
-// plateScale is nominal physical scale on tangent plane (radians / TPA physical units)
-bool pmAstromReadWCS (pmFPA *fpa, pmChip *chip, const psMetadata *header, double plateScale)
+// pixelScale is microns per pixel
+bool pmAstromReadWCS (pmFPA *fpa, pmChip *chip, const psMetadata *header, double pixelScale)
 {
     pmAstromWCS *wcs = pmAstromWCSfromHeader (header);
@@ -40 +40 @@
     }
 
-    bool status = pmAstromWCStoFPA (fpa, chip, wcs, plateScale);
+    bool status = pmAstromWCStoFPA (fpa, chip, wcs, pixelScale);
 
     psFree (wcs);
@@ -67 +67 @@
 
 // interpret chip header WCS as bilevel chip components
-// plateScale is nominal physical scale on tangent plane (radians / TPA physical units)
 bool pmAstromReadBilevelChip (pmChip *chip, const psMetadata *header)
 {
@@ -82 +81 @@
 
 // convert toFPA / toSky components to traditional WCS
-// plateScale is nominal physical scale on tangent plane (microns / arcsecond)
 bool pmAstromReadBilevelMosaic (pmFPA *fpa, const psMetadata *header)
 {
@@ -359 +357 @@
 
 // interpret header WCS (only handles traditional WCS for the moment)
-// plateScale is nominal physical scale on tangent plane (radians / TPA physical units)
-bool pmAstromWCStoFPA (pmFPA *fpa, pmChip *chip, const pmAstromWCS *wcs, double plateScale)
+// pixelScale is the pixel size in microns per pixel
+bool pmAstromWCStoFPA (pmFPA *fpa, pmChip *chip, const pmAstromWCS *wcs, double pixelScale)
 {
     psPlaneTransform *toFPA;
@@ -367 +365 @@
      * wcs->trans, which will convert X,Y in pixels to L,M in degrees.  we also have the following 
      * elements defined:
-     * type (CTYPE)
+     * type (projection type)
      * crval1,2 (in RA,DEC degrees)
      * crpix1,2 
      * cdelt1,2 (in degrees / pixel)
-     * plateScale (radians / physical TPA units)
+     * pixelScale (microns / pixel)
      * 
      * now we convert wcs->trans to toFPA, which is different from wcs->trans in 3 important ways:
-     * 1) the output is in pixels (not degrees): divide by cdelt1,2 raised to an appropriate power
+     * 1) the output is in microns (not degrees): divide by cdelt1,2
      * 2) X,Y are applied directly, without an applied Xo,Yo offset
      * 3) there is an allowed Lo,Mo term ([0][0] coefficients)
      */
 
-    // convert wcs->trans to a matrix which yields L,M in pixels
-    double cdelt1 = hypot (wcs->trans->x->coeff[1][0], wcs->trans->x->coeff[0][1]);
-    double cdelt2 = hypot (wcs->trans->y->coeff[1][0], wcs->trans->y->coeff[0][1]);
-    for (int i = 0; i <= wcs->trans->x->nX; i++) {
-        for (int j = 0; j <= wcs->trans->x->nX; j++) {
-            wcs->trans->x->coeff[i][j] /= cdelt1;
-            wcs->trans->y->coeff[i][j] /= cdelt2;
-        }
-    }
-
     // create transformation with 0,0 reference pixel
     toFPA = psPlaneTransformSetCenter (NULL, wcs->trans, -wcs->crpix1, -wcs->crpix2);
 
-    // scale from FPA to TPA (microns / pixel)
-    double pdelt1 = cdelt1*PM_RAD_DEG / plateScale;
-    double pdelt2 = cdelt2*PM_RAD_DEG / plateScale;
+    // modify scale of toFPA to yield L,M in microns
+    double cdelt1 = hypot (toFPA->x->coeff[1][0], toFPA->x->coeff[0][1]);
+    double cdelt2 = hypot (toFPA->y->coeff[1][0], toFPA->y->coeff[0][1]);
+    for (int i = 0; i <= toFPA->x->nX; i++) {
+        for (int j = 0; j <= toFPA->x->nX; j++) {
+            toFPA->x->coeff[i][j] *= pixelScale/cdelt1;
+            toFPA->y->coeff[i][j] *= pixelScale/cdelt2;
+        }
+    }
+
+    // scale from FPA to TPA (degrees / micron)
+    double pdelt1 = cdelt1 / pixelScale;
+    double pdelt2 = cdelt2 / pixelScale;
     float rX = 1.0;
     float rY = 1.0;
 
-    // projection from TPA to SKY
-    psProjection *toSky = psProjectionAlloc (wcs->toSky->R, wcs->toSky->D, plateScale, plateScale, wcs->toSky->type);
+    // projection from TPA ("linear" degrees) to SKY (radians)
+    psProjection *toSky = psProjectionAlloc (wcs->toSky->R, wcs->toSky->D, PM_RAD_DEG, PM_RAD_DEG, wcs->toSky->type);
 
     if (fpa->toSky == NULL) {
@@ -558 +556 @@
 
     // cdelt1,2 convert from pixels->degrees
-    double cdelt1 = fpa->toTPA->x->coeff[1][0]*fpa->toSky->Xs*PM_DEG_RAD;
-    double cdelt2 = fpa->toTPA->y->coeff[0][1]*fpa->toSky->Ys*PM_DEG_RAD;
+    double cdelt1 = fpa->toTPA->x->coeff[1][0];
+    double cdelt2 = fpa->toTPA->y->coeff[0][1];
     wcs->cdelt1 = cdelt1;
     wcs->cdelt2 = cdelt2;