Changeset 13015

Timestamp:

Apr 24, 2007, 4:40:57 PM (19 years ago)

Author:

Paul Price

Message:

Adding WCS, so that psastro can be used.

File:

• trunk/ppSim/src/ppSimLoop.c (modified) (2 diffs)

trunk/ppSim/src/ppSimLoop.c

@@ -283 +283 @@
         M_PI / 3600.0 / 180.0; // Plate scale (radians/pixel)
     float zp = psMetadataLookupF32(NULL, config->arguments, "ZEROPOINT"); // Photometric zero point
+    float ra0 = psMetadataLookupF32(NULL, config->arguments, "RA"); // Boresight RA (radians)
+    float dec0 = psMetadataLookupF32(NULL, config->arguments, "DEC"); // Boresight Dec (radians)
+    float pa = psMetadataLookupF32(NULL, config->arguments, "PA"); // Position angle (radians)
 
     // Add catalogue stars
     if (type == PPSIM_TYPE_OBJECT) {
-        float ra0 = psMetadataLookupF32(NULL, config->arguments, "RA"); // Boresight RA (radians)
-        float dec0 = psMetadataLookupF32(NULL, config->arguments, "DEC"); // Boresight Dec (radians)
-        float pa = psMetadataLookupF32(NULL, config->arguments, "PA"); // Position angle (radians)
-
         // Read catalogue stars using psastro
         psMetadata *astroRecipe = psMetadataLookupPtr(NULL, config->recipes, PSASTRO_RECIPE);
@@ -651 +650 @@
                              "Binning in y", binning);
 
+            // Set up WCS in header
+            pmHDU *hdu = pmHDUGetLowest(fpa, chip, cell);
+            if (!hdu->header) {
+                hdu->header = psMetadataAlloc();
+            }
+            pmAstromWCS *wcs = pmAstromWCSAlloc(1, 1); // WCS structure
+            wcs->toSky = psProjectionAlloc(ra0, dec0, scale * xParityChip, scale * yParityChip, PS_PROJ_TAN);
+            wcs->cdelt1 = scale * PM_DEG_RAD * xParityChip;
+            wcs->cdelt2 = scale * PM_DEG_RAD * yParityChip;
+            wcs->crpix1 = fpa2cell(0.0, x0Cell, xParityCell, binning, x0Chip, xParityChip);
+            wcs->crpix2 = fpa2cell(0.0, y0Cell, yParityCell, binning, y0Chip, yParityChip);
+            wcs->trans->x->coeff[1][0] = cos(pa) * wcs->cdelt1;
+            wcs->trans->x->coeff[0][1] = -sin(pa) * wcs->cdelt1;
+            wcs->trans->y->coeff[1][0] = sin(pa) * wcs->cdelt2;
+            wcs->trans->y->coeff[0][1] = cos(pa) * wcs->cdelt2;
+            // These aren't used by pmAstromWCStoHeader, but set them anyway
+            wcs->trans->x->coeff[0][0] = ra0;
+            wcs->trans->y->coeff[0][0] = dec0;
+            wcs->trans->x->coeff[1][1] = 0.0;
+            wcs->trans->y->coeff[1][1] = 0.0;
+
+            pmAstromWCStoHeader(hdu->header, wcs);
+            psFree(wcs);
+
             if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) {
                 psError(PS_ERR_IO, false, "Unable to write file.");