IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Nov 23, 2005, 1:54:30 PM (21 years ago)
Author:
gusciora
Message:

I added tests for, and a workaround from a possible bug with the psProject()
psDeproject() code, so that the various psAstrometry transforms to and from
the sky work.

File:
1 edited

Legend:

Unmodified
Added
Removed

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

    r5553 r5587  
    88*  @author GLG, MHPCC
    99*
    10 *  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
    11 *  @date $Date: 2005-11-19 01:11:03 $
     10*  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
     11*  @date $Date: 2005-11-23 23:54:30 $
    1212*
    1313*  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    1414*/
     15
     16
     17
     18
     19
     20
     21
     22
     23
     24
     25
     26
     27
     28
     29
     30
     31
     32
     33
     34
     35
     36
     37
     38
     39
     40
     41
     42
     43
     44
     45
    1546
    1647/******************************************************************************/
    … …  
    2253#include "psDB.h"
    2354#include "pmAstrometry.h"
    24 
    25 void PS_PRINT_POLY2D(psPolynomial2D *poly)
    26 {
    27     for (psS32 x = 0 ; x < poly->nX+1 ; x++) {
    28         for (psS32 y = 0 ; y < poly->nY+1 ; y++) {
    29             printf("    (%.2f x^%d y^%d)\n", poly->coeff[x][y], x, y);
    30         }
    31     }
    32 }
    33 
    34 void PS_PRINT_PLANE_TRANSFORM(psPlaneTransform *pt)
    35 {
    36     printf("---------------------- Plane Transform ----------------------\n");
    37     printf("x:\n");
    38     PS_PRINT_POLY2D(pt->x);
    39     printf("y:\n");
    40     PS_PRINT_POLY2D(pt->y);
    41 }
    4255
    4356/*****************************************************************************
    … …  
    7083        psFree(readout->mask);
    7184        psFree(readout->weight);
     85        //
     86        // XXX: Not sure if this is the right way to do things.  Currently the psListAdd()
     87        // increase the memory reference counter to the list data.  So, we
     88        // iterate through the list, and decrement the reference counters.
     89        //
     90        if (1) {
     91            if ((readout->bias != NULL) && (readout->bias->head != NULL)) {
     92                psListElem *tmpElem = (psListElem *) readout->bias->head;
     93                while (NULL != tmpElem) {
     94                    psMemDecrRefCounter((psImage *) tmpElem->data);
     95                    tmpElem = tmpElem->next;
     96                }
     97            }
     98        }
    7299        psFree(readout->bias);
    73100        psFree(readout->analysis);
    … …  
    84111        psFree(cell->toSky);
    85112        psFree(cell->concepts);
    86         psFree(cell->camera);
    87113        psFree(cell->analysis);
    88114        //
    … …  
    195221    }
    196222    tmpCell->camera = cameradata;
    197     psMemIncrRefCounter((psMetadata *) cameradata);
    198223    tmpCell->analysis = NULL;
    199224    tmpCell->readouts = psArrayAlloc(0);
    … …  
    247272    tmpFPA->analysis = NULL;
    248273    tmpFPA->camera = camera;
    249     psMemIncrRefCounter((psPtr) camera);
    250274    tmpFPA->chips = psArrayAlloc(0);
    251275    tmpFPA->private = NULL;
    … …  
    338362
    339363    // Convert to those chip coordinates.
    340     pmCoordFPAToChip(&chipCoord, fpaCoord, tmpChip);
     364    psPlane *rc = pmCoordFPAToChip(&chipCoord, fpaCoord, tmpChip);
     365    if (rc == NULL) {
     366        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine Chip coords.\n");
     367        return(NULL);
     368    }
    341369
    342370    // Determine which cell contains those chip coordinates.
    343371    outCell = pmCellInChip(&chipCoord, tmpChip);
     372    if (outCell == NULL) {
     373        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the cell.\n");
     374        return(NULL);
     375    }
    344376
    345377    return (outCell);
    … …  
    379411    }
    380412
    381     // XXX: Print warning here?
     413    psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the chip.\n");
    382414    return (NULL);
    383415}
    … …  
    414446            chipToCell = p_psPlaneTransformLinearInvert(tmpCell->toChip);
    415447        } else {
    416             // XXX: Generate warning message.
    417             // XXX: Can we use the following function to derive a transform?
    418             // chipToCell = psPlaneTransformInvert(NULL, tmpCell->toChip, NULL, -1);
    419         }
    420 
    421         PS_ASSERT_PTR_NON_NULL(chipToCell, NULL);
     448            psLogMsg(__func__, PS_LOG_WARN, "WARNING: non-linear cell->chip transforms are not yet implemented.\n");
     449            //chipToCell = psPlaneTransformInvert(NULL, tmpCell->toChip, NULL, -1);
     450            chipToCell = NULL;
     451        }
     452        if (chipToCell == NULL) {
     453            psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not invert the Cell->toChip transform.\n");
     454            return(NULL);
     455        }
    422456        psArray* readouts = tmpCell->readouts;
    423457
    … …  
    442476    }
    443477
     478    psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the cell.\n");
    444479    return (NULL);
    445480}
    446481
     482
     483psPlane* pmCoordCellToFPA(
     484    psPlane* fpaCoord,
     485    const psPlane* cellCoord,
     486    const pmCell* cell)
     487{
     488    PS_ASSERT_PTR_NON_NULL(cellCoord, NULL);
     489    PS_ASSERT_PTR_NON_NULL(cell, NULL);
     490
     491    psPlane *rc = psPlaneTransformApply(fpaCoord, cell->toFPA, cellCoord);
     492    if (rc == NULL) {
     493        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform cell coords to FPA coords.\n");
     494    }
     495    return(rc);
     496}
     497
     498
     499psPlane* pmCoordChipToFPA(
     500    psPlane* outCoord,
     501    const psPlane* inCoord,
     502    const pmChip* chip)
     503{
     504    PS_ASSERT_PTR_NON_NULL(inCoord, NULL);
     505    PS_ASSERT_PTR_NON_NULL(chip, NULL);
     506
     507    psPlane *rc = psPlaneTransformApply(outCoord, chip->toFPA, inCoord);
     508    if (rc == NULL) {
     509        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform chip coords to FPA coords.\n");
     510    }
     511    return(rc);
     512}
     513
     514
     515psPlane* pmCoordFPAToChip(
     516    psPlane* chipCoord,
     517    const psPlane* fpaCoord,
     518    const pmChip* chip)
     519{
     520    PS_ASSERT_PTR_NON_NULL(fpaCoord, NULL);
     521    PS_ASSERT_PTR_NON_NULL(chip, NULL);
     522    PS_ASSERT_PTR_NON_NULL(chip->fromFPA, NULL);
     523
     524    psPlane *rc = psPlaneTransformApply(chipCoord, chip->fromFPA, fpaCoord);
     525    if (rc == NULL) {
     526        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform FPA coords to Chip coords.\n");
     527    }
     528    return(rc);
     529}
    447530
    448531psPlane* pmCoordCellToChip(
    … …  
    454537    PS_ASSERT_PTR_NON_NULL(cell, NULL);
    455538
    456     return (psPlaneTransformApply(outCoord, cell->toChip, inCoord));
    457 }
    458 
    459 
    460 psPlane* pmCoordChipToFPA(
    461     psPlane* outCoord,
    462     const psPlane* inCoord,
    463     const pmChip* chip)
    464 {
    465     PS_ASSERT_PTR_NON_NULL(inCoord, NULL);
    466     PS_ASSERT_PTR_NON_NULL(chip, NULL);
    467 
    468     return (psPlaneTransformApply(outCoord, chip->toFPA, inCoord));
    469 }
    470 
     539    psPlane *rc = psPlaneTransformApply(outCoord, cell->toChip, inCoord);
     540    if (rc == NULL) {
     541        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform Cell coords to Chip coords.\n");
     542    }
     543    return(rc);
     544}
     545
     546psPlane* pmCoordChipToCell(
     547    psPlane* cellCoord,
     548    const psPlane* chipCoord,
     549    const pmCell* cell)
     550{
     551    PS_ASSERT_PTR_NON_NULL(chipCoord, NULL);
     552    PS_ASSERT_PTR_NON_NULL(cell, NULL);
     553    PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
     554
     555    pmCell *tmpCell = pmCellInChip(chipCoord, cell->parent);
     556    if (tmpCell == NULL) {
     557        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the proper cell.\n");
     558        return(NULL);
     559    }
     560
     561    psPlaneTransform *tmpChipToCell = NULL;
     562    PS_ASSERT_PTR_NON_NULL(tmpCell->toChip, NULL);
     563    if (true ==  p_psIsProjectionLinear(tmpCell->toChip)) {
     564        tmpChipToCell = p_psPlaneTransformLinearInvert(tmpCell->toChip);
     565    } else {
     566        psLogMsg(__func__, PS_LOG_WARN, "WARNING: non-linear cell->chip transforms are not yet implemented.\n");
     567        // XXX: tmpChipToCell = psPlaneTransformInvert(NULL, tmpCell->toChip, NULL, -1);
     568        tmpChipToCell = NULL;
     569    }
     570    if (tmpChipToCell == NULL) {
     571        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not invert the Cell->toChip transform.\n");
     572        return(NULL);
     573    }
     574
     575    psPlane *rc = psPlaneTransformApply(cellCoord, tmpChipToCell, chipCoord);
     576    if (rc == NULL) {
     577        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform Chip coords to Cell coords.\n");
     578    }
     579    psFree(tmpChipToCell);
     580    return(rc);
     581}
    471582
    472583psPlane* pmCoordFPAToTP(
    … …  
    480591    PS_ASSERT_PTR_NON_NULL(fpa, NULL);
    481592
    482     return(psPlaneDistortApply(outCoord, fpa->toTangentPlane, inCoord, color, magnitude));
     593    psPlane *rc = psPlaneDistortApply(outCoord, fpa->toTangentPlane, inCoord, color, magnitude);
     594    if (rc == NULL) {
     595        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform FPA coords to tangent plane coords.\n");
     596    }
     597    return(rc);
     598}
     599
     600psPlane* pmCoordTPToFPA(
     601    psPlane* fpaCoord,
     602    const psPlane* tpCoord,
     603    double color,
     604    double magnitude,
     605    const pmFPA* fpa)
     606{
     607    PS_ASSERT_PTR_NON_NULL(tpCoord, NULL);
     608    PS_ASSERT_PTR_NON_NULL(fpa, NULL);
     609    PS_ASSERT_PTR_NON_NULL(fpa->fromTangentPlane, NULL);
     610
     611    psPlane *rc = psPlaneDistortApply(fpaCoord, fpa->fromTangentPlane, tpCoord, color, magnitude);
     612    if (rc == NULL) {
     613        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform tangent plane coords to FPA coords.\n");
     614    }
     615    return(rc);
     616}
     617
     618
     619/*****************************************************************************
     620XXXDeproject(outSphere, coord, projection): This private routine is a wrapper
     621for p_psDeproject().  The reason: p_psDeproject() and p_psProject() combined
     622do not seem to produce the original coordinates when they even though they
     623should.  XXXDeproject() simply negates the ->r and ->d members of the output
     624psSphere if the input ->y is larger than 0.0.  I don't know why it works.
     625 
     626I'm guessing the p_psProject() and p_psDeproject() functions have bugs.
     627 *****************************************************************************/
     628psSphere* XXXDeproject(
     629    psSphere *outSphere,
     630    const psPlane* coord,
     631    const psProjection* projection)
     632{
     633    psSphere *rc = p_psDeproject(outSphere, coord, projection);
     634
     635    if (coord->y >= 0.0) {
     636        rc->d = -rc->d;
     637        rc->r = -rc->r;
     638    }
     639
     640    return(rc);
    483641}
    484642
    … …  
    496654    PS_ASSERT_PTR_NON_NULL(projection, NULL);
    497655
    498     return(p_psDeproject(outSphere, tpCoord, projection));
     656    psSphere *rc = XXXDeproject(outSphere, tpCoord, projection);
     657    if (rc == NULL) {
     658        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform tangent plane coords to sky coords.\n");
     659    }
     660    return(rc);
    499661}
    500662/*****************************************************************************
    … …  
    511673    PS_ASSERT_PTR_NON_NULL(projection, NULL);
    512674
    513     return(p_psProject(tpCoord, in, projection));
    514 }
    515 
    516 
    517 psPlane* pmCoordCellToFPA(
    518     psPlane* fpaCoord,
    519     const psPlane* cellCoord,
    520     const pmCell* cell)
    521 {
    522     PS_ASSERT_PTR_NON_NULL(cellCoord, NULL);
    523     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    524 
    525     return (psPlaneTransformApply(fpaCoord, cell->toFPA, cellCoord));
     675    psPlane *rc = p_psProject(tpCoord, in, projection);
     676    if (rc == NULL) {
     677        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform sky to tangent plane coords.\n");
     678    }
     679    return(rc);
    526680}
    527681
    … …  
    545699    psPlane fpaCoord;
    546700    psPlane tpCoord;
     701    psPlane *rc;
    547702    pmFPA* parFPA = (cell->parent)->parent;
    548703
    549704    // Convert the input cell coordinates to FPA coordinates.
    550     psPlaneTransformApply(&fpaCoord, cell->toFPA, cellCoord);
     705    rc = psPlaneTransformApply(&fpaCoord, cell->toFPA, cellCoord);
     706    if (rc == NULL) {
     707        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform cell coords to FPA coords.\n");
     708        return(NULL);
     709    }
    551710
    552711    // Convert the FPA coordinates to tangent plane Coordinates.
    553     psPlaneDistortApply(&tpCoord, parFPA->toTangentPlane, &fpaCoord, color, magnitude);
     712    rc = psPlaneDistortApply(&tpCoord, parFPA->toTangentPlane, &fpaCoord, color, magnitude);
     713    if (rc == NULL) {
     714        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform FPA coords to tangent plane coords.\n");
     715        return(NULL);
     716    }
    554717
    555718    // Convert the tangent plane Coordinates to sky coordinates.
    556     return(pmCoordTPToSky(NULL, &tpCoord, parFPA->projection));
    557 }
    558 
     719    psSphere *rc2 = pmCoordTPToSky(skyCoord, &tpCoord, parFPA->projection);
     720    if (rc2 == NULL) {
     721        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform cell coords to sky coords.\n");
     722    }
     723
     724    return(rc2);
     725}
     726
     727psPlane* pmCoordSkyToCell(
     728    psPlane* cellCoord,
     729    const psSphere* skyCoord,
     730    float color,
     731    float magnitude,
     732    const pmCell* cell)
     733{
     734    PS_ASSERT_PTR_NON_NULL(skyCoord, NULL);
     735    PS_ASSERT_PTR_NON_NULL(cell, NULL);
     736    PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
     737    PS_ASSERT_PTR_NON_NULL(cell->parent->parent, NULL);
     738    pmChip *parChip = cell->parent;
     739    pmFPA *parFPA = parChip->parent;
     740    psPlane tpCoord;
     741    psPlane fpaCoord;
     742    psPlane chipCoord;
     743    psPlane *rc;
     744
     745    // Convert the skyCoords to tangent plane coords.
     746    rc = pmCoordSkyToTP(&tpCoord, skyCoord, parFPA->projection);
     747    if (rc == NULL) {
     748        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine tangent plane coords.\n");
     749        return(NULL);
     750    }
     751
     752    // Convert the tangent plane coords to FPA coords.
     753    rc = pmCoordTPToFPA(&fpaCoord, &tpCoord, color, magnitude, parFPA);
     754    if (rc == NULL) {
     755        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine FPA coords.\n");
     756        return(NULL);
     757    }
     758
     759    // Convert the FPA coords to chip coords.
     760    rc = pmCoordFPAToChip(&chipCoord, &fpaCoord, parChip);
     761    if (rc == NULL) {
     762        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine chip coords.\n");
     763        return(NULL);
     764    }
     765
     766    // Convert the chip coords to cell coords.
     767    rc = pmCoordChipToCell(cellCoord, &chipCoord, cell);
     768    if (rc == NULL) {
     769        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine cell coords.\n");
     770        return(NULL);
     771    }
     772
     773    return (cellCoord);
     774}
    559775
    560776// XXX: This has not been tested.
    … …  
    569785    PS_ASSERT_PTR_NON_NULL(cell->toSky, NULL);
    570786    psPlane outPlane;
    571     psPlaneTransformApply(&outPlane, cell->toSky, cellCoord);
     787    psPlane *rc;
     788    rc = psPlaneTransformApply(&outPlane, cell->toSky, cellCoord);
     789    if (rc == NULL) {
     790        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform cell coords to sky coords.\n");
     791        return(NULL);
     792    }
     793
    572794    psSphere *out = outSphere;
    573795    if (out == NULL) {
    … …  
    577799    out->d = outPlane.x;
    578800
    579     return(NULL);
     801    return(out);
    580802}
    581803
    … …  
    599821    skyPlane.x = skyCoord->d;
    600822
    601     return(psPlaneTransformApply(cellCoord, cell->toSky, &skyPlane));
    602 }
    603 
    604 
    605 psPlane* pmCoordTPToFPA(
    606     psPlane* fpaCoord,
    607     const psPlane* tpCoord,
    608     double color,
    609     double magnitude,
    610     const pmFPA* fpa)
    611 {
    612     PS_ASSERT_PTR_NON_NULL(tpCoord, NULL);
    613     PS_ASSERT_PTR_NON_NULL(fpa, NULL);
    614     PS_ASSERT_PTR_NON_NULL(fpa->fromTangentPlane, NULL);
    615 
    616     return (psPlaneDistortApply(fpaCoord, fpa->fromTangentPlane, tpCoord, color, magnitude));
    617 }
    618 
    619 
    620 psPlane* pmCoordFPAToChip(
    621     psPlane* chipCoord,
    622     const psPlane* fpaCoord,
    623     const pmChip* chip)
    624 {
    625     PS_ASSERT_PTR_NON_NULL(fpaCoord, NULL);
    626     PS_ASSERT_PTR_NON_NULL(chip, NULL);
    627     PS_ASSERT_PTR_NON_NULL(chip->fromFPA, NULL);
    628 
    629     chipCoord = psPlaneTransformApply(chipCoord, chip->fromFPA, fpaCoord);
    630     return(chipCoord);
    631 }
    632 
    633 psPlane* pmCoordChipToCell(
    634     psPlane* cellCoord,
    635     const psPlane* chipCoord,
    636     const pmCell* cell)
    637 {
    638     PS_ASSERT_PTR_NON_NULL(chipCoord, NULL);
    639     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    640     PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
    641 
    642     pmCell *tmpCell = pmCellInChip(chipCoord, cell->parent);
    643     PS_ASSERT_PTR_NON_NULL(tmpCell->toChip, NULL);
    644     psPlaneTransform *tmpChipToCell = p_psPlaneTransformLinearInvert(tmpCell->toChip);
    645     PS_ASSERT_PTR_NON_NULL(tmpChipToCell, NULL);
    646     cellCoord = psPlaneTransformApply(cellCoord, tmpChipToCell, chipCoord);
    647     psFree(tmpChipToCell);
     823    psPlane *rc = psPlaneTransformApply(cellCoord, cell->toSky, &skyPlane);
     824    if (rc == NULL) {
     825        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform sky to cell coords.\n");
     826    }
    648827    return(cellCoord);
    649828}
    650 
    651 
    652 psPlane* pmCoordSkyToCell(
    653     psPlane* cellCoord,
    654     const psSphere* skyCoord,
    655     float color,
    656     float magnitude,
    657     const pmCell* cell)
    658 {
    659     PS_ASSERT_PTR_NON_NULL(skyCoord, NULL);
    660     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    661     PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
    662     PS_ASSERT_PTR_NON_NULL(cell->parent->parent, NULL);
    663     pmChip *parChip = cell->parent;
    664     pmFPA *parFPA = parChip->parent;
    665 
    666     // Convert the skyCoords to tangent plane coords.
    667     psPlane *tpCoord = pmCoordSkyToTP(NULL, skyCoord, parFPA->projection);
    668 
    669     // Convert the tangent plane coords to FPA coords.
    670     psPlane *fpaCoord = pmCoordTPToFPA(NULL, tpCoord, color, magnitude, parFPA);
    671 
    672     // Convert the FPA coords to chip coords.
    673     psPlane *chipCoord = pmCoordFPAToChip(NULL, fpaCoord, parChip);
    674 
    675     // Convert the chip coords to cell coords.
    676     cellCoord = pmCoordChipToCell(cellCoord, chipCoord, cell);
    677 
    678     psFree(tpCoord);
    679     psFree(fpaCoord);
    680     psFree(chipCoord);
    681 
    682     return (cellCoord);
    683 }
    684 
Note: See TracChangeset for help on using the changeset viewer.