IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 5718


Ignore:
Timestamp:
Dec 6, 2005, 9:55:31 PM (20 years ago)
Author:
eugene
Message:

finshed deblender, dropped local copies of psModule code

Location:
trunk/psphot
Files:
1 added
10 deleted
11 edited

Legend:

Unmodified
Added
Removed

  • trunk/psphot/Makefile

    r5672 r5718  
    3939$(SRC)/psphotMagnitudes.$(ARCH).o  \
    4040$(SRC)/psphotImageBackground.$(ARCH).o  \
     41$(SRC)/psphotBasicDeblend.$(ARCH).o  \
     42$(SRC)/psphotModelGroupInit.$(ARCH).o  \
     43$(SRC)/pmSourceContour.$(ARCH).o \
    4144$(SRC)/psLine.$(ARCH).o            \
    4245$(SRC)/psModulesUtils.$(ARCH).o    \

  • trunk/psphot/src/models/pmModel_WAUSS.c

    r4954 r5718  
    5959    if (flux <= 0) return (1.0);
    6060    if (params->data.F32[1] <= 0) return (1.0);
    61     if (flux >= params->data.F32[1] <= 0) return (1.0);
    6261
    6362    psF64 sigma  = sqrt(2.0) * hypot (1.0 / params->data.F32[4], 1.0 / params->data.F32[5]);
    … …  
    6665}
    6766
    68 bool psModelGuess_WAUSS (psModel *model, psSource *source) {
     67bool psModelGuess_WAUSS (psModel *model, pmSource *source) {
    6968
    7069    psVector *params = model->params;

  • trunk/psphot/src/pmSourceContour.c

    r5672 r5718  
    100100new implementation of source contour function
    101101*****************************************************************************/
    102 psArray *pmSourceContour_EAM (psImage *image, int x, int y, float threshold) {
     102psArray *pmSourceContour_EAM (psImage *image, int xc, int yc, float threshold) {
    103103
    104104    int xR, yR, x0, x1, x0s, x1s;
    … …  
    106106    psTrace(__func__, 3, "---- %s() begin ----\n", __func__);
    107107    PS_ASSERT_PTR_NON_NULL(image, false);
     108
     109    int x = xc - image->col0;
     110    int y = yc - image->row0;
    108111
    109112    // the requested point must be within the contour
    … …  
    141144            x0 = findContourPos (image, threshold, xR, yR, RIGHT, x1);
    142145            if (x0 == x1) {
    143                 fprintf (stderr, "top: %d (%d - %d)\n", yR, xR, x1);
     146                // fprintf (stderr, "top: %d (%d - %d)\n", yR, xR, x1);
    144147                goto pt1;
    145148            }
    … …  
    150153            x1 = findContourNeg (image, threshold, xR, yR, RIGHT);
    151154        }
    152         fprintf (stderr, "pos: %d (%d - %d)\n", yR, x0, x1);
     155        // fprintf (stderr, "pos: %d (%d - %d)\n", yR, x0, x1);
    153156
    154157        xVec->data.F32[Npt + 0] = image->col0 + x0;
    … …  
    175178            x0 = findContourPos (image, threshold, xR, yR, RIGHT, x1);
    176179            if (x0 == x1) {
    177                 fprintf (stderr, "top: %d (%d - %d)\n", yR, xR, x1);
     180                // fprintf (stderr, "top: %d (%d - %d)\n", yR, xR, x1);
    178181                goto pt2;
    179182            }
    … …  
    184187            x1 = findContourNeg (image, threshold, xR, yR, RIGHT);
    185188        }
    186         fprintf (stderr, "neg: %d (%d - %d)\n", yR, x0, x1);
     189        // fprintf (stderr, "neg: %d (%d - %d)\n", yR, x0, x1);
    187190
    188191        xVec->data.F32[Npt + 0] = image->col0 + x0;
    … …  
    202205    yVec->n = Npt;
    203206
    204     fprintf (stderr, "done\n");
     207    // fprintf (stderr, "done\n");
    205208    psArray *tmpArray = psArrayAlloc(2);
    206209    tmpArray->data[0] = (psPtr *) xVec;

  • trunk/psphot/src/psModulesUtils.h

    r5125 r5718  
    55# include <strings.h>  // for strcasecmp
    66# include <pslib.h>
     7# include <pmObjects.h>
     8# include <pmModelGroup.h>
    79# include "psLibUtils.h"
    8 # include "pmObjects_EAM.h"
    9 # include "pmModelGroup.h"
    1010
    1111// psModule extra utilities

  • trunk/psphot/src/psphot.c

    r5672 r5718  
    4747    if (!strcasecmp (breakPt, "CLASS")) exit (0);
    4848
     49    psphotBasicDeblend (sources, config, sky);
     50    if (!strcasecmp (breakPt, "DEBLEND")) exit (0);
     51
    4952    // source analysis is done in S/N order (brightest first)
    5053    sources = psArraySort (sources, psphotSortBySN);

  • trunk/psphot/src/psphot.h

    r5672 r5718  
    2121
    2222// top-level psphot functions
    23 psMetadata  *psphotArguments (int *argc, char **argv);
    24 eamReadout  *psphotSetup (psMetadata *config);
    25 psStats     *psphotImageStats (eamReadout *imdata, psMetadata *config);
    26 psArray     *psphotSourceStats (eamReadout *imdata, psMetadata *config, psArray *allpeaks);
    27 pmPSF       *psphotChoosePSF (psMetadata *config, psArray *sources, psStats *sky);
    28 bool         psphotApplyPSF (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
    29 bool         psphotFixedPSF (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
    30 bool         psphotFitGalaxies (eamReadout *imdata, psMetadata *config, psArray *sources, psStats *skyStats);
    31 void         psphotOutput (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
     23psMetadata     *psphotArguments (int *argc, char **argv);
     24eamReadout     *psphotSetup (psMetadata *config);
     25psStats        *psphotImageStats (eamReadout *imdata, psMetadata *config);
     26psArray        *psphotSourceStats (eamReadout *imdata, psMetadata *config, psArray *allpeaks);
     27pmPSF          *psphotChoosePSF (psMetadata *config, psArray *sources, psStats *sky);
     28bool            psphotApplyPSF (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
     29bool            psphotFixedPSF (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
     30bool            psphotFitGalaxies (eamReadout *imdata, psMetadata *config, psArray *sources, psStats *skyStats);
     31void            psphotOutput (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
    3232
    33 bool         psphotMarkPSF (pmSource *source, float shapeNsigma, float minSN, float maxChi, float SATURATE);
    34 bool         psphotSubtractPSF (pmSource *source);
    35 int          psphotSortBySN (const void **a, const void **b);
    36 int          psphotSortByY (const void **a, const void **b);
    37 int          psphotSaveImage (psMetadata *header, psImage *image, char *filename);
    38 bool         psphotDefinePixels (pmSource *mySource, const eamReadout *imdata, psF32 x, psF32 y, psF32 Radius);
     33bool            psphotMarkPSF (pmSource *source, float shapeNsigma, float minSN, float maxChi, float SATURATE);
     34bool            psphotSubtractPSF (pmSource *source);
     35int             psphotSortBySN (const void **a, const void **b);
     36int             psphotSortByY (const void **a, const void **b);
     37int             psphotSaveImage (psMetadata *header, psImage *image, char *filename);
     38bool            psphotDefinePixels (pmSource *mySource, const eamReadout *imdata, psF32 x, psF32 y, psF32 Radius);
     39void            psphotModelGroupInit (void);
    3940
    40 bool         pmSourceFitFixed (pmSource *source, pmModel *model);
    41 psArray     *pmPeaksSigmaLimit (eamReadout *imdata, psMetadata *config, psStats *sky);
    42 pmModel     *pmSourceMagnitudes (pmSource *source, pmPSF *psf, float apRadius);
    43 pmModel     *pmSourceSelectModel (pmSource *source);
     41bool            pmSourceFitFixed (pmSource *source, pmModel *model);
     42psArray        *pmPeaksSigmaLimit (eamReadout *imdata, psMetadata *config, psStats *sky);
     43pmModel        *pmSourceMagnitudes (pmSource *source, pmPSF *psf, float apRadius);
     44pmModel        *pmSourceSelectModel (pmSource *source);
    4445
    4546// eamReadout functions
    46 eamReadout *eamReadoutAlloc (psImage *image, psImage *noise, psImage *mask, psMetadata *header);
     47eamReadout     *eamReadoutAlloc (psImage *image, psImage *noise, psImage *mask, psMetadata *header);
    4748
    4849// output functions
    49 bool         pmSourcesWriteText (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
    50 bool         pmSourcesWriteOBJ  (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
    51 bool         pmSourcesWriteCMP  (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
    52 bool         pmSourcesWriteCMF  (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
    53 bool         pmSourcesWriteSX   (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
     50bool            pmSourcesWriteText (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
     51bool            pmSourcesWriteOBJ  (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
     52bool            pmSourcesWriteCMP  (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
     53bool            pmSourcesWriteCMF  (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
     54bool            pmSourcesWriteSX   (eamReadout *imdata, psMetadata *config, char *filename, psArray *sources, pmPSF *psf, psStats *sky);
    5455
    55 bool         pmPeaksWriteText (psArray *sources, char *filename);
    56 bool         pmMomentsWriteText (psArray *sources, char *filename);
    57 bool         pmModelWritePSFs (psArray *sources, psMetadata *config, char *filename, pmPSF *psf);
    58 bool         pmModelWriteFLTs (psArray *sources, char *filename);
    59 bool         pmModelWriteNULLs (psArray *sources, char *filename);
     56bool            pmPeaksWriteText (psArray *sources, char *filename);
     57bool            pmMomentsWriteText (psArray *sources, char *filename);
     58bool            pmModelWritePSFs (psArray *sources, psMetadata *config, char *filename, pmPSF *psf);
     59bool            pmModelWriteFLTs (psArray *sources, char *filename);
     60bool            pmModelWriteNULLs (psArray *sources, char *filename);
    6061
    61 int          pmSourcesDophotType (pmSource *source);
     62int             pmSourcesDophotType (pmSource *source);
    6263
    6364// psphotModelTest functions
    64 psMetadata *modelTestArguments (int *argc, char **argv);
    65 bool modelTestFitSource (eamReadout *imdata, psMetadata *config);
    66 bool psphotEnsemblePSF (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
    67 float psphotCrossProduct (pmSource *Mi, pmSource *Mj);
     65psMetadata     *modelTestArguments (int *argc, char **argv);
     66bool            modelTestFitSource (eamReadout *imdata, psMetadata *config);
     67bool            psphotEnsemblePSF (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
     68float           psphotCrossProduct (pmSource *Mi, pmSource *Mj);
    6869psPolynomial2D *psphotImageBackground (eamReadout *imdata, psMetadata *config, psStats *sky);
    69 bool psphotReapplyPSF (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
     70bool            psphotReapplyPSF (eamReadout *imdata, psMetadata *config, psArray *sources, pmPSF *psf, psStats *sky);
    7071
    71 pmModel *pmModelCopy (pmModel *model);
    72 psArray *pmSourceContour_EAM (psImage *image, int x, int y, float threshold);
    73 psMetadata *psphotTestArguments (int *argc, char **argv);
     72pmModel        *pmModelCopy (pmModel *model);
     73psArray        *pmSourceContour_EAM (psImage *image, int x, int y, float threshold);
     74psMetadata     *psphotTestArguments (int *argc, char **argv);
     75bool            psphotBasicDeblend (psArray *sources, psMetadata *config, psStats *sky);

  • trunk/psphot/src/psphotBasicDeblend.c

    r5704 r5718  
    11# include "psphot.h"
    22
    3 // try PSF models and select best option
     3bool psphotBasicDeblend (psArray *sources, psMetadata *config, psStats *sky) {
    44
    5 pmPSF *psphotBasicDeblend (psMetadata *config, psArray *sources, psStats *skystats)
    6 {
     5    int N;
     6    bool status;
     7    float threshold;
     8    pmSource *source, *testSource;
    79
    8     // source analysis is done in S/N order (brightest first)
     10    FILE *f = fopen ("deblend.dat", "w");
     11    if (f == NULL) psAbort ("psphot", "can't open deblend.dat output file");
     12    int Nblend = 0;
     13
     14    psTimerStart ("psphot");
     15
     16    // this should be added to the PM_SOURCE flags:
     17    int PM_SOURCE_BLEND = PM_SOURCE_OTHER + 1;
     18
     19    float FRACTION = psMetadataLookupF32 (&status, config, "DEBLEND_PEAK_FRACTION");
     20    float NSIGMA   = psMetadataLookupF32 (&status, config, "DEBLEND_SKY_NSIGMA");
     21    float minThreshold = NSIGMA*sky->sampleStdev;
     22    fprintf (stderr, "min threshold: %f\n", minThreshold);
     23
     24    // we need sources spatially-sorted to find overlaps
    925    sources = psArraySort (sources, psphotSortByY);
    1026
     27    // source analysis is done in peak order (brightest first)
     28    // we use an index for this so the spatial sorting is kept
    1129    psVector *SN = psVectorAlloc (sources->n, PS_DATA_F32);
    1230    for (int i = 0; i < SN->n; i++) {
    13       source = sources->data[i];
    14       SN->data.F32[i] = source->moments->SN;
     31        source = sources->data[i];
     32        SN->data.F32[i] = source->peak->counts;
    1533    }
    1634    psVector *index = psVectorSortIndex (NULL, SN);
    17    
     35    // this results in an index of increasing SN
     36
     37    // temporary array for overlapping objects we find
    1838    psArray *overlap = psArrayAlloc (100);
    1939
    20     for (int i = 0; i < sources->n; i++) {
    21       N = index->data.U32[i];
    22       source = sources->data[N];
    23       overlap->n = 0;
     40    // XXX make sure this results in decreasing, not increasing, SN
     41    for (int i = sources->n - 1; i >= 0; i--) {
     42        N = index->data.U32[i];
     43        source = sources->data[N];
    2444
    25       // search backwards for overlapping sources
    26       for (int j = N - 1; j >= 0; j--) {
    27         testSource = sources->data[j];
    28         if (testSource->moments->x < source->pixels->col0) continue;
    29         if (testSource->moments->x >= source->pixels->col0 + source->pixels->numCols) continue;
    30         if (testSource->moments->y < source->pixels->row0) break;
    31         if (testSource->moments->y >= source->pixels->row0 + source->pixels->numRows) {
    32           fprintf (stderr, "warning: invalid condition\n");
    33           continue;
     45        if (source->type == PM_SOURCE_BLEND) continue;
     46
     47        overlap->n = 0;
     48
     49        // search backwards for overlapping sources
     50        for (int j = N - 1; j >= 0; j--) {
     51            testSource = sources->data[j];
     52            if (testSource->peak->x <  source->pixels->col0) continue;
     53            if (testSource->peak->x >= source->pixels->col0 + source->pixels->numCols) continue;
     54            if (testSource->peak->y <  source->pixels->row0) break;
     55            if (testSource->peak->y >= source->pixels->row0 + source->pixels->numRows) {
     56                fprintf (stderr, "warning: invalid condition\n");
     57                continue;
     58            }
     59            psArrayAdd (overlap, 100, testSource);
    3460        }
    35         psArrayAdd (overlap, 100, testSource);
    36       }
    3761
    38       // search forwards for overlapping sources
    39       for (int j = N + 1; j < sources->n; j++) {
    40         testSource = sources->data[j];
    41         if (testSource->moments->x < source->pixels->col0) continue;
    42         if (testSource->moments->x >= source->pixels->col0 + source->pixels->numCols) continue;
    43         if (testSource->moments->y < source->pixels->row0) {
    44           fprintf (stderr, "warning: invalid condition\n");
    45           continue;
     62        // search forwards for overlapping sources
     63        for (int j = N + 1; j < sources->n; j++) {
     64            testSource = sources->data[j];
     65            if (testSource->peak->x <  source->pixels->col0) continue;
     66            if (testSource->peak->x >= source->pixels->col0 + source->pixels->numCols) continue;
     67            if (testSource->peak->y <  source->pixels->row0) {
     68                fprintf (stderr, "warning: invalid condition\n");
     69                continue;
     70            }
     71            if (testSource->peak->y >= source->pixels->row0 + source->pixels->numRows) break;
     72            psArrayAdd (overlap, 100, testSource);
    4673        }
    47         if (testSource->moments->y >= source->pixels->row0 + source->pixels->numRows) break;
    48         psArrayAdd (overlap, 100, testSource);
    49       }
    5074
    51       // generate source contour (1/4 peak counts)
    52       if (overlap->n > 0) {
    53         // XXX EAM : make the 1/4 user-defined.
    54         threshold = 0.25 * (source->moments->Peak - source->moments->Sky) + source->moments->Sky);
    55         psArray *contour = pmSourceContour_EAM (source->pixels, source->moments->x, source->moments->x, threshold);
     75        // generate source contour (1/4 peak counts)
     76        if (overlap->n > 0) {
     77            // XXX EAM : make the 1/4 user-defined.
     78            // XXX keep threshold from dropping too low (N*sky.sigma)
     79            threshold = FRACTION * (source->peak->counts - source->moments->Sky) + source->moments->Sky;
     80            threshold = PS_MAX (threshold, minThreshold);
    5681
    57         // XXX EAM : the order is wrong: need to check j, j+1 contour points
    58         psVector *xv = contour->data[0];
    59         psVector *yv = contour->data[1];
    60         for (int j = 0; j < xv->n; j++) {
    61           for (int k = 0; k < overlap->n; k++) {
    62             testSource = overlap->data[k];
    63             if (fabs(yv->data.F32[j] - testSource->moments->y) > 0.5) continue;
    64             if (xv->data.F32[j] > testSource->moments->x) continue;
    65             if (xv->data.F32[j] < testSource->moments->x) continue;
    66             // XXX EAM : mark source with flag (blended)
    67           }
    68         } 
     82            // XXX EAM : should the contour input coordinate be in parent or subimage coords? parent, for now
     83            psArray *contour = pmSourceContour_EAM (source->pixels, source->peak->x, source->peak->y, threshold);
     84            if (contour == NULL) {
     85                fprintf (stderr, "below threshold? invalid peak?\n");
     86                continue;
     87            }
     88
     89            // the source contour consists of two vectors, xv and yv.  the contour is
     90            // a series of coordinate pairs, (xv[i],yv[i]) & (xv[i+1],yv[i+1]).  both
     91            // coordinate pairs have the same yv[] value, with xv[i] corresponding to
     92            // the left boundary and xv[i+1] corresponding to the right boundary
     93
     94            psVector *xv = contour->data[0];
     95            psVector *yv = contour->data[1];
     96            for (int k = 0; k < overlap->n; k++) {
     97                testSource = overlap->data[k];
     98                if (testSource->peak->counts > source->peak->counts) continue;
     99                for (int j = 0; j < xv->n; j+=2) {
     100                    if (fabs(yv->data.F32[j] - testSource->peak->y) > 0.5) continue;
     101                    if (xv->data.F32[j+0] > testSource->peak->x) break;
     102                    if (xv->data.F32[j+1] < testSource->peak->x) break;
     103
     104                    int xp0 = source->moments->x - source->pixels->col0;
     105                    int xp1 = source->peak->x - source->pixels->col0;
     106                    int xp2 = testSource->moments->x - source->pixels->col0;
     107                    int xp3 = testSource->peak->x - source->pixels->col0;
     108
     109                    int yp0 = source->moments->y - source->pixels->row0;
     110                    int yp1 = source->peak->y - source->pixels->row0;
     111                    int yp2 = testSource->moments->y - testSource->pixels->row0;
     112                    int yp3 = testSource->peak->y - testSource->pixels->row0;
     113
     114                    fprintf (f, "%d %d (%f, %f) :  %d %d (%f, %f)  vs %f\n",
     115                             source->peak->x, source->peak->y,
     116                             source->pixels->data.F32[yp0][xp0], source->pixels->data.F32[yp1][xp1],
     117                             testSource->peak->x, testSource->peak->y,
     118                             testSource->pixels->data.F32[yp2][xp2], testSource->pixels->data.F32[yp3][xp3], threshold
     119                             );
     120
     121                    testSource->type = PM_SOURCE_BLEND;
     122                    Nblend ++;
     123                    j = xv->n;
     124                }
     125            } 
     126        }
    69127    }
     128    fclose (f);
     129    psLogMsg ("psphot.deblend", 3, "identified %d blended objects (%f)\n", Nblend, psTimerMark ("psphot"));
     130    return true;
     131}
    70132
    71     psLogMsg ("psphot.pspsf", 3, "selected psf model %s, ApResid: %f +/- %f\n", modelName, psf->ApResid, psf->dApResid);
    72133
    73     return (psf);
    74 }
     134

  • trunk/psphot/src/psphotEnsemblePSF.c

    r5704 r5718  
    1313    sources = psArraySort (sources, psphotSortByY);
    1414
     15    // this should be added to the PM_SOURCE flags:
     16    int PM_SOURCE_BLEND = PM_SOURCE_OTHER + 1;
     17
    1518    float OUTER_RADIUS     = psMetadataLookupF32 (&status, config, "SKY_OUTER_RADIUS");
    1619    float PSF_FIT_NSIGMA   = psMetadataLookupF32 (&status, config, "PSF_FIT_NSIGMA");
    … …  
    2326    // pre-calculate all model pixels
    2427    psArray  *models = psArrayAlloc (sources->n);
     28    psVector *index  = psVectorAlloc (sources->n, PS_TYPE_U32);
     29    models->n = 0;
     30    index->n = 0;
     31
    2532    for (int i = 0; i < sources->n; i++) {
    2633        pmSource *inSource = sources->data[i];
     34
     35        // skip non-astronomical objects (very likely defects)
     36        // XXX EAM : should we try these anyway?
     37        if (inSource->type == PM_SOURCE_BLEND) continue;
     38        if (inSource->type == PM_SOURCE_DEFECT) continue;
     39        if (inSource->type == PM_SOURCE_SATURATED) continue;
     40
    2741        pmSource *otSource = pmSourceAlloc ();
    2842
    … …  
    7286        // save source in list
    7387        otSource->modelPSF = model;
    74         models->data[i] = otSource;
     88        index->data.U32[models->n] = i;
     89        psArrayAdd (models, 100, otSource);
    7590    }
    7691    psLogMsg ("psphot.emsemble", 4, "built models: %f (%d objects)\n", psTimerMark ("psphot"), sources->n);
    … …  
    7994
    8095    // fill out the sparse matrix
    81     psSparse *sparse = psSparseAlloc (sources->n, 100);
    82     for (int i = 0; i < sources->n; i++) {
    83         pmSource *Fi = sources->data[i];
     96    psSparse *sparse = psSparseAlloc (models->n, 100);
     97    for (int i = 0; i < models->n; i++) {
     98        int N = index->data.U32[i];
     99        pmSource *Fi = sources->data[N];
    84100        pmSource *Mi = models->data[i];
    85101
    … …  
    93109
    94110        // loop over all other stars following this one
    95         for (int j = i + 1; j < sources->n; j++) {
     111        for (int j = i + 1; j < models->n; j++) {
    96112            pmSource *Mj = models->data[j];
    97113
    98             // XXX double check that this is working!  should not break here
    99             check me here;
    100             if (Mi->pixels->col0 + Mi->pixels->numCols < Mj->pixels->col0) break;
     114            // skip over disjoint source images, break after last possible overlap
     115            if (Mi->pixels->row0 + Mi->pixels->numRows < Mj->pixels->row0) break;
     116            if (Mj->pixels->row0 + Mj->pixels->numRows < Mi->pixels->row0) continue;
     117            if (Mi->pixels->col0 + Mi->pixels->numCols < Mj->pixels->col0) continue;
    101118            if (Mj->pixels->col0 + Mj->pixels->numCols < Mi->pixels->col0) continue;
    102             if (Mi->pixels->row0 + Mi->pixels->numRows < Mj->pixels->row0) continue;
    103             if (Mj->pixels->row0 + Mj->pixels->numRows < Mi->pixels->row0) continue;
    104119           
    105120            // got an overlap; calculate cross-product and add to output array
    … …  
    115130
    116131    // adjust models, set sources and subtract
    117     for (int i = 0; i < sources->n; i++) {
    118         pmSource *Fi = sources->data[i];
     132    for (int i = 0; i < models->n; i++) {
     133        int N = index->data.U32[i];
     134        pmSource *Fi = sources->data[N];
    119135        pmSource *Mi = models->data[i];
    120136

  • trunk/psphot/src/psphotImageBackground.c

    r5672 r5718  
    5858    for (int j = 0; j < image->numRows; j++) {
    5959        for (int i = 0; i < image->numCols; i++) {
    60             if (mask->data.U8[j][i]) {
    61                 image->data.F32[j][i] = 0;
    62             } else {
     60            if (!mask->data.U8[j][i]) {
    6361                Sky = psPolynomial2DEval (skyModel, i, j);
    6462                image->data.F32[j][i] -= Sky;
    … …  
    6664        }
    6765    }
     66
    6867    psLogMsg ("psphot", 5, "back: %f sec (fit model)\n", psTimerMark ("psphot"));
    6968

  • trunk/psphot/src/psphotModelGroupInit.c

    r5672 r5718  
    11# include "psphot.h"
    22
    3 # include "models/pmModel_WAUSS.c"
     3# include "models/pmModel_TAUSS.c"
    44
    55static pmModelGroup userModels[] = {
    6     {"PS_MODEL_WAUSS", 9, pmModelFunc_WAUSS,  pmModelFlux_WAUSS,  pmModelRadius_WAUSS,  pmModelLimits_WAUSS,  pmModelGuess_WAUSS, pmModelFromPSF_WAUSS, pmModelFitStatus_WAUSS},
     6    {"PS_MODEL_TAUSS", 9, pmModelFunc_TAUSS,  pmModelFlux_TAUSS,  pmModelRadius_TAUSS,  pmModelLimits_TAUSS,  pmModelGuess_TAUSS, pmModelFromPSF_TAUSS, pmModelFitStatus_TAUSS},
    77};
    88
    … …  
    1717        pmModelGroupAdd (&userModels[i]);
    1818    }
    19     return true;
     19    return;
    2020}

  • trunk/psphot/src/psphotReapplyPSF.c

    r5672 r5718  
    1414    sources = psArraySort (sources, psphotSortBySN);
    1515   
    16     // we may set this differently here from the value used to mark likely saturated stars
     16     // this should be added to the PM_SOURCE flags:
     17    int PM_SOURCE_BLEND = PM_SOURCE_OTHER + 1;
     18
     19   // we may set this differently here from the value used to mark likely saturated stars
    1720    float SATURATION       = psMetadataLookupF32 (&status, config, "SATURATION");
    1821    float PSF_MIN_SN       = psMetadataLookupF32 (&status, config, "PSF_MIN_SN");
    … …  
    2629        // skip non-astronomical objects (very likely defects)
    2730        // XXX EAM : should we try these anyway?
     31        if (source->type == PM_SOURCE_BLEND) continue;
    2832        if (source->type == PM_SOURCE_DEFECT) continue;
    2933        if (source->type == PM_SOURCE_SATURATED) continue;
Note: See TracChangeset for help on using the changeset viewer.