IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 26964


Ignore:
Timestamp:
Feb 16, 2010, 3:36:17 PM (16 years ago)
Author:
Paul Price
Message:

Make zero point a 'concept' so that it is easily passed from stage to stage. Checked that this results in calibrated magnitudes out of ppSub.

Location:
trunk
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • trunk/ippconfig/gpc1/format_mef.config

    r23624 r26964  
    179179        FPA.AZ          STR     AZ
    180180        FPA.TEMP        STR     DETTEM
     181        FPA.ZP          STR     MAG_ZP
    181182        CHIP.ID         STR     DETECTOR
    182183        CELL.XBIN       STR     CCDSUM

  • trunk/ippconfig/gpc1/format_relphot.config

    r23200 r26964  
    116116        FPA.AZ          STR     AZ
    117117        FPA.TEMP        STR     DETTEM
     118        FPA.ZP          STR     MAG_ZP
    118119        CHIP.ID         STR     DETECTOR
    119120        CELL.XBIN       STR     CCDSUM

  • trunk/psModules/src/concepts/pmConcepts.c

    r26002 r26964  
    292292        conceptRegisterF32("FPA.FOCUS", "Telescope focus", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA);
    293293        conceptRegisterF32("FPA.AIRMASS", "Airmass at boresight", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA);
    294         // XXX p_pmConceptParse_FPA_FILTER -> p_pmConceptParse_FPA_FILTERID (and Format as well)?
     294        // XXX p_pmConceptParse_FPA_FILTER -> p_pmConceptParse_FPA_FILTERID (and Format as well)?
    295295        conceptRegisterStr("FPA.FILTERID", "Filter used (parsed, abstract name)", p_pmConceptParse_FPA_FILTER, p_pmConceptFormat_FPA_FILTER, NULL, false, PM_FPA_LEVEL_FPA);
    296296        conceptRegisterStr("FPA.FILTER", "Filter used (instrument name)", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA);
    … …  
    334334        conceptRegisterF32("FPA.TELTEMP.EXTRA", "Telescope Temperatures: extra", p_pmConceptParse_TELTEMPS, NULL, NULL, false, PM_FPA_LEVEL_FPA);
    335335        conceptRegisterF32("FPA.PON.TIME", "Power On Time", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA);
    336         conceptRegisterS32("FPA.BURNTOOL.APPLIED", "[T=applied] Burn streaks applied to image data", p_pmConceptParse_BTOOLAPP,p_pmConceptFormat_BTOOLAPP,NULL,false,PM_FPA_LEVEL_FPA);
    337         //      conceptRegisterBool("FPA.BURNTOOL.APPLIED", "Status of burntool processing", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA);
     336        conceptRegisterS32("FPA.BURNTOOL.APPLIED", "[T=applied] Burn streaks applied to image data", p_pmConceptParse_BTOOLAPP,p_pmConceptFormat_BTOOLAPP,NULL,false,PM_FPA_LEVEL_FPA);
     337        //      conceptRegisterBool("FPA.BURNTOOL.APPLIED", "Status of burntool processing", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA);
    338338        conceptRegisterF32("FPA.EXPOSURE", "Exposure time (sec)", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA);
     339        conceptRegisterF32("FPA.ZP", "Magnitude zero point", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA);
    339340    }
    340341    if (!conceptsChip) {

  • trunk/psModules/src/concepts/pmConceptsAverage.c

    r22706 r26964  
    3636
    3737    double time      = 0.0;             // Time of observation
     38    double zp        = 0.0;             // Zero point
    3839    psTimeType timeSys = 0;             // Time system
    3940    char *filter     = NULL;            // Filter
    … …  
    6364        psTimeConvert(fpaTime, PS_TIME_TAI);
    6465        time       += psTimeToMJD(fpaTime);
     66
     67        zp += psMetadataLookupF32(NULL, fpa->concepts, "FPA.ZP");
     68
    6569        if (num == 1) {
    6670            timeSys = psMetadataLookupS32(NULL, fpa->concepts, "FPA.TIMESYS");
    … …  
    8589
    8690    time /= (double)num;
     91    zp /= (double)num;
    8792
    8893    MD_UPDATE(target->concepts, "FPA.TIMESYS", S32, timeSys);
    … …  
    9297    MD_UPDATE_STR(target->concepts, "FPA.INSTRUMENT", instrument);
    9398    MD_UPDATE_STR(target->concepts, "FPA.DETECTOR", detector);
     99    MD_UPDATE(target->concepts, "FPA.ZP", F32, zp);
    94100
    95101    // FPA.TIME needs special care

  • trunk/psModules/src/objects/pmSourceIO_CMF_PS1_V1.c

    r25754 r26964  
    6969    float magOffset = NAN;
    7070    float exptime   = psMetadataLookupF32 (&status1, fpa->concepts, "FPA.EXPOSURE");
    71     float zeropt    = psMetadataLookupF32 (&status2, imageHeader, "ZPT_OBS");
    72     float zeroptErr = psMetadataLookupF32 (&status2, imageHeader, "ZPT_ERR");
     71    float zeropt    = psMetadataLookupF32(&status2, fpa->concepts, "FPA.ZP");
     72    if (!isfinite(zeropt)) {
     73        zeropt    = psMetadataLookupF32 (&status2, imageHeader, "ZPT_OBS");
     74    }
    7375    if (status1 && status2 && (exptime > 0.0)) {
    7476        magOffset = zeropt + 2.5*log10(exptime);
    7577    }
     78    float zeroptErr = psMetadataLookupF32 (&status2, imageHeader, "ZPT_ERR");
    7679
    7780    // if the sequence is defined, write these in seq order; otherwise
    … …  
    477480        }
    478481
    479 # endif 
     482# endif
    480483        psArrayAdd (table, 100, row);
    481484        psFree (row);

  • trunk/psModules/src/objects/pmSourceIO_CMF_PS1_V2.c

    r26893 r26964  
    6969    float magOffset = NAN;
    7070    float exptime   = psMetadataLookupF32 (&status1, fpa->concepts, "FPA.EXPOSURE");
    71     float zeropt    = psMetadataLookupF32 (&status2, imageHeader, "ZPT_OBS");
    72     float zeroptErr = psMetadataLookupF32 (&status2, imageHeader, "ZPT_ERR");
     71    float zeropt    = psMetadataLookupF32(&status2, fpa->concepts, "FPA.ZP");
     72    if (!isfinite(zeropt)) {
     73        zeropt    = psMetadataLookupF32 (&status2, imageHeader, "ZPT_OBS");
     74    }
    7375    if (status1 && status2 && (exptime > 0.0)) {
    7476        magOffset = zeropt + 2.5*log10(exptime);
    7577    }
     78    float zeroptErr = psMetadataLookupF32 (&status2, imageHeader, "ZPT_ERR");
    7679
    7780    // if the sequence is defined, write these in seq order; otherwise

  • trunk/psastro/src/psastroZeroPoint.c

    r26558 r26964  
    11/** @file psastroMosaicZeroPoint.c
    22 *
    3  *  @brief 
     3 *  @brief
    44 *
    55 *  @ingroup libpsastro
    … …  
    2424    psMetadata *recipe  = psMetadataLookupPtr (NULL, config->recipes, PSASTRO_RECIPE);
    2525    if (!recipe) {
    26         psError(PSASTRO_ERR_CONFIG, false, "Can't find PSASTRO recipe!\n");
    27         return false;
     26        psError(PSASTRO_ERR_CONFIG, false, "Can't find PSASTRO recipe!\n");
     27        return false;
    2828    }
    2929
    … …  
    3434    pmFPAfile *input = psMetadataLookupPtr (NULL, config->files, "PSASTRO.INPUT");
    3535    if (!input) {
    36         psError(PSASTRO_ERR_CONFIG, false, "Can't find input data!\n");
    37         return false;
     36        psError(PSASTRO_ERR_CONFIG, false, "Can't find input data!\n");
     37        return false;
    3838    }
    3939    pmFPA *fpa = input->fpa;
    … …  
    5555    // really error-out here?  or just skip?
    5656    if (!psastroZeroPointFromRecipe (&zeropt, &exptime, NULL, fpa, recipe)) {
    57         psLogMsg ("psastro", PS_LOG_INFO, "failed to load zeropt data from recipe");
    58         return false;
     57        psLogMsg ("psastro", PS_LOG_INFO, "failed to load zeropt data from recipe");
     58        return false;
    5959    }
    6060
    … …  
    6767        psTrace ("psastro", 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process);
    6868        if (!chip->process || !chip->file_exists) { continue; }
    69         if (!chip->toFPA) { continue; }
    70 
    71         while ((cell = pmFPAviewNextCell (view, fpa, 1)) != NULL) {
     69        if (!chip->toFPA) { continue; }
     70
     71        while ((cell = pmFPAviewNextCell (view, fpa, 1)) != NULL) {
    7272            psTrace ("psastro", 4, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process);
    7373            if (!cell->process || !cell->file_exists) { continue; }
    7474
    75             // process each of the readouts
    76             // XXX there can only be one readout per chip, right?
    77             while ((readout = pmFPAviewNextReadout (view, fpa, 1)) != NULL) {
    78                 if (! readout->data_exists) { continue; }
    79 
    80                 // calculate dMag for the matched stars
    81                 dMag = psastroZeroPointReadoutAccum (dMag, readout, exptime);
    82 
    83                 if (!byExposure) {
    84                     // calculate dMag for the matched stars just for this readout (well, chip)
    85                     psMetadata *header = psMetadataLookupMetadata (&status, readout->analysis, "PSASTRO.HEADER");
    86                     psastroZeroPointAnalysis (header, dMag, zeropt, recipe);
    87                     psFree (dMag);
    88                     dMag = NULL;
    89                 }
    90             }
    91         }
     75            // process each of the readouts
     76            // XXX there can only be one readout per chip, right?
     77            while ((readout = pmFPAviewNextReadout (view, fpa, 1)) != NULL) {
     78                if (! readout->data_exists) { continue; }
     79
     80                // calculate dMag for the matched stars
     81                dMag = psastroZeroPointReadoutAccum (dMag, readout, exptime);
     82
     83                if (!byExposure) {
     84                    // calculate dMag for the matched stars just for this readout (well, chip)
     85                    psMetadata *header = psMetadataLookupMetadata (&status, readout->analysis, "PSASTRO.HEADER");
     86                    psastroZeroPointAnalysis (header, dMag, zeropt, recipe);
     87                    psFree (dMag);
     88                    dMag = NULL;
     89                }
     90            }
     91        }
    9292    }
    9393
    9494    if (byExposure) {
    95         psMetadata *header = psMetadataLookupMetadata (&status, fpa->analysis, "PSASTRO.HEADER");
    96         if (!header) {
    97           header = psMetadataAlloc ();
    98           psMetadataAddMetadata (fpa->analysis, PS_LIST_TAIL, "PSASTRO.HEADER",  PS_META_REPLACE, "psastro header stats", header);
    99           psFree (header);
    100         }
    101         psastroZeroPointAnalysis (header, dMag, zeropt, recipe);
    102         psFree (dMag);
    103         dMag = NULL;
    104     }
     95        psMetadata *header = psMetadataLookupMetadata (&status, fpa->analysis, "PSASTRO.HEADER");
     96        if (!header) {
     97          header = psMetadataAlloc ();
     98          psMetadataAddMetadata (fpa->analysis, PS_LIST_TAIL, "PSASTRO.HEADER",  PS_META_REPLACE, "psastro header stats", header);
     99          psFree (header);
     100        }
     101        psastroZeroPointAnalysis (header, dMag, zeropt, recipe);
     102        psFree (dMag);
     103        dMag = NULL;
     104    }
     105
     106    psMetadataItem *item = psMetadataLookup(fpa->concepts, "FPA.ZP");
     107    item->data.F32 = zeropt;
    105108
    106109    psFree (view);
    … …  
    127130
    128131    if (!dMag) {
    129         dMag  = psVectorAllocEmpty (100, PS_TYPE_F32);
     132        dMag  = psVectorAllocEmpty (100, PS_TYPE_F32);
    130133    }
    131134
    … …  
    164167    psStats *stats = NULL;
    165168    if (useMean) {
    166         // stats structure and mask for use in measuring the clipping statistic
    167         // this analysis has too few data points to use the robust median method
    168         stats = psStatsAlloc (PS_STAT_CLIPPED_MEAN | PS_STAT_CLIPPED_STDEV);
    169         if (!psVectorStats (stats, dMag, NULL, NULL, 0)) {
    170             psError(PS_ERR_UNKNOWN, false, "failure to measure zero point by mean");
    171             return false;
    172         }
     169        // stats structure and mask for use in measuring the clipping statistic
     170        // this analysis has too few data points to use the robust median method
     171        stats = psStatsAlloc (PS_STAT_CLIPPED_MEAN | PS_STAT_CLIPPED_STDEV);
     172        if (!psVectorStats (stats, dMag, NULL, NULL, 0)) {
     173            psError(PS_ERR_UNKNOWN, false, "failure to measure zero point by mean");
     174            return false;
     175        }
    173176    } else {
    174         stats = psastroStatsPercentile (dMag, recipe);
    175         if (!stats) {
    176             psError(PS_ERR_UNKNOWN, false, "failure to measure zero point by edge");
    177             return false;
    178         }
     177        stats = psastroStatsPercentile (dMag, recipe);
     178        if (!stats) {
     179            psError(PS_ERR_UNKNOWN, false, "failure to measure zero point by edge");
     180            return false;
     181        }
    179182    }
    180183    fprintf (stderr, "zero point %f +/- %f using %ld stars; transparency %f\n", stats->clippedMean, stats->clippedStdev, dMag->n, zeropt - stats->clippedMean);
    … …  
    243246    // search for the 'blue' edge of the dMag distribution:
    244247    // the distribution is not a normal population, but instead has a broad range with fairly hard edges.
    245     // construct a histogram and look for the 
     248    // construct a histogram and look for the
    246249
    247250    bool status;
    … …  
    256259    float min = NAN, max = NAN;         // Mimimum and maximum values
    257260
    258     psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS); 
     261    psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS);
    259262
    260263    // Get the minimum and maximum values
    … …  
    267270    // If all data points have the same value, then we set the appropriate members of stats and return.
    268271    if (fabs(max - min) <= FLT_EPSILON) {
    269         stats->clippedMean = min;
    270         stats->clippedStdev = NAN;
    271         stats->results |= PS_STAT_CLIPPED_MEAN;
    272         stats->results |= PS_STAT_CLIPPED_STDEV;
    273         return stats;
    274     }
    275    
    276     // Define the histogram bin size. 
     272        stats->clippedMean = min;
     273        stats->clippedStdev = NAN;
     274        stats->results |= PS_STAT_CLIPPED_MEAN;
     275        stats->results |= PS_STAT_CLIPPED_STDEV;
     276        return stats;
     277    }
     278
     279    // Define the histogram bin size.
    277280    float binSize = MAG_RESOLUTION;
    278281    long numBins = (max - min) / binSize; // Number of bins
    … …  
    293296    float S2 = 0.0;
    294297    for (int i = 0; i < edgeSample; i++) {
    295        
    296         // generate the subset vector
    297         for (long i = 0; i < nSubset; i++) {
    298             double frnd = psRandomUniform(rng);
    299             int entry = PS_MIN(myVector->n - 1, PS_MAX(0, myVector->n * frnd));
    300 
    301             subset->data.F32[i] = myVector->data.F32[entry];
    302         }
    303 
    304         float value = psastroStatsPercentileValue (histogram, cumulative, subset, edgeFraction);
    305 
    306         Sum += value;
    307         S2 += value*value;
     298
     299        // generate the subset vector
     300        for (long i = 0; i < nSubset; i++) {
     301            double frnd = psRandomUniform(rng);
     302            int entry = PS_MIN(myVector->n - 1, PS_MAX(0, myVector->n * frnd));
     303
     304            subset->data.F32[i] = myVector->data.F32[entry];
     305        }
     306
     307        float value = psastroStatsPercentileValue (histogram, cumulative, subset, edgeFraction);
     308
     309        Sum += value;
     310        S2 += value*value;
    308311    }
    309312    psTrace("psastro", 6, "subset stats: Sum: %f, S2: %f, Npts: %d\n", Sum, S2, edgeSample);
    … …  
    345348    psVectorInit (histogram->nums, 0);
    346349    if (!psVectorHistogram(histogram, myVector, NULL, NULL, 0)) {
    347         // if psVectorHistogram returns false, we have a programming error
    348         psAbort ("Unable to generate histogram for psastroZeroPointAnalysis");
     350        // if psVectorHistogram returns false, we have a programming error
     351        psAbort ("Unable to generate histogram for psastroZeroPointAnalysis");
    349352    }
    350353    if (psTraceGetLevel("psastro") >= 8) {
    351         PS_VECTOR_PRINT_F32(histogram->bounds);
    352         PS_VECTOR_PRINT_F32(histogram->nums);
     354        PS_VECTOR_PRINT_F32(histogram->bounds);
     355        PS_VECTOR_PRINT_F32(histogram->nums);
    353356    }
    354357
    … …  
    357360    cumulative->nums->data.F32[0] = histogram->nums->data.F32[0];
    358361    for (long i = 1; i < histogram->nums->n; i++) {
    359         cumulative->nums->data.F32[i] = cumulative->nums->data.F32[i-1] + histogram->nums->data.F32[i];
    360         cumulative->bounds->data.F32[i-1] = histogram->bounds->data.F32[i];
     362        cumulative->nums->data.F32[i] = cumulative->nums->data.F32[i-1] + histogram->nums->data.F32[i];
     363        cumulative->bounds->data.F32[i-1] = histogram->bounds->data.F32[i];
    361364    }
    362365    if (psTraceGetLevel("psastro") >= 8) {
    363         PS_VECTOR_PRINT_F32(cumulative->bounds);
    364         PS_VECTOR_PRINT_F32(cumulative->nums);
     366        PS_VECTOR_PRINT_F32(cumulative->bounds);
     367        PS_VECTOR_PRINT_F32(cumulative->nums);
    365368    }
    366369
    … …  
    403406    if (item->type != PS_DATA_METADATA_MULTI) ESCAPE ("PHOTCODE.DATA not a multi");
    404407
    405     // PHOTCODE.DATA is a multi of metadata items 
     408    // PHOTCODE.DATA is a multi of metadata items
    406409    psListIterator *iter = psListIteratorAlloc(item->data.list, PS_LIST_HEAD, false);
    407410
    408411    psMetadataItem *refItem = NULL;
    409412    while ((refItem = psListGetAndIncrement (iter))) {
    410         if (refItem->type != PS_DATA_METADATA) ESCAPE ("PHOTCODE.DATA entry is not a metadata folder");
    411    
    412         char *refFilter = psMetadataLookupStr (&status, refItem->data.md, "FILTER");
    413         if (!status) {
    414             // psLogMsg ("psastro", PS_LOG_INFO, "a PHOTCODE.DATA recipe folder is missing FILTER");
    415             continue;
    416         }
    417 
    418         // does this entry match the current filter?
    419         if (strcmp (refFilter, filter)) continue;
    420 
    421         psLogMsg ("psastro", PS_LOG_DETAIL, "PHOTCODE.DATA found for filter %s", filter);
    422 
    423         *zeropt = psMetadataLookupF32 (&status, refItem->data.md, "ZEROPT");
    424         if (!status) {
    425             psLogMsg ("psastro", PS_LOG_INFO, "a PHOTCODE.DATA recipe folder is missing ZEROPT");
    426             continue;
    427         }
    428         if (ghostMaxMag) {
    429             *ghostMaxMag = psMetadataLookupF32 (&status, refItem->data.md, "GHOST_MAX_MAG");
    430             if (!status) {
    431                 psLogMsg ("psastro", PS_LOG_INFO, "a PHOTCODE.DATA recipe folder is missing GHOST_MAX_MAG");
    432                 continue;
    433             }
    434         }
    435         psFree (iter);
    436         return true;
     413        if (refItem->type != PS_DATA_METADATA) ESCAPE ("PHOTCODE.DATA entry is not a metadata folder");
     414
     415        char *refFilter = psMetadataLookupStr (&status, refItem->data.md, "FILTER");
     416        if (!status) {
     417            // psLogMsg ("psastro", PS_LOG_INFO, "a PHOTCODE.DATA recipe folder is missing FILTER");
     418            continue;
     419        }
     420
     421        // does this entry match the current filter?
     422        if (strcmp (refFilter, filter)) continue;
     423
     424        psLogMsg ("psastro", PS_LOG_DETAIL, "PHOTCODE.DATA found for filter %s", filter);
     425
     426        *zeropt = psMetadataLookupF32 (&status, refItem->data.md, "ZEROPT");
     427        if (!status) {
     428            psLogMsg ("psastro", PS_LOG_INFO, "a PHOTCODE.DATA recipe folder is missing ZEROPT");
     429            continue;
     430        }
     431        if (ghostMaxMag) {
     432            *ghostMaxMag = psMetadataLookupF32 (&status, refItem->data.md, "GHOST_MAX_MAG");
     433            if (!status) {
     434                psLogMsg ("psastro", PS_LOG_INFO, "a PHOTCODE.DATA recipe folder is missing GHOST_MAX_MAG");
     435                continue;
     436            }
     437        }
     438        psFree (iter);
     439        return true;
    437440    }
    438441    psFree (iter);

  • trunk/pswarp/src/pswarpLoop.c

    r26896 r26964  
    367367    }
    368368
     369    // Update ZP from the astrometry
     370    {
     371        psMetadataItem *item = psMetadataLookup(outFPA->concepts, "FPA.ZP");
     372        item->data.F32 = psMetadataLookupF32(NULL, astrom->fpa->concepts, "FPA.ZP");
     373    }
     374
    369375    pmHDU *hdu = outFPA->hdu;           ///< HDU for the output warped image
    370376
Note: See TracChangeset for help on using the changeset viewer.