IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Feb 17, 2006, 7:13:42 AM (20 years ago)
Author:
magnier
Message:

bulk merge of eam_rel9_p0 onto this branch

File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/rel10_ifa/psModules/src/astrom/pmAstrometry.c

    r6205 r6448  
    1010* XXX: We should review the extent of the warning messages on these functions
    1111* when the transformations are not successful.
    12 * 
     12*
    1313* XXX: Should we implement non-linear cell->chip transforms?
    14 * 
    15 *  @version $Revision: 1.12 $ $Name: not supported by cvs2svn $
    16 *  @date $Date: 2006-01-26 21:10:50 $
     14*
     15*  @version $Revision: 1.12.4.1 $ $Name: not supported by cvs2svn $
     16*  @date $Date: 2006-02-17 17:13:41 $
    1717*
    1818*  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    2525#include <math.h>
    2626#include "pslib.h"
     27
    2728#include "pmAstrometry.h"
     29#include "pmConcepts.h"
    2830
    2931/*****************************************************************************
    … …  
    6062        psFree(readout->mask);
    6163        psFree(readout->weight);
    62         //
    63         // XXX: Not sure if this is the right way to do things.  Currently the psListAdd()
    64         // increase the memory reference counter to the list data.  So, we
    65         // iterate through the list, and decrement the reference counters.
    66         //
    67         if (1) {
    68             if ((readout->bias != NULL) && (readout->bias->head != NULL)) {
    69                 psListElem *tmpElem = (psListElem *) readout->bias->head;
    70                 while (NULL != tmpElem) {
    71                     psMemDecrRefCounter((psImage *) tmpElem->data);
    72                     tmpElem = tmpElem->next;
    73                 }
    74             }
    75         }
    76         psFree(readout->bias);
    7764        psFree(readout->analysis);
     65        #if 0
     66
    7867        psFree(readout->parent);
     68        #endif
     69
     70        readout->parent = NULL;
    7971    }
    8072}
    … …  
    8880        psFree(cell->concepts);
    8981        psFree(cell->analysis);
     82        psFree(cell->camera);
    9083        //
    9184        // Set the parent to NULL in all cell->readouts before psFree(cell->readouts)
    9285        // in order to avoid memory reference counter problems.
    9386        //
     87        #if 0
     88
    9489        for (psS32 i = 0 ; i < cell->readouts->n ; i++) {
    9590            pmReadout *tmpReadout = (pmReadout *) cell->readouts->data[i];
    … …  
    9994            }
    10095        }
     96        psFree(cell->parent);
     97        #endif
     98
     99        cell->parent = NULL;
     100
    101101        psFree(cell->readouts);
    102         psFree(cell->parent);
    103102        psFree(cell->hdu);
     103
    104104    }
    105105}
    … …  
    116116        // in order to avoid memory reference counter problems.
    117117        //
     118        #if 0
     119
    118120        for (psS32 i = 0 ; i < chip->cells->n ; i++) {
    119121            pmCell *tmpCell = (pmCell *) chip->cells->data[i];
    … …  
    123125            }
    124126        }
     127        psFree(chip->parent);
     128        #endif
     129
     130        chip->parent = NULL;
    125131        psFree(chip->cells);
    126         psFree(chip->parent);
    127132        psFree(chip->hdu);
    128133    }
    … …  
    143148        // in order to avoid memory reference counter problems.
    144149        //
     150        #if 0
     151
    145152        for (psS32 i = 0 ; i < fpa->chips->n ; i++) {
    146153            pmChip *tmpChip = (pmChip *) fpa->chips->data[i];
    … …  
    150157            }
    151158        }
     159        #endif
    152160        psFree(fpa->chips);
    153161        psFree(fpa->hdu);
    … …  
    156164}
    157165
     166void p_pmHDUFree(p_pmHDU *hdu)
     167{
     168    if (hdu) {
     169        psFree(hdu->extname);
     170        psFree(hdu->header);
     171        psFree(hdu->images);
     172        psFree(hdu->masks);
     173        psFree(hdu->weights);
     174    }
     175}
     176
    158177// XXX: Verify these default values for row0, col0, rowBins, colBins
     178// PAP: These values may disappear in the future in favour of values in parent->concepts?
    159179pmReadout *pmReadoutAlloc(pmCell *cell)
    160180{
    161181    pmReadout *tmpReadout = (pmReadout *) psAlloc(sizeof(pmReadout));
    162182
    163     tmpReadout->col0 = -1;
    164     tmpReadout->row0 = -1;
    165     tmpReadout->colBins = -1;
    166     tmpReadout->rowBins = -1;
     183    tmpReadout->col0 = 0;
     184    tmpReadout->row0 = 0;
     185    tmpReadout->colBins = 0;
     186    tmpReadout->rowBins = 0;
    167187    tmpReadout->image = NULL;
    168188    tmpReadout->mask = NULL;
    169189    tmpReadout->weight = NULL;
    170     tmpReadout->bias = NULL;
    171190    tmpReadout->analysis = psMetadataAlloc();
    172191    tmpReadout->parent = cell;
    … …  
    179198
    180199// XXX: Verify these default values for row0, col0.
     200// PAP: These values may disappear in the future in favour of values in the "concepts"?
    181201pmCell *pmCellAlloc(
    182202    pmChip *chip,
    183     psMetadata *cameradata,
    184     psString name)
     203    psMetadata *cameraData,
     204    const char *name)
    185205{
    186206    pmCell *tmpCell = (pmCell *) psAlloc(sizeof(pmCell));
    187207
    188     tmpCell->col0 = -1;
    189     tmpCell->row0 = -1;
     208    tmpCell->col0 = 0;
     209    tmpCell->row0 = 0;
    190210    tmpCell->toChip = NULL;
    191211    tmpCell->toFPA = NULL;
    … …  
    196216        psLogMsg(__func__, PS_LOG_WARN, "WARNING: Could not add CELL.NAME to metadata.\n");
    197217    }
    198     tmpCell->camera = cameradata;
    199     tmpCell->analysis = NULL;
     218    tmpCell->camera = psMemIncrRefCounter(cameraData);
     219    tmpCell->analysis = psMetadataAlloc();
    200220    tmpCell->readouts = psArrayAlloc(0);
    201221    tmpCell->parent = chip;
    … …  
    203223        chip->cells = psArrayAdd(chip->cells, 1, (psPtr) tmpCell);
    204224    }
    205     tmpCell->valid = false;
     225    tmpCell->process = true;            // All cells are processed by default
     226    tmpCell->exists = false;            // Not yet read in
    206227    tmpCell->hdu = NULL;
     228
     229    pmConceptsBlankCell(tmpCell);
    207230
    208231    psMemSetDeallocator(tmpCell, (psFreeFunc) cellFree);
    … …  
    211234
    212235// XXX: Verify these default values for row0, col0.
     236// PAP: row0, col0 may disappear in the future in favour of storing values in the "concepts".
    213237pmChip *pmChipAlloc(
    214238    pmFPA *fpa,
    215     psString name)
     239    const char *name)
    216240{
    217241    pmChip *tmpChip = (pmChip *) psAlloc(sizeof(pmChip));
    218242
    219     tmpChip->col0 = -1;
    220     tmpChip->row0 = -1;
     243    tmpChip->col0 = 0;
     244    tmpChip->row0 = 0;
    221245    tmpChip->toFPA = NULL;
    222246    tmpChip->fromFPA = NULL;
    … …  
    226250        psLogMsg(__func__, PS_LOG_WARN, "WARNING: Could not add CHIP.NAME to metadata.\n");
    227251    }
    228     tmpChip->analysis = NULL;
     252    tmpChip->analysis = psMetadataAlloc();
    229253    tmpChip->cells = psArrayAlloc(0);
    230254    tmpChip->parent = fpa;
    … …  
    232256        fpa->chips = psArrayAdd(fpa->chips, 1, (psPtr) tmpChip);
    233257    }
    234     tmpChip->valid = false;
     258    tmpChip->process = true;            // Work on all chips, by default
     259    tmpChip->exists = false;            // Not read in yet
    235260    tmpChip->hdu = NULL;
     261
     262    pmConceptsBlankChip(tmpChip);
    236263
    237264    psMemSetDeallocator(tmpChip, (psFreeFunc) chipFree);
    … …  
    248275    tmpFPA->concepts = psMetadataAlloc();
    249276    tmpFPA->analysis = NULL;
    250     tmpFPA->camera = camera;
     277    tmpFPA->camera = psMemIncrRefCounter((psPtr)camera);
    251278    tmpFPA->chips = psArrayAlloc(0);
    252279    tmpFPA->hdu = NULL;
    253280    tmpFPA->phu = NULL;
    254281
     282    pmConceptsBlankFPA(tmpFPA);
     283
    255284    psMemSetDeallocator(tmpFPA, (psFreeFunc) FPAFree);
    256285    return(tmpFPA);
     286}
     287
     288p_pmHDU *p_pmHDUAlloc(const char *extname)
     289{
     290    p_pmHDU *hdu = psAlloc(sizeof(p_pmHDU));
     291    psMemSetDeallocator(hdu, (psFreeFunc)p_pmHDUFree);
     292
     293    hdu->extname = psStringCopy(extname);
     294    hdu->header = NULL;
     295    hdu->images = NULL;
     296    hdu->masks = NULL;
     297    hdu->weights = NULL;
     298
     299    return hdu;
    257300}
    258301
    … …  
    317360
    318361
    319 /*****************************************************************************/
    320 /* FUNCTION IMPLEMENTATION - PUBLIC                                          */
    321 /*****************************************************************************/
    322 
    323 pmCell* pmCellInFPA(
    324     const psPlane* fpaCoord,
    325     const pmFPA* FPA)
    326 {
    327     PS_ASSERT_PTR_NON_NULL(fpaCoord, NULL);
    328     PS_ASSERT_PTR_NON_NULL(FPA, NULL);
    329 
    330     pmChip* tmpChip = NULL;
    331     psPlane chipCoord;
    332     pmCell* outCell = NULL;
    333 
    334     // Determine which chip contains the fpaCoords.
    335     tmpChip = pmChipInFPA(fpaCoord, FPA);
    336     if (tmpChip == NULL) {
    337         return(NULL);
    338     }
    339 
    340     // Convert to those chip coordinates.
    341     psPlane *rc = pmCoordFPAToChip(&chipCoord, fpaCoord, tmpChip);
    342     if (rc == NULL) {
    343         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine Chip coords.\n");
    344         return(NULL);
    345     }
    346 
    347     // Determine which cell contains those chip coordinates.
    348     outCell = pmCellInChip(&chipCoord, tmpChip);
    349     if (outCell == NULL) {
    350         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the cell.\n");
    351         return(NULL);
    352     }
    353 
    354     return (outCell);
    355 }
    356 
    357 pmChip* pmChipInFPA(
    358     const psPlane* fpaCoord,
    359     const pmFPA* FPA)
    360 {
    361     PS_ASSERT_PTR_NON_NULL(fpaCoord, NULL);
    362     PS_ASSERT_PTR_NON_NULL(FPA, NULL);
    363     PS_ASSERT_PTR_NON_NULL(FPA->chips, NULL);
    364 
    365     psArray* chips = FPA->chips;
    366     psS32 nChips = chips->n;
    367     psPlane chipCoord;
    368     pmCell *tmpCell = NULL;
    369 
    370     //
    371     // Loop through every chip in this FPA.  Convert the original FPA
    372     // coordinates to chip coordinates for that chip.  Then, determine if any
    373     // cells in that chip contain those chip coordinates.
    374     // XXX: Depending on the number of chips, and their topology, there may be
    375     // a much more efficient way of doing this.
    376     //
    377     for (psS32 i = 0; i < nChips; i++) {
    378         pmChip* tmpChip = chips->data[i];
    379         PS_ASSERT_PTR_NON_NULL(tmpChip, NULL);
    380         PS_ASSERT_PTR_NON_NULL(tmpChip->fromFPA, NULL);
    381 
    382         psPlaneTransformApply(&chipCoord, tmpChip->fromFPA, fpaCoord);
    383 
    384         tmpCell = pmCellInChip(&chipCoord, tmpChip);
    385         if (tmpCell != NULL) {
    386             return(tmpChip);
    387         }
    388     }
    389 
    390     psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the chip.\n");
    391     return (NULL);
    392 }
    393 
    394 
    395 pmCell* pmCellInChip(
    396     const psPlane* chipCoord,
    397     const pmChip* chip)
    398 {
    399     PS_ASSERT_PTR_NON_NULL(chipCoord, NULL);
    400     PS_ASSERT_PTR_NON_NULL(chip, NULL);
    401 
    402     psPlane cellCoord;
    403     psArray* cells;
    404 
    405     cells = chip->cells;
    406     if (cells == NULL) {
    407         return NULL;
    408     }
    409 
    410     //
    411     // We loop over each cell in the chip.  We transform the chipCoord into
    412     // a cellCoord for that cell and determine if that cellCoord is valid.
    413     // If so, then we return that cell.
    414     // XXX: Depending on the number of cells, and their topology, there may be
    415     // a much more efficient way of doing this.
    416     //
    417     for (psS32 i = 0; i < cells->n; i++) {
    418         pmCell* tmpCell = (pmCell* ) cells->data[i];
    419         PS_ASSERT_PTR_NON_NULL(tmpCell, NULL);
    420 
    421         psPlaneTransform *chipToCell = NULL;
    422         if (true ==  p_psIsProjectionLinear(tmpCell->toChip)) {
    423             chipToCell = p_psPlaneTransformLinearInvert(tmpCell->toChip);
    424         } else {
    425             psLogMsg(__func__, PS_LOG_WARN, "WARNING: non-linear cell->chip transforms are not yet implemented.\n");
    426             //chipToCell = psPlaneTransformInvert(NULL, tmpCell->toChip, NULL, -1);
    427             chipToCell = NULL;
    428         }
    429         if (chipToCell == NULL) {
    430             psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not invert the Cell->toChip transform.\n");
    431             return(NULL);
    432         }
    433         psArray* readouts = tmpCell->readouts;
    434 
    435         if (readouts != NULL) {
    436             for (psS32 j = 0; j < readouts->n; j++) {
    437                 pmReadout* tmpReadout = readouts->data[j];
    438                 PS_ASSERT_READOUT_NON_NULL(tmpReadout, NULL);
    439 
    440                 psPlaneTransformApply(&cellCoord,
    441                                       chipToCell,
    442                                       chipCoord);
    443 
    444                 if (checkValidImageCoords(cellCoord.x,
    445                                           cellCoord.y,
    446                                           tmpReadout->image)) {
    447                     psFree(chipToCell);
    448                     return (tmpCell);
    449                 }
    450             }
    451         }
    452         psFree(chipToCell);
    453     }
    454 
    455     //psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the cell.\n");
    456     return (NULL);
    457 }
    458 
    459 
    460 psPlane* pmCoordCellToFPA(
    461     psPlane* fpaCoord,
    462     const psPlane* cellCoord,
    463     const pmCell* cell)
    464 {
    465     PS_ASSERT_PTR_NON_NULL(cellCoord, NULL);
    466     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    467 
    468     psPlane *rc = psPlaneTransformApply(fpaCoord, cell->toFPA, cellCoord);
    469     if (rc == NULL) {
    470         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform cell coords to FPA coords.\n");
    471     }
    472     return(rc);
    473 }
    474 
    475 
    476 psPlane* pmCoordChipToFPA(
    477     psPlane* outCoord,
    478     const psPlane* inCoord,
    479     const pmChip* chip)
    480 {
    481     PS_ASSERT_PTR_NON_NULL(inCoord, NULL);
    482     PS_ASSERT_PTR_NON_NULL(chip, NULL);
    483 
    484     psPlane *rc = psPlaneTransformApply(outCoord, chip->toFPA, inCoord);
    485     if (rc == NULL) {
    486         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform chip coords to FPA coords.\n");
    487     }
    488     return(rc);
    489 }
    490 
    491 
    492 psPlane* pmCoordFPAToChip(
    493     psPlane* chipCoord,
    494     const psPlane* fpaCoord,
    495     const pmChip* chip)
    496 {
    497     PS_ASSERT_PTR_NON_NULL(fpaCoord, NULL);
    498     PS_ASSERT_PTR_NON_NULL(chip, NULL);
    499     PS_ASSERT_PTR_NON_NULL(chip->fromFPA, NULL);
    500 
    501     psPlane *rc = psPlaneTransformApply(chipCoord, chip->fromFPA, fpaCoord);
    502     if (rc == NULL) {
    503         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform FPA coords to Chip coords.\n");
    504     }
    505     return(rc);
    506 }
    507 
    508 psPlane* pmCoordCellToChip(
    509     psPlane* outCoord,
    510     const psPlane* inCoord,
    511     const pmCell* cell)
    512 {
    513     PS_ASSERT_PTR_NON_NULL(inCoord, NULL);
    514     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    515 
    516     psPlane *rc = psPlaneTransformApply(outCoord, cell->toChip, inCoord);
    517     if (rc == NULL) {
    518         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform Cell coords to Chip coords.\n");
    519     }
    520     return(rc);
    521 }
    522 
    523 psPlane* pmCoordChipToCell(
    524     psPlane* cellCoord,
    525     const psPlane* chipCoord,
    526     const pmCell* cell)
    527 {
    528     PS_ASSERT_PTR_NON_NULL(chipCoord, NULL);
    529     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    530     PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
    531 
    532     pmCell *tmpCell = pmCellInChip(chipCoord, cell->parent);
    533     if (tmpCell == NULL) {
    534         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the proper cell.\n");
    535         return(NULL);
    536     }
    537 
    538     psPlaneTransform *tmpChipToCell = NULL;
    539     PS_ASSERT_PTR_NON_NULL(tmpCell->toChip, NULL);
    540     if (true ==  p_psIsProjectionLinear(tmpCell->toChip)) {
    541         tmpChipToCell = p_psPlaneTransformLinearInvert(tmpCell->toChip);
    542     } else {
    543         psLogMsg(__func__, PS_LOG_WARN, "WARNING: non-linear cell->chip transforms are not yet implemented.\n");
    544         // XXX: tmpChipToCell = psPlaneTransformInvert(NULL, tmpCell->toChip, NULL, -1);
    545         tmpChipToCell = NULL;
    546     }
    547     if (tmpChipToCell == NULL) {
    548         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not invert the Cell->toChip transform.\n");
    549         return(NULL);
    550     }
    551 
    552     psPlane *rc = psPlaneTransformApply(cellCoord, tmpChipToCell, chipCoord);
    553     if (rc == NULL) {
    554         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform Chip coords to Cell coords.\n");
    555     }
    556     psFree(tmpChipToCell);
    557     return(rc);
    558 }
    559 
    560 psPlane* pmCoordFPAToTP(
    561     psPlane* outCoord,
    562     const psPlane* inCoord,
    563     double color,
    564     double magnitude,
    565     const pmFPA* fpa)
    566 {
    567     PS_ASSERT_PTR_NON_NULL(inCoord, NULL);
    568     PS_ASSERT_PTR_NON_NULL(fpa, NULL);
    569 
    570     psPlane *rc = psPlaneDistortApply(outCoord, fpa->toTangentPlane, inCoord, color, magnitude);
    571     if (rc == NULL) {
    572         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform FPA coords to tangent plane coords.\n");
    573     }
    574     return(rc);
    575 }
    576 
    577 psPlane* pmCoordTPToFPA(
    578     psPlane* fpaCoord,
    579     const psPlane* tpCoord,
    580     double color,
    581     double magnitude,
    582     const pmFPA* fpa)
    583 {
    584     PS_ASSERT_PTR_NON_NULL(tpCoord, NULL);
    585     PS_ASSERT_PTR_NON_NULL(fpa, NULL);
    586     PS_ASSERT_PTR_NON_NULL(fpa->fromTangentPlane, NULL);
    587 
    588     psPlane *rc = psPlaneDistortApply(fpaCoord, fpa->fromTangentPlane, tpCoord, color, magnitude);
    589     if (rc == NULL) {
    590         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform tangent plane coords to FPA coords.\n");
    591     }
    592     return(rc);
    593 }
    594 
    595 
    596 /*****************************************************************************
    597 XXXDeproject(outSphere, coord, projection): This private routine is a wrapper
    598 for p_psDeproject().  The reason: p_psDeproject() and p_psProject() combined
    599 do not seem to produce the original coordinates when they even though they
    600 should.  XXXDeproject() simply negates the ->r and ->d members of the output
    601 psSphere if the input ->y is larger than 0.0.  I don't know why it works.
    602  
    603 I'm guessing the p_psProject() and p_psDeproject() functions have bugs.
    604  
    605 XXX: It appears that p_psProject() and p_psDeproject() have been fixed.
    606 Remove this.
    607  *****************************************************************************/
    608 psSphere* XXXDeproject(
    609     psSphere *outSphere,
    610     const psPlane* coord,
    611     const psProjection* projection)
    612 {
    613     psSphere *rc = p_psDeproject(outSphere, coord, projection);
    614 
    615     if (coord->y >= 0.0) {
    616         rc->d = -rc->d;
    617         rc->r = -rc->r;
    618     }
    619 
    620     return(rc);
    621 }
    622 
    623 /*****************************************************************************
    624   *****************************************************************************/
    625 psSphere* pmCoordTPToSky(
    626     psSphere* outSphere,
    627     const psPlane* tpCoord,
    628     const psProjection *projection)
    629 {
    630     PS_ASSERT_PTR_NON_NULL(tpCoord, NULL);
    631     PS_ASSERT_PTR_NON_NULL(projection, NULL);
    632 
    633     //    psSphere *rc = XXXDeproject(outSphere, tpCoord, projection);
    634     psSphere *rc = p_psDeproject(outSphere, tpCoord, projection);
    635     if (rc == NULL) {
    636         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform tangent plane coords to sky coords.\n");
    637     }
    638     return(rc);
    639 }
    640 
    641362/*****************************************************************************
    642363 *****************************************************************************/
    643 psPlane* pmCoordSkyToTP(
    644     psPlane* tpCoord,
    645     const psSphere* in,
    646     const psProjection *projection)
    647 {
    648     PS_ASSERT_PTR_NON_NULL(in, NULL);
    649     PS_ASSERT_PTR_NON_NULL(projection, NULL);
    650 
    651     psPlane *rc = p_psProject(tpCoord, in, projection);
    652     if (rc == NULL) {
    653         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform sky to tangent plane coords.\n");
    654     }
    655     return(rc);
    656 }
    657 
    658 /*****************************************************************************
    659  *****************************************************************************/
    660 psSphere* pmCoordCellToSky(
    661     psSphere* skyCoord,
    662     const psPlane* cellCoord,
    663     double color,
    664     double magnitude,
    665     const pmCell* cell)
    666 {
    667     PS_ASSERT_PTR_NON_NULL(cellCoord, NULL);
    668     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    669     PS_ASSERT_PTR_NON_NULL(cell->toFPA, NULL);
    670     PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
    671     PS_ASSERT_PTR_NON_NULL(cell->parent->parent, NULL);
    672     PS_ASSERT_PTR_NON_NULL(cell->parent->parent->toTangentPlane, NULL);
    673     PS_ASSERT_PTR_NON_NULL(cell->parent->parent->projection, NULL);
    674     psPlane fpaCoord;
    675     psPlane tpCoord;
    676     psPlane *rc;
    677     pmFPA* parFPA = (cell->parent)->parent;
    678 
    679     // Convert the input cell coordinates to FPA coordinates.
    680     rc = psPlaneTransformApply(&fpaCoord, cell->toFPA, cellCoord);
    681     if (rc == NULL) {
    682         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform cell coords to FPA coords.\n");
    683         return(NULL);
    684     }
    685 
    686     // Convert the FPA coordinates to tangent plane Coordinates.
    687     rc = psPlaneDistortApply(&tpCoord, parFPA->toTangentPlane, &fpaCoord, color, magnitude);
    688     if (rc == NULL) {
    689         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform FPA coords to tangent plane coords.\n");
    690         return(NULL);
    691     }
    692 
    693     // Convert the tangent plane Coordinates to sky coordinates.
    694     psSphere *rc2 = pmCoordTPToSky(skyCoord, &tpCoord, parFPA->projection);
    695     if (rc2 == NULL) {
    696         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform cell coords to sky coords.\n");
    697     }
    698 
    699     return(rc2);
    700 }
    701 
    702 /*****************************************************************************
    703  *****************************************************************************/
    704 psPlane* pmCoordSkyToCell(
    705     psPlane* cellCoord,
    706     const psSphere* skyCoord,
    707     float color,
    708     float magnitude,
    709     const pmCell* cell)
    710 {
    711     PS_ASSERT_PTR_NON_NULL(skyCoord, NULL);
    712     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    713     PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
    714     PS_ASSERT_PTR_NON_NULL(cell->parent->parent, NULL);
    715     pmChip *parChip = cell->parent;
    716     pmFPA *parFPA = parChip->parent;
    717     psPlane tpCoord;
    718     psPlane fpaCoord;
    719     psPlane chipCoord;
    720     psPlane *rc;
    721 
    722     // Convert the skyCoords to tangent plane coords.
    723     rc = pmCoordSkyToTP(&tpCoord, skyCoord, parFPA->projection);
    724     if (rc == NULL) {
    725         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine tangent plane coords.\n");
    726         return(NULL);
    727     }
    728 
    729     // Convert the tangent plane coords to FPA coords.
    730     rc = pmCoordTPToFPA(&fpaCoord, &tpCoord, color, magnitude, parFPA);
    731     if (rc == NULL) {
    732         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine FPA coords.\n");
    733         return(NULL);
    734     }
    735 
    736     // Convert the FPA coords to chip coords.
    737     rc = pmCoordFPAToChip(&chipCoord, &fpaCoord, parChip);
    738     if (rc == NULL) {
    739         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine chip coords.\n");
    740         return(NULL);
    741     }
    742 
    743     // Convert the chip coords to cell coords.
    744     rc = pmCoordChipToCell(cellCoord, &chipCoord, cell);
    745     if (rc == NULL) {
    746         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine cell coords.\n");
    747         return(NULL);
    748     }
    749 
    750     return (cellCoord);
    751 }
    752 
    753 /*****************************************************************************
    754  *****************************************************************************/
    755 psSphere* pmCoordCellToSkyQuick(
    756     psSphere* outSphere,
    757     const psPlane* cellCoord,
    758     const pmCell* cell)
    759 {
    760     PS_ASSERT_PTR_NON_NULL(cellCoord, NULL);
    761     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    762     PS_ASSERT_PTR_NON_NULL(cell->toSky, NULL);
    763     psPlane outPlane;
    764     psPlane *rc;
    765     rc = psPlaneTransformApply(&outPlane, cell->toSky, cellCoord);
    766     if (rc == NULL) {
    767         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform cell coords to sky coords.\n");
    768         return(NULL);
    769     }
    770 
    771     psSphere *out = outSphere;
    772     if (out == NULL) {
    773         out = psSphereAlloc();
    774     }
    775     out->r = outPlane.y;
    776     out->d = outPlane.x;
    777 
    778     return(out);
    779 }
    780 
    781 /*****************************************************************************
    782  *****************************************************************************/
    783 psPlane* pmCoordSkyToCellQuick(
    784     psPlane* cellCoord,
    785     const psSphere* skyCoord,
    786     const pmCell* cell)
    787 {
    788     PS_ASSERT_PTR_NON_NULL(skyCoord, NULL);
    789     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    790     PS_ASSERT_PTR_NON_NULL(cell->toSky, NULL);
    791     psPlane skyPlane;
    792     skyPlane.y = skyCoord->r;
    793     skyPlane.x = skyCoord->d;
    794 
    795     psPlane *rc = psPlaneTransformApply(cellCoord, cell->toSky, &skyPlane);
    796     if (rc == NULL) {
    797         psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform sky to cell coords.\n");
    798     }
    799     return(cellCoord);
     364
     365// Set cells within a chip to be processed or not
     366static bool setCellsProcess(const pmChip *chip, // Chip of interest
     367                            bool process  // Process this chip?
     368                           )
     369{
     370    PS_ASSERT_PTR_NON_NULL(chip, false);
     371
     372    psArray *cells = chip->cells;       // Component cells
     373    if (! cells) {
     374        return false;
     375    }
     376    for (int i = 0; i < cells->n; i++) {
     377        pmCell *tmpCell = cells->data[i]; // Cell of interest
     378        if (tmpCell) {
     379            tmpCell->process = process;
     380        }
     381    }
     382
     383    return true;
    800384}
    801385
    … …  
    807391{
    808392    PS_ASSERT_PTR_NON_NULL(fpa, false);
    809     if ((fpa->chips == NULL) || (chipNum >= fpa->chips->n)) {
     393
     394    psArray *chips = fpa->chips;        // Component chips
     395    if ((chips == NULL) || (chipNum >= chips->n)) {
    810396        return(false);
    811397    }
    812     psBool rc = true;
    813 
    814     for (psS32 i = 0 ; i < fpa->chips->n ; i++) {
    815         pmChip *tmpChip = (pmChip *) fpa->chips->data[i];
     398
     399    for (int i = 0 ; i < chips->n ; i++) {
     400        pmChip *tmpChip = (pmChip *) chips->data[i];
    816401        if (tmpChip == NULL) {
    817             rc = false;
     402            continue;
     403        }
     404        if (i == chipNum) {
     405            tmpChip->process = true;
     406            setCellsProcess(tmpChip, true);
    818407        } else {
    819             if (i == chipNum) {
    820                 tmpChip->valid = true;
    821             } else {
    822                 tmpChip->valid = false;
    823             }
    824         }
    825     }
    826 
    827     return(rc);
     408            tmpChip->process = false;
     409            setCellsProcess(tmpChip, false);
     410        }
     411
     412    }
     413
     414    return true;
    828415}
    829416
    … …  
    831418/*****************************************************************************
    832419XXX: The SDRS is ambiguous on a few things:
    833     Whether or not the other chips should be set valid=true.
    834     Should we return the number of chip valid=true before or after they're set,
     420    Whether or not the other chips should be set process=true. [PAP: No]
     421    Should we return the number of chip process=true before or after they're set, [PAP: After]
    835422 *****************************************************************************/
    836423/**
    837  * 
    838  * pmFPAExcludeChip shall set valid to false only for the specified chip
     424 *
     425 * pmFPAExcludeChip shall set process to false only for the specified chip
    839426 * number (chipNum). In the event that the specified chip number does not exist
    840427 * within the fpa, the function shall generate a warning, and perform no action.
    841  * The function shall return the number of chips within the fpa that have valid
     428 * The function shall return the number of chips within the fpa that have process
    842429 * set to true.
    843  * 
     430 *
    844431 */
    845432int pmFPAExcludeChip(
    … …  
    849436    PS_ASSERT_PTR_NON_NULL(fpa, false);
    850437
    851     if (fpa->chips == NULL) {
     438    psArray *chips = fpa->chips;        // Component chips
     439    if (chips == NULL) {
    852440        psLogMsg(__func__, PS_LOG_WARN, "WARNING: fpa->chips == NULL\n");
    853441        return(0);
    854442    }
    855     if ((chipNum >= fpa->chips->n) || (NULL == (pmChip *) fpa->chips->data[chipNum])) {
     443    if ((chipNum >= chips->n) || (NULL == (pmChip *) chips->data[chipNum])) {
    856444        psLogMsg(__func__, PS_LOG_WARN, "WARNING: the specified chip (%d) does not exist.\n", chipNum);
    857445        return(0);
    858446    }
    859447
    860     psS32 numChips = 0;
    861     for (psS32 i = 0 ; i < fpa->chips->n ; i++) {
    862         pmChip *tmpChip = (pmChip *) fpa->chips->data[i];
     448    int numChips = 0;                   // Number of chips to be processed
     449    for (int i = 0 ; i < chips->n ; i++) {
     450        pmChip *tmpChip = (pmChip *) chips->data[i]; // Chip of interest
    863451        if (tmpChip != NULL) {
    864452            if (i == chipNum) {
    865                 tmpChip->valid = false;
    866             } else {
    867                 tmpChip->valid = true;
     453                tmpChip->process = false;
     454                setCellsProcess(tmpChip, false); // Wipe out the cell as well
     455            } else if (tmpChip->process) {
    868456                numChips++;
    869457            }
    … …  
    874462}
    875463
     464
     465bool pmCellSetWeights(pmCell *cell // Cell for which to set weights
     466                     )
     467{
     468    float gain = psMetadataLookupF32(NULL, cell->concepts, "CELL.GAIN"); // Cell gain
     469    float readnoise = psMetadataLookupF32(NULL, cell->concepts, "CELL.READNOISE"); // Cell read noise
     470    psRegion *trimsec = psMetadataLookupPtr(NULL, cell->concepts, "CELL.TRIMSEC"); // Trim section
     471
     472    p_pmHDU *hdu = cell->hdu;           // The data unit, containing the weight and mask originals
     473    if (!hdu) {
     474        pmChip *chip = cell->parent;    // The parent chip
     475        if (chip->hdu) {
     476            hdu = chip->hdu;
     477        } else {
     478            pmFPA *fpa = chip->parent;  // The parent FPA
     479            hdu = fpa->hdu;
     480            if (!hdu) {
     481                psError(PS_ERR_UNKNOWN, true, "Unable to find the HDU in the hierarchy!\n");
     482                return false;
     483            }
     484        }
     485    }
     486
     487    psArray *pixels = hdu->images;      // Array of images
     488    psArray *weights = hdu->weights;    // Array of weight images
     489    psArray *masks = hdu->masks;        // Array of mask images
     490    // Generate the weights and masks if required
     491    if (! weights) {
     492        weights = psArrayAlloc(pixels->n);
     493        for (int i = 0; i < pixels->n; i++) {
     494            psImage *image = pixels->data[i];
     495            weights->data[i] = psImageAlloc(image->numCols, image->numRows, PS_TYPE_F32);
     496            psImageInit(weights->data[i], 0.0);
     497        }
     498        hdu->weights = weights;
     499    }
     500    if (! masks) {
     501        masks = psArrayAlloc(pixels->n);
     502        for (int i = 0; i < pixels->n; i++) {
     503            psImage *image = pixels->data[i];
     504            masks->data[i] = psImageAlloc(image->numCols, image->numRows, PS_TYPE_U8);
     505            psImageInit(masks->data[i], 0);
     506        }
     507        hdu->masks = masks;
     508    }
     509
     510    // Set the pixels
     511    psArray *readouts = cell->readouts; // Array of readouts
     512    for (int i = 0; i < readouts->n; i++) {
     513        pmReadout *readout = readouts->data[i]; // The readout of interest
     514
     515        if (! readout->weight) {
     516            readout->weight = weights->data[i];
     517        }
     518        if (! readout->mask) {
     519            readout->mask = masks->data[i];
     520        }
     521
     522        // Mask is already set to 0
     523
     524        // Set weight image to the variance = g*f + rn^2
     525        psImage *image = psImageSubset(readout->image, *trimsec); // The pixels
     526        psImage *weight = psImageSubset(readout->weight, *trimsec); // The weight map
     527        psBinaryOp(weight, image, "*", psScalarAlloc(gain, PS_TYPE_F32));
     528        psBinaryOp(weight, weight, "+", psScalarAlloc(readnoise*readnoise, PS_TYPE_F32));
     529    }
     530
     531    return true;
     532}
     533
     534
Note: See TracChangeset for help on using the changeset viewer.