IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 6753


Ignore:
Timestamp:
Mar 31, 2006, 4:58:31 PM (20 years ago)
Author:
eugene
Message:

fixed up low-level issues with sky levels, sky bias

Location:
trunk/psphot/src
Files:
10 edited

Legend:

Unmodified
Added
Removed
  • trunk/psphot/src/psphot.h

    r6727 r6753  
    44# include <pslib.h>
    55# include <psmodules.h>
    6 
    7 # if (0)
    8 # include <psLine.h>
    9 # include <psSparse.h>
    10 # include <psEllipse.h>
    11 # include <psAdditionals.h>
    12 # include <pmConfig.h>
    13 # include <pmHDU.h>
    14 # include <pmHDUUtils.h>
    15 # include <pmFPA.h>
    16 # include <pmFPAview.h>
    17 # include <pmFPAfile.h>
    18 # include <pmFPARead.h>
    19 # include <pmFPAConstruct.h>
    20 # include <pmMaskBadPixels.h>
    21 # include <pmConcepts.h>
    22 
    23 # include <pmPeaks.h>
    24 # include <pmMoments.h>
    25 # include <pmModel.h>
    26 # include <pmSource.h>
    27 # include <pmSourceIO.h>
    28 # include <pmSourceSky.h>
    29 # include <pmSourceFitModel.h>
    30 # include <pmSourceFitSet.h>
    31 # include <pmSourceContour.h>
    32 # include <pmGrowthCurve.h>
    33 # include <pmPSF.h>
    34 # include <pmPSFtry.h>
    35 # include <pmModelGroup.h>
    36 # include <pmSourcePhotometry.h>
    37 # endif
    386
    397# define PSPHOT_RECIPE "PSPHOT"
  • trunk/psphot/src/psphotApResid.c

    r6715 r6753  
    1717    // measure the aperture loss as a function of radius for PSF
    1818    float REF_RADIUS = psMetadataLookupF32 (&status, config, "PSF_REF_RADIUS");
    19     psf->growth = pmGrowthCurveAlloc (3.0, REF_RADIUS, 0.1);
     19    psf->growth = pmGrowthCurveAlloc (3.0, 50.0, REF_RADIUS);
    2020    psphotGrowthCurve (readout, psf);
    2121   
     
    4343        if (source->moments->SN < 2*FIT_SN_LIM) continue;
    4444
    45         // get uncorrected magnitudes in scaled apertures
    46         model = pmSourceMagnitudes (source, NULL, 0);
     45        // get growth-corrected, apTrend-uncorrected magnitudes in scaled apertures
     46        model = pmSourceMagnitudes (source, psf, true, false);
    4747        if (model == NULL) continue;
    4848
     49        apResid->data.F64[Npsf] = source->apMag - source->psfMag;
     50
     51        xPos->data.F64[Npsf]    = model->params->data.F32[2];
     52        yPos->data.F64[Npsf]    = model->params->data.F32[3];
     53
     54        flux->data.F64[Npsf]    = pow(10.0, -0.4*source->psfMag);
     55        r2rflux->data.F64[Npsf] = PS_SQR(model->radiusTMP) / flux->data.F64[Npsf];
     56
    4957        mask->data.U8[Npsf] = 0;
    50         xPos->data.F64[Npsf] = model->params->data.F32[2];
    51         yPos->data.F64[Npsf] = model->params->data.F32[3];
    52 
    53         flux->data.F64[Npsf] = pow(10.0, -0.4*source->psfMag);
    54         r2rflux->data.F64[Npsf] = PS_SQR(model->radius) / flux->data.F64[Npsf];
    55        
    56         apResid->data.F64[Npsf] = source->apMag + pmGrowthCurveCorrect (psf->growth, model->radius) - source->psfMag;
    5758
    5859        // XXX sanity clip?
     
    8384    stats->max = 3.0;
    8485
     86    // no correction
     87    if (!strcasecmp (ApTrendOption, "NONE")) {
     88        // remove ApTrend fit from pmPSFtry
     89        psf->ApTrend->coeff[0][0][0][0] = 0;
     90    }
     91
    8592    // constant only
    8693    if (!strcasecmp (ApTrendOption, "CONSTANT")) {
     
    202209    psVector *resid = (psVector *) psBinaryOp (NULL, (void *) apResid, "-", (void *) fit);
    203210
    204 # if (0)
     211# if (1)
    205212    FILE *fout = fopen ("resid.dat", "w");
    206213    for (int i = 0; i < resid->n; i++) {
     
    241248
    242249    psLogMsg ("psphot.apresid", 3, "measure full-frame aperture residual: %f sec\n", psTimerMark ("psphot"));
    243     psLogMsg ("psphot.apresid", 4, "aperture residual: %f +/- %f : %f bias\n", psf->ApResid, psf->dApResid, psf->skyBias);
     250    psLogMsg ("psphot.apresid", 4, "aperture residual: %f +/- %f : %f bias, %f skysat\n", psf->ApResid, psf->dApResid, psf->skyBias, psf->skySat);
    244251    psLogMsg ("psphot.apresid", 4, "apresid trends: %f %f %f %f %f\n",
    245252              1e3*psf->ApTrend->coeff[1][0][0][0],
  • trunk/psphot/src/psphotBlendFit.c

    r6715 r6753  
    3131        if (source->type == PM_SOURCE_TYPE_DEFECT) continue;
    3232        if (source->type == PM_SOURCE_TYPE_SATURATED) continue;
     33
     34        // skip DBL second sources (ie, added by psphotFitBlob)
     35        if (source->mode &  PM_SOURCE_MODE_PAIR) continue;
    3336
    3437        // limit selection to some SN limit
  • trunk/psphot/src/psphotEnsemblePSF.c

    r6715 r6753  
    118118        // fill in the model pixel values
    119119        psImageInit (flux, 0.0);
    120         psImageKeepCircle (mask, x, y, model->radius, "OR", PM_SOURCE_MASK_MARKED);
     120        psImageKeepCircle (mask, x, y, model->radiusTMP, "OR", PM_SOURCE_MASK_MARKED);
    121121        pmModelAdd (flux, mask, model, false, false);
    122122
     
    209209        y = model->params->data.F32[3];
    210210
    211         psImageKeepCircle (Fi->mask, x, y, model->radius, "OR", PM_SOURCE_MASK_MARKED);
     211        psImageKeepCircle (Fi->mask, x, y, model->radiusTMP, "OR", PM_SOURCE_MASK_MARKED);
    212212        pmSourceChisq (model, Fi->pixels, Fi->mask, Fi->weight);
    213         psImageKeepCircle (Fi->mask, x, y, model->radius, "AND", ~PM_SOURCE_MASK_MARKED);
     213        psImageKeepCircle (Fi->mask, x, y, model->radiusTMP, "AND", ~PM_SOURCE_MASK_MARKED);
    214214    }
    215215
  • trunk/psphot/src/psphotGrowthCurve.c

    r6715 r6753  
    5252
    5353        psImageKeepCircle (mask, xc, yc, radius, "AND", ~PM_SOURCE_MASK_MARKED);
    54 
    5554        psf->growth->apMag->data.F32[i] = apMag;
    5655    }
  • trunk/psphot/src/psphotImageMedian.c

    r6715 r6753  
    211211    psLogMsg ("psphot", 3, "build resampled image: %f sec\n", psTimerMark ("psphot"));
    212212
     213    psphotSaveImage (NULL, background, "back.fits");
     214
    213215    // back-sub image pixels, from output background file (don't create if not requested)
    214216    psImage *backSub = pmFPAfileReadoutImage (config->files, view, "PSPHOT.BACKSUB", 0, 0, PS_TYPE_F32);
  • trunk/psphot/src/psphotMagnitudes.c

    r6715 r6753  
    33bool psphotMagnitudes (psArray *sources, psMetadata *config, pmPSF *psf) {
    44
    5     bool status;
    6 
    75    psTimerStart ("psphot");
    86
    9     float RADIUS = psMetadataLookupF32 (&status, config, "AP_RADIUS");
    107    for (int i = 0; i < sources->n; i++) {
    118        pmSource *source = (pmSource *) sources->data[i];
    12         pmSourceMagnitudes (source, psf, RADIUS);
     9        pmSourceMagnitudes (source, psf, true, true);
    1310    }   
    1411
     
    1613    return true;
    1714}
    18 
    19 // XXX should this function use RADIUS or should it
  • trunk/psphot/src/psphotRadiusChecks.c

    r6715 r6753  
    11# include "psphot.h"
     2# define RADIUS_TYPE float
    23
    34static float PSF_FIT_NSIGMA;
     
    2324
    2425    // set the fit radius based on the object flux limit and the model
    25     model->radius = modelRadiusPSF (model->params, PSF_FIT_NSIGMA*moments->dSky) + PSF_FIT_PADDING;
    26     if (isnan(model->radius)) psAbort ("apply_psf_model", "error in radius");
     26    model->radiusTMP = (RADIUS_TYPE) (modelRadiusPSF (model->params, PSF_FIT_NSIGMA*moments->dSky) + PSF_FIT_PADDING);
     27    if (isnan(model->radiusTMP)) psAbort ("apply_psf_model", "error in radius");
    2728       
    28     if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
    29         model->radius *= 2;
     29    if (source->mode & PM_SOURCE_MODE_SATSTAR) {
     30        model->radiusTMP *= 2;
    3031    }
    3132
    32     bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radius);
     33    bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radiusTMP);
    3334    return status;
    3435}
     
    4041
    4142    // set the fit radius based on the object flux limit and the model
    42     model->radius = modelRadiusPSF (model->params, PSF_FIT_NSIGMA*moments->dSky) + dR + PSF_FIT_PADDING;
    43     if (isnan(model->radius)) psAbort ("apply_psf_model", "error in radius");
     43    model->radiusTMP = (RADIUS_TYPE) (modelRadiusPSF (model->params, PSF_FIT_NSIGMA*moments->dSky) + dR + PSF_FIT_PADDING);
     44    if (isnan(model->radiusTMP)) psAbort ("apply_psf_model", "error in radius");
    4445       
    4546    if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
    46         model->radius *= 2;
     47        model->radiusTMP *= 2;
    4748    }
    4849
    49     bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radius);
     50    bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radiusTMP);
    5051    return status;
    5152}
     
    7374
    7475    // set the fit radius based on the object flux limit and the model
    75     model->radius = modelRadiusEXT (model->params, EXT_FIT_NSIGMA*moments->dSky) + EXT_FIT_PADDING;
    76     if (isnan(model->radius)) psAbort (__func__, "error in radius");
     76    model->radiusTMP = (RADIUS_TYPE) (modelRadiusEXT (model->params, EXT_FIT_NSIGMA*moments->dSky) + EXT_FIT_PADDING);
     77    if (isnan(model->radiusTMP)) psAbort (__func__, "error in radius");
    7778
    7879    // redefine the pixels if needed
    79     bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radius);
     80    bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radiusTMP);
    8081    return status;
    8182}
  • trunk/psphot/src/psphotReadout.c

    r6727 r6753  
    5656
    5757    // linear PSF fit to remaining peaks
    58     psphotEnsemblePSF (readout, sources, recipe, psf, TRUE);
     58    // psphotEnsemblePSF (readout, sources, recipe, psf, TRUE);
    5959
    6060    // measure aperture photometry corrections
     
    8383}
    8484
    85     # if (1)
    86     fprintf (stderr, "making fake sources\n");
    87     sources = psphotFakeSources ();
    88     psMetadata *tmp = psMetadataAlloc ();
    89     status = psMetadataAdd (readout->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_DATA_ARRAY,    "psphot sources", sources);
    90     status = psMetadataAdd (readout->analysis, PS_LIST_TAIL, "PSPHOT.HEADER",  PS_DATA_METADATA, "header stats", tmp);
    91     psFree (sources);
    92     psFree (tmp);
    93     return true;
    94     # endif
     85# if (0)
     86fprintf (stderr, "making fake sources\n");
     87sources = psphotFakeSources ();
     88psMetadata *tmp = psMetadataAlloc ();
     89status = psMetadataAdd (readout->analysis, PS_LIST_TAIL, "PSPHOT.SOURCES", PS_DATA_ARRAY,    "psphot sources", sources);
     90status = psMetadataAdd (readout->analysis, PS_LIST_TAIL, "PSPHOT.HEADER",  PS_DATA_METADATA, "header stats", tmp);
     91psFree (sources);
     92psFree (tmp);
     93return true;
     94# endif
  • trunk/psphot/src/psphotSourceFits.c

    r6715 r6753  
    1717
    1818    // fit PSF model (set/unset the pixel mask)
    19     psImageKeepCircle (source->mask, x, y, PSF->radius, "OR", PM_SOURCE_MASK_MARKED);
     19    psImageKeepCircle (source->mask, x, y, PSF->radiusTMP, "OR", PM_SOURCE_MASK_MARKED);
    2020    pmSourceFitModel (source, PSF, true);
    21     psImageKeepCircle (source->mask, x, y, PSF->radius, "AND", ~PM_SOURCE_MASK_MARKED);
     21    psImageKeepCircle (source->mask, x, y, PSF->radiusTMP, "AND", ~PM_SOURCE_MASK_MARKED);
    2222
    2323    // does the PSF model succeed?
     
    189189
    190190    // fit EXT (not PSF) model (set/unset the pixel mask)
    191     psImageKeepCircle (source->mask, x, y, PSF->radius, "OR", PM_SOURCE_MASK_MARKED);
     191    psImageKeepCircle (source->mask, x, y, PSF->radiusTMP, "OR", PM_SOURCE_MASK_MARKED);
    192192    pmSourceFitSet (source, modelSet, true);
    193     psImageKeepCircle (source->mask, x, y, PSF->radius, "AND", ~PM_SOURCE_MASK_MARKED);
     193    psImageKeepCircle (source->mask, x, y, PSF->radiusTMP, "AND", ~PM_SOURCE_MASK_MARKED);
    194194
    195195    return (modelSet);
     
    209209
    210210    // fit EXT (not PSF) model (set/unset the pixel mask)
    211     psImageKeepCircle (source->mask, x, y, EXT->radius, "OR", PM_SOURCE_MASK_MARKED);
     211    psImageKeepCircle (source->mask, x, y, EXT->radiusTMP, "OR", PM_SOURCE_MASK_MARKED);
    212212    pmSourceFitModel (source, EXT, false);
    213     psImageKeepCircle (source->mask, x, y, EXT->radius, "AND", ~PM_SOURCE_MASK_MARKED);
     213    psImageKeepCircle (source->mask, x, y, EXT->radiusTMP, "AND", ~PM_SOURCE_MASK_MARKED);
    214214
    215215    return (EXT);
     
    272272
    273273    // fit PSF model (set/unset the pixel mask)
    274     psImageKeepCircle (source->mask, x, y, PSF->radius, "OR", PM_SOURCE_MASK_MARKED);
     274    psImageKeepCircle (source->mask, x, y, PSF->radiusTMP, "OR", PM_SOURCE_MASK_MARKED);
    275275    pmSourceFitSet (source, modelSet, true);
    276     psImageKeepCircle (source->mask, x, y, PSF->radius, "AND", ~PM_SOURCE_MASK_MARKED);
     276    psImageKeepCircle (source->mask, x, y, PSF->radiusTMP, "AND", ~PM_SOURCE_MASK_MARKED);
    277277
    278278    // evaluate the blend objects, subtract if good, free otherwise
Note: See TracChangeset for help on using the changeset viewer.