Changeset 36208 for branches/eam_branches/ipp-20130904/psModules/src/objects/pmSource.c

Timestamp:

Oct 13, 2013, 1:31:44 PM (13 years ago)

Author:

eugene

Message:

strip out old pmSourceMoments code (long ago moved to pmSourceMoments.c)

File:

• branches/eam_branches/ipp-20130904/psModules/src/objects/pmSource.c (modified) (2 diffs)

branches/eam_branches/ipp-20130904/psModules/src/objects/pmSource.c

@@ -692 +692 @@
             // why do we recalculate moments here?
             // we already attempt to do this in psphotSourceStats
-            // pmSourceMoments (source, INNER_RADIUS);
             Nsatstar ++;
             continue;
@@ -806 +805 @@
     return true;
 }
-
-/******************************************************************************
-pmSourceMoments(source, radius): this function takes a subImage defined in the
-pmSource data structure, along with the peak location, and determines the
-various moments associated with that peak.
-
-Requires the following to have been created:
-    pmSource
-    pmSource->peak
-    pmSource->pixels
-    pmSource->variance
-    pmSource->mask
-
-XXX: The peak calculations are done in image coords, not subImage coords.
-
-XXX EAM : this version clips input pixels on S/N
-XXX EAM : this version returns false for several reasons
-*****************************************************************************/
-# define VALID_RADIUS(X,Y,RAD2) (((RAD2) >= (PS_SQR(X) + PS_SQR(Y))) ? 1 : 0)
-
-/*** this been moved to pmSourceMoments.c ***/
-# if (0)
-bool pmSourceMoments(pmSource *source,
-                     psF32 radius)
-{
-    psTrace("psModules.objects", 10, "---- begin ----\n");
-    PS_ASSERT_PTR_NON_NULL(source, false);
-    PS_ASSERT_PTR_NON_NULL(source->peak, false);
-    PS_ASSERT_PTR_NON_NULL(source->pixels, false);
-    PS_ASSERT_FLOAT_LARGER_THAN(radius, 0.0, false);
-
-    //
-    // XXX: Verify the setting for sky if source->moments == NULL.
-    //
-    psF32 sky = 0.0;
-    if (source->moments == NULL) {
-        source->moments = pmMomentsAlloc();
-    } else {
-        sky = source->moments->Sky;
-    }
-
-    //
-    // Sum = SUM (z - sky)
-    // X1  = SUM (x - xc)*(z - sky)
-    // X2  = SUM (x - xc)^2 * (z - sky)
-    // XY  = SUM (x - xc)*(y - yc)*(z - sky)
-    //
-    psF32 peakPixel = -PS_MAX_F32;
-    psS32 numPixels = 0;
-    psF32 Sum = 0.0;
-    psF32 Var = 0.0;
-    psF32 X1 = 0.0;
-    psF32 Y1 = 0.0;
-    psF32 X2 = 0.0;
-    psF32 Y2 = 0.0;
-    psF32 XY = 0.0;
-    psF32 x  = 0;
-    psF32 y  = 0;
-    psF32 R2 = PS_SQR(radius);
-
-    psF32 xPeak = source->peak->x;
-    psF32 yPeak = source->peak->y;
-    psF32 xOff = source->pixels->col0 - source->peak->x;
-    psF32 yOff = source->pixels->row0 - source->peak->y;
-
-    // XXX why do I get different results for these two methods of finding Sx?
-    // XXX Sx, Sy would be better measured if we clip pixels close to sky
-    // XXX Sx, Sy can still be imaginary, so we probably need to keep Sx^2?
-    // We loop through all pixels in this subimage (source->pixels), and for each
-    // pixel that is not masked, AND within the radius of the peak pixel, we
-    // proceed with the moments calculation.  need to do two loops for a
-    // numerically stable result.  first loop: get the sums.
-    // XXX EAM : mask == 0 is valid
-
-    for (psS32 row = 0; row < source->pixels->numRows ; row++) {
-
-        psF32 *vPix = source->pixels->data.F32[row];
-        psF32 *vWgt = source->variance->data.F32[row];
-        psImageMaskType *vMsk = (source->maskObj == NULL) ? NULL : source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA[row];
-
-        for (psS32 col = 0; col < source->pixels->numCols ; col++, vPix++, vWgt++) {
-            if (vMsk) {
-                if (*vMsk) {
-                    vMsk++;
-                    psTrace("psModules.objects", 10, "Ignoring pixel %d,%d due to mask: %d\n",
-                            col, row, (int)*vMsk);
-                    continue;
-                }
-                vMsk++;
-            }
-            if (isnan(*vPix)) continue;
-
-            psF32 xDiff = col + xOff;
-            psF32 yDiff = row + yOff;
-
-            // radius is just a function of (xDiff, yDiff)
-            if (!VALID_RADIUS(xDiff, yDiff, R2)) {
-#if 1
-                psTrace("psModules.objects", 10, "Ignoring pixel %d,%d due to position: %f %f\n",
-                        col, row, xDiff, yDiff);
-#endif
-                continue;
-            }
-
-            psF32 pDiff = *vPix - sky;
-            psF32 wDiff = *vWgt;
-
-            // XXX EAM : check for valid S/N in pixel
-            // XXX EAM : should this limit be user-defined?
-#if 1
-            if (PS_SQR(pDiff) < wDiff) {
-                psTrace("psModules.objects", 10, "Ignoring pixel %d,%d due to insignificance: %f, %f\n",
-                        col, row, pDiff, wDiff);
-                continue;
-            }
-#endif
-
-            Var += wDiff;
-            Sum += pDiff;
-
-            psF32 xWght = xDiff * pDiff;
-            psF32 yWght = yDiff * pDiff;
-
-            X1  += xWght;
-            Y1  += yWght;
-
-            XY  += xDiff * yWght;
-            X2  += xDiff * xWght;
-            Y2  += yDiff * yWght;
-
-            peakPixel = PS_MAX (*vPix, peakPixel);
-            numPixels++;
-        }
-    }
-
-    // if we have less than (1/4) of the possible pixels, force a retry
-    // XXX EAM - the limit is a bit arbitrary.  make it user defined?
-    if ((numPixels < 0.75*R2) || (Sum <= 0)) {
-        psTrace ("psModules.objects", 3, "insufficient valid pixels (%d vs %d; %f) for source\n",
-                 numPixels, (int)(0.75*R2), Sum);
-        psTrace("psModules.objects", 10, "---- end (false) ----\n");
-        return (false);
-    }
-
-    psTrace ("psModules.objects", 4, "sky: %f  Sum: %f  X1: %f  Y1: %f  X2: %f  Y2: %f  XY: %f  Npix: %d\n",
-             sky, Sum, X1, Y1, X2, Y2, XY, numPixels);
-
-    //
-    // first moment X  = X1/Sum + xc
-    // second moment X = sqrt (X2/Sum - (X1/Sum)^2)
-    // Sxy             = XY / Sum
-    //
-    x = X1/Sum;
-    y = Y1/Sum;
-    if ((fabs(x) > radius) || (fabs(y) > radius)) {
-        psTrace ("psModules.objects", 3, "large centroid swing; invalid peak %d, %d\n",
-                 source->peak->x, source->peak->y);
-        psTrace("psModules.objects", 10, "---- end(false)  ----\n");
-        return (false);
-    }
-
-    source->moments->Mx = x + xPeak;
-    source->moments->My = y + yPeak;
-
-    // XXX EAM : Sxy needs to have x*y subtracted
-    source->moments->Mxy = XY/Sum - x*y;
-    source->moments->Sum = Sum;
-    source->moments->SN  = Sum / sqrt(Var);
-    source->moments->Peak = peakPixel;
-    source->moments->nPixels = numPixels;
-
-    // XXX EAM : these values can be negative, so we need to limit the range
-    // XXX EAM : make the use of this consistent: should this be the second moment or sqrt?
-    // source->moments->Mxx = sqrt(PS_MAX(X2/Sum - PS_SQR(x), 0));
-    // source->moments->Myy = sqrt(PS_MAX(Y2/Sum - PS_SQR(y), 0));
-    source->moments->Mxx = PS_MAX(X2/Sum - PS_SQR(x), 0);
-    source->moments->Myy = PS_MAX(Y2/Sum - PS_SQR(y), 0);
-
-    psTrace ("psModules.objects", 4,
-             "sky: %f  Sum: %f  Mx: %f  My: %f  Mxx: %f  Myy: %f  Mxy: %f\n",
-             sky, Sum, source->moments->Mx, source->moments->My,
-             source->moments->Mxx, source->moments->Myy, source->moments->Mxy);
-
-    psTrace("psModules.objects", 10, "---- end ----\n");
-    return(true);
-}
-# endif
 
 // construct a realization of the source model