IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Mar 3, 2011, 3:11:09 PM (15 years ago)
Author:
eugene
Message:

fix the calculation of chisq values for constant errors (re-calculated, do not attempt to correct raw value)

Location:
branches/eam_branches/ipp-20110213/psphot
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • branches/eam_branches/ipp-20110213/psphot/doc/stack.txt

    r30779 r30784  
    4949     1652559         2000000    pmSubtractionMatch.c:189
    5050
    51  
    52      2132413         4000000    pmFPAfileDefine.c:1275
     51  chisq:
    5352     2133004         4000000    pmFPACopy.c:71
    5453     2133010         4000000    pmFPACopy.c:71
    5554     2133007         2000000    pmFPACopy.c:71
     55
     56  backmdl: 
     57     2132413         4000000    pmFPAfileDefine.c:1275
     58
     59  ???
    5660     8614548         2000000    psBinaryOp.c:502
    5761     8617313         2000000    psBinaryOp.c:502
     62
     63  pmSource (modelFlx):
     64     2775 * 6k = 16.9M
     65
     66  pmSource (maskObj):
     67     2775 * 3k =  8.4M
     68
     69  (span + footprints account for ~5 - 10 M, rest in little things)
    5870
    597120101221
  • branches/eam_branches/ipp-20110213/psphot/src/psphotFitSourcesLinear.c

    r30764 r30784  
    312312    for (int i = 0; i < fitSources->n; i++) {
    313313        pmSource *source = fitSources->data[i];
    314         if (source->mode & PM_SOURCE_MODE_NONLINEAR_FIT) continue;
    315314        pmModel *model = pmSourceGetModel (NULL, source);
    316         pmSourceChisq (model, source->pixels, source->maskObj, source->variance, maskVal, covarFactor, 1);
     315        if (!(source->mode & PM_SOURCE_MODE_NONLINEAR_FIT)) {
     316            model->nPar = 1; // LINEAR-only sources have 1 parameter; NONLINEAR sources have their original value
     317        }
     318        pmSourceChisq (model, source->pixels, source->maskObj, source->variance, maskVal);
    317319    }
    318320    psLogMsg ("psphot.ensemble", PS_LOG_MINUTIA, "get chisqs: %f sec (%d elements)\n", psTimerMark ("psphot.linear"), sparse->Nelem);
  • branches/eam_branches/ipp-20110213/psphot/src/psphotFitSourcesLinearStack.c

    r30764 r30784  
    163163    for (int i = 0; i < fitSources->n; i++) {
    164164        pmSource *source = fitSources->data[i];
    165         if (source->mode & PM_SOURCE_MODE_NONLINEAR_FIT) continue;
    166165        pmModel *model = pmSourceGetModel (NULL, source);
    167         pmSourceChisq (model, source->pixels, source->maskObj, source->variance, maskVal, COVAR_FACTOR, 1);
     166        if (!(source->mode & PM_SOURCE_MODE_NONLINEAR_FIT)) {
     167            model->nPar = 1; // LINEAR-only sources have 1 parameter; NONLINEAR sources have their original value
     168        }
     169        pmSourceChisq (model, source->pixels, source->maskObj, source->variance, maskVal);
    168170    }
    169171    psLogMsg ("psphot.ensemble", PS_LOG_MINUTIA, "get chisqs: %f sec (%d elements)\n", psTimerMark ("psphot.linear"), sparse->Nelem);
  • branches/eam_branches/ipp-20110213/psphot/src/psphotReadout.c

    r30749 r30784  
    99}
    1010
     11// for now, let's store the detections on the readout->analysis for each readout
     12bool psphotDumpChisqs (pmConfig *config, const pmFPAview *view, const char *filerule)
     13{
     14    static int npass = 0;
     15    char filename[64];
     16
     17    bool status = true;
     18
     19    int num = psphotFileruleCount(config, filerule);
     20
     21    snprintf (filename, 64, "chisq.%02d.dat", npass);
     22    FILE *f = fopen (filename, "w");
     23
     24    // loop over the available readouts
     25    for (int i = 0; i < num; i++) {
     26
     27        // find the currently selected readout
     28        pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, i); // File of interest
     29        psAssert (file, "missing file?");
     30
     31        pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
     32        psAssert (readout, "missing readout?");
     33
     34        pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
     35        psAssert (detections, "missing detections?");
     36
     37        psArray *sources = detections->allSources;
     38        psAssert (sources, "missing sources?");
     39
     40        for (int i = 0; i < sources->n; i++) {
     41            pmSource *source = sources->data[i];
     42            if (!source) continue;
     43
     44            pmModel *model = pmSourceGetModel (NULL, source);
     45            if (!model) continue;
     46       
     47            if (source->mode & PM_SOURCE_MODE_NONLINEAR_FIT) {
     48                fprintf (f, "%f %f %f %d %d %f  1 NONLINEAR\n", model->mag, model->params->data.F32[1], model->chisq, model->nDOF, model->nPix, model->chisqNorm);
     49            } else {
     50                fprintf (f, "%f %f %f %d %d %f  0 LINEAR\n", model->mag, model->params->data.F32[1], model->chisq, model->nDOF, model->nPix, model->chisqNorm);
     51            }
     52        }
     53    }
     54    fclose (f);
     55    npass ++;
     56
     57    return true;
     58}
     59
    1160bool psphotReadout(pmConfig *config, const pmFPAview *view, const char *filerule) {
    1261
     
    147196    // linear PSF fit to source peaks, subtract the models from the image (in PSF mask)
    148197    psphotFitSourcesLinear (config, view, filerule, false); // pass 1 (detections->allSources)
     198    psphotDumpChisqs (config, view, filerule);
    149199
    150200    // identify CRs and extended sources (only unmeasured sources are measured)
     
    157207    // replace model flux, adjust mask as needed, fit, subtract the models (full stamp)
    158208    psphotBlendFit (config, view, filerule); // pass 1 (detections->allSources)
     209    psphotDumpChisqs (config, view, filerule);
    159210
    160211    // replace all sources
     
    164215    // NOTE : apply to ALL sources (extended + psf)
    165216    psphotFitSourcesLinear (config, view, filerule, true); // pass 2 (detections->allSources)
     217    psphotDumpChisqs (config, view, filerule);
    166218
    167219    // if we only do one pass, skip to extended source analysis
     
    209261        // NOTE: apply to ALL sources
    210262        psphotFitSourcesLinear (config, view, filerule, true); // pass 3 (detections->allSources)
     263        psphotDumpChisqs (config, view, filerule);
    211264    }
    212265
  • branches/eam_branches/ipp-20110213/psphot/src/psphotSourceFits.c

    r30624 r30784  
    101101    if (!isfinite(PSF->params->data.F32[PM_PAR_I0])) psAbort("nan in fit");
    102102
    103     // correct model chisq for flux trend
    104     double chiTrend = psPolynomial1DEval (psf->ChiTrend, PSF->params->data.F32[PM_PAR_I0]);
    105     PSF->chisqNorm = PSF->chisq / chiTrend;
    106 
    107103    // evaluate the blend objects, subtract if good, free otherwise
    108104    for (int i = 1; i < modelSet->n; i++) {
     
    111107
    112108        if (!isfinite(model->params->data.F32[PM_PAR_I0])) psAbort("nan in fit");
    113 
    114         // correct model chisq for flux trend
    115         chiTrend = psPolynomial1DEval (psf->ChiTrend, model->params->data.F32[PM_PAR_I0]);
    116         model->chisqNorm = model->chisq / chiTrend;
    117109
    118110        // if this one failed, skip it
     
    159151bool psphotFitPSF (pmReadout *readout, pmSource *source, pmPSF *psf, pmSourceFitOptions *fitOptions, psImageMaskType maskVal, psImageMaskType markVal) {
    160152
    161     double chiTrend;
    162153    pmSourceFitOptions options = *fitOptions;
    163154
     
    182173    // clear the circular mask
    183174    psImageMaskPixels (source->maskObj, "AND", PS_NOT_IMAGE_MASK(markVal));
    184 
    185     // correct model chisq for flux trend
    186     chiTrend = psPolynomial1DEval (psf->ChiTrend, PSF->params->data.F32[PM_PAR_I0]);
    187     PSF->chisqNorm = PSF->chisq / chiTrend;
    188175
    189176    // does the PSF model succeed?
     
    225212    bool okEXT, okDBL;
    226213    float chiEXT, chiDBL;
    227     double chiTrend;
    228214    pmModel *ONE = NULL;
    229215    pmSource *tmpSrc = NULL;
     
    271257
    272258        // correct first model chisqs for flux trend
    273         chiDBL = NAN;
    274259        ONE = DBL->data[0];
    275260        if (ONE) {
    276261            if (!isfinite(ONE->params->data.F32[PM_PAR_I0])) psAbort("nan in fit");
    277             chiTrend = psPolynomial1DEval (psf->ChiTrend, ONE->params->data.F32[1]);
    278             ONE->chisqNorm = ONE->chisq / chiTrend;
    279             chiDBL = ONE->chisq / ONE->nDOF; // save chisq for double-star/galaxy comparison
     262            chiDBL = ONE->chisqNorm; // save chisq for double-star/galaxy comparison
    280263            ONE->fitRadius = radius;
    281264        }
     
    285268        if (ONE) {
    286269            if (!isfinite(ONE->params->data.F32[PM_PAR_I0])) psAbort("nan in fit");
    287             chiTrend = psPolynomial1DEval (psf->ChiTrend, ONE->params->data.F32[1]);
    288             ONE->chisqNorm = ONE->chisq / chiTrend;
    289270            ONE->fitRadius = radius;
    290271        }
     
    298279
    299280        okEXT = psphotEvalEXT (tmpSrc, EXT);
    300         chiEXT = EXT ? EXT->chisq / EXT->nDOF : NAN;
     281        chiEXT = EXT ? EXT->chisqNorm : NAN;
    301282    }
    302283
Note: See TracChangeset for help on using the changeset viewer.