IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Mar 19, 2007, 12:48:11 PM (19 years ago)
Author:
Paul Price
Message:

Avoiding warning due to "unused variable".

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/psLib/src/math/psStats.c

    r12458 r12495  
    1313 * use ->min and ->max (PS_STAT_USE_RANGE)
    1414 *
    15  *  @version $Revision: 1.203 $ $Name: not supported by cvs2svn $
    16  *  @date $Date: 2007-03-16 00:33:55 $
     15 *  @version $Revision: 1.204 $ $Name: not supported by cvs2svn $
     16 *  @date $Date: 2007-03-19 22:48:11 $
    1717 *
    1818 *  Copyright 2006 IfA, University of Hawaii
     
    6767(0.5 * (HISTOGRAM->bounds->data.F32[(BIN_NUM)] + HISTOGRAM->bounds->data.F32[(BIN_NUM)+1]))
    6868
    69 // set the bin closest to the corresponding value.  if USE_END is +/- 1, 
     69// set the bin closest to the corresponding value.  if USE_END is +/- 1,
    7070// out-of-range saturates on lower/upper bin REGARDLESS of actual value
    7171#define PS_BIN_FOR_VALUE(RESULT, VECTOR, VALUE, USE_END) { \
    72         psVectorBinaryDisectResult result; \
    73         psScalar tmpScalar; \
    74         tmpScalar.type.type = PS_TYPE_F32; \
    75         tmpScalar.data.F32 = (VALUE); \
    76         RESULT = psVectorBinaryDisect (&result, VECTOR, &tmpScalar); \
    77         switch (result) { \
    78           case PS_BINARY_DISECT_PASS: \
    79             break; \
    80           case PS_BINARY_DISECT_OUTSIDE_RANGE: \
    81             psTrace ("psLib.math", 6, "selected bin outside range"); \
    82             if (USE_END == -1) { RESULT = 0; } \
    83             if (USE_END == +1) { RESULT = VECTOR->n - 1; } \
    84             break; \
    85           case PS_BINARY_DISECT_INVALID_INPUT: \
    86           case PS_BINARY_DISECT_INVALID_TYPE: \
    87             psAbort ("programming error"); \
    88             break; \
     72        psVectorBinaryDisectResult result; \
     73        psScalar tmpScalar; \
     74        tmpScalar.type.type = PS_TYPE_F32; \
     75        tmpScalar.data.F32 = (VALUE); \
     76        RESULT = psVectorBinaryDisect (&result, VECTOR, &tmpScalar); \
     77        switch (result) { \
     78          case PS_BINARY_DISECT_PASS: \
     79            break; \
     80          case PS_BINARY_DISECT_OUTSIDE_RANGE: \
     81            psTrace ("psLib.math", 6, "selected bin outside range"); \
     82            if (USE_END == -1) { RESULT = 0; } \
     83            if (USE_END == +1) { RESULT = VECTOR->n - 1; } \
     84            break; \
     85          case PS_BINARY_DISECT_INVALID_INPUT: \
     86          case PS_BINARY_DISECT_INVALID_TYPE: \
     87            psAbort ("programming error"); \
     88            break; \
    8989        } }
    9090
     
    9292        float dX, dY, Xo, Yo, Xt; \
    9393        if (BIN == BOUNDS->n - 1) { \
    94             dX = 0.5*(BOUNDS->data.F32[BIN+1] - BOUNDS->data.F32[BIN-1]); \
    95             dY = VECTOR->data.F32[BIN] - VECTOR->data.F32[BIN-1]; \
    96             Xo = 0.5*(BOUNDS->data.F32[BIN+1] + BOUNDS->data.F32[BIN]); \
    97             Yo = VECTOR->data.F32[BIN]; \
     94            dX = 0.5*(BOUNDS->data.F32[BIN+1] - BOUNDS->data.F32[BIN-1]); \
     95            dY = VECTOR->data.F32[BIN] - VECTOR->data.F32[BIN-1]; \
     96            Xo = 0.5*(BOUNDS->data.F32[BIN+1] + BOUNDS->data.F32[BIN]); \
     97            Yo = VECTOR->data.F32[BIN]; \
    9898        } else { \
    99             dX = 0.5*(BOUNDS->data.F32[BIN+2] - BOUNDS->data.F32[BIN]); \
    100             dY = VECTOR->data.F32[BIN+1] - VECTOR->data.F32[BIN]; \
    101             Xo = 0.5*(BOUNDS->data.F32[BIN+1] + BOUNDS->data.F32[BIN]); \
    102             Yo = VECTOR->data.F32[BIN]; \
    103         } \
    104         if (dY != 0) { \
    105             Xt = (VALUE - Yo)*dX/dY + Xo; \
    106         } else { \
    107             Xt = Xo; \
    108         } \
    109         Xt = PS_MIN (BOUNDS->data.F32[BIN+1], PS_MAX(BOUNDS->data.F32[BIN], Xt)); \
     99            dX = 0.5*(BOUNDS->data.F32[BIN+2] - BOUNDS->data.F32[BIN]); \
     100            dY = VECTOR->data.F32[BIN+1] - VECTOR->data.F32[BIN]; \
     101            Xo = 0.5*(BOUNDS->data.F32[BIN+1] + BOUNDS->data.F32[BIN]); \
     102            Yo = VECTOR->data.F32[BIN]; \
     103        } \
     104        if (dY != 0) { \
     105            Xt = (VALUE - Yo)*dX/dY + Xo; \
     106        } else { \
     107            Xt = Xo; \
     108        } \
     109        Xt = PS_MIN (BOUNDS->data.F32[BIN+1], PS_MAX(BOUNDS->data.F32[BIN], Xt)); \
    110110        psTrace("psLib.math", 6, "(Xo, Yo, dX, dY, Xt, Yt) is (%.2f %.2f %.2f %.2f %.2f %.2f)\n", \
    111111                Xo, Yo, dX, dY, Xt, VALUE); \
     
    166166*****************************************************************************/
    167167    static bool vectorSampleMean(const psVector* myVector,
    168                                 const psVector* errors,
    169                                 const psVector* maskVector,
    170                                 psMaskType maskVal,
    171                                 psStats* stats)
     168                                const psVector* errors,
     169                                const psVector* maskVector,
     170                                psMaskType maskVal,
     171                                psStats* stats)
    172172{
    173173    psTrace("psLib.math", 4, "---- %s() begin ----\n", __func__);
     
    234234*****************************************************************************/
    235235    static long vectorMinMax(const psVector* myVector,
    236                              const psVector* maskVector,
    237                              psMaskType maskVal,
    238                              psStats* stats
    239         )
     236                             const psVector* maskVector,
     237                             psMaskType maskVal,
     238                             psStats* stats
     239        )
    240240{
    241241    psTrace("psLib.math", 4, "---- %s() begin ----\n", __func__);
     
    377377*****************************************************************************/
    378378    static bool vectorSampleStdev(const psVector* myVector,
    379                                   const psVector* errors,
    380                                   const psVector* maskVector,
    381                                   psMaskType maskVal,
    382                                   psStats* stats)
     379                                  const psVector* errors,
     380                                  const psVector* maskVector,
     381                                  psMaskType maskVal,
     382                                  psStats* stats)
    383383{
    384384    psTrace("psLib.math", 4, "---- %s() begin ----\n", __func__);
     
    536536            for (long j = 0; j < myVector->n; j++) {
    537537                if (!tmpMask->data.U8[j] &&
    538                     fabsf(myVector->data.F32[j] - clippedMean) > stats->clipSigma * errors->data.F32[j]) {
     538                    fabsf(myVector->data.F32[j] - clippedMean) > stats->clipSigma * errors->data.F32[j]) {
    539539                    tmpMask->data.U8[j] = 0xff;
    540540                    psTrace("psLib.math", 10, "Clipped %ld: %f +/- %f\n", j,
     
    547547            for (long j = 0; j < myVector->n; j++) {
    548548                if (!tmpMask->data.U8[j] &&
    549                     fabsf(myVector->data.F32[j] - clippedMean) > (stats->clipSigma * clippedStdev)) {
     549                    fabsf(myVector->data.F32[j] - clippedMean) > (stats->clipSigma * clippedStdev)) {
    550550                    tmpMask->data.U8[j] = 0xff;
    551551                    psTrace("psLib.math", 10, "Clipped %ld: %f\n", j, myVector->data.F32[j]);
     
    667667            // Data range calculation failed
    668668            psError(PS_ERR_UNKNOWN, false, "Failed to calculate the min/max of the input vector.\n");
    669             goto escape;
     669            goto escape;
    670670        }
    671671        psTrace("psLib.math", 6, "Data min/max is (%.2f, %.2f)\n", min, max);
     
    731731        psTrace("psLib.math", 6, "Total data points is %ld\n", totalDataPoints);
    732732
    733         PS_BIN_FOR_VALUE(binMedian, cumulative->nums, totalDataPoints/2.0, 0);
     733        PS_BIN_FOR_VALUE(binMedian, cumulative->nums, totalDataPoints/2.0, 0);
    734734        psTrace("psLib.math", 6, "The median bin is %ld (%.2f to %.2f)\n", binMedian,
    735735                cumulative->bounds->data.F32[binMedian], cumulative->bounds->data.F32[binMedian+1]);
     
    737737        // ADD step 3: Interpolate to the exact 50% position: this is the robust histogram median.
    738738        stats->robustMedian = fitQuadraticSearchForYThenReturnX(cumulative->bounds, cumulative->nums,
    739                                                                 binMedian, totalDataPoints/2.0);
     739                                                                binMedian, totalDataPoints/2.0);
    740740        if (isnan(stats->robustMedian)) {
    741741            psError(PS_ERR_UNKNOWN, false, "Failed to fit a quadratic and calculate the 50-percent position.\n");
    742             goto escape;
     742            goto escape;
    743743        }
    744744        psTrace("psLib.math", 6, "Current robust median is %f\n", stats->robustMedian);
    745745
    746746        // ADD step 4: Find the bins which contains the 15.8655% (-1 sigma) and 84.1345% (+1 sigma) data points
    747         // also find the 30.8538% (-0.5 sigma) and 69.1462% (+0.5 sigma) points
    748         PS_BIN_FOR_VALUE(binLo, cumulative->nums, totalDataPoints * 0.158655f, -1);
    749         PS_BIN_FOR_VALUE(binHi, cumulative->nums, totalDataPoints * 0.841345f, +1);
    750 
    751         PS_BIN_FOR_VALUE(binL2, cumulative->nums, totalDataPoints * 0.308538f, -1);
    752         PS_BIN_FOR_VALUE(binH2, cumulative->nums, totalDataPoints * 0.691462f, +1);
     747        // also find the 30.8538% (-0.5 sigma) and 69.1462% (+0.5 sigma) points
     748        PS_BIN_FOR_VALUE(binLo, cumulative->nums, totalDataPoints * 0.158655f, -1);
     749        PS_BIN_FOR_VALUE(binHi, cumulative->nums, totalDataPoints * 0.841345f, +1);
     750
     751        PS_BIN_FOR_VALUE(binL2, cumulative->nums, totalDataPoints * 0.308538f, -1);
     752        PS_BIN_FOR_VALUE(binH2, cumulative->nums, totalDataPoints * 0.691462f, +1);
    753753
    754754        psTrace("psLib.math", 6, "The 15.8655%% and 84.1345%% data point bins are (%ld, %ld).\n",
     
    759759        if ((binLo < 0) || (binHi < 0)) {
    760760            psError(PS_ERR_UNKNOWN, false, "Failed to calculate the 15.8655%% and 84.1345%% data points.\n");
    761             goto escape;
     761            goto escape;
    762762        }
    763763
     
    769769
    770770        // find the +0.5 and -0.5 sigma points with linear interpolation.  binLo and binHi are the bins
    771         // containing the value of interest.  constrain the answer to be within the current bin.
    772         // (extrapolation should not be needed and will result in errors)
     771        // containing the value of interest.  constrain the answer to be within the current bin.
     772        // (extrapolation should not be needed and will result in errors)
    773773        float binLoF32, binHiF32, binL2F32, binH2F32;
    774         PS_BIN_INTERPOLATE (binLoF32, cumulative->nums, cumulative->bounds, binL2, totalDataPoints * 0.158655f);
    775         PS_BIN_INTERPOLATE (binHiF32, cumulative->nums, cumulative->bounds, binH2, totalDataPoints * 0.841345f);
    776         PS_BIN_INTERPOLATE (binL2F32, cumulative->nums, cumulative->bounds, binL2, totalDataPoints * 0.308538f);
    777         PS_BIN_INTERPOLATE (binH2F32, cumulative->nums, cumulative->bounds, binH2, totalDataPoints * 0.691462f);
     774        PS_BIN_INTERPOLATE (binLoF32, cumulative->nums, cumulative->bounds, binL2, totalDataPoints * 0.158655f);
     775        PS_BIN_INTERPOLATE (binHiF32, cumulative->nums, cumulative->bounds, binH2, totalDataPoints * 0.841345f);
     776        PS_BIN_INTERPOLATE (binL2F32, cumulative->nums, cumulative->bounds, binL2, totalDataPoints * 0.308538f);
     777        PS_BIN_INTERPOLATE (binH2F32, cumulative->nums, cumulative->bounds, binH2, totalDataPoints * 0.691462f);
    778778
    779779        // report +/- 1 sigma points
     
    810810            // Free the histograms; they will be recreated on the next iteration, with new bounds
    811811            psFree(histogram);
    812             histogram = NULL;
     812            histogram = NULL;
    813813
    814814            psFree(cumulative);
    815             cumulative = NULL;
     815            cumulative = NULL;
    816816        } else {
    817817            // We've got the bin size correct now
     
    823823    // XXX test lines while studying algorithm errors
    824824    // fprintf (stderr, "robust stats test %7.1f +/- %7.1f : %4ld %4ld %4ld %4ld %4ld  : %f %f %f\n",
    825     // stats->robustMedian, stats->robustStdev, 
     825    // stats->robustMedian, stats->robustStdev,
    826826    // binLo, binL2, binMedian, binH2, binHi, binSize, max, min);
    827827
     
    835835    // positions.
    836836    psF32 binLo25F32 = fitQuadraticSearchForYThenReturnX(cumulative->bounds, cumulative->nums,
    837                                                         binLo25, totalDataPoints * 0.25f);
     837                                                        binLo25, totalDataPoints * 0.25f);
    838838    psF32 binHi25F32 = fitQuadraticSearchForYThenReturnX(cumulative->bounds, cumulative->nums,
    839                                                         binHi25, totalDataPoints * 0.75f);
     839                                                        binHi25, totalDataPoints * 0.75f);
    840840    if (isnan(binLo25F32) || isnan(binHi25F32)) {
    841841        psError(PS_ERR_UNKNOWN, false,
    842842                "could not determine the robustUQ: fitQuadraticSearchForYThenReturnX() returned a NAN.\n");
    843         goto escape;
     843        goto escape;
    844844    }
    845845
     
    851851    for (long i = 0 ; i < myVector->n ; i++) {
    852852        if (!mask->data.U8[i] &&
    853             (binLo25F32 <= myVector->data.F32[i]) && (binHi25F32 >= myVector->data.F32[i])) {
     853            (binLo25F32 <= myVector->data.F32[i]) && (binHi25F32 >= myVector->data.F32[i])) {
    854854            N50++;
    855855        }
     
    986986        }
    987987
    988         // select the min and max bins, saturating on the lower and upper end-points
    989         long binMin, binMax;
    990         PS_BIN_FOR_VALUE (binMin, histogram->bounds, guessMean - minFitSigma*guessStdev, -1);
    991         PS_BIN_FOR_VALUE (binMax, histogram->bounds, guessMean + maxFitSigma*guessStdev, +1);
     988        // select the min and max bins, saturating on the lower and upper end-points
     989        long binMin, binMax;
     990        PS_BIN_FOR_VALUE (binMin, histogram->bounds, guessMean - minFitSigma*guessStdev, -1);
     991        PS_BIN_FOR_VALUE (binMax, histogram->bounds, guessMean + maxFitSigma*guessStdev, +1);
    992992
    993993        // Generate the variables that will be used in the Gaussian fitting
     
    11641164        }
    11651165
    1166         // select the min and max bins, saturating on the lower and upper end-points
    1167         long binMin, binMax;
    1168         PS_BIN_FOR_VALUE (binMin, histogram->bounds, guessMean - minFitSigma*guessStdev, -1);
    1169         PS_BIN_FOR_VALUE (binMax, histogram->bounds, guessMean + maxFitSigma*guessStdev, +1);
     1166        // select the min and max bins, saturating on the lower and upper end-points
     1167        long binMin, binMax;
     1168        PS_BIN_FOR_VALUE (binMin, histogram->bounds, guessMean - minFitSigma*guessStdev, -1);
     1169        PS_BIN_FOR_VALUE (binMax, histogram->bounds, guessMean + maxFitSigma*guessStdev, +1);
    11701170
    11711171        // Generate the variables that will be used in the Gaussian fitting
     
    13571357        }
    13581358
    1359         // select the min and max bins, saturating on the lower and upper end-points
    1360         long binMin, binMax;
    1361         PS_BIN_FOR_VALUE (binMin, histogram->bounds, guessMean - minFitSigma*guessStdev, -1);
    1362         PS_BIN_FOR_VALUE (binMax, histogram->bounds, guessMean + maxFitSigma*guessStdev, +1);
     1359        // select the min and max bins, saturating on the lower and upper end-points
     1360        long binMin, binMax;
     1361        PS_BIN_FOR_VALUE (binMin, histogram->bounds, guessMean - minFitSigma*guessStdev, -1);
     1362        PS_BIN_FOR_VALUE (binMax, histogram->bounds, guessMean + maxFitSigma*guessStdev, +1);
    13631363
    13641364        // search for mode (peak of histogram within range mean-2sigma - mean+2sigma
     
    14941494
    14951495            // calculate upper mean & stdev from parabolic fit -- use this as the result
     1496#ifndef PS_NO_TRACE
    14961497            float upperStdev = sqrt(-0.5/poly->coeff[2]);
    14971498            float upperMean = poly->coeff[1]*PS_SQR(upperStdev);
     
    14991500            psTrace("psLib.math", 6, "The upper mean  is %f.\n", upperMean);
    15001501            psTrace("psLib.math", 6, "The upper stdev is %f.\n", upperStdev);
     1502#endif
    15011503
    15021504            psFree (poly);
     
    20432045*****************************************************************************/
    20442046static psF32 fitQuadraticSearchForYThenReturnX(const psVector *xVec,
    2045                                                psVector *yVec,
    2046                                                psS32 binNum,
    2047                                                psF32 yVal
     2047                                               psVector *yVec,
     2048                                               psS32 binNum,
     2049                                               psF32 yVal
    20482050    )
    20492051{
     
    20842086        //
    20852087        if (! (((y->data.F64[0] <= yVal) && (yVal <= y->data.F64[2])) ||
    2086                ((y->data.F64[2] <= yVal) && (yVal <= y->data.F64[0]))) ) {
     2088               ((y->data.F64[2] <= yVal) && (yVal <= y->data.F64[0]))) ) {
    20872089            psError(PS_ERR_BAD_PARAMETER_VALUE, true,
    20882090                    _("Specified yVal, %g, is not within y-range, %g to %g."),
     
    20942096        //
    20952097        if (((y->data.F64[0] < y->data.F64[1]) && !(y->data.F64[1] <= y->data.F64[2])) ||
    2096             ((y->data.F64[0] > y->data.F64[1]) && !(y->data.F64[1] >= y->data.F64[2]))) {
     2098            ((y->data.F64[0] > y->data.F64[1]) && !(y->data.F64[1] >= y->data.F64[2]))) {
    20972099            psError(PS_ERR_UNKNOWN, true,
    20982100                    "This routine must be called with monotonically increasing or decreasing data points.\n");
Note: See TracChangeset for help on using the changeset viewer.