IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Oct 2, 2009, 3:12:47 PM (17 years ago)
Author:
eugene
Message:

check in changes from gene@development branch : extensive changes to moments calculation, psf model generation, aperture residuals; improvements to extended source analysis

Location:
trunk/psphot/src
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • trunk/psphot/src

    • Property svn:ignore

      •  

        old new  
        1818psphotVersionDefinitions.h
        1919psphotMomentsStudy
         20psphotPetrosianStudy

  • trunk/psphot/src/psphotRadiusChecks.c

    r21183 r25755  
    77                                        // and a per-object radius is calculated)
    88
     9static float PSF_APERTURE = 0;  // radius to use in PSF aperture mags
     10
     11
    912bool psphotInitRadiusPSF(const psMetadata *recipe, const pmModelType type) {
    1013
    … …  
    1316    PSF_FIT_NSIGMA  = psMetadataLookupF32(&status, recipe, "PSF_FIT_NSIGMA");
    1417    PSF_FIT_PADDING = psMetadataLookupF32(&status, recipe, "PSF_FIT_PADDING");
    15     PSF_FIT_RADIUS =  psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS");
     18    PSF_FIT_RADIUS  =  psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS");
     19    PSF_APERTURE    =  psMetadataLookupF32(&status, recipe, "PSF_APERTURE");
    1620
    1721    return true;
    … …  
    3438            radiusFit = model->modelRadius(model->params, 1.0);
    3539        }
     40        model->fitRadius = (RADIUS_TYPE)(radiusFit + PSF_FIT_PADDING);
     41    } else {
     42        model->fitRadius = radiusFit;
    3643    }
    37     model->radiusFit = (RADIUS_TYPE)(radiusFit + PSF_FIT_PADDING);
    38 
    39     if (isnan(model->radiusFit)) psAbort("error in radius");
     44    if (isnan(model->fitRadius)) psAbort("error in radius");
    4045       
    4146    if (source->mode & PM_SOURCE_MODE_SATSTAR) {
    42         model->radiusFit *= 2;
     47        model->fitRadius *= 2;
    4348    }
    4449
    45     bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);
     50    // radius used to measure aperture photometry
     51    source->apRadius = PSF_APERTURE;
     52
     53    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius);
    4654
    4755    // set the mask to flag the excluded pixels
    48     psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
     56    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius, "OR", markVal);
    4957    return status;
    5058}
    … …  
    5866
    5967    // set the fit radius based on the object flux limit and the model
    60     model->radiusFit = (RADIUS_TYPE) (model->modelRadius (model->params, PSF_FIT_NSIGMA*moments->dSky) + dR + PSF_FIT_PADDING);
    61     if (isnan(model->radiusFit)) psAbort("error in radius");
    62        
     68    float radiusFit = PSF_FIT_RADIUS;
     69    if (radiusFit <= 0) {               // use fixed radius
     70        if (moments == NULL) {
     71            radiusFit = model->modelRadius(model->params, PSF_FIT_NSIGMA*moments->dSky);
     72        } else {
     73            radiusFit = model->modelRadius(model->params, 1.0);
     74        }
     75        model->fitRadius = (RADIUS_TYPE)(radiusFit + PSF_FIT_PADDING);
     76    } else {
     77        model->fitRadius = radiusFit;
     78    }
     79    if (isnan(model->fitRadius)) psAbort("error in radius");
     80
     81    // above sets a radius for a single star, bump by blend separation
     82    model->fitRadius += dR;
     83
    6384    if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
    64         model->radiusFit *= 2;
     85        model->fitRadius *= 2;
    6586    }
    6687
    67     bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);
     88    // radius used to measure aperture photometry
     89    source->apRadius = PSF_APERTURE;
     90
     91    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius);
    6892
    6993    // set the mask to flag the excluded pixels
    70     psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
     94    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius, "OR", markVal);
    7195    return status;
    7296}
    … …  
    7498static float EXT_FIT_NSIGMA;
    7599static float EXT_FIT_PADDING;
     100static float EXT_FIT_MAX_RADIUS;
    76101
    77102bool psphotInitRadiusEXT (psMetadata *recipe, pmModelType type) {
    … …  
    79104    bool status;
    80105
    81     EXT_FIT_NSIGMA   = psMetadataLookupF32 (&status, recipe, "EXT_FIT_NSIGMA");
    82     EXT_FIT_PADDING  = psMetadataLookupF32 (&status, recipe, "EXT_FIT_PADDING");
     106    EXT_FIT_NSIGMA     = psMetadataLookupF32 (&status, recipe, "EXT_FIT_NSIGMA");
     107    EXT_FIT_PADDING    = psMetadataLookupF32 (&status, recipe, "EXT_FIT_PADDING");
     108    EXT_FIT_MAX_RADIUS = psMetadataLookupF32 (&status, recipe, "EXT_FIT_MAX_RADIUS");
    83109
    84110    return true;
    … …  
    86112
    87113// call this function whenever you (re)-define the EXT model
     114float psphotSetRadiusEXT (pmReadout *readout, pmSource *source, psImageMaskType markVal) {
     115
     116    psAssert (source, "source not defined??");
     117    psAssert (source->peak, "peak not defined??");
     118
     119    pmPeak *peak = source->peak;
     120
     121    // set the radius based on the footprint:
     122    if (!peak->footprint) goto escape;
     123    pmFootprint *footprint = peak->footprint;
     124    if (!footprint->spans) goto escape;
     125    if (footprint->spans->n < 1) goto escape;
     126
     127    // find the max radius
     128    float radius = 0.0;
     129    for (int j = 0; j < footprint->spans->n; j++) {
     130        pmSpan *span = footprint->spans->data[j];
     131
     132        float dY  = span->y  - peak->yf;
     133        float dX0 = span->x0 - peak->xf;
     134        float dX1 = span->x1 - peak->xf;
     135
     136        radius = PS_MAX (radius, hypot(dY, dX0));
     137        radius = PS_MAX (radius, hypot(dY, dX1));
     138    }
     139
     140    radius += EXT_FIT_PADDING;
     141    if (isnan(radius)) psAbort("error in radius");
     142
     143    radius = PS_MIN (radius, EXT_FIT_MAX_RADIUS);
     144
     145    // redefine the pixels if needed
     146    pmSourceRedefinePixels (source, readout, peak->xf, peak->yf, radius);
     147
     148    // set the mask to flag the excluded pixels
     149    psImageKeepCircle (source->maskObj, peak->xf, peak->yf, radius, "OR", markVal);
     150    return radius;
     151
     152escape:
     153    return NAN;
     154    // bool result = psphotCheckRadiusEXT (readout, source, model, markVal);
     155    // return result;
     156}
     157
     158// alternative EXT radius based on model guess (for use without footprints)
    88159bool psphotCheckRadiusEXT (pmReadout *readout, pmSource *source, pmModel *model, psImageMaskType markVal) {
     160
     161    psAbort ("do not use this function");
    89162
    90163    psF32 *PAR = model->params->data.F32;
    … …  
    96169    float rawRadius = model->modelRadius (model->params, EXT_FIT_NSIGMA*moments->dSky);
    97170
    98     model->radiusFit = rawRadius + EXT_FIT_PADDING;
    99     if (isnan(model->radiusFit)) psAbort("error in radius");
     171    model->fitRadius = rawRadius + EXT_FIT_PADDING;
     172    if (isnan(model->fitRadius)) psAbort("error in radius");
    100173
    101174    // redefine the pixels if needed
    102     bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit);
     175    bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius);
    103176
    104177    // set the mask to flag the excluded pixels
    105     psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", markVal);
     178    psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->fitRadius, "OR", markVal);
    106179    return status;
    107180}
Note: See TracChangeset for help on using the changeset viewer.