IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 10873


Ignore:
Timestamp:
Jan 1, 2007, 11:05:46 AM (19 years ago)
Author:
magnier
Message:

changed def of TPA to same units as FPA (microns)

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/psModules/src/astrom/pmAstrometryWCS.c

    r10859 r10873  
    77 *  @author EAM, IfA
    88 *
    9  *  @version $Revision: 1.14 $ $Name: not supported by cvs2svn $
    10  *  @date $Date: 2006-12-30 03:27:19 $
     9 *  @version $Revision: 1.15 $ $Name: not supported by cvs2svn $
     10 *  @date $Date: 2007-01-01 21:05:46 $
    1111 *
    1212 *  Copyright 2006 Institute for Astronomy, University of Hawaii
    … …  
    4848// convert toFPA / toSky components to pmAstromWCS
    4949// tolerance is convergence for inversion of non-linear terms in pixels
    50 bool pmAstromWriteWCS (psMetadata *header, const pmFPA *fpa, const pmChip *chip, float tol)
     50bool pmAstromWriteWCS (psMetadata *header, const pmFPA *fpa, const pmChip *chip, double tol)
    5151{
    5252    // XXX require chip->toFPA->x->nX == chip->toFPA->x->nY
    … …  
    9595
    9696// convert chip->toFPA components to bilevel WCS
    97 bool pmAstromWriteBilevelChip (psMetadata *header, const pmChip *chip, float tol)
     97bool pmAstromWriteBilevelChip (psMetadata *header, const pmChip *chip, double tol)
    9898{
    9999    pmAstromWCS *wcs = pmAstromWCSBilevelChipFromFPA (chip, tol);
    … …  
    107107
    108108// convert fpa->toTPA, fpa->toSky components to bilevel WCS
    109 bool pmAstromWriteBilevelMosaic (psMetadata *header, const pmFPA *fpa, float tol)
     109bool pmAstromWriteBilevelMosaic (psMetadata *header, const pmFPA *fpa, double tol)
    110110{
    111111    pmAstromWCS *wcs = pmAstromWCSBilevelMosaicFromFPA (fpa, tol);
    … …  
    198198    // XXX add check for CROTA2
    199199    int fitOrder = psMetadataLookupS32 (&isPoly, header, "NPLYTERM");
    200     psMetadataLookupF32 (&pcKeys, header, "PC001001");
    201     psMetadataLookupF32 (&cdKeys, header, "CD1_1");
     200    psMetadataLookupF64 (&pcKeys, header, "PC001001");
     201    psMetadataLookupF64 (&cdKeys, header, "CD1_1");
    202202
    203203    if (cdKeys && pcKeys) {
    … …  
    231231    // and then define a transformation from degrees to degrees
    232232
    233     wcs->crval1 = psMetadataLookupF32 (&status, header, "CRVAL1");
    234     wcs->crval2 = psMetadataLookupF32 (&status, header, "CRVAL2");
    235     wcs->crpix1 = psMetadataLookupF32 (&status, header, "CRPIX1");
    236     wcs->crpix2 = psMetadataLookupF32 (&status, header, "CRPIX2");
     233    wcs->crval1 = psMetadataLookupF64 (&status, header, "CRVAL1");
     234    wcs->crval2 = psMetadataLookupF64 (&status, header, "CRVAL2");
     235    wcs->crpix1 = psMetadataLookupF64 (&status, header, "CRPIX1");
     236    wcs->crpix2 = psMetadataLookupF64 (&status, header, "CRPIX2");
    237237    wcs->toSky = psProjectionAlloc (wcs->crval1*PM_RAD_DEG, wcs->crval2*PM_RAD_DEG, PM_RAD_DEG, PM_RAD_DEG, type);
    238238
    239239    // test the CDELTi varient
    240240    if (pcKeys) {
    241         wcs->cdelt1 = psMetadataLookupF32 (&status, header, "CDELT1");
    242         wcs->cdelt2 = psMetadataLookupF32 (&status, header, "CDELT2");
     241        wcs->cdelt1 = psMetadataLookupF64 (&status, header, "CDELT1");
     242        wcs->cdelt2 = psMetadataLookupF64 (&status, header, "CDELT2");
    243243
    244244        // test the CROTAi varient:
    245245        // XXX double check lambda..
    246         double rotate = psMetadataLookupF32 (&status, header, "CROTA2");
     246        double rotate = psMetadataLookupF64 (&status, header, "CROTA2");
    247247        if (status) {
    248248            wcs->trans->x->coeff[1][0] = +wcs->cdelt1 * cos(rotate*PM_RAD_DEG); // == PC1_1
    … …  
    254254
    255255        // test the PC00i00j varient:
    256         wcs->trans->x->coeff[1][0] = wcs->cdelt1 * psMetadataLookupF32 (&status, header, "PC001001"); // == PC1_1
    257         wcs->trans->x->coeff[0][1] = wcs->cdelt2 * psMetadataLookupF32 (&status, header, "PC001002"); // == PC1_2
    258         wcs->trans->y->coeff[1][0] = wcs->cdelt1 * psMetadataLookupF32 (&status, header, "PC002001"); // == PC2_1
    259         wcs->trans->y->coeff[0][1] = wcs->cdelt2 * psMetadataLookupF32 (&status, header, "PC002002"); // == PC2_2
     256        wcs->trans->x->coeff[1][0] = wcs->cdelt1 * psMetadataLookupF64 (&status, header, "PC001001"); // == PC1_1
     257        wcs->trans->x->coeff[0][1] = wcs->cdelt2 * psMetadataLookupF64 (&status, header, "PC001002"); // == PC1_2
     258        wcs->trans->y->coeff[1][0] = wcs->cdelt1 * psMetadataLookupF64 (&status, header, "PC002001"); // == PC2_1
     259        wcs->trans->y->coeff[0][1] = wcs->cdelt2 * psMetadataLookupF64 (&status, header, "PC002002"); // == PC2_2
    260260
    261261        if (isPoly) {
    … …  
    272272                    }
    273273                    sprintf (name, "PCA1X%1dY%1d", i, j);
    274                     wcs->trans->x->coeff[i][j] = pow(wcs->cdelt1, i) * pow(wcs->cdelt2, j) * psMetadataLookupF32 (&status, header, name);
     274                    wcs->trans->x->coeff[i][j] = pow(wcs->cdelt1, i) * pow(wcs->cdelt2, j) * psMetadataLookupF64 (&status, header, name);
    275275                    sprintf (name, "PCA2X%1dY%1d", i, j);
    276                     wcs->trans->y->coeff[i][j] = pow(wcs->cdelt1, i) * pow(wcs->cdelt2, j) * psMetadataLookupF32 (&status, header, name);
     276                    wcs->trans->y->coeff[i][j] = pow(wcs->cdelt1, i) * pow(wcs->cdelt2, j) * psMetadataLookupF64 (&status, header, name);
    277277                }
    278278            }
    … …  
    285285        wcs->cdelt1 = 1.0;
    286286        wcs->cdelt2 = 1.0;
    287         wcs->trans->x->coeff[1][0] = psMetadataLookupF32 (&status, header, "CD1_1"); // == PC1_1
    288         wcs->trans->x->coeff[0][1] = psMetadataLookupF32 (&status, header, "CD1_2"); // == PC1_2
    289         wcs->trans->y->coeff[1][0] = psMetadataLookupF32 (&status, header, "CD2_1"); // == PC2_1
    290         wcs->trans->y->coeff[0][1] = psMetadataLookupF32 (&status, header, "CD2_2"); // == PC2_2
     287        wcs->trans->x->coeff[1][0] = psMetadataLookupF64 (&status, header, "CD1_1"); // == PC1_1
     288        wcs->trans->x->coeff[0][1] = psMetadataLookupF64 (&status, header, "CD1_2"); // == PC1_2
     289        wcs->trans->y->coeff[1][0] = psMetadataLookupF64 (&status, header, "CD2_1"); // == PC2_1
     290        wcs->trans->y->coeff[0][1] = psMetadataLookupF64 (&status, header, "CD2_2"); // == PC2_2
    291291        return wcs;
    292292    }
    … …  
    328328
    329329    // test the PC00i00j varient:
    330     psMetadataAddF32 (header, PS_LIST_TAIL, "PC001001", PS_META_REPLACE, "", wcs->trans->x->coeff[1][0] / cdelt1); // == PC1_1
    331     psMetadataAddF32 (header, PS_LIST_TAIL, "PC001002", PS_META_REPLACE, "", wcs->trans->x->coeff[0][1] / cdelt2); // == PC1_2
    332     psMetadataAddF32 (header, PS_LIST_TAIL, "PC002001", PS_META_REPLACE, "", wcs->trans->y->coeff[1][0] / cdelt1); // == PC2_1
    333     psMetadataAddF32 (header, PS_LIST_TAIL, "PC002002", PS_META_REPLACE, "", wcs->trans->y->coeff[0][1] / cdelt2); // == PC2_2
     330    psMetadataAddF64 (header, PS_LIST_TAIL, "PC001001", PS_META_REPLACE, "", wcs->trans->x->coeff[1][0] / cdelt1); // == PC1_1
     331    psMetadataAddF64 (header, PS_LIST_TAIL, "PC001002", PS_META_REPLACE, "", wcs->trans->x->coeff[0][1] / cdelt2); // == PC1_2
     332    psMetadataAddF64 (header, PS_LIST_TAIL, "PC002001", PS_META_REPLACE, "", wcs->trans->y->coeff[1][0] / cdelt1); // == PC2_1
     333    psMetadataAddF64 (header, PS_LIST_TAIL, "PC002002", PS_META_REPLACE, "", wcs->trans->y->coeff[0][1] / cdelt2); // == PC2_2
    334334
    335335    // Elixir-style polynomial terms
    … …  
    345345                    continue;
    346346                sprintf (name, "PCA1X%1dY%1d", i, j);
    347                 psMetadataAddF32 (header, PS_LIST_TAIL, name, PS_META_REPLACE, "", wcs->trans->x->coeff[i][j] / pow(wcs->cdelt1, i) / pow(wcs->cdelt2, j));
     347                psMetadataAddF64 (header, PS_LIST_TAIL, name, PS_META_REPLACE, "", wcs->trans->x->coeff[i][j] / pow(wcs->cdelt1, i) / pow(wcs->cdelt2, j));
    348348                sprintf (name, "PCA2X%1dY%1d", i, j);
    349                 psMetadataAddF32 (header, PS_LIST_TAIL, name, PS_META_REPLACE, "", wcs->trans->y->coeff[i][j] / pow(wcs->cdelt1, i) / pow(wcs->cdelt2, j));
     349                psMetadataAddF64 (header, PS_LIST_TAIL, name, PS_META_REPLACE, "", wcs->trans->y->coeff[i][j] / pow(wcs->cdelt1, i) / pow(wcs->cdelt2, j));
    350350            }
    351351        }
    … …  
    390390    }
    391391
    392     // scale from FPA to TPA (degrees / micron)
     392    // scale from TPA (microns) to degrees
    393393    double pdelt1 = cdelt1 / pixelScale;
    394394    double pdelt2 = cdelt2 / pixelScale;
    395     float rX = 1.0;
    396     float rY = 1.0;
    397 
    398     // projection from TPA ("linear" degrees) to SKY (radians)
    399     psProjection *toSky = psProjectionAlloc (wcs->toSky->R, wcs->toSky->D, PM_RAD_DEG, PM_RAD_DEG, wcs->toSky->type);
     395
     396    // projection from TPA (linear microns) to SKY (radians)
     397    psProjection *toSky = psProjectionAlloc (wcs->toSky->R, wcs->toSky->D, PM_RAD_DEG*pdelt1, PM_RAD_DEG*pdelt2, wcs->toSky->type);
    400398
    401399    if (fpa->toSky == NULL) {
    402400        fpa->toTPA = psPlaneTransformIdentity (1);
    403401        fpa->fromTPA = psPlaneTransformIdentity (1);
    404         fpa->toTPA->x->coeff[1][0] = pdelt1;
    405         fpa->toTPA->y->coeff[0][1] = pdelt2;
    406         fpa->fromTPA->x->coeff[1][0] = 1.0 / pdelt1;
    407         fpa->fromTPA->y->coeff[0][1] = 1.0 / pdelt2;
    408402        fpa->toSky = toSky;
    409403    } else {
    … …  
    420414        // convert from pixels on this chip to pixels on reference chip
    421415        // rX has units of refpixels / pixel
    422         rX = pdelt1 / fpa->toTPA->x->coeff[1][0];
    423         rY = pdelt2 / fpa->toTPA->y->coeff[0][1];
     416        double rX = toSky->Xs / fpa->toSky->Xs;
     417        double rY = toSky->Ys / fpa->toSky->Ys;
    424418        for (int i = 0; i <= fpa->toTPA->x->nX; i++) {
    425419            for (int j = 0; j <= fpa->toTPA->x->nY; j++) {
    … …  
    527521// convert toFPA / toSky components to pmAstromWCS
    528522// tolerance is allowed error in center solution in pixels
    529 pmAstromWCS *pmAstromWCSfromFPA (const pmFPA *fpa, const pmChip *chip, float tol)
     523pmAstromWCS *pmAstromWCSfromFPA (const pmFPA *fpa, const pmChip *chip, double tol)
    530524{
    531525
    … …  
    549543    psPlane *center = psPlaneTransformGetCenter (chip->toFPA, tol);
    550544
     545    // create wcs transform from toFPA, result converts pixels to microns
    551546    // adjust wcs transform to use center as reference coordinate
    552547    psPlaneTransformSetCenter (wcs->trans, chip->toFPA, center->x, center->y);
    … …  
    557552    psFree (center);
    558553
    559     // cdelt1,2 convert from pixels->degrees
    560     double cdelt1 = fpa->toTPA->x->coeff[1][0];
    561     double cdelt2 = fpa->toTPA->y->coeff[0][1];
    562     wcs->cdelt1 = cdelt1;
    563     wcs->cdelt2 = cdelt2;
    564 
    565     // convert wcs->trans to a matrix which yields L,M in pixels
     554    // pdelt1,2 converts from microns->degrees
     555    double pdelt1 = fpa->toSky->Xs * PM_DEG_RAD;
     556    double pdelt2 = fpa->toSky->Ys * PM_DEG_RAD;
     557
     558    // cdelt1,2 converts from pixels->degrees
     559    wcs->cdelt1 = pdelt1 * hypot (chip->toFPA->x->coeff[1][0], chip->toFPA->x->coeff[0][1]);
     560    wcs->cdelt2 = pdelt2 * hypot (chip->toFPA->y->coeff[1][0], chip->toFPA->y->coeff[0][1]);
     561
     562    // convert wcs->trans to a matrix which yields L,M in degrees
    566563    for (int i = 0; i <= wcs->trans->x->nX; i++) {
    567564        for (int j = 0; j <= wcs->trans->x->nX; j++) {
    568             wcs->trans->x->coeff[i][j] *= cdelt1;
    569             wcs->trans->y->coeff[i][j] *= cdelt2;
     565            wcs->trans->x->coeff[i][j] *= pdelt1;
     566            wcs->trans->y->coeff[i][j] *= pdelt2;
    570567        }
    571568    }
    … …  
    580577
    581578// convert the chip-level toFPA to a wcs polynomial transformation
    582 pmAstromWCS *pmAstromWCSBilevelChipFromFPA (const pmChip *chip, float tol)
     579pmAstromWCS *pmAstromWCSBilevelChipFromFPA (const pmChip *chip, double tol)
    583580{
    584581    // XXX require chip->toFPA->x->nX == chip->toFPA->x->nY
    … …  
    616613
    617614// convert the fpa-level toTPA, toSky to a wcs polynomial transformation
    618 pmAstromWCS *pmAstromWCSBilevelMosaicFromFPA (const pmFPA *fpa, float tol)
     615pmAstromWCS *pmAstromWCSBilevelMosaicFromFPA (const pmFPA *fpa, double tol)
    619616{
    620617    // XXX require fpa->toTPA->x->nX == fpa->toTPA->x->nY
Note: See TracChangeset for help on using the changeset viewer.