IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 5565


Ignore:
Timestamp:
Nov 21, 2005, 5:33:37 PM (20 years ago)
Author:
eugene
Message:

reorg

Location:
trunk/psastro
Files:
1 added
11 edited

Legend:

Unmodified
Added
Removed

  • trunk/psastro/Makefile

    r5560 r5565  
    2727$(SRC)/psastroArguments.$(ARCH).o  \
    2828$(SRC)/psastroUtils.$(ARCH).o      \
     29$(SRC)/psModUtils.$(ARCH).o        \
    2930$(SRC)/psLibUtils.$(ARCH).o        \
    3031$(SRC)/pmAstrom.$(ARCH).o          \
    … …  
    3435$(SRC)/psastro-mktest.$(ARCH).o    \
    3536$(SRC)/psastroArguments.$(ARCH).o  \
     37$(SRC)/psastroIO.$(ARCH).o         \
     38$(SRC)/psModUtils.$(ARCH).o        \
    3639$(SRC)/psLibUtils.$(ARCH).o        \
    3740$(SRC)/pmAstrom.$(ARCH).o
    38 
    39 # $(SRC)/psModUtils.$(ARCH).o      \
    4041
    4142psastro: $(BIN)/psastro.$(ARCH)
    … …  
    4849
    4950INSTALL = psastro psastro-mktest
     51
     52all: psastro.install psastro-mktest.install
    5053
    5154# dependancy rules for binary code #########################

  • trunk/psastro/src/pmAstrom.c

    r5560 r5565  
    1313}
    1414
    15 psPlane *psCoordReadoutToCell (psPlane *cellpix, psPlane *readpix, pmReadout *readout) {
    16 
    17     if (cellpix == NULL) {
    18         cellpix = psPlaneAlloc ();
    19     }
    20 
    21     cellpix->x = readpix->x*readout->colBins + readout->col0;
    22     cellpix->y = readpix->y*readout->rowBins + readout->row0;
    23 
    24     return (cellpix);
    25 }
    26 
    27 psPlane *psCoordCellToReadout (psPlane *readpix, psPlane *cellpix, pmReadout *readout) {
    28 
    29     if (readpix == NULL) {
    30         readpix = psPlaneAlloc ();
    31     }
    32 
    33     readpix->x = (cellpix->x - readout->col0) / readout->colBins;
    34     readpix->y = (cellpix->y - readout->row0) / readout->rowBins;
    35 
    36     return (readpix);
    37 }
    38 
    39 psPlane* psCoordChipToCell_EAM(psPlane* cellCoord,
    40                            const psPlane* chipCoord,
    41                            const pmCell* cell)
     15// sort by mag (descending)
     16int pmAstromObjSortByMag (const void **a, const void **b)
    4217{
    43     PS_ASSERT_PTR_NON_NULL(chipCoord, NULL);
    44     PS_ASSERT_PTR_NON_NULL(cell, NULL);
    45     PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
    46 
    47     // XXX EAM : why was this being done?
    48     // pmCell *tmpCell = pmCellInChip(chipCoord, cell->parent);
    49     PS_ASSERT_PTR_NON_NULL(cell->toChip, NULL);
    50     psPlaneTransform *tmpChipToCell = p_psPlaneTransformLinearInvert(cell->toChip);
    51     PS_ASSERT_PTR_NON_NULL(tmpChipToCell, NULL);
    52     cellCoord = psPlaneTransformApply(cellCoord, tmpChipToCell, chipCoord);
    53     psFree(tmpChipToCell);
    54     return(cellCoord);
    55 }
    56 
    57 bool psastroProjectFPA (pmFPA *fpa, char *starlist, bool toSky) {
    58 
    59     bool status;
    60 
    61     for (int i = 0; i < fpa->chips->n; i++) {
    62         pmChip *chip = fpa->chips->data[i];
    63         for (int j = 0; j < chip->cells->n; j++) {
    64             pmCell *cell = chip->cells->data[j];
    65             for (int k = 0; k < cell->readouts->n; k++) {
    66                 pmReadout *readout = cell->readouts->data[k];
    67                 psArray *stars = psMetadataLookupPtr (&status, readout->analysis, starlist);
    68                 if (toSky) {
    69                     psastroProjectRawstars (stars, readout);
    70                 } else {
    71                     psastroProjectRefstars (stars, readout);
    72                 }
    73             }
    74         }
    75     }
    76     return true;
    77 }
    78  
    79 bool psastroProjectRawstars (psArray *stars, pmReadout *readout) {
    80 
    81     pmCell *cell = readout->parent;
    82     pmChip *chip = cell->parent;
    83     pmFPA  *fpa  = chip->parent;
    84 
    85     for (int i = 0; i < stars->n; i++) {
    86         pmAstromObj *star = stars->data[i];
    87         psCoordReadoutToCell (&star->cell, &star->pix, readout);
    88         psCoordCellToChip (&star->chip, &star->cell, cell);
    89         psCoordChipToFPA (&star->FP, &star->chip, chip);
    90         psCoordFPAToTP (&star->TP, &star->FP, 0.0, 0.0, fpa);
    91         psCoordTPToSky (&star->sky, &star->TP, fpa->projection);
    92     }
    93     return true;
    94 }
    95  
    96 bool psastroProjectRefstars (psArray *stars, pmReadout *readout) {
    97 
    98     pmCell *cell = readout->parent;
    99     pmChip *chip = cell->parent;
    100     pmFPA  *fpa  = chip->parent;
    101 
    102     for (int i = 0; i < stars->n; i++) {
    103         pmAstromObj *star = stars->data[i];
    104         psCoordSkyToTP (&star->TP, &star->sky, fpa->projection);
    105         psCoordTPToFPA (&star->FP, &star->TP, 0.0, 0.0, fpa);
    106         psCoordFPAToChip (&star->chip, &star->FP, chip);
    107         psCoordChipToCell_EAM (&star->cell, &star->chip, cell);
    108         psCoordCellToReadout (&star->pix, &star->cell, readout);
    109     }
    110     return true;
    111 }
    112  
     18    pmAstromObj *A = *(pmAstromObj **)a;
     19    pmAstromObj *B = *(pmAstromObj **)b;
     20
     21    psF32 diff = A->Mag - B->Mag;
     22    if (diff > FLT_EPSILON) return (-1);
     23    if (diff < FLT_EPSILON) return (+1);
     24    return (0);
     25}
     26
    11327psArray *pmAstromRadiusMatch (psArray *st1, psArray *st2, psMetadata *config) {
    11428

  • trunk/psastro/src/pmAstrom.h

    r5560 r5565  
    2929} pmAstromStats;
    3030
     31/* in pmAstrom.c */
    3132pmAstromObj      *pmAstromObjAlloc (void);
    3233pmAstromObj      *pmAstromObjCopy (pmAstromObj *old);
    … …  
    3536psArray          *pmAstromRotateObj (psArray *old, psPlane center, double angle);
    3637psPlaneTransform *pmAstromMatchFit (psPlaneTransform *map, psArray *st1, psArray *st2, psArray *match, psMetadata *config);
     38int               pmAstromObjSortByFPX (const void **a, const void **b);
     39int               pmAstromObjSortByMag (const void **a, const void **b);
    3740
    38 int               pmAstromObjSortByFPX (const void **a, const void **b);
    39 
    40 bool              psastroProjectFPA (pmFPA *fpa, char *starlist, bool toSky);
    41 bool              psastroProjectRawstars (psArray *stars, pmReadout *readout);
    42 bool              psastroProjectRefstars (psArray *stars, pmReadout *readout);
    43 
    44 psPlane* psCoordChipToCell_EAM(
    45     psPlane* out,                      ///< a plane struct to recycle. If NULL, a new struct is created
    46     const psPlane* in,                 ///< the Chip coordinate
    47     const pmCell* cell                 ///< the cell of interest
    48 );
    49 
    50 psPlaneTransform *p_psPlaneTransformLinearInvert_EAM(psPlaneTransform *transform);
     41/* in pmAstromGrid.c */
     42pmAstromStats     pmAstromGridAngle (psArray *st1, psArray *st2, psMetadata *config);
     43psPlaneTransform *pmAstromGridApply (psPlaneTransform *map, pmAstromStats stat);
     44pmAstromStats     pmAstromGridMatch (psArray *st1, psArray *st2, psMetadata *config);

  • trunk/psastro/src/psLibUtils.c

    r5562 r5565  
    171171    return(out);
    172172}
     173
     174// returns the rotation term, forcing positive parity
     175double psPlaneTransformGetRotation (psPlaneTransform *map) {
     176
     177    if (map->x->nX < 1) return 0;
     178    if (map->x->nY < 1) return 0;
     179
     180    if (map->y->nX < 1) return 0;
     181    if (map->y->nY < 1) return 0;
     182   
     183    double pc1_1 = map->x->coeff[1][0];
     184    double pc1_2 = map->x->coeff[0][1];
     185    double pc2_1 = map->y->coeff[1][0];
     186    double pc2_2 = map->y->coeff[0][1];
     187
     188    double px = SIGN (pc1_1);
     189    double py = SIGN (pc2_2);
     190
     191    // both x and y terms imply an angle. take the average
     192    double t1 = -atan2 (px*pc1_2, px*pc1_1);
     193    double t2 = +atan2 (py*pc2_1, py*pc2_2);
     194   
     195    // careful near -pi,+pi boundary...
     196    if (t1 - t2 > M_PI/2) t2 += 2*M_PI;
     197    if (t2 - t1 > M_PI/2) t1 += 2*M_PI;
     198
     199    double theta = 0.5*(t1 + t2);
     200    while (theta < M_PI) theta += 2*M_PI;
     201    while (theta > M_PI) theta -= 2*M_PI;
     202   
     203    return (theta);
     204}

  • trunk/psastro/src/psLibUtils.h

    r5562 r5565  
    11# ifndef PS_LIB_UTILS
    22# define PS_LIB_UTILS
    3 
    4 // structure to carry a dynamic string
    5 typedef struct {
    6     int NLINE;
    7     int Nline;
    8     char *line;
    9 } psLine;
    103
    114# define psMemCopy(A)(psMemIncrRefCounter((A)))
    … …  
    136# define RAD_DEG  0.017453292519943
    147
     8# define SIGN(X)  (((X) == 0) ? 0 : ((fabs((double)(X))) / (X)))
     9
    1510// psLib extra utilities
    1611psS32        psLogArguments (int *argc, char **argv);   // added to SDRS (part of psArgumentVerbosity)
    1712psS32        psTraceArguments (int *argc, char **argv); // added to SDRS (part of psArgumentVerbosity)
    1813
    19 // psLine functions -- keep out for now?
    20 psLine      *psLineAlloc (int Nline);
    21 bool         psLineInit (psLine *line);
    22 bool         psLineAdd (psLine *line, char *format, ...);
     14psPolynomial2D *psPolynomial2DCopy (psPolynomial2D *input);
     15psPolynomial4D *psPolynomial4DCopy (psPolynomial4D *input);
     16psPlaneDistort *psPlaneDistortCopy (psPlaneDistort *input);
     17psPlaneTransform *psPlaneTransformCopy (psPlaneTransform *input);
     18psProjection *psProjectionCopy (psProjection *input);
     19psPlaneDistort *psPlaneDistortInvert(psPlaneDistort *distort);
     20double psPlaneTransformGetRotation (psPlaneTransform *map);
    2321
    2422# endif

  • trunk/psastro/src/psModUtils.c

    r5560 r5565  
    11# include "psastro.h"
    22
    3 // these functions temporarily replace existing psModules code
     3psPlane *psCoordReadoutToCell (psPlane *cellpix, psPlane *readpix, pmReadout *readout) {
    44
    5 // template sample alloc/free pair:
    6 # if (0)
    7 static void psFooFree (psFoo *foo) {
     5    if (cellpix == NULL) {
     6        cellpix = psPlaneAlloc ();
     7    }
    88
    9   if (foo == NULL) return;
    10   return;
    11 }
    12 psFoo *psFooAlloc (int i1, int i2) {
     9    cellpix->x = readpix->x*readout->colBins + readout->col0;
     10    cellpix->y = readpix->y*readout->rowBins + readout->row0;
    1311
    14   psFoo *foo = psAlloc (sizeof(psFoo));
    15   psMemSetDeallocator(foo, (psFreeFunc) psFooFree);
    16 
    17   return (foo);
    18 }
    19 # endif
    20 
    21 // pmFPA
    22 static void pmFPAFree (pmFPA *fpa) {
    23 
    24   if (fpa == NULL) return;
    25 
    26   psFree (fpa->toSky);
    27   psFree (fpa->toTPA);
    28   psFree (fpa->fromTPA);
    29   psFree (fpa->chips);
    30 
    31   return;
     12    return (cellpix);
    3213}
    3314
    34 pmFPA *pmFPAAlloc (void) {
     15psPlane *psCoordCellToReadout (psPlane *readpix, psPlane *cellpix, pmReadout *readout) {
    3516
    36   pmFPA *fpa = psAlloc (sizeof(pmFPA));
    37   psMemSetDeallocator(fpa, (psFreeFunc) pmFPAFree);
     17    if (readpix == NULL) {
     18        readpix = psPlaneAlloc ();
     19    }
    3820
    39   fpa->toSky   = NULL;
    40   fpa->toTPA   = NULL;
    41   fpa->fromTPA = NULL;
    42   fpa->chips   = NULL;
     21    readpix->x = (cellpix->x - readout->col0) / readout->colBins;
     22    readpix->y = (cellpix->y - readout->row0) / readout->rowBins;
    4323
    44   return (fpa);
     24    return (readpix);
    4525}
    4626
    47 // pmChip
    48 static void pmChipFree (pmChip *chip) {
     27psPlane* psCoordChipToCell_EAM(psPlane* cellCoord,
     28                               const psPlane* chipCoord,
     29                               const pmCell* cell)
     30{
     31    PS_ASSERT_PTR_NON_NULL(chipCoord, NULL);
     32    PS_ASSERT_PTR_NON_NULL(cell, NULL);
     33    PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
    4934
    50   if (chip == NULL) return;
    51 
    52   psFree (chip->toFPA);
    53   psFree (chip->fromFPA);
    54   psFree (chip->cells);
    55 
    56   return;
     35    // XXX EAM : why was this being done?
     36    // pmCell *tmpCell = pmCellInChip(chipCoord, cell->parent);
     37    PS_ASSERT_PTR_NON_NULL(cell->toChip, NULL);
     38    psPlaneTransform *tmpChipToCell = p_psPlaneTransformLinearInvert(cell->toChip);
     39    PS_ASSERT_PTR_NON_NULL(tmpChipToCell, NULL);
     40    cellCoord = psPlaneTransformApply(cellCoord, tmpChipToCell, chipCoord);
     41    psFree(tmpChipToCell);
     42    return(cellCoord);
    5743}
    5844
    59 pmChip *pmChipAlloc (void) {
     45bool psastroProjectFPA (pmFPA *fpa, char *starlist, bool toSky) {
    6046
    61   pmChip *chip = psAlloc (sizeof(pmChip));
    62   psMemSetDeallocator(chip, (psFreeFunc) pmChipFree);
     47    bool status;
    6348
    64   chip->toFPA   = NULL;
    65   chip->fromFPA = NULL;
    66   chip->cells   = NULL;
     49    for (int i = 0; i < fpa->chips->n; i++) {
     50        pmChip *chip = fpa->chips->data[i];
     51        for (int j = 0; j < chip->cells->n; j++) {
     52            pmCell *cell = chip->cells->data[j];
     53            for (int k = 0; k < cell->readouts->n; k++) {
     54                pmReadout *readout = cell->readouts->data[k];
     55                psArray *stars = psMetadataLookupPtr (&status, readout->analysis, starlist);
     56                if (toSky) {
     57                    psastroProjectRawstars (stars, readout);
     58                } else {
     59                    psastroProjectRefstars (stars, readout);
     60                }
     61            }
     62        }
     63    }
     64    return true;
     65}
     66 
     67bool psastroProjectRawstars (psArray *stars, pmReadout *readout) {
    6768
    68   return (chip);
     69    pmCell *cell = readout->parent;
     70    pmChip *chip = cell->parent;
     71    pmFPA  *fpa  = chip->parent;
     72
     73    for (int i = 0; i < stars->n; i++) {
     74        pmAstromObj *star = stars->data[i];
     75        psCoordReadoutToCell (&star->cell, &star->pix, readout);
     76        psCoordCellToChip (&star->chip, &star->cell, cell);
     77        psCoordChipToFPA (&star->FP, &star->chip, chip);
     78        psCoordFPAToTP (&star->TP, &star->FP, 0.0, 0.0, fpa);
     79        psCoordTPToSky (&star->sky, &star->TP, fpa->projection);
     80    }
     81    return true;
     82}
     83 
     84bool psastroProjectRefstars (psArray *stars, pmReadout *readout) {
     85
     86    pmCell *cell = readout->parent;
     87    pmChip *chip = cell->parent;
     88    pmFPA  *fpa  = chip->parent;
     89
     90    for (int i = 0; i < stars->n; i++) {
     91        pmAstromObj *star = stars->data[i];
     92        psCoordSkyToTP (&star->TP, &star->sky, fpa->projection);
     93        psCoordTPToFPA (&star->FP, &star->TP, 0.0, 0.0, fpa);
     94        psCoordFPAToChip (&star->chip, &star->FP, chip);
     95        psCoordChipToCell_EAM (&star->cell, &star->chip, cell);
     96        psCoordCellToReadout (&star->pix, &star->cell, readout);
     97    }
     98    return true;
     99}
     100 
     101pmFPA *pmFPACopyAstrom (pmFPA *inFPA) {
     102
     103    pmFPA *fpa = pmFPAAlloc (NULL, NULL);
     104
     105    // copy FPA astrometry data
     106    fpa->toSky   = psProjectionCopy (inFPA->toSky);
     107    fpa->toTPA   = psPlaneDistortCopy (inFPA->toTPA);
     108    fpa->fromTPA = psPlaneDistortCopy (inFPA->fromTPA);
     109
     110    psArrayRealloc (fpa->chips, inFPA->chips->n);
     111    for (int i = 0; i < inFPA->chips->n; i++) {
     112        pmChip *inChip = inFPA->chips->data[i];
     113        pmChip *chip = pmChipAlloc (fpa);
     114
     115        // copy Chip astrometry data
     116        chip->toFPA = psPlaneTransformCopy (inChip->toFPA);
     117        chip->fromFPA = psPlaneTransformCopy (inChip->fromFPA);
     118
     119        psArrayRealloc (chip->cells, inChip->cells->n);
     120        for (int j = 0; j < inChip->cells->n; j++) {
     121            pmCell *inCell = inChip->cells->data[j];
     122            pmCell *cell = pmCellAlloc (chip);
     123
     124            // cell.header is a view on inCell.header
     125            cell->header = psMemCopy (inCell->header);
     126
     127            // copy Cell astrometry data
     128            cell->toChip = psPlaneTransformCopy (inCell->toChip);
     129
     130            psArrayRealloc (cell->readouts, inCell->readouts->n);
     131            for (int k = 0; k < inCell->readouts->n; k++) {
     132                pmReadout *inReadout = inCell->readouts->data[k];
     133                pmReadout *readout = pmReadoutAlloc (cell);
     134
     135                // copy Readout data
     136                *readout = *inReadout;
     137
     138                cell->readouts->data[k] = readout;
     139            }
     140            chip->cells->data[j] = cell;
     141        }
     142        fpa->chips->data[i] = chip;
     143    }
     144    return (fpa);
    69145}
    70146
    71 // pmCell
    72 static void pmCellFree (pmCell *cell) {
    73 
    74   if (cell == NULL) return;
    75 
    76   psFree (cell->toChip);
    77   psFree (cell->readouts);
    78 
    79   return;
    80 }
    81 
    82 pmCell *pmCellAlloc (void) {
    83 
    84   pmCell *cell = psAlloc (sizeof(pmCell));
    85   psMemSetDeallocator(cell, (psFreeFunc) pmCellFree);
    86 
    87   cell->toChip   = NULL;
    88   cell->readouts = NULL;
    89 
    90   return (cell);
    91 }
    92 
    93 // pmReadout
    94 static void pmReadoutFree (pmReadout *readout) {
    95 
    96   if (readout == NULL) return;
    97 
    98   psFree (readout->stars);
    99 
    100   return;
    101 }
    102 
    103 pmReadout *pmReadoutAlloc (void) {
    104 
    105   pmReadout *readout = psAlloc (sizeof(pmReadout));
    106   psMemSetDeallocator(readout, (psFreeFunc) pmReadoutFree);
    107 
    108   readout->col0 = 0;
    109   readout->row0 = 0;
    110   readout->colBins = 1;
    111   readout->rowBins = 1;
    112   readout->stars = NULL;
    113 
    114   return (readout);
    115 }
    116 

  • trunk/psastro/src/psastro-mktest.c

    r5562 r5565  
    11# include "psastro.h"
    2 
    3 bool testWriteCMP (psMetadata *header, char *filename, psArray *sources);
    4 bool testWriteRef (char *filename, psArray *sources);
    5 bool testWriteRaw (char *filename, psArray *sources);
    62
    73int main (int argc, char **argv) {
    … …  
    173169    exit (0);
    174170}
    175 
    176 // elixir-style pseudo FITS table (header + ascii list)
    177 bool testWriteCMP (psMetadata *header, char *filename, psArray *sources) {
    178 
    179     int i;
    180     psMetadataItem *mdi;
    181 
    182     psMetadataAdd (header, PS_LIST_TAIL, "NSTARS", PS_DATA_S32    | PS_META_REPLACE, "foo", sources->n);
    183 
    184     // create file, write-out header
    185     unlink (filename);
    186     psFits *fits = psFitsAlloc (filename);
    187 
    188     // set NAXIS to 0 : CFITSIO requires isolated header to have NAXIS = 0
    189     mdi = psMetadataLookup (header, "NAXIS");
    190     mdi->data.S32 = 0;
    191     mdi->type = PS_DATA_S32;
    192 
    193     psFitsWriteHeader (header, fits);
    194     psFree (fits);
    195 
    196     // re-open, add data to end of file
    197     FILE *f = fopen (filename, "a+");
    198     if (f == NULL) {
    199         psLogMsg ("WriteSourceOBJ", 3, "can't open output file for output %s\n", filename);
    200         return false;
    201     }
    202     fseek (f, 0, SEEK_END);
    203 
    204     for (i = 0; i < sources->n; i++) {
    205        
    206         pmAstromObj *star = sources->data[i];
    207 
    208         fprintf (f, "%6.1f %6.1f %6.3f %03d %2d %3.1f %6.3f %6.3f %6.2f %6.2f %5.1f\n",
    209                  star->pix.x, star->pix.y, star->Mag, (int)(star->dMag*1000), 1, 2.5, 11.0, 12.0, 2.0, 0.5, 15.0);
    210     }
    211     fclose (f);
    212     return true;
    213 }
    214 
    215 // elixir-style pseudo FITS table (header + ascii list)
    216 bool testWriteRef (char *filename, psArray *sources) {
    217 
    218     int i;
    219 
    220     // re-open, add data to end of file
    221     FILE *f = fopen (filename, "w");
    222     if (f == NULL) {
    223         psLogMsg ("WriteSourceOBJ", 3, "can't open output file for output %s\n", filename);
    224         return false;
    225     }
    226 
    227     for (i = 0; i < sources->n; i++) {
    228        
    229         pmAstromObj *star = sources->data[i];
    230 
    231         fprintf (f, "%10.7f %10.7f %6.3f %6.3f\n", star->sky.r, star->sky.d+10*RAD_DEG/3600, 10.0, 0.05);
    232     }
    233     fclose (f);
    234     return true;
    235 }
    236 
    237 // elixir-style pseudo FITS table (header + ascii list)
    238 bool testWriteRaw (char *filename, psArray *sources) {
    239 
    240     int i;
    241 
    242     // re-open, add data to end of file
    243     FILE *f = fopen (filename, "w");
    244     if (f == NULL) {
    245         psLogMsg ("WriteSourceOBJ", 3, "can't open output file for output %s\n", filename);
    246         return false;
    247     }
    248 
    249     for (i = 0; i < sources->n; i++) {
    250        
    251         pmAstromObj *star = sources->data[i];
    252 
    253         fprintf (f, "%8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %10.6f %10.6f   %8.2f %8.2f\n",
    254                  star->pix.x, star->pix.y,
    255                  star->cell.x, star->cell.y,
    256                  star->chip.x, star->chip.y,
    257                  star->FP.x, star->FP.y,
    258                  star->TP.x, star->TP.y,
    259                  star->sky.r*DEG_RAD, star->sky.d*DEG_RAD,
    260                  star->Mag, star->dMag);
    261     }
    262     fclose (f);
    263     return true;
    264 }

  • trunk/psastro/src/psastro.h

    r5560 r5565  
    1111# define toSky projection
    1212
    13 # if (0)
    14 typedef struct
    15 {
    16     psProjection *toSky;                ///< Projection from tangent plane to sky
    17     psPlaneDistort *toTPA;              ///< Transformation from focal plane to tangent plane
    18     psPlaneDistort *fromTPA;            ///< Transformation from focal plane to tangent plane
    19     psArray *chips;                     ///< The chips
    20 }
    21 pmFPA;
    22 
    23 typedef struct
    24 {
    25     // Astrometric transformations
    26     psPlaneTransform *toFPA;            ///< Transformation from chip to FPA coordinates
    27     psPlaneTransform *fromFPA;          ///< Transformation from FPA to chip coordinates
    28     psArray *cells;                     ///< The cells
    29     pmFPA *fpa;
    30 }
    31 pmChip;
    32 
    33 typedef struct
    34 {
    35     // Astrometric transformations
    36     psPlaneTransform *toChip;           ///< Transformations from cell to chip coordinates
    37     psArray *readouts;                  ///< The readouts (referred to by number)
    38     psMetadata *header;                 ///< header always corresponds to a cell
    39     pmChip *chip;
    40 }
    41 pmCell;
    42 
    43 typedef struct
    44 {
    45     // Position on the cell
    46     int col0;                           ///< Offset from the left of chip.
    47     int row0;                           ///< Offset from the bottom of chip.
    48     int colBins;                        ///< Amount of binning in x-dimension
    49     int rowBins;                        ///< Amount of binning in y-dimension
    50     psArray *stars;                     ///< sources detected / measured on readout
    51     pmCell *cell;
    52 }
    53 pmReadout;
    54 # endif
    55 
    56 pmAstromStats     pmAstromGridAngle (psArray *st1, psArray *st2, psMetadata *config);
    57 psPlaneTransform *pmAstromGridApply (psPlaneTransform *map, pmAstromStats stat);
    58 pmAstromStats     pmAstromGridMatch (psArray *st1, psArray *st2, psMetadata *config);
    59 int               pmAstromObjSortByMag (const void **a, const void **b);
    6013bool              pmCellInterpretWCS (pmCell *cell, psMetadata *header) ;
    6114pmFPA            *pmFPACopyAstrom (pmFPA *inFPA);
    … …  
    7023psMetadata *testArguments (int *argc, char **argv);
    7124
    72 # if (0)
    73 pmFPA            *pmFPAAlloc (void);
    74 pmChip           *pmChipAlloc (void);
    75 pmCell           *pmCellAlloc (void);
    76 pmReadout        *pmReadoutAlloc (void);
    77 # endif
    78 
    79 # define SIGN(X)  (((X) == 0) ? 0 : ((fabs((double)(X))) / (X)))
    80 
    81 psPolynomial2D *psPolynomial2DCopy (psPolynomial2D *input);
    82 psPolynomial4D *psPolynomial4DCopy (psPolynomial4D *input);
    83 psPlaneDistort *psPlaneDistortCopy (psPlaneDistort *input);
    84 psPlaneTransform *psPlaneTransformCopy (psPlaneTransform *input);
    85 psProjection *psProjectionCopy (psProjection *input);
    8625psPlane *psCoordReadoutToCell (psPlane *cellpix, psPlane *readpix, pmReadout *readout);
    8726psPlane *psCoordCellToReadout (psPlane *readpix, psPlane *cellpix, pmReadout *readout);
    88 double psPlaneTransformGetRotation (psPlaneTransform *map);
     27psPlane* psCoordChipToCell_EAM(psPlane* cellCoord, const psPlane* chipCoord, const pmCell* cell);
    8928
    90 psPlaneDistort *psPlaneDistortInvert(psPlaneDistort *distort);
    91 psPlaneTransform *psPlaneTransformLinearInvert(psPlaneTransform *transform);
     29bool testWriteCMP (psMetadata *header, char *filename, psArray *sources);
     30bool testWriteRef (char *filename, psArray *sources);
     31bool testWriteRaw (char *filename, psArray *sources);
     32
     33bool              psastroProjectFPA (pmFPA *fpa, char *starlist, bool toSky);
     34bool              psastroProjectRawstars (psArray *stars, pmReadout *readout);
     35bool              psastroProjectRefstars (psArray *stars, pmReadout *readout);
     36
     37psPlane* psCoordChipToCell_EAM(
     38    psPlane* out,                      ///< a plane struct to recycle. If NULL, a new struct is created
     39    const psPlane* in,                 ///< the Chip coordinate
     40    const pmCell* cell                 ///< the cell of interest
     41);
     42
     43psPlaneTransform *p_psPlaneTransformLinearInvert_EAM(psPlaneTransform *transform);
     44pmFPA *pmFPACopyAstrom (pmFPA *inFPA);

  • trunk/psastro/src/psastroArguments.c

    r5560 r5565  
    4646
    4747static void usage_test (void) {
    48     fprintf (stderr, "USAGE: psastro-mktest (config) (cmpfile) (reflist) (rawfile) (rawfile)\n");
     48    fprintf (stderr, "USAGE: psastro-mktest (config) (cmpfile) (catfile) (rawfile) (reffile)\n");
    4949    exit (2);
    5050}

  • trunk/psastro/src/psastroIO.c

    r5562 r5565  
    124124}
    125125
     126// elixir-style pseudo FITS table (header + ascii list)
     127bool testWriteCMP (psMetadata *header, char *filename, psArray *sources) {
     128
     129    int i;
     130    psMetadataItem *mdi;
     131
     132    psMetadataAdd (header, PS_LIST_TAIL, "NSTARS", PS_DATA_S32    | PS_META_REPLACE, "foo", sources->n);
     133
     134    // create file, write-out header
     135    unlink (filename);
     136    psFits *fits = psFitsAlloc (filename);
     137
     138    // set NAXIS to 0 : CFITSIO requires isolated header to have NAXIS = 0
     139    mdi = psMetadataLookup (header, "NAXIS");
     140    mdi->data.S32 = 0;
     141    mdi->type = PS_DATA_S32;
     142
     143    psFitsWriteHeader (header, fits);
     144    psFree (fits);
     145
     146    // re-open, add data to end of file
     147    FILE *f = fopen (filename, "a+");
     148    if (f == NULL) {
     149        psLogMsg ("WriteSourceOBJ", 3, "can't open output file for output %s\n", filename);
     150        return false;
     151    }
     152    fseek (f, 0, SEEK_END);
     153
     154    for (i = 0; i < sources->n; i++) {
     155       
     156        pmAstromObj *star = sources->data[i];
     157
     158        fprintf (f, "%6.1f %6.1f %6.3f %03d %2d %3.1f %6.3f %6.3f %6.2f %6.2f %5.1f\n",
     159                 star->pix.x, star->pix.y, star->Mag, (int)(star->dMag*1000), 1, 2.5, 11.0, 12.0, 2.0, 0.5, 15.0);
     160    }
     161    fclose (f);
     162    return true;
     163}
     164
     165// write out raw objects
     166bool testWriteRaw (char *filename, psArray *sources) {
     167
     168    int i;
     169
     170    // re-open, add data to end of file
     171    FILE *f = fopen (filename, "w");
     172    if (f == NULL) {
     173        psLogMsg ("WriteSourceOBJ", 3, "can't open output file for output %s\n", filename);
     174        return false;
     175    }
     176
     177    for (i = 0; i < sources->n; i++) {
     178       
     179        pmAstromObj *star = sources->data[i];
     180
     181        fprintf (f, "%8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %10.6f %10.6f   %8.2f %8.2f\n",
     182                 star->pix.x, star->pix.y,
     183                 star->cell.x, star->cell.y,
     184                 star->chip.x, star->chip.y,
     185                 star->FP.x, star->FP.y,
     186                 star->TP.x, star->TP.y,
     187                 star->sky.r*DEG_RAD, star->sky.d*DEG_RAD,
     188                 star->Mag, star->dMag);
     189    }
     190    fclose (f);
     191    return true;
     192}
     193
     194// elixir-style pseudo FITS table (header + ascii list)
     195bool testWriteRef (char *filename, psArray *sources) {
     196
     197    int i;
     198
     199    // re-open, add data to end of file
     200    FILE *f = fopen (filename, "w");
     201    if (f == NULL) {
     202        psLogMsg ("WriteSourceOBJ", 3, "can't open output file for output %s\n", filename);
     203        return false;
     204    }
     205
     206    for (i = 0; i < sources->n; i++) {
     207       
     208        pmAstromObj *star = sources->data[i];
     209
     210        fprintf (f, "%10.7f %10.7f %6.3f %6.3f\n", star->sky.r, star->sky.d+10*RAD_DEG/3600, 10.0, 0.05);
     211    }
     212    fclose (f);
     213    return true;
     214}
     215
     216void psastroDumpStars (psArray *sources) {
     217
     218    for (int i = 0; i < sources->n; i++) {
     219        pmAstromObj *star = sources->data[i];
     220
     221        fprintf (stderr, "%8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %10.6f %10.6f   %8.2f %8.2f\n",
     222                 star->pix.x, star->pix.y,
     223                 star->cell.x, star->cell.y,
     224                 star->chip.x, star->chip.y,
     225                 star->FP.x, star->FP.y,
     226                 star->TP.x, star->TP.y,
     227                 star->sky.r*DEG_RAD, star->sky.d*DEG_RAD,
     228                 star->Mag, star->dMag);
     229    }
     230}

  • trunk/psastro/src/psastroUtils.c

    r5560 r5565  
    11# include "psastro.h"
    22# define RENORM 0
    3 
    4 bool testWriteRaw (char *filename, psArray *sources);
    5 
    6 void psastroDumpStars (psArray *sources) {
    7 
    8     for (int i = 0; i < sources->n; i++) {
    9         pmAstromObj *star = sources->data[i];
    10 
    11         fprintf (stderr, "%8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %10.6f %10.6f   %8.2f %8.2f\n",
    12                  star->pix.x, star->pix.y,
    13                  star->cell.x, star->cell.y,
    14                  star->chip.x, star->chip.y,
    15                  star->FP.x, star->FP.y,
    16                  star->TP.x, star->TP.y,
    17                  star->sky.r*DEG_RAD, star->sky.d*DEG_RAD,
    18                  star->Mag, star->dMag);
    19     }
    20 }
    21 
    22 // sort by mag (descending)
    23 int pmAstromObjSortByMag (const void **a, const void **b)
    24 {
    25     pmAstromObj *A = *(pmAstromObj **)a;
    26     pmAstromObj *B = *(pmAstromObj **)b;
    27 
    28     psF32 diff = A->Mag - B->Mag;
    29     if (diff > FLT_EPSILON) return (-1);
    30     if (diff < FLT_EPSILON) return (+1);
    31     return (0);
    32 }
    333
    344bool psastroSelectBrightStars (pmFPA *fpa, psMetadata *config) {
    … …  
    6232}
    6333
    64 pmFPA *pmFPACopyAstrom (pmFPA *inFPA) {
    65 
    66     pmFPA *fpa = pmFPAAlloc (NULL, NULL);
    67 
    68     // copy FPA astrometry data
    69     fpa->toSky   = psProjectionCopy (inFPA->toSky);
    70     fpa->toTPA   = psPlaneDistortCopy (inFPA->toTPA);
    71     fpa->fromTPA = psPlaneDistortCopy (inFPA->fromTPA);
    72 
    73     psArrayRealloc (fpa->chips, inFPA->chips->n);
    74     for (int i = 0; i < inFPA->chips->n; i++) {
    75         pmChip *inChip = inFPA->chips->data[i];
    76         pmChip *chip = pmChipAlloc (fpa);
    77 
    78         // copy Chip astrometry data
    79         chip->toFPA = psPlaneTransformCopy (inChip->toFPA);
    80         chip->fromFPA = psPlaneTransformCopy (inChip->fromFPA);
    81 
    82         psArrayRealloc (chip->cells, inChip->cells->n);
    83         for (int j = 0; j < inChip->cells->n; j++) {
    84             pmCell *inCell = inChip->cells->data[j];
    85             pmCell *cell = pmCellAlloc (chip);
    86 
    87             // cell.header is a view on inCell.header
    88             cell->header = psMemCopy (inCell->header);
    89 
    90             // copy Cell astrometry data
    91             cell->toChip = psPlaneTransformCopy (inCell->toChip);
    92 
    93             psArrayRealloc (cell->readouts, inCell->readouts->n);
    94             for (int k = 0; k < inCell->readouts->n; k++) {
    95                 pmReadout *inReadout = inCell->readouts->data[k];
    96                 pmReadout *readout = pmReadoutAlloc (cell);
    97 
    98                 // copy Readout data
    99                 *readout = *inReadout;
    100 
    101                 cell->readouts->data[k] = readout;
    102             }
    103             chip->cells->data[j] = cell;
    104         }
    105         fpa->chips->data[i] = chip;
    106     }
    107     return (fpa);
    108 }
    109 
    11034// measure per-chip astrometry terms
    11135bool psastroChipAstrom (pmFPA *fpa, psMetadata *config) {
    … …  
    15377                psastroProjectRefstars (refstars, readout);
    15478
    155                 testWriteRaw ("ref.inp", refstars);
    156                 testWriteRaw ("raw.inp", rawstars);
     79                // testWriteRaw ("ref.inp", refstars);
     80                // testWriteRaw ("raw.inp", rawstars);
    15781
    15882                // fprintf (stderr, "rawstars:\n");
    … …  
    17296                psastroProjectRefstars (refstars, readout);
    17397
    174                 testWriteRaw ("ref.dat", refstars);
    175                 testWriteRaw ("raw.dat", rawstars);
     98                // testWriteRaw ("ref.dat", refstars);
     99                // testWriteRaw ("raw.dat", rawstars);
    176100   
    177101                // use small radius to match stars
    … …  
    284208    return true;
    285209}
    286 
    287 // returns the rotation term, forcing positive parity
    288 double psPlaneTransformGetRotation (psPlaneTransform *map) {
    289 
    290     if (map->x->nX < 1) return 0;
    291     if (map->x->nY < 1) return 0;
    292 
    293     if (map->y->nX < 1) return 0;
    294     if (map->y->nY < 1) return 0;
    295    
    296     double pc1_1 = map->x->coeff[1][0];
    297     double pc1_2 = map->x->coeff[0][1];
    298     double pc2_1 = map->y->coeff[1][0];
    299     double pc2_2 = map->y->coeff[0][1];
    300 
    301     double px = SIGN (pc1_1);
    302     double py = SIGN (pc2_2);
    303 
    304     // both x and y terms imply an angle. take the average
    305     double t1 = -atan2 (px*pc1_2, px*pc1_1);
    306     double t2 = +atan2 (py*pc2_1, py*pc2_2);
    307    
    308     // careful near -pi,+pi boundary...
    309     if (t1 - t2 > M_PI/2) t2 += 2*M_PI;
    310     if (t2 - t1 > M_PI/2) t1 += 2*M_PI;
    311 
    312     double theta = 0.5*(t1 + t2);
    313     while (theta < M_PI) theta += 2*M_PI;
    314     while (theta > M_PI) theta -= 2*M_PI;
    315    
    316     return (theta);
    317 }
    318 
    319 // elixir-style pseudo FITS table (header + ascii list)
    320 bool testWriteRaw (char *filename, psArray *sources) {
    321 
    322     int i;
    323 
    324     // re-open, add data to end of file
    325     FILE *f = fopen (filename, "w");
    326     if (f == NULL) {
    327         psLogMsg ("WriteSourceOBJ", 3, "can't open output file for output %s\n", filename);
    328         return false;
    329     }
    330 
    331     for (i = 0; i < sources->n; i++) {
    332        
    333         pmAstromObj *star = sources->data[i];
    334 
    335         fprintf (f, "%8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %8.2f %8.2f   %10.6f %10.6f   %8.2f %8.2f\n",
    336                  star->pix.x, star->pix.y,
    337                  star->cell.x, star->cell.y,
    338                  star->chip.x, star->chip.y,
    339                  star->FP.x, star->FP.y,
    340                  star->TP.x, star->TP.y,
    341                  star->sky.r*DEG_RAD, star->sky.d*DEG_RAD,
    342                  star->Mag, star->dMag);
    343     }
    344     fclose (f);
    345     return true;
    346 }
Note: See TracChangeset for help on using the changeset viewer.