IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Jun 19, 2007, 4:22:26 PM (19 years ago)
Author:
Paul Price
Message:

Extensive changes to APIs to allow use of a nominated value to mask
against (the maskVal). Previously, the mask values were either
hard-coded (e.g., PM_MASK_SAT) or taken as anything non-zero. The
code is tested under psphot (which has similar changes) and does not
crash, but neither is it successful in marking all bad pixels (EAM
will investigate). For this reason, I have left the "gutter" pixels
(cell gaps) set to 0 instead of NAN in pmFPAMosaic.

File:
1 edited

Legend:

Unmodified
Added
Removed

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

    r13803 r13898  
    55 *  @author EAM, IfA
    66 *
    7  *  @version $Revision: 1.42 $ $Name: not supported by cvs2svn $
    8  *  @date $Date: 2007-06-13 23:41:51 $
     7 *  @version $Revision: 1.43 $ $Name: not supported by cvs2svn $
     8 *  @date $Date: 2007-06-20 02:22:26 $
    99 *
    1010 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    9696
    9797// generate a pmPSFtry with a copy of the test PSF sources
    98 pmPSFtry *pmPSFtryModel (psArray *sources, char *modelName, float RADIUS, bool poissonErrors, psPolynomial2D *psfTrendMask, bool applyWeights)
     98pmPSFtry *pmPSFtryModel (psArray *sources, char *modelName, float RADIUS, bool poissonErrors, psPolynomial2D *psfTrendMask, bool applyWeights, psMaskType maskVal, psMaskType mark)
    9999{
    100100    bool status;
    … …  
    116116        if (source->modelEXT == NULL) {
    117117            psError(PS_ERR_UNKNOWN, false, "failed to build model");
    118             psFree(psfTry);
     118            psFree(psfTry);
    119119            return NULL;
    120120        }
    121121
    122122        // set object mask to define valid pixels
    123         psImageKeepCircle (source->maskObj, source->peak->x, source->peak->y, RADIUS, "OR", PM_MASK_MARK);
     123        psImageKeepCircle (source->maskObj, source->peak->x, source->peak->y, RADIUS, "OR", mark);
    124124
    125125        // fit model as EXT, not PSF
    126         status = pmSourceFitModel (source, source->modelEXT, PM_SOURCE_FIT_EXT);
     126        status = pmSourceFitModel (source, source->modelEXT, PM_SOURCE_FIT_EXT, maskVal);
    127127
    128128        // exclude the poor fits
    … …  
    139139    if (!pmPSFFromPSFtry (psfTry, applyWeights)) {
    140140        psError(PS_ERR_UNKNOWN, false, "failed to construct a psf model from collection of sources");
    141         psFree(psfTry);
     141        psFree(psfTry);
    142142        return NULL;
    143143    }
    … …  
    154154        // set shape for this model based on PSF
    155155        source->modelPSF = pmModelFromPSF (source->modelEXT, psfTry->psf);
    156         if (source->modelPSF == NULL) {
     156        if (source->modelPSF == NULL) {
    157157            psfTry->mask->data.U8[i] = PSFTRY_MASK_BAD_MODEL;
    158             abort();
    159             continue;
    160         }
     158            abort();
     159            continue;
     160        }
    161161        source->modelPSF->radiusFit = RADIUS;
    162162
    163163        // set object mask to define valid pixels
    164         psImageKeepCircle (source->maskObj, source->peak->x, source->peak->y, RADIUS, "OR", PM_MASK_MARK);
     164        psImageKeepCircle (source->maskObj, source->peak->x, source->peak->y, RADIUS, "OR", mark);
    165165
    166166        // fit the PSF model to the source
    167         status = pmSourceFitModel (source, source->modelPSF, PM_SOURCE_FIT_PSF);
     167        status = pmSourceFitModel (source, source->modelPSF, PM_SOURCE_FIT_PSF, maskVal);
    168168
    169169        // skip poor fits
    170170        if (!status) {
    171171            psfTry->mask->data.U8[i] = PSFTRY_MASK_PSF_FAIL;
    172             psTrace ("psModules.pmPSFtry", 4, "dropping %d (%d,%d) : failed PSF fit\n", i, source->peak->x, source->peak->y);
     172            psTrace ("psModules.pmPSFtry", 4, "dropping %d (%d,%d) : failed PSF fit\n", i, source->peak->x, source->peak->y);
    173173            continue;
    174174        }
    175175
    176176        // XXX : use a different aperture radius from the fit radius?
    177         if (!pmSourceMagnitudes (source, psfTry->psf, PM_SOURCE_PHOT_INTERP)) {
     177        if (!pmSourceMagnitudes (source, psfTry->psf, PM_SOURCE_PHOT_INTERP, maskVal, mark)) {
    178178            psfTry->mask->data.U8[i] = PSFTRY_MASK_BAD_PHOT;
    179             psTrace ("psModules.pmPSFtry", 4, "dropping %d (%d,%d) : poor photometry\n", i, source->peak->x, source->peak->y);
     179            psTrace ("psModules.pmPSFtry", 4, "dropping %d (%d,%d) : poor photometry\n", i, source->peak->x, source->peak->y);
    180180            continue;
    181181        }
    … …  
    203203            flux->data.F64[i] = 0.0;
    204204            chisq->data.F64[i] = 0.0;
    205             mask->data.U8[i] = 1;
     205            mask->data.U8[i] = 0xff;
    206206        } else {
    207207            flux->data.F64[i] = source->modelPSF->params->data.F32[PM_PAR_I0];
    … …  
    217217
    218218    // linear clipped fit of chisq trend vs flux
    219     bool result = psVectorClipFitPolynomial1D(psfTry->psf->ChiTrend, stats, mask, 1, chisq, NULL, flux);
     219    bool result = psVectorClipFitPolynomial1D(psfTry->psf->ChiTrend, stats, mask, 0xff, chisq, NULL, flux);
    220220    psLogMsg ("pmPSFtry", 4, "chisq vs flux fit: %f +/- %f\n", stats->robustMedian, stats->robustStdev);
    221221
    … …  
    289289        psError(PS_ERR_UNKNOWN, false, "Failed to fit clipped poly");
    290290
    291         psFree(poly);
    292         psFree(r2rflux);
    293         psFree(stats);
     291        psFree(poly);
     292        psFree(r2rflux);
     293        psFree(stats);
    294294
    295295        return false;
    … …  
    373373        y->data.F64[i] = source->modelEXT->params->data.F32[PM_PAR_YPOS];
    374374
    375         // weight by the error on the source flux
    376         if (dz) {
    377             dz->data.F64[i] = source->modelEXT->dparams->data.F32[PM_PAR_I0];
    378         }
     375        // weight by the error on the source flux
     376        if (dz) {
     377            dz->data.F64[i] = source->modelEXT->dparams->data.F32[PM_PAR_I0];
     378        }
    379379    }
    380380
    … …  
    399399    psVector *e2   = psVectorAlloc (psfTry->sources->n, PS_TYPE_F64);
    400400    for (int i = 0; i < psfTry->sources->n; i++) {
    401         pmSource *source = psfTry->sources->data[i];
    402         if (source->modelEXT == NULL) continue;
    403 
    404         psEllipsePol pol = pmPSF_ModelToFit (source->modelEXT->params->data.F32);
    405 
    406         e0->data.F64[i] = pol.e0;
    407         e1->data.F64[i] = pol.e1;
    408         e2->data.F64[i] = pol.e2;
     401        pmSource *source = psfTry->sources->data[i];
     402        if (source->modelEXT == NULL) continue;
     403
     404        psEllipsePol pol = pmPSF_ModelToFit (source->modelEXT->params->data.F32);
     405
     406        e0->data.F64[i] = pol.e0;
     407        e1->data.F64[i] = pol.e1;
     408        e2->data.F64[i] = pol.e2;
    409409    }
    410410
    … …  
    413413    for (int i = 0; i < stats->clipIter; i++) {
    414414        psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E0], stats, psfTry->mask, 0xff, e0, dz, x, y);
    415         psTrace ("psModules.pmPSFtry", 4, "clipped E0 : keeping %ld of %ld\n", stats->clippedNvalues, e0->n);
     415        psTrace ("psModules.pmPSFtry", 4, "clipped E0 : keeping %ld of %ld\n", stats->clippedNvalues, e0->n);
    416416        psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E1], stats, psfTry->mask, 0xff, e1, dz, x, y);
    417         psTrace ("psModules.pmPSFtry", 4, "clipped E1 : keeping %ld of %ld\n", stats->clippedNvalues, e1->n);
     417        psTrace ("psModules.pmPSFtry", 4, "clipped E1 : keeping %ld of %ld\n", stats->clippedNvalues, e1->n);
    418418        psVectorClipFitPolynomial2D (psf->params_NEW->data[PM_PAR_E2], stats, psfTry->mask, 0xff, e2, dz, x, y);
    419         psTrace ("psModules.pmPSFtry", 4, "clipped E2 : keeping %ld of %ld\n", stats->clippedNvalues, e2->n);
     419        psTrace ("psModules.pmPSFtry", 4, "clipped E2 : keeping %ld of %ld\n", stats->clippedNvalues, e2->n);
    420420    }
    421421
    422422    // XXX temporary dump of the psf parameters
    423423    if (psTraceGetLevel("psModules.objects") >= 4) {
    424         FILE *f = fopen ("pol.dat", "w");
    425         for (int i = 0; i < e0->n; i++) {
    426             fprintf (f, "%f %f  :  %f %f %f  : %d\n",
    427                      x->data.F64[i], y->data.F64[i],
    428                      e0->data.F64[i], e1->data.F64[i], e2->data.F64[i], psfTry->mask->data.U8[i]);
    429         }
    430         fclose (f);
     424        FILE *f = fopen ("pol.dat", "w");
     425        for (int i = 0; i < e0->n; i++) {
     426            fprintf (f, "%f %f  :  %f %f %f  : %d\n",
     427                     x->data.F64[i], y->data.F64[i],
     428                     e0->data.F64[i], e1->data.F64[i], e2->data.F64[i], psfTry->mask->data.U8[i]);
     429        }
     430        fclose (f);
    431431    }
    432432
    … …  
    437437    // skip the unfitted parameters (X, Y, Io, Sky) and the shape parameters (SXX, SYY, SXY)
    438438    for (int i = 0; i < psf->params_NEW->n; i++) {
    439         switch (i) {
    440           case PM_PAR_SKY:
    441           case PM_PAR_I0:
    442           case PM_PAR_XPOS:
    443           case PM_PAR_YPOS:           
    444           case PM_PAR_SXX:           
    445           case PM_PAR_SYY:           
    446           case PM_PAR_SXY:           
    447             continue;
    448           default:
    449             break;
    450         }
     439        switch (i) {
     440          case PM_PAR_SKY:
     441          case PM_PAR_I0:
     442          case PM_PAR_XPOS:
     443          case PM_PAR_YPOS:
     444          case PM_PAR_SXX:
     445          case PM_PAR_SYY:
     446          case PM_PAR_SXY:
     447            continue;
     448          default:
     449            break;
     450        }
    451451
    452452        // select the per-object fitted data for this parameter
    … …  
    463463        if (!psVectorClipFitPolynomial2D(psf->params_NEW->data[i], stats, psfTry->mask, 0xff, z, NULL, x, y)) {
    464464            psError(PS_ERR_UNKNOWN, false, "failed to build psf model for parameter %d", i);
    465             psFree(stats);
    466             psFree(x);
    467             psFree(y);
    468             psFree(z);
    469             psFree(dz);
     465            psFree(stats);
     466            psFree(x);
     467            psFree(y);
     468            psFree(z);
     469            psFree(dz);
    470470            return false;
    471471        }
Note: See TracChangeset for help on using the changeset viewer.