IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 8046


Ignore:
Timestamp:
Aug 1, 2006, 4:15:32 PM (20 years ago)
Author:
Paul Price
Message:

Initialising values to avoid compiler error.

Location:
trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/psastro/src/psastroAstromGuess.c

    r7332 r8046  
    44// XXX WCS keywords in the headers corresponding to the chips.
    55// XXX this function is both converting the header WCS astrometry terms to the fpa terms
    6 //     and applying the astrometry to the detected objects. 
     6//     and applying the astrometry to the detected objects.
    77bool psastroAstromGuess (pmConfig *config) {
    88
     
    1010    bool status = false;
    1111    bool isMosaic = false;
    12     double RAmin, RAmax, RAminSky, RAmaxSky;
    13     double DECmin, DECmax;
     12    double RAmin = NAN, RAmax = NAN, RAminSky = NAN, RAmaxSky = NAN;
     13    double DECmin = NAN, DECmax = NAN;
    1414
    1515    pmChip *chip = NULL;
     
    2323    psMetadata *recipe  = psMetadataLookupPtr (NULL, config->recipes, "PSASTRO");
    2424    if (!recipe) {
    25         psErrorStackPrint(stderr, "Can't find PSASTRO recipe!\n");
    26         exit(EXIT_FAILURE);
     25        psErrorStackPrint(stderr, "Can't find PSASTRO recipe!\n");
     26        exit(EXIT_FAILURE);
    2727    }
    2828
     
    3030    pmFPAfile *input = psMetadataLookupPtr (NULL, config->files, "PSASTRO.INPUT");
    3131    if (!input) {
    32         psErrorStackPrint(stderr, "Can't find input data!\n");
    33         exit(EXIT_FAILURE);
     32        psErrorStackPrint(stderr, "Can't find input data!\n");
     33        exit(EXIT_FAILURE);
    3434    }
    3535
     
    4444        psTrace (__func__, 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process);
    4545        if (!chip->process || !chip->file_exists) { continue; }
    46        
     46
    4747        // read WCS data from the corresponding header
    4848        pmHDU *hdu = pmFPAviewThisHDU (view, fpa);
    4949
    50         pmAstromReadWCS (fpa, chip, hdu->header, plateScale, isMosaic);
    51         if (newFPA) {
    52             newFPA = false;
    53             RAminSky = fpa->projection->R - M_PI;
    54             RAmaxSky = fpa->projection->R + M_PI;
    55             RAmin = RAmax = fpa->projection->R;
    56             DECmin = DECmax = fpa->projection->D;
    57         }
     50        pmAstromReadWCS (fpa, chip, hdu->header, plateScale, isMosaic);
     51        if (newFPA) {
     52            newFPA = false;
     53            RAminSky = fpa->projection->R - M_PI;
     54            RAmaxSky = fpa->projection->R + M_PI;
     55            RAmin = RAmax = fpa->projection->R;
     56            DECmin = DECmax = fpa->projection->D;
     57        }
    5858
    59         // apply the new WCS guess data to all of the data in the readouts
    60         // XXX should this go into a different function? this would separate WCS interpretation from application
    61         while ((cell = pmFPAviewNextCell (view, fpa, 1)) != NULL) {
     59        // apply the new WCS guess data to all of the data in the readouts
     60        // XXX should this go into a different function? this would separate WCS interpretation from application
     61        while ((cell = pmFPAviewNextCell (view, fpa, 1)) != NULL) {
    6262            psTrace (__func__, 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process);
    6363            if (!cell->process || !cell->file_exists) { continue; }
    6464
    65             // process each of the readouts
    66             while ((readout = pmFPAviewNextReadout (view, fpa, 1)) != NULL) {
    67                 if (! readout->data_exists) { continue; }
     65            // process each of the readouts
     66            while ((readout = pmFPAviewNextReadout (view, fpa, 1)) != NULL) {
     67                if (! readout->data_exists) { continue; }
    6868
    69                 psArray *rawstars = psMetadataLookupPtr (NULL, readout->analysis, "PSASTRO.RAWSTARS");
    70                 if (rawstars == NULL) { continue; }
     69                psArray *rawstars = psMetadataLookupPtr (NULL, readout->analysis, "PSASTRO.RAWSTARS");
     70                if (rawstars == NULL) { continue; }
    7171
    72                 for (int i = 0; i < rawstars->n; i++) {
    73                     pmAstromObj *raw = rawstars->data[i];
    74        
    75                     psPlaneTransformApply (raw->FP, chip->toFPA, raw->chip);
    76                     psPlaneDistortApply (raw->TP, fpa->toTangentPlane, raw->FP, 0.0, 0.0);
    77                     p_psDeproject (raw->sky, raw->TP, fpa->projection);
     72                for (int i = 0; i < rawstars->n; i++) {
     73                    pmAstromObj *raw = rawstars->data[i];
    7874
    79                     if (i < 0) {
    80                         fprintf (stderr, "up: %f,%f -> %f,%f -> %f,%f -> %f,%f\n",
    81                                  raw->chip->x, raw->chip->y,
    82                                  raw->FP->x, raw->FP->y,
    83                                  raw->TP->x, raw->TP->y,
    84                                  raw->sky->r, raw->sky->d);
     75                    psPlaneTransformApply (raw->FP, chip->toFPA, raw->chip);
     76                    psPlaneDistortApply (raw->TP, fpa->toTangentPlane, raw->FP, 0.0, 0.0);
     77                    p_psDeproject (raw->sky, raw->TP, fpa->projection);
    8578
    86                         psPlane *fp = psPlaneAlloc();
    87                         psPlane *tp = psPlaneAlloc();
    88                         psPlane *ch = psPlaneAlloc();
     79                    if (i < 0) {
     80                        fprintf (stderr, "up: %f,%f -> %f,%f -> %f,%f -> %f,%f\n",
     81                                 raw->chip->x, raw->chip->y,
     82                                 raw->FP->x, raw->FP->y,
     83                                 raw->TP->x, raw->TP->y,
     84                                 raw->sky->r, raw->sky->d);
    8985
    90                         p_psProject (tp, raw->sky, fpa->projection);
    91                         psPlaneDistortApply (fp, fpa->fromTangentPlane, tp, 0.0, 0.0);
    92                         psPlaneTransformApply (ch, chip->fromFPA, fp);
     86                        psPlane *fp = psPlaneAlloc();
     87                        psPlane *tp = psPlaneAlloc();
     88                        psPlane *ch = psPlaneAlloc();
    9389
    94                         fprintf (stderr, "dn: %f,%f <- %f,%f <- %f,%f <- %f,%f\n",
    95                                  ch->x, ch->y,
    96                                  fp->x, fp->y,
    97                                  tp->x, tp->y,
    98                                  raw->sky->r, raw->sky->d);
    99                     }
     90                        p_psProject (tp, raw->sky, fpa->projection);
     91                        psPlaneDistortApply (fp, fpa->fromTangentPlane, tp, 0.0, 0.0);
     92                        psPlaneTransformApply (ch, chip->fromFPA, fp);
    10093
    101                     // rationalize ra to sky range centered on boresite
    102                     while (raw->sky->r < RAminSky) raw->sky->r += M_PI;
    103                     while (raw->sky->r > RAmaxSky) raw->sky->r -= M_PI;
     94                        fprintf (stderr, "dn: %f,%f <- %f,%f <- %f,%f <- %f,%f\n",
     95                                 ch->x, ch->y,
     96                                 fp->x, fp->y,
     97                                 tp->x, tp->y,
     98                                 raw->sky->r, raw->sky->d);
     99                    }
    104100
    105                     RAmin = PS_MIN (raw->sky->r, RAmin);
    106                     RAmax = PS_MAX (raw->sky->r, RAmax);
     101                    // rationalize ra to sky range centered on boresite
     102                    while (raw->sky->r < RAminSky) raw->sky->r += M_PI;
     103                    while (raw->sky->r > RAmaxSky) raw->sky->r -= M_PI;
    107104
    108                     DECmin = PS_MIN (raw->sky->d, DECmin);
    109                     DECmax = PS_MAX (raw->sky->d, DECmax);
    110                 }
    111             }
    112         }
     105                    RAmin = PS_MIN (raw->sky->r, RAmin);
     106                    RAmax = PS_MAX (raw->sky->r, RAmax);
     107
     108                    DECmin = PS_MIN (raw->sky->d, DECmin);
     109                    DECmax = PS_MAX (raw->sky->d, DECmax);
     110                }
     111            }
     112        }
    113113    }
    114114
     
    119119    psMetadataAdd (recipe, PS_LIST_TAIL, "DEC_MIN", PS_DATA_F32 | PS_META_REPLACE, "", DECmin);
    120120    psMetadataAdd (recipe, PS_LIST_TAIL, "DEC_MAX", PS_DATA_F32 | PS_META_REPLACE, "", DECmax);
    121    
     121
    122122    psFree (view);
    123123    return true;
  • trunk/psphot/src/psphotChoosePSF.c

    r7758 r8046  
    22
    33// try PSF models and select best option
    4 pmPSF *psphotChoosePSF (pmReadout *readout, psArray *sources, psMetadata *recipe) { 
     4pmPSF *psphotChoosePSF (pmReadout *readout, psArray *sources, psMetadata *recipe) {
    55
    6     bool            status;
    7     char           *modelName;
    8     pmPSF          *psf = NULL;
    9     pmPSFtry       *try = NULL;
    10     psArray        *stars = NULL;
     6    bool            status;
     7    char           *modelName;
     8    pmPSF          *psf = NULL;
     9    pmPSFtry       *try = NULL;
     10    psArray        *stars = NULL;
    1111
    1212    psTimerStart ("psphot");
     
    3232    // select the candidate PSF stars (pointers to original sources)
    3333    for (int i = 0; (i < sources->n) && (stars->n < NSTARS); i++) {
    34         pmSource *source = sources->data[i];
    35         if (source->mode & PM_SOURCE_MODE_PSFSTAR) psArrayAdd (stars, 200, source);
     34        pmSource *source = sources->data[i];
     35        if (source->mode & PM_SOURCE_MODE_PSFSTAR) psArrayAdd (stars, 200, source);
    3636    }
    3737    psLogMsg ("psphot.pspsf", 4, "selected candidate %d PSF objects\n", stars->n);
     
    4747
    4848    if (mdi->type == PS_DATA_STRING) {
    49         list = psListAlloc(NULL);
    50         psListAdd (list, PS_LIST_HEAD, mdi);
     49        list = psListAlloc(NULL);
     50        psListAdd (list, PS_LIST_HEAD, mdi);
    5151    } else {
    52         if (mdi->type != PS_DATA_METADATA_MULTI) psAbort ("psphotChoosePSF", "missing PSF_MODEL selection");
    53         list = psMemIncrRefCounter(mdi->data.list);
     52        if (mdi->type != PS_DATA_METADATA_MULTI) psAbort ("psphotChoosePSF", "missing PSF_MODEL selection");
     53        list = psMemIncrRefCounter(mdi->data.list);
    5454    }
    5555
     
    5959
    6060    // try each model option listed in config
    61     psListIterator *iter = psListIteratorAlloc (list, PS_LIST_HEAD, FALSE); 
    62     for (int i = 0; i < models->n; i++) { 
    63         psMetadataItem *item = psListGetAndIncrement (iter);
    64         modelName = item->data.V;
    65         models->data[i] = pmPSFtryModel (stars, modelName, RADIUS, POISSON_ERRORS);
     61    psListIterator *iter = psListIteratorAlloc (list, PS_LIST_HEAD, FALSE);
     62    for (int i = 0; i < models->n; i++) {
     63        psMetadataItem *item = psListGetAndIncrement (iter);
     64        modelName = item->data.V;
     65        models->data[i] = pmPSFtryModel (stars, modelName, RADIUS, POISSON_ERRORS);
    6666    }
    6767    psFree (iter);
     
    7575    float bestM = try->psf->dApResid;
    7676    for (int i = 1; i < models->n; i++) {
    77         try = models->data[i];
    78         float M = try->psf->dApResid;
    79         if (M < bestM) {
    80             bestM = M;
    81             bestN = i;
    82         }
     77        try = models->data[i];
     78        float M = try->psf->dApResid;
     79        if (M < bestM) {
     80            bestM = M;
     81            bestN = i;
     82        }
    8383    }
    84    
     84
    8585    // use the best model:
    8686    try = models->data[bestN];
     
    8888    // unset the PSFSTAR flag for stars not used for PSF model
    8989    for (int i = 0; i < try->sources->n; i++) {
    90         pmSource *source = try->sources->data[i];
    91         if (try->mask->data.U8[i] & PSFTRY_MASK_EXT_FAIL) {
    92             source->mode &= ~PM_SOURCE_MODE_PSFSTAR;
    93         }
     90        pmSource *source = try->sources->data[i];
     91        if (try->mask->data.U8[i] & PSFTRY_MASK_EXT_FAIL) {
     92            source->mode &= ~PM_SOURCE_MODE_PSFSTAR;
     93        }
    9494    }
    9595
     
    130130    // get the model full-width at half-max
    131131    psF64 FWHM_X = 2*modelRadiusFunc (modelPSF->params, 0.5);
    132    
     132
    133133    shape.sx  = modelPSF->params->data.F32[4];
    134134    shape.sy  = modelPSF->params->data.F32[5];
     
    141141    psMetadataAdd (recipe, PS_LIST_TAIL, "FWHM_Y", PS_DATA_F32 | PS_META_REPLACE, "PSF FWHM Minor axis", FWHM_Y);
    142142    psMetadataAdd (recipe, PS_LIST_TAIL, "ANGLE",  PS_DATA_F32 | PS_META_REPLACE, "PSF angle",           axes.theta);
    143    
     143
    144144    psFree (modelEXT);
    145145    psFree (modelPSF);
     
    150150bool psphotMomentsStats (pmReadout *readout, psMetadata *recipe, psArray *sources) {
    151151
    152     // without the PSF model, we can only rely on the source->moments 
     152    // without the PSF model, we can only rely on the source->moments
    153153    // as a measure of the FWHM values
    154154
    155     double FWHM_X, FWHM_Y, FWHM_T;
    156     int FWHM_N;
     155    double FWHM_X = 0.0, FWHM_Y = 0.0, FWHM_T = 0.0;
     156    int FWHM_N = 0;
    157157
    158158    psEllipseMoments moments;
    159159    psEllipseAxes axes;
    160160
    161     FWHM_N = 0;
    162     FWHM_X = FWHM_Y = 0.0;
     161    for (int i = 0; i < sources->n; i++) {
     162        pmSource *source = sources->data[i];
     163        if (!(source->mode & PM_SOURCE_MODE_PSFSTAR)) continue;
    163164
    164     for (int i = 0; i < sources->n; i++) {
    165         pmSource *source = sources->data[i];
    166         if (!(source->mode & PM_SOURCE_MODE_PSFSTAR)) continue;
    167        
    168         // moments->Sx,Sy are roughly sigma_x,y
    169         moments.x2 = PS_SQR(source->moments->Sx);
    170         moments.y2 = PS_SQR(source->moments->Sy);
    171         moments.xy = source->moments->Sxy;
     165        // moments->Sx,Sy are roughly sigma_x,y
     166        moments.x2 = PS_SQR(source->moments->Sx);
     167        moments.y2 = PS_SQR(source->moments->Sy);
     168        moments.xy = source->moments->Sxy;
    172169
    173         axes = psEllipseMomentsToAxes (moments);
     170        axes = psEllipseMomentsToAxes (moments);
    174171
    175         FWHM_X += axes.major * 2.35;
    176         FWHM_Y += axes.minor * 2.35;
    177         FWHM_T += axes.theta;
    178         FWHM_N ++;
     172        FWHM_X += axes.major * 2.35;
     173        FWHM_Y += axes.minor * 2.35;
     174        FWHM_T += axes.theta;
     175        FWHM_N ++;
    179176    }
    180177
     
    186183    psMetadataAdd (recipe, PS_LIST_TAIL, "FWHM_Y", PS_DATA_F32 | PS_META_REPLACE, "PSF FWHM Minor axis", FWHM_Y);
    187184    psMetadataAdd (recipe, PS_LIST_TAIL, "ANGLE",  PS_DATA_F32 | PS_META_REPLACE, "PSF angle",           FWHM_T);
    188    
     185
    189186    return true;
    190187}
Note: See TracChangeset for help on using the changeset viewer.