IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 36747 for trunk/psphot/src


Ignore:
Timestamp:
May 13, 2014, 2:40:41 PM (12 years ago)
Author:
bills
Message:

change psphotFullForce galaxy shapes analysis to use different radius
trial values based on the signal to nosie of the source as measured
by psphotStack

Location:
trunk/psphot/src
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/psphot/src/psphot.h

    r36718 r36747  
    432432
    433433typedef struct {
     434    float   Q;
     435    float   NSigma;
     436    float   clampSN;
     437#ifdef notdef
    434438    float fRmajorMin;
    435439    float fRmajorMax;
    … …  
    438442    float fRminorMax;
    439443    float fRminorDel;
     444    float snMinAllTrials;
     445    float reduceTrialsFactor;
     446    bool  reducedTrials;
     447#endif
    440448} psphotGalaxyShapeOptions;
    441449
    … …  
    539547bool psphotGalaxyShapeSource (pmPCMdata *pcm, pmSource *source, pmSourceGalaxyFits *galaxyFits, psImageMaskType maskVal, int psfSize, bool saveResults);
    540548psphotGalaxyShapeOptions *psphotGalaxyShapeOptionsAlloc();
     549void psphotGalaxyShapeOptionsSet(pmSource *source, pmSourceGalaxyFits *galaxyFits, psphotGalaxyShapeOptions *defaultOptions);
    541550
    542551bool psphotRadialProfileWings (pmConfig *config, const pmFPAview *view, const char *filerule);

  • trunk/psphot/src/psphotFullForceSummaryReadout.c

    r36731 r36747  
    33
    44typedef struct {
     5    psF32   Q;
     6    psF32   NSigma;
     7    psF32   clampSN;
     8#ifdef notdef
     9    // these are per object now
    510    int     numTrials;
    6     psF64   fRMajorMin;
    7     psF64   fRMajorMax;
    8     psF64   fRMajorDel;
    9     psF64   fRMinorMin;
    10     psF64   fRMinorMax;
    11     psF64   fRMinorDel;
    12     psVector    *fRMajor;
    13     psVector    *fRMinor;
     11    psF64   fRmajorMin;
     12    psF64   fRmajorMax;
     13    psF64   fRmajorDel;
     14    psF64   fRminorMin;
     15    psF64   fRminorMax;
     16    psF64   fRminorDel;
     17    psVector    *fRmajor;
     18    psVector    *fRminor;
     19#endif
    1420    psArray *zeroPt;    // zero points for each input used with exptime to scale flux
    1521    psArray *exptime;   // exposure times for each input
    1622    psArray *cffTables; // one for each model type index is (model_type + 1) entry 0 is no extended model (star)
    17 } galaxyShapeOptions;
    18 
    19 
    20 
    21 static pmSource *psphotFullForceSummarizeObject(pmConfig *config, pmPhotObj *obj, psVector *fluxScaleFactor, galaxyShapeOptions *options);
     23} galaxyShapeSummaryOptions;
     24
     25
     26
     27static pmSource *psphotFullForceSummarizeObject(pmConfig *config, pmPhotObj *obj, psVector *fluxScaleFactor, galaxyShapeSummaryOptions *options);
    2228static pmPhotObj *findObjectForSource(psArray **pObjects, pmSource *source);
    2329
    24 static bool setOptions(galaxyShapeOptions *options, pmReadout *readout, psMetadata *recipe, bool saveVectors);
    25 static bool checkOptions(galaxyShapeOptions *options, pmReadout *readout, psMetadata *recipe);
     30static bool setOptions(galaxyShapeSummaryOptions *options, pmReadout *readout, psMetadata *recipe, bool saveVectors);
     31static bool checkOptions(galaxyShapeSummaryOptions *options, pmReadout *readout, psMetadata *recipe);
    2632
    2733
    … …  
    5460
    5561    psVector *fluxScaleFactor = psVectorAlloc(num, PS_TYPE_F32);
    56     galaxyShapeOptions options;
     62    galaxyShapeSummaryOptions options;
    5763    psArray *objects = NULL;
    5864    for (int index = 0; index < num; index++) {
    … …  
    174180}
    175181
    176 static pmSource *psphotFullForceSummarizeObject(pmConfig *config, pmPhotObj *obj, psVector *fluxScaleFactor, galaxyShapeOptions *options) {
     182static pmSource *psphotFullForceSummarizeObject(pmConfig *config, pmPhotObj *obj, psVector *fluxScaleFactor, galaxyShapeSummaryOptions *options) {
    177183
    178184    pmSource *outSrc = NULL;
    … …  
    220226                outSrc->extpars =  psMemIncrRefCounter(source->extpars);
    221227            }
    222             outSrc->modelEXT = psFree(outSrc->modelEXT);
    223228        }
    224229        if (nModels && !outSrc->modelFits && source->modelFits && source->modelFits->n == nModels) {
    … …  
    264269        psVector *sumInvSig2 = NULL;
    265270        psVector *numerator   = NULL;
     271        psVector *fRmajorVec = NULL;
     272        psVector *fRminorVec = NULL;
    266273        psF32   totalNPix = 0;
    267274        long    vectorLength = 0;
     275        psF32   fRmajorMin = 0;
     276        psF32   fRmajorMax = 0;
     277        psF32   fRmajorDel = 0;
     278        psF32   fRminorMin = 0;
     279        psF32   fRminorMax = 0;
     280        psF32   fRminorDel = 0;
    268281        for (int i=0; i < obj->sources->n; i++) {
    269282            pmSource *source = obj->sources->data[i];
    … …  
    282295                galaxyFits = NULL;
    283296            }
    284             // end of proposed function
    285297
    286298            if (galaxyFits && isfinite(galaxyFits->nPix) && galaxyFits->chisq->n) {
    … …  
    294306                    numerator   = psVectorAlloc(vectorLength, PS_TYPE_F32);
    295307                    psVectorInit(numerator, 0.0);
     308
     309                    // Initialize the fractional radius vectors
     310                    // save these to initialize the output galaxy fits
     311                    fRmajorMin = galaxyFits->fRmajorMin;
     312                    fRmajorMax = galaxyFits->fRmajorMax;
     313                    fRmajorDel = galaxyFits->fRmajorDel;
     314                    fRminorMin = galaxyFits->fRminorMin;
     315                    fRminorMax = galaxyFits->fRminorMax;
     316                    fRminorDel = galaxyFits->fRminorDel;
     317                    fRmajorVec = psVectorAlloc(vectorLength, PS_TYPE_F32);
     318                    fRminorVec = psVectorAlloc(vectorLength, PS_TYPE_F32);
     319                    int k = 0;
     320                    for (float fRmajor = fRmajorMin; fRmajor < fRmajorMax + 0.5*fRmajorDel; fRmajor += fRmajorDel) {
     321                        for (float fRminor = fRminorMin; fRminor < fRminorMax + 0.5*fRminorDel; fRminor += fRminorDel) {
     322                            fRminorVec->data.F32[k] = fRminor;
     323                            fRmajorVec->data.F32[k] = fRmajor;
     324                            k++;
     325                        }
     326                    }
     327                    psAssert(k == vectorLength, "mismatched vectors");
    296328                }
    297329
    298330                // Die if the lengths of the vectors is not the same in all sources
     331#ifdef notdef
    299332                psAssert(vectorLength == options->numTrials, "length of chisq vector %ld does not match options %d",
    300333                    vectorLength, options->numTrials);
     334#endif
    301335                psAssert(galaxyFits->chisq->n == vectorLength, "length of chisq vectors do not match %ld %ld",
    302336                             galaxyFits->chisq->n, vectorLength);
    … …  
    339373                galaxyFits->chisq = psVectorRecycle(galaxyFits->chisq, vectorLength, PS_TYPE_F32);
    340374
     375            galaxyFits->fRmajorMin = fRmajorMin;
     376            galaxyFits->fRmajorMax = fRmajorMax;
     377            galaxyFits->fRmajorDel = fRmajorDel;
     378            galaxyFits->fRminorMin = fRminorMin;
     379            galaxyFits->fRminorMax = fRminorMax;
     380            galaxyFits->fRminorDel = fRminorDel;
    341381
    342382            // fill the summary galaxyFits vectors and find the trial with the minimum value for chisq
    … …  
    369409                bool fitMajor = true;
    370410                bool fitMinor = true;
    371                 psF64 fRMajorBest = options->fRMajor->data.F32[min_k];
    372                 if ((fabs(fRMajorBest - options->fRMajorMin) < options->fRMajorDel) ||
    373                     (fabs(fRMajorBest - options->fRMajorMax) < options->fRMajorDel)) {
     411                psF64 fRmajorBest = fRmajorVec->data.F32[min_k];
     412                if ((fabs(fRmajorBest - galaxyFits->fRmajorMin) < galaxyFits->fRmajorDel) ||
     413                    (fabs(fRmajorBest - galaxyFits->fRmajorMax) < galaxyFits->fRmajorDel)) {
    374414                    fitMajor = false;
    375415                }
    376                 psF64 fRMinorBest = options->fRMinor->data.F32[min_k];
    377                 if ((fabs(fRMinorBest - options->fRMinorMin) < options->fRMinorDel) ||
    378                     (fabs(fRMinorBest - options->fRMinorMax) < options->fRMinorDel)) {
     416                psF64 fRminorBest = fRminorVec->data.F32[min_k];
     417                if ((fabs(fRminorBest - galaxyFits->fRminorMin) < galaxyFits->fRminorDel) ||
     418                    (fabs(fRminorBest - galaxyFits->fRminorMax) < galaxyFits->fRminorDel)) {
    379419                    fitMinor = false;
    380420                }
    … …  
    396436                    psVector *dFlux = psVectorAllocEmpty(NUM_TRIALS_INIT, PS_TYPE_F64);
    397437
    398                     // XXX: use a recipe parameter
    399                     psF64 maxDeltaMaj = 2.2 * options->fRMajorDel;
    400                     psF64 maxDeltaMin = 2.2 * options->fRMinorDel;
     438                    // XXX: use a recipe parameter instead of 2.2
     439                    psF64 maxDeltaMaj = 2.2 * galaxyFits->fRmajorDel;
     440                    psF64 maxDeltaMin = 2.2 * galaxyFits->fRminorDel;
    401441
    402442                    psF64 majorMin = NAN;
    … …  
    408448                    for (int k = 0; k < vectorLength; k++) {
    409449
    410                         if (fabs(options->fRMajor->data.F32[k] - fRMajorBest) < maxDeltaMaj &&
    411                             fabs(options->fRMinor->data.F32[k] - fRMinorBest) < maxDeltaMin) {
     450                        if (fabs(fRmajorVec->data.F32[k] - fRmajorBest) < maxDeltaMaj &&
     451                            fabs(fRminorVec->data.F32[k] - fRminorBest) < maxDeltaMin) {
    412452
    413453                            if (isfinite(chisqVec->data.F32[k]) &&
    … …  
    417457                                // compute major and minor radius vectors from nominal and trial fractions
    418458                                // also find the ranges in the vectors
    419                                 psF64 thisMajor = R_MAJOR * options->fRMajor->data.F32[k];
     459                                psF64 thisMajor = R_MAJOR * fRmajorVec->data.F32[k];
    420460                                if (!isfinite(majorMin) || thisMajor < majorMin) {
    421461                                    majorMin = thisMajor;
    … …  
    426466                                psVectorAppend(major, thisMajor);
    427467
    428                                 psF64 thisMinor = R_MINOR * options->fRMinor->data.F32[k];
     468                                psF64 thisMinor = R_MINOR * fRminorVec->data.F32[k];
    429469                                if (!isfinite(minorMin) || thisMinor < minorMin) {
    430470                                    minorMin = thisMinor;
    … …  
    449489                        }
    450490                    }
     491
    451492
    452493                    // see if we ever get too few good values (haven't seen this happen)
    … …  
    521562#ifdef PRINTVALS
    522563                        fprintf (stderr, "%4d %3d %3ld | %6.3f %6.3f %6.3f %4.2f |  %6.3f %6.3f %6.3f %4.2f | %7.4f %7.4f %7.1f %7.1f\n",
    523                             outSrc->seq, min_k, major->n, R_MAJOR, fRMajorBest*R_MAJOR, major0, major0/R_MAJOR, R_MINOR, fRMinorBest*R_MINOR, minor0, minor0/R_MINOR, minChisq, chisq0, fluxBest, flux0);
     564                            outSrc->seq, min_k, major->n, R_MAJOR, fRmajorBest*R_MAJOR, major0, major0/R_MAJOR, R_MINOR, fRminorBest*R_MINOR, minor0, minor0/R_MINOR, minChisq, chisq0, fluxBest, flux0);
    524565#endif
    525566
    … …  
    546587#ifdef PRINTVALS
    547588                        fprintf(stderr, "%4d %3d %3ld | %6.3f %6.3f %6.3f %6.3f bad fit\n",
    548                             outSrc->seq, min_k, major->n, R_MAJOR, fRMajorBest, R_MINOR, fRMinorBest);
     589                            outSrc->seq, min_k, major->n, R_MAJOR, fRmajorBest, R_MINOR, fRminorBest);
    549590#endif
    550591                        psErrorClear();
    … …  
    563604                    // No fit
    564605                    // Set the values based on the trial with the best chisq
    565                     axes.major = R_MAJOR * fRMajorBest;
    566                     axes.minor = R_MINOR * fRMinorBest;
     606                    axes.major = R_MAJOR * fRmajorBest;
     607                    axes.minor = R_MINOR * fRminorBest;
    567608                    model->chisq = minChisq;
    568609                    flux0 = fluxBest;
    … …  
    570611#ifdef PRINTVALS
    571612                    fprintf(stderr, "%4d %3d %3ld | %6.3f %6.3f %6.3f %6.3f skip fit\n",
    572                             outSrc->seq, min_k, 0L, R_MAJOR, fRMajorBest, R_MINOR, fRMinorBest);
     613                            outSrc->seq, min_k, 0L, R_MAJOR, fRmajorBest, R_MINOR, fRminorBest);
    573614#endif
    574615                }
    … …  
    583624            }
    584625        }
     626        psFree(fRminorVec);
     627        psFree(fRmajorVec);
    585628    }
    586629
    … …  
    595638    }
    596639
    597 static bool setOptions(galaxyShapeOptions *opt, pmReadout *readout, psMetadata *recipe, bool makeVectors) {
     640static bool setOptions(galaxyShapeSummaryOptions *opt, pmReadout *readout, psMetadata *recipe, bool makeVectors) {
    598641    bool status;
    599642    bool useAnalysis;   // fall back to recipe if we dont' find values in analysis. Probably should no longer do this
    600643
     644#ifdef notdef
    601645    psMetadataLookupF32(&useAnalysis, readout->analysis, "GALAXY_SHAPES_FR_MAJOR_MIN");
     646#endif
     647    psMetadataLookupF32(&useAnalysis, readout->analysis, "GALAXY_SHAPES_Q");
    602648    psMetadata *md = useAnalysis ? readout->analysis : recipe;
    603649
    604     GETVAL(fRMajorMin, "GALAXY_SHAPES_FR_MAJOR_MIN");
    605     GETVAL(fRMajorMax, "GALAXY_SHAPES_FR_MAJOR_MAX");
    606     GETVAL(fRMajorDel, "GALAXY_SHAPES_FR_MAJOR_DEL");
    607     GETVAL(fRMinorMin, "GALAXY_SHAPES_FR_MINOR_MIN");
    608     GETVAL(fRMinorMax, "GALAXY_SHAPES_FR_MINOR_MAX");
    609     GETVAL(fRMinorDel, "GALAXY_SHAPES_FR_MINOR_DEL");
    610 
    611     opt->numTrials = ceil((opt->fRMajorMax - opt->fRMajorMin + 0.5*opt->fRMajorDel) / opt->fRMajorDel) *
    612                          ceil((opt->fRMinorMax - opt->fRMinorMin + 0.5*opt->fRMinorDel) / opt->fRMinorDel) ;
     650    GETVAL(Q, "GALAXY_SHAPES_Q");
     651    GETVAL(NSigma, "GALAXY_SHAPES_NSIGMA");
     652    GETVAL(clampSN, "GALAXY_SHAPES_CLAMP_SN");
     653
     654#ifdef notdef
     655    // these are per object now
     656    GETVAL(fRmajorMin, "GALAXY_SHAPES_FR_MAJOR_MIN");
     657    GETVAL(fRmajorMax, "GALAXY_SHAPES_FR_MAJOR_MAX");
     658    GETVAL(fRmajorDel, "GALAXY_SHAPES_FR_MAJOR_DEL");
     659    GETVAL(fRminorMin, "GALAXY_SHAPES_FR_MINOR_MIN");
     660    GETVAL(fRminorMax, "GALAXY_SHAPES_FR_MINOR_MAX");
     661    GETVAL(fRminorDel, "GALAXY_SHAPES_FR_MINOR_DEL");
     662
     663    opt->numTrials = ceil((opt->fRmajorMax - opt->fRmajorMin + 0.5*opt->fRmajorDel) / opt->fRmajorDel) *
     664                         ceil((opt->fRminorMax - opt->fRminorMin + 0.5*opt->fRminorDel) / opt->fRminorDel) ;
    613665
    614666    if (makeVectors) {
    615         opt->fRMinor = psVectorAlloc(opt->numTrials, PS_TYPE_F32);
    616         opt->fRMajor = psVectorAlloc(opt->numTrials, PS_TYPE_F32);
     667        opt->fRminor = psVectorAlloc(opt->numTrials, PS_TYPE_F32);
     668        opt->fRmajor = psVectorAlloc(opt->numTrials, PS_TYPE_F32);
    617669        int i = 0;
    618         for (float fRMajor = opt->fRMajorMin; fRMajor < opt->fRMajorMax + 0.5*opt->fRMajorDel;
    619                 fRMajor += opt->fRMajorDel) {
    620             for (float fRMinor = opt->fRMinorMin; fRMinor < opt->fRMinorMax + 0.5*opt->fRMinorDel;
    621                     fRMinor += opt->fRMinorDel) {
    622                 opt->fRMinor->data.F32[i] = fRMinor;
    623                 opt->fRMajor->data.F32[i] = fRMajor;
     670        for (float fRmajor = opt->fRmajorMin; fRmajor < opt->fRmajorMax + 0.5*opt->fRmajorDel;
     671                fRmajor += opt->fRmajorDel) {
     672            for (float fRminor = opt->fRminorMin; fRminor < opt->fRminorMax + 0.5*opt->fRminorDel;
     673                    fRminor += opt->fRminorDel) {
     674                opt->fRminor->data.F32[i] = fRminor;
     675                opt->fRmajor->data.F32[i] = fRmajor;
    624676                i++;
    625677            }
    … …  
    627679        psAssert(i == opt->numTrials, "Something's wrong with my loop got %d entries expected %d", i, opt->numTrials);
    628680    } else {
    629         opt->fRMinor = NULL;
    630         opt->fRMajor = NULL;
    631     }
     681        opt->fRminor = NULL;
     682        opt->fRmajor = NULL;
     683    }
     684#endif
    632685       
    633686    return true;
    … …  
    640693    }
    641694
    642 static bool checkOptions(galaxyShapeOptions *options, pmReadout *readout, psMetadata *recipe) {
    643     galaxyShapeOptions thisReadoutsOptions;
     695static bool checkOptions(galaxyShapeSummaryOptions *options, pmReadout *readout, psMetadata *recipe) {
     696    galaxyShapeSummaryOptions thisReadoutsOptions;
    644697   
    645698    if (!setOptions(&thisReadoutsOptions, readout, recipe, false)) {
    … …  
    647700        return false;
    648701    }
     702    CHECKVAL(options, thisReadoutsOptions, Q, "mismatched Q");
     703    CHECKVAL(options, thisReadoutsOptions, NSigma, "mismatched NSIGMA");
     704    CHECKVAL(options, thisReadoutsOptions, clampSN, "mismatched cleampSN");
     705
     706#ifdef notdef
     707    // these are per object now
    649708    CHECKVAL(options, thisReadoutsOptions, numTrials, "mismatched number of trials")
    650     CHECKVAL(options, thisReadoutsOptions, fRMajorMin, "mismatched fRMajorMin")
    651     CHECKVAL(options, thisReadoutsOptions, fRMajorMax, "mismatched fRMajorMax")
    652     CHECKVAL(options, thisReadoutsOptions, fRMajorDel, "mismatched fRMajorDel")
    653     CHECKVAL(options, thisReadoutsOptions, fRMinorMin, "mismatched fRMinorMin")
    654     CHECKVAL(options, thisReadoutsOptions, fRMinorMax, "mismatched fRMinorMax")
    655     CHECKVAL(options, thisReadoutsOptions, fRMinorDel, "mismatched fRMinorDel")
     709    CHECKVAL(options, thisReadoutsOptions, fRmajorMin, "mismatched fRmajorMin")
     710    CHECKVAL(options, thisReadoutsOptions, fRmajorMax, "mismatched fRmajorMax")
     711    CHECKVAL(options, thisReadoutsOptions, fRmajorDel, "mismatched fRmajorDel")
     712    CHECKVAL(options, thisReadoutsOptions, fRminorMin, "mismatched fRminorMin")
     713    CHECKVAL(options, thisReadoutsOptions, fRminorMax, "mismatched fRminorMax")
     714    CHECKVAL(options, thisReadoutsOptions, fRminorDel, "mismatched fRminorDel")
     715#endif
    656716
    657717    return true;

  • trunk/psphot/src/psphotGalaxyShape.c

    r36718 r36747  
    8080    assert (markVal);
    8181
    82     // maskVal is used to test for rejected pixels, and must include markVal
    8382    maskVal |= markVal;
    8483
    85     // bit-mask to mark pixels not used in analysis
    8684    psphotGalaxyShapeOptions *opt = psphotGalaxyShapeOptionsAlloc();
     85    opt->Q = psMetadataLookupF32(&status, recipe, "GALAXY_SHAPES_Q");
     86    psAssert (status, "missing GALAXY_SHAPES_Q");
     87    opt->NSigma = psMetadataLookupF32(&status, recipe, "GALAXY_SHAPES_NSIGMA");
     88    psAssert (status, "missing GALAXY_SHAPES_NSIGMA");
     89    opt->clampSN = psMetadataLookupF32(&status, recipe, "GALAXY_SHAPES_CLAMP_SN");
     90    psAssert (status, "missing GALAXY_SHAPES_NSIGMA");
     91
     92#ifdef notdef
    8793    opt->fRmajorMin = psMetadataLookupF32(&status, recipe, "GALAXY_SHAPES_FR_MAJOR_MIN"); if (!status) opt->fRmajorMin = 0.5;
    8894    opt->fRmajorMax = psMetadataLookupF32(&status, recipe, "GALAXY_SHAPES_FR_MAJOR_MAX"); if (!status) opt->fRmajorMax = 2.0;
    … …  
    9197    opt->fRminorMax = psMetadataLookupF32(&status, recipe, "GALAXY_SHAPES_FR_MINOR_MAX"); if (!status) opt->fRminorMax = 2.0;
    9298    opt->fRminorDel = psMetadataLookupF32(&status, recipe, "GALAXY_SHAPES_FR_MINOR_DEL"); if (!status) opt->fRminorDel = 0.1;
     99#endif
    93100
    94101    // what fraction of the PSF is used? (radius in pixels : 2 -> 5x5 box)
    … …  
    238245bool psphotGalaxyShapeGrid (pmSource *source, pmSourceFitOptions *fitOptions, psphotGalaxyShapeOptions *opt, psImageMaskType maskVal, int psfSize) {
    239246
     247
    240248    for (int iModel = 0 ; iModel < source->modelFits->n; iModel++) {
    241249        pmModel *model = source->modelFits->data[iModel];
    … …  
    252260
    253261        // double check that the guess is carried along...
    254 
    255         // I have some source guess (e0, e1, e2)
    256         psEllipseAxes guessAxes = pmPSF_ModelToAxes (PAR, modelType);
    257 
    258262        if (!source->galaxyFits) {
    259263            source->galaxyFits = psArrayAllocEmpty(1);
    260             // source->galaxyFits = pmSourceGalaxyFitsAlloc();
    261264        }
    262265        pmSourceGalaxyFits *galaxyFits = pmSourceGalaxyFitsAlloc();
    … …  
    264267        psFree(galaxyFits);
    265268
     269        // set the options for this source
     270        psphotGalaxyShapeOptionsSet(source, galaxyFits, opt);
     271
     272        // I have some source guess (e0, e1, e2)
     273        psEllipseAxes guessAxes = pmPSF_ModelToAxes (PAR, modelType);
     274
    266275        float fRmajorBest = NAN;
    267276        float fRminorBest = NAN;
    268277        float chisqBest = NAN;
    269         for (float fRmajor = opt->fRmajorMin; fRmajor < opt->fRmajorMax + 0.5*opt->fRmajorDel; fRmajor += opt->fRmajorDel) {
    270             for (float fRminor = opt->fRminorMin; fRminor < opt->fRminorMax + 0.5*opt->fRminorDel; fRminor += opt->fRminorDel) {
     278        for (float fRmajor = galaxyFits->fRmajorMin; fRmajor < galaxyFits->fRmajorMax + 0.5*galaxyFits->fRmajorDel; fRmajor += galaxyFits->fRmajorDel) {
     279            for (float fRminor = galaxyFits->fRminorMin; fRminor < galaxyFits->fRminorMax + 0.5*galaxyFits->fRminorDel; fRminor += galaxyFits->fRminorDel) {
    271280     
    272281                psEllipseAxes testAxes = guessAxes;
    … …  
    405414    psMemSetDeallocator(opt, (psFreeFunc) psphotGalaxyShapeOptionsFree);
    406415   
     416#ifdef notdef
    407417    opt->fRmajorMin = 0.5;
    408418    opt->fRmajorMax = 2.0;
    … …  
    411421    opt->fRminorMax = 2.0;
    412422    opt->fRminorDel = 0.1;
     423    opt->snMinAllTrials = 20;      // for sources with cff SN < 20 ...
     424    opt->reduceTrialsFactor = 2; // cut the number of trials in half (double the deltas)
     425#endif
    413426
    414427    return opt;
    415428}
    416429
     430void psphotGalaxyShapeOptionsSet(pmSource *source, pmSourceGalaxyFits *galaxyFits, psphotGalaxyShapeOptions *options)
     431{
     432    // XXX: put these in recipe
     433    // doesn't make sense to use too small of a spacing
     434    float clampedSN = source->extSN < options->clampSN ? source->extSN : options->clampSN;
     435
     436    float f_del = options->Q / clampedSN;
     437    float f_min = 1 - options->NSigma * f_del;
     438    float f_max = 1 + options->NSigma * f_del;
     439
     440    galaxyFits->fRmajorMin = f_min;
     441    galaxyFits->fRmajorMax = f_max;
     442    galaxyFits->fRmajorDel = f_del;
     443    galaxyFits->fRminorMin = f_min;
     444    galaxyFits->fRminorMax = f_max;
     445    galaxyFits->fRminorDel = f_del;
     446}
     447
Note: See TracChangeset for help on using the changeset viewer.