IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 24876


Ignore:
Timestamp:
Jul 21, 2009, 10:46:14 AM (17 years ago)
Author:
eugene
Message:

drop deprecated unthreaded code (had been ifdefed away); replace pmModelFromPSF calls with pmModelFromPSFforXY; drop wildGuess, replace with better logic on moments vs peaks

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/psphot/src/psphotGuessModels.c

    r21519 r24876  
    66// 2) loop over the sources once and associate them with their cell
    77// 3) define the threaded function to work with sources for a given cell
    8 
    9 // A guess for when the moments aren't available
    10 static pmModel *wildGuess(pmSource *source, // Source for which to guess
    11                           pmPSF *psf    // The point-spread function
    12     )
    13 {
    14     pmModel *model = pmModelAlloc(psf->type);
    15     psF32 *PAR = model->params->data.F32;
    16     PAR[PM_PAR_SKY]  = 0;
    17     // XXX get this from the image pixels
    18     PAR[PM_PAR_I0]   = source->peak->flux;
    19     PAR[PM_PAR_XPOS] = source->peak->xf;
    20     PAR[PM_PAR_YPOS] = source->peak->yf;
    21     return model;
    22 }
    238
    249// construct an initial PSF model for each object
     
    8166            }
    8267            psFree(job);
    83 
    84 # if (0)               
    85                 if (!psphotGuessModel_Unthreaded (readout, cells->data[j], psf, maskVal, markVal)) {
    86                     psError(PS_ERR_UNKNOWN, false, "Unable to guess model.");
    87                     return false;
    88                 }
    89 # endif
    9068        }
    9169
     
    152130        nSrc ++;
    153131       
    154         // XXX if a source is faint, it will not have moments measured.
    155         // it must be modelled as a PSF.  In this case, we need to use
    156         // the peak centroid to get the coordinates and get the peak flux
    157         // from the image?
    158         pmModel *modelEXT;
    159         if (!source->moments) {
    160             modelEXT = wildGuess(source, psf);
     132        // the guess central intensity comes from the peak:
     133        float Io = source->peak->flux;
     134
     135        // We have two options to get a guess for the object position: the position from the
     136        // peak and the position from the moments.  Use the peak position if (a) there are no
     137        // moments and (b) the sources is not saturated
     138
     139        bool useMoments = false;
     140        useMoments = (source->mode & PM_SOURCE_MODE_SATSTAR);  // we only want to try if SATSTAR is set, but..
     141        useMoments = (useMoments && source->moments);          // can't if there are no moments
     142        useMoments = (useMoments && source->moments->nPixels); // can't if the moments were not measured
     143        useMoments = (useMoments && !(source->mode && PM_SOURCE_MODE_MOMENTS_FAILURE)); // can't if the moments failed...
     144
     145        float Xo, Yo;
     146        if (useMoments) {
     147            Xo = source->moments->Mx;
     148            Yo = source->moments->My;
    161149        } else {
    162             // use the source moments, etc to guess basic model parameters
    163             modelEXT = pmSourceModelGuess (source, psf->type); // ALLOC X5
    164             if (!modelEXT) {
    165                 modelEXT = wildGuess(source, psf);
    166             }
    167             // these valuse are set in pmSourceModelGuess, should this rule be in there as well?
    168             if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
    169                 modelEXT->params->data.F32[PM_PAR_XPOS] = source->moments->Mx;
    170                 modelEXT->params->data.F32[PM_PAR_YPOS] = source->moments->My;
    171             } else {
    172                 modelEXT->params->data.F32[PM_PAR_XPOS] = source->peak->xf;
    173                 modelEXT->params->data.F32[PM_PAR_YPOS] = source->peak->yf;
    174             }
     150            Xo = source->peak->xf;
     151            Yo = source->peak->yf;
    175152        }
    176153
    177         // set PSF parameters for this model (apply 2D shape model)
    178         pmModel *modelPSF = pmModelFromPSF (modelEXT, psf); // ALLOC X5
     154        // set PSF parameters for this model (apply 2D shape model to coordinates Xo, Yo)
     155        pmModel *modelPSF = pmModelFromPSFforXY(psf, Xo, Yo, Io);
     156
    179157        if (modelPSF == NULL) {
    180             psWarning ("Failed to determine PSF model at r,c = (%d,%d); trying centre of image",
    181                     source->peak->y, source->peak->x);
     158            psWarning ("Failed to determine PSF model at (%f,%f); trying image center", Xo, Yo);
    182159
    183             // Try the center of the image
    184             modelEXT->params->data.F32[PM_PAR_XPOS] = 0.5*readout->image->numCols;
    185             modelEXT->params->data.F32[PM_PAR_YPOS] = 0.5*readout->image->numRows;
    186             modelPSF = pmModelFromPSF (modelEXT, psf);
     160            float Xc = 0.5*readout->image->numCols;
     161            float Yc = 0.5*readout->image->numRows;
     162            pmModel *modelPSF = pmModelFromPSFforXY(psf, Xc, Yc, Io);
    187163            if (modelPSF == NULL) {
    188164                psError(PSPHOT_ERR_PSF, false, "Failed to determine PSF model at center of image");
    189                 psFree(modelEXT);
    190165                return false;
    191166            }
     167
     168            // Now set the object position at the expected location:
     169            modelPSF->params->data.F32[PM_PAR_XPOS] = Xo;
     170            modelPSF->params->data.F32[PM_PAR_YPOS] = Yo;
    192171            source->mode |= PM_SOURCE_MODE_BADPSF;
    193172        }
    194         psFree (modelEXT); // FREE (x3)
    195173
    196         // XXX need to define the guess flux?
    197174        // set the fit radius based on the object flux limit and the model
    198175        // this function affects the mask pixels
     
    209186    return true;
    210187}
    211 
    212 # if (0)
    213 // construct models only for sources in the specified region
    214 bool psphotGuessModel_Unthreaded (pmReadout *readout, psArray *sources, pmPSF *psf, psImageMaskType maskVal, psImageMaskType markVal) {
    215 
    216     int nSrc = 0;
    217 
    218     for (int i = 0; i < sources->n; i++) {
    219         pmSource *source = sources->data[i];
    220 
    221         // XXXX this is just for a test: use this to mark sources for which the model is measured
    222         // check later that all are used.
    223         source->mode |= PM_SOURCE_MODE_EXT_LIMIT;
    224 
    225         // skip non-astronomical objects (very likely defects)
    226         if (source->type == PM_SOURCE_TYPE_DEFECT) continue;
    227         if (source->type == PM_SOURCE_TYPE_SATURATED) continue;
    228         if (!source->peak) continue;
    229 
    230         nSrc ++;
    231        
    232         // XXX if a source is faint, it will not have moments measured.
    233         // it must be modelled as a PSF.  In this case, we need to use
    234         // the peak centroid to get the coordinates and get the peak flux
    235         // from the image?
    236         pmModel *modelEXT;
    237         if (!source->moments) {
    238             modelEXT = wildGuess(source, psf);
    239         } else {
    240             // use the source moments, etc to guess basic model parameters
    241             modelEXT = pmSourceModelGuess (source, psf->type); // ALLOC
    242             if (!modelEXT) {
    243                 modelEXT = wildGuess(source, psf);
    244             }
    245             // these valuse are set in pmSourceModelGuess, should this rule be in there as well?
    246             if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
    247                 modelEXT->params->data.F32[PM_PAR_XPOS] = source->moments->Mx;
    248                 modelEXT->params->data.F32[PM_PAR_YPOS] = source->moments->My;
    249             } else {
    250                 modelEXT->params->data.F32[PM_PAR_XPOS] = source->peak->xf;
    251                 modelEXT->params->data.F32[PM_PAR_YPOS] = source->peak->yf;
    252             }
    253         }
    254 
    255         // set PSF parameters for this model (apply 2D shape model)
    256         pmModel *modelPSF = pmModelFromPSF (modelEXT, psf); // ALLOC
    257         if (modelPSF == NULL) {
    258             psError(PSPHOT_ERR_PSF, false,
    259                     "Failed to determine PSF model at r,c = (%d,%d); trying centre of image",
    260                     source->peak->y, source->peak->x);
    261             //
    262             // Try the centre of the image
    263             //
    264             modelEXT->params->data.F32[PM_PAR_XPOS] = 0.5*readout->image->numCols;
    265             modelEXT->params->data.F32[PM_PAR_YPOS] = 0.5*readout->image->numRows;
    266             modelPSF = pmModelFromPSF (modelEXT, psf);
    267             if (modelPSF == NULL) {
    268                 psError(PSPHOT_ERR_PSF, false,
    269                         "Failed to determine PSF model at centre of image");
    270                 psFree(modelEXT);
    271                 return false;
    272             }
    273 
    274             source->mode |= PM_SOURCE_MODE_BADPSF;
    275         }
    276         psFree (modelEXT);
    277 
    278         // XXX need to define the guess flux?
    279         // set the fit radius based on the object flux limit and the model
    280         // this function affects the mask pixels
    281         psphotCheckRadiusPSF (readout, source, modelPSF, markVal);
    282 
    283         // set the source PSF model
    284         source->modelPSF = modelPSF;
    285         source->modelPSF->residuals = psf->residuals;
    286 
    287         pmSourceCacheModel (source, maskVal);
    288 
    289     }
    290 
    291     return true;
    292 }
    293 # endif
Note: See TracChangeset for help on using the changeset viewer.