IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Oct 24, 2007, 9:08:24 AM (19 years ago)
Author:
Paul Price
Message:

Cleaned up code. Some required code was hidden in !PS_NO_TRACE
blocks, which prevented optimized compilation.

File:
1 edited

Legend:

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

    r15126 r15370  
    1313 * use ->min and ->max (PS_STAT_USE_RANGE)
    1414 *
    15  *  @version $Revision: 1.218 $ $Name: not supported by cvs2svn $
    16  *  @date $Date: 2007-09-29 23:58:55 $
     15 *  @version $Revision: 1.219 $ $Name: not supported by cvs2svn $
     16 *  @date $Date: 2007-10-24 19:08:24 $
    1717 *
    1818 *  Copyright 2006 IfA, University of Hawaii
     
    152152// static prototypes:
    153153static psF32 minimizeLMChi2Gauss1D(psVector *deriv, const psVector *params, const psVector *coords);
    154 static psF32 fitQuadraticSearchForYThenReturnX(const psVector *xVec, psVector *yVec, psS32 binNum, psF32 yVal);
     154static psF32 fitQuadraticSearchForYThenReturnX(const psVector *xVec, psVector *yVec, psS32 binNum,
     155                                               psF32 yVal);
    155156
    156157/******************************************************************************
     
    309310
    310311    // Allocate temporary vectors for the data.
    311     psVector *outVector = psVectorAllocEmpty(inVector->n, PS_TYPE_F32); // Vector containing the unmasked values
     312    psVector *outVector = psVectorAllocEmpty(inVector->n, PS_TYPE_F32); // Vector containing unmasked values
    312313    psF32 *output = outVector->data.F32; // Dereference the vector
    313314
     
    401402    // If the mean is NAN, then generate a warning and set the stdev to NAN.
    402403    if (isnan(stats->sampleMean)) {
    403         psTrace(TRACE, 5, "WARNING: vectorSampleStdev(): sample mean is NAN.  Setting stats->sampleStdev = NAN.\n");
     404        psTrace(TRACE, 5, "WARNING: vectorSampleStdev(): sample mean is NAN.  "
     405                "Setting stats->sampleStdev = NAN.\n");
    404406        stats->sampleStdev = NAN;
    405407        return true;
     
    443445        // Assume that the user knows what he's doing when he masks out everything --> no error.
    444446        stats->sampleStdev = NAN;
    445         psTrace(TRACE, 5, "WARNING: vectorSampleStdev(): no valid psVector elements (%ld).  Setting stats->sampleStdev = NAN.\n", count);
     447        psTrace(TRACE, 5, "WARNING: vectorSampleStdev(): no valid psVector elements (%ld).  "
     448                "Setting stats->sampleStdev = NAN.\n", count);
    446449        return true;
    447450    }
    448451    if (count == 1) {
    449452        stats->sampleStdev = 0.0;
    450         psTrace(TRACE, 5, "WARNING: vectorSampleStdev(): only one valid psVector elements (%ld).  Setting stats->sampleStdev = 0.0.\n", count);
     453        psTrace(TRACE, 5, "WARNING: vectorSampleStdev(): only one valid psVector elements (%ld).  "
     454                "Setting stats->sampleStdev = 0.0.\n", count);
    451455        return true;
    452456    }
     
    623627        // a) Exclude all values x_i for which |x_i - x| > K * stdev
    624628        if (errors) {
    625             // XXXX if we convert errors to variance, this should square the other terms (A*A faster than sqrt(A))
     629            // XXXX if we convert errors to variance, this should square the other terms (A*A faster than
     630            // sqrt(A))
    626631            for (long j = 0; j < myVector->n; j++) {
    627632                if (!tmpMask->data.U8[j] &&
     
    749754    // Iterate to get the best bin size
    750755    for (int iterate = 1; iterate > 0; iterate++) {
    751         psTrace(TRACE, 6, "-------------------- Iterating on Bin size.  Iteration number %d --------------------\n", iterate);
     756        psTrace(TRACE, 6,
     757                "-------------------- Iterating on Bin size.  Iteration number %d --------------------\n",
     758                iterate);
    752759
    753760        // Get the minimum and maximum values
     
    832839                                                                binMedian, totalDataPoints/2.0);
    833840        if (isnan(stats->robustMedian)) {
    834             psError(PS_ERR_UNKNOWN, false, "Failed to fit a quadratic and calculate the 50-percent position.\n");
     841            psError(PS_ERR_UNKNOWN, false,
     842                    "Failed to fit a quadratic and calculate the 50-percent position.\n");
    835843            goto escape;
    836844        }
    837845        psTrace(TRACE, 6, "Current robust median is %f\n", stats->robustMedian);
    838846
    839         // ADD step 4: Find the bins which contains the 15.8655% (-1 sigma) and 84.1345% (+1 sigma) data points
    840         // also find the 30.8538% (-0.5 sigma) and 69.1462% (+0.5 sigma) points
     847        // ADD step 4: Find the bins which contains the 15.8655% (-1 sigma) and 84.1345% (+1 sigma) data
     848        // points also find the 30.8538% (-0.5 sigma) and 69.1462% (+0.5 sigma) points
    841849        PS_BIN_FOR_VALUE(binLo, cumulative->nums, totalDataPoints * 0.158655f, 0);
    842850        PS_BIN_FOR_VALUE(binHi, cumulative->nums, totalDataPoints * 0.841345f, 0);
     
    858866        }
    859867
    860         // ADD step 4b: Interpolate Sigma (linearly) to find these two positions exactly: these are the 1sigma positions.
     868        // ADD step 4b: Interpolate Sigma (linearly) to find these two positions exactly: these are the 1sigma
     869        // positions.
    861870        psTrace(TRACE, 6, "binLo is %ld.  Nums at that bin and the next are (%.2f, %.2f)\n",
    862871                binLo, cumulative->nums->data.F32[binLo], cumulative->nums->data.F32[binLo+1]);
     
    868877        // (extrapolation should not be needed and will result in errors)
    869878        float binLoF32, binHiF32, binL2F32, binH2F32, binL4F32, binH4F32;
    870         PS_BIN_INTERPOLATE (binLoF32, cumulative->nums, cumulative->bounds, binLo, totalDataPoints * 0.158655f);
    871         PS_BIN_INTERPOLATE (binHiF32, cumulative->nums, cumulative->bounds, binHi, totalDataPoints * 0.841345f);
    872         PS_BIN_INTERPOLATE (binL2F32, cumulative->nums, cumulative->bounds, binL2, totalDataPoints * 0.308538f);
    873         PS_BIN_INTERPOLATE (binH2F32, cumulative->nums, cumulative->bounds, binH2, totalDataPoints * 0.691462f);
    874         PS_BIN_INTERPOLATE (binL4F32, cumulative->nums, cumulative->bounds, binL4, totalDataPoints * 0.022481f);
    875         PS_BIN_INTERPOLATE (binH4F32, cumulative->nums, cumulative->bounds, binH4, totalDataPoints * 0.977519f);
     879        PS_BIN_INTERPOLATE (binLoF32, cumulative->nums, cumulative->bounds, binLo,
     880                            totalDataPoints * 0.158655f);
     881        PS_BIN_INTERPOLATE (binHiF32, cumulative->nums, cumulative->bounds, binHi,
     882                            totalDataPoints * 0.841345f);
     883        PS_BIN_INTERPOLATE (binL2F32, cumulative->nums, cumulative->bounds, binL2,
     884                            totalDataPoints * 0.308538f);
     885        PS_BIN_INTERPOLATE (binH2F32, cumulative->nums, cumulative->bounds, binH2,
     886                            totalDataPoints * 0.691462f);
     887        PS_BIN_INTERPOLATE (binL4F32, cumulative->nums, cumulative->bounds, binL4,
     888                            totalDataPoints * 0.022481f);
     889        PS_BIN_INTERPOLATE (binH4F32, cumulative->nums, cumulative->bounds, binH4,
     890                            totalDataPoints * 0.977519f);
    876891
    877892        // report +/- 1 sigma points
     
    891906        float sigma4 = (binH4F32 - binL4F32) / 4.0;
    892907
    893         // take the smallest of the three: if we have a clump with wide outliers, sigma2 and
    894         // sigma4 will be biased high; if we have a bi-modal distribution, sigma1 and sigma2
    895         // will be biased high.
    896         sigma = PS_MIN (sigma1, PS_MIN (sigma2, sigma4));
     908        // take the smallest of the three: if we have a clump with wide outliers, sigma2 and
     909        // sigma4 will be biased high; if we have a bi-modal distribution, sigma1 and sigma2
     910        // will be biased high.
     911        sigma = PS_MIN (sigma1, PS_MIN (sigma2, sigma4));
    897912
    898913        psTrace(TRACE, 6, "The 2x sigma is %f.\n", sigma1);
     
    13281343
    13291344        if (poly->coeff[2] >= 0.0) {
    1330             psTrace(TRACE, 6, "Parabolic fit results: %f + %f x + %f x^2\n", poly->coeff[0], poly->coeff[1], poly->coeff[2]);
     1345            psTrace(TRACE, 6, "Parabolic fit results: %f + %f x + %f x^2\n",
     1346                    poly->coeff[0], poly->coeff[1], poly->coeff[2]);
    13311347            psError(PS_ERR_UNKNOWN, false, "fit did not converge\n");
    13321348            psFree(x);
     
    13441360            done = true;
    13451361        } else {
    1346             psTrace(TRACE, 6, "Parabolic fit results: %f + %f x + %f x^2\n", poly->coeff[0], poly->coeff[1], poly->coeff[2]);
     1362            psTrace(TRACE, 6, "Parabolic fit results: %f + %f x + %f x^2\n",
     1363                    poly->coeff[0], poly->coeff[1], poly->coeff[2]);
    13471364            psTrace(TRACE, 6, "The new guess mean  is %f.\n", guessMean);
    13481365            psTrace(TRACE, 6, "The new guess stdev is %f.\n", guessStdev);
     
    14901507                valPeak = histogram->nums->data.F32[binPeak];
    14911508            }
    1492             psTrace (TRACE, 6, "(%f = %.0f) ", histogram->bounds->data.F32[i], histogram->nums->data.F32[i]);
    1493         }
    1494         psTrace (TRACE, 6, "\n");
     1509            psTrace (TRACE, 6, "(%f = %.0f) ", histogram->bounds->data.F32[i], histogram->nums->data.F32[i]);
     1510        }
     1511        psTrace (TRACE, 6, "\n");
    14951512
    14961513        // assume a reasonably well-defined gaussian-like population; run from peak out until val < 0.25*peak
    14971514
    14981515        psTrace(TRACE, 6, "The clipped numBins is %ld\n", binMax - binMin);
    1499         psTrace(TRACE, 6, "The clipped min is %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binMin), binMin);
    1500         psTrace(TRACE, 6, "The clipped max is %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binMax - 1), binMax - 1);
    1501         psTrace(TRACE, 6, "The clipped peak is %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binPeak), binPeak);
     1516        psTrace(TRACE, 6, "The clipped min is %f (%ld)\n",
     1517                PS_BIN_MIDPOINT(histogram, binMin), binMin);
     1518        psTrace(TRACE, 6, "The clipped max is %f (%ld)\n",
     1519                PS_BIN_MIDPOINT(histogram, binMax - 1), binMax - 1);
     1520        psTrace(TRACE, 6, "The clipped peak is %f (%ld)\n",
     1521                PS_BIN_MIDPOINT(histogram, binPeak), binPeak);
    15021522        psTrace(TRACE, 6, "The clipped peak value is %f\n", histogram->nums->data.F32[binPeak]);
    15031523
     
    15131533                }
    15141534            }
    1515             psTrace(TRACE, 6, "Lower bound for lower half: %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binS), binS);
    1516             psTrace(TRACE, 6, "Upper bound for lower half: %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binE), binE);
     1535            psTrace(TRACE, 6, "Lower bound for lower half: %f (%ld)\n",
     1536                    PS_BIN_MIDPOINT(histogram, binS), binS);
     1537            psTrace(TRACE, 6, "Upper bound for lower half: %f (%ld)\n",
     1538                    PS_BIN_MIDPOINT(histogram, binE), binE);
    15171539
    15181540            psVector *y = psVectorAllocEmpty(binE - binS, PS_TYPE_F32); // Vector of coordinates
     
    15231545                    continue;
    15241546                x->data.F32[j] = PS_BIN_MIDPOINT(histogram, i);
    1525                 y->data.F32[j] = log(histogram->nums->data.F32[i]); // note this is the natural log: expected distribution is A exp(-(x-xo)^2/2sigma^2)
     1547                // note this is the natural log: expected distribution is A exp(-(x-xo)^2/2sigma^2)
     1548                y->data.F32[j] = log(histogram->nums->data.F32[i]);
    15261549                j++;
    15271550            }
     
    15441567
    15451568            if (poly->coeff[2] >= 0.0) {
    1546                 psTrace(TRACE, 6, "Failed parabolic fit: %f + %f x + %f x^2\n", poly->coeff[0], poly->coeff[1], poly->coeff[2]);
     1569                psTrace(TRACE, 6, "Failed parabolic fit: %f + %f x + %f x^2\n",
     1570                        poly->coeff[0], poly->coeff[1], poly->coeff[2]);
    15471571                psFree(poly);
    15481572                psFree(histogram);
    15491573                psFree(statsMinMax);
    15501574
    1551                 // sometimes, the guessStdev is much too large.  in this case, the entire real population tends to be found
    1552                 // in a single bin.  make one attempt to recover by dropping the guessStdev down by a jump and trying again
     1575                // sometimes, the guessStdev is much too large.  in this case, the entire real population
     1576                // tends to be found in a single bin.  make one attempt to recover by dropping the guessStdev
     1577                // down by a jump and trying again
    15531578                if (iteration == 0) {
    15541579                    guessStdev = 0.25*guessStdev;
     
    15681593                done = true;
    15691594            }
    1570             psTrace(TRACE, 6, "Parabolic Lower fit results: %f + %f x + %f x^2\n", poly->coeff[0], poly->coeff[1], poly->coeff[2]);
     1595            psTrace(TRACE, 6, "Parabolic Lower fit results: %f + %f x + %f x^2\n",
     1596                    poly->coeff[0], poly->coeff[1], poly->coeff[2]);
    15711597            psTrace(TRACE, 6, "The lower mean  is %f.\n", guessMean);
    15721598            psTrace(TRACE, 6, "The lower stdev is %f.\n", guessStdev);
     
    15751601        }
    15761602
    1577         // for test, measure the same result for the upper section
     1603        // for test, measure the same result for the upper section
    15781604        {
    15791605            // fit the upper half of the distribution
     
    15871613                }
    15881614            }
    1589             psTrace(TRACE, 6, "Lower bound for upper half: %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binS), binS);
    1590             psTrace(TRACE, 6, "Upper bound for upper half: %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binE), binE);
     1615            psTrace(TRACE, 6, "Lower bound for upper half: %f (%ld)\n",
     1616                    PS_BIN_MIDPOINT(histogram, binS), binS);
     1617            psTrace(TRACE, 6, "Upper bound for upper half: %f (%ld)\n",
     1618                    PS_BIN_MIDPOINT(histogram, binE), binE);
    15911619
    15921620            psVector *y = psVectorAllocEmpty(binE - binS, PS_TYPE_F32); // Vector of coordinates
     
    15971625                    continue;
    15981626                x->data.F32[j] = PS_BIN_MIDPOINT(histogram, i);
    1599                 y->data.F32[j] = log(histogram->nums->data.F32[i]); // note this is the natural log: expected distribution is A exp(-(x-xo)^2/2sigma^2)
     1627                // note this is the natural log: expected distribution is A exp(-(x-xo)^2/2sigma^2)
     1628                y->data.F32[j] = log(histogram->nums->data.F32[i]);
    16001629                j++;
    16011630            }
     
    16181647
    16191648            // calculate upper mean & stdev from parabolic fit -- ignore this value
    1620 #ifndef PS_NO_TRACE
    16211649            float upperStdev = sqrt(-0.5/poly->coeff[2]);
    16221650            float upperMean = poly->coeff[1]*PS_SQR(upperStdev);
    1623             psTrace(TRACE, 6, "Parabolic Upper fit results: %f + %f x + %f x^2\n", poly->coeff[0], poly->coeff[1], poly->coeff[2]);
     1651#ifndef PS_NO_TRACE
     1652            psTrace(TRACE, 6, "Parabolic Upper fit results: %f + %f x + %f x^2\n",
     1653                    poly->coeff[0], poly->coeff[1], poly->coeff[2]);
    16241654            psTrace(TRACE, 6, "The upper mean  is %f.\n", upperMean);
    16251655            psTrace(TRACE, 6, "The upper stdev is %f.\n", upperStdev);
    16261656#endif
    16271657
    1628             // if the resulting value is outside of the range binMin - binMax, use the upper value
    1629             if (done && (guessMean > PS_BIN_MIDPOINT(histogram, binMax - 1))) {
    1630                 guessMean = upperMean;
    1631                 guessStdev = upperStdev;
    1632             }
     1658            // if the resulting value is outside of the range binMin - binMax, use the upper value
     1659            if (done && (guessMean > PS_BIN_MIDPOINT(histogram, binMax - 1))) {
     1660                guessMean = upperMean;
     1661                guessStdev = upperStdev;
     1662            }
    16331663
    16341664            psFree (poly);
     
    17181748        // Calculate the number of bins.
    17191749        // XXX can we calculate the binMin, binMax **before** building this histogram?
    1720         // if the range is too absurd, adjust numBins & binSize
     1750        // if the range is too absurd, adjust numBins & binSize
    17211751        long numBins = PS_MAX (50, PS_MIN (10000, (max - min) / binSize));
    1722         binSize = (max - min) / (float) numBins;
     1752        binSize = (max - min) / (float) numBins;
    17231753        psTrace(TRACE, 6, "The new min/max values are (%f, %f).\n", min, max);
    17241754        psTrace(TRACE, 6, "The new bin size is %f.\n", binSize);
     
    17741804                valPeak = histogram->nums->data.F32[binPeak];
    17751805            }
    1776             psTrace (TRACE, 6, "(%f = %.0f) ", histogram->bounds->data.F32[i], histogram->nums->data.F32[i]);
    1777         }
    1778         psTrace (TRACE, 6, "\n");
     1806            psTrace (TRACE, 6, "(%f = %.0f) ", histogram->bounds->data.F32[i], histogram->nums->data.F32[i]);
     1807        }
     1808        psTrace (TRACE, 6, "\n");
    17791809
    17801810        // assume a reasonably well-defined gaussian-like population; run from peak out until val < 0.25*peak
     
    17821812        psTrace(TRACE, 6, "The clipped numBins is %ld\n", binMax - binMin);
    17831813        psTrace(TRACE, 6, "The clipped min is %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binMin), binMin);
    1784         psTrace(TRACE, 6, "The clipped max is %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binMax - 1), binMax - 1);
     1814        psTrace(TRACE, 6, "The clipped max is %f (%ld)\n",
     1815                PS_BIN_MIDPOINT(histogram, binMax - 1), binMax - 1);
    17851816        psTrace(TRACE, 6, "The clipped peak is %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binPeak), binPeak);
    17861817        psTrace(TRACE, 6, "The clipped peak value is %f\n", histogram->nums->data.F32[binPeak]);
     
    18041835                }
    18051836            }
    1806             psTrace(TRACE, 6, "Lower bound for lower half: %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binS), binS);
    1807             psTrace(TRACE, 6, "Upper bound for lower half: %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binE), binE);
     1837            psTrace(TRACE, 6, "Lower bound for lower half: %f (%ld)\n",
     1838                    PS_BIN_MIDPOINT(histogram, binS), binS);
     1839            psTrace(TRACE, 6, "Upper bound for lower half: %f (%ld)\n",
     1840                    PS_BIN_MIDPOINT(histogram, binE), binE);
    18081841
    18091842            psVector *y = psVectorAllocEmpty(binE - binS, PS_TYPE_F32); // Vector of coordinates
     
    18141847                    continue;
    18151848                x->data.F32[j] = PS_BIN_MIDPOINT(histogram, i);
    1816                 y->data.F32[j] = log(histogram->nums->data.F32[i]); // note this is the natural log: expected distribution is A exp(-(x-xo)^2/2sigma^2)
     1849                // note this is the natural log: expected distribution is A exp(-(x-xo)^2/2sigma^2)
     1850                y->data.F32[j] = log(histogram->nums->data.F32[i]);
    18171851                j++;
    18181852            }
     
    18351869
    18361870            if (poly->coeff[2] >= 0.0) {
    1837                 psTrace(TRACE, 6, "Failed parabolic fit: %f + %f x + %f x^2\n", poly->coeff[0], poly->coeff[1], poly->coeff[2]);
     1871                psTrace(TRACE, 6, "Failed parabolic fit: %f + %f x + %f x^2\n",
     1872                        poly->coeff[0], poly->coeff[1], poly->coeff[2]);
    18381873                psFree(poly);
    18391874                psFree(histogram);
    18401875                psFree(statsMinMax);
    18411876
    1842                 // sometimes, the guessStdev is much too large.  in this case, the entire real population tends to be found
    1843                 // in a single bin.  make one attempt to recover by dropping the guessStdev down by a jump and trying again
     1877                // sometimes, the guessStdev is much too large.  in this case, the entire real population
     1878                // tends to be found in a single bin.  make one attempt to recover by dropping the guessStdev
     1879                // down by a jump and trying again
    18441880                if (iteration == 0) {
    18451881                    guessStdev = 0.25*guessStdev;
     
    18591895                done = true;
    18601896            }
    1861             psTrace(TRACE, 6, "Parabolic Lower fit results: %f + %f x + %f x^2\n", poly->coeff[0], poly->coeff[1], poly->coeff[2]);
     1897            psTrace(TRACE, 6, "Parabolic Lower fit results: %f + %f x + %f x^2\n",
     1898                    poly->coeff[0], poly->coeff[1], poly->coeff[2]);
    18621899            psTrace(TRACE, 6, "The lower mean  is %f.\n", guessMean);
    18631900            psTrace(TRACE, 6, "The lower stdev is %f.\n", guessStdev);
     
    18661903        }
    18671904
    1868         // if we converge on a solution outside the range binMin - binMax, use a more conservative range
    1869         float minValue = PS_BIN_MIDPOINT(histogram, binMin);
    1870         float maxValue = PS_BIN_MIDPOINT(histogram, binMax - 1);
    1871 
    1872         if (done && ((guessMean < minValue) || (guessMean > maxValue))) {
     1905        // if we converge on a solution outside the range binMin - binMax, use a more conservative range
     1906        float minValue = PS_BIN_MIDPOINT(histogram, binMin);
     1907        float maxValue = PS_BIN_MIDPOINT(histogram, binMax - 1);
     1908
     1909        if (done && ((guessMean < minValue) || (guessMean > maxValue))) {
    18731910            psTrace(TRACE, 6, "Inconsistent result, re-trying the fit\n");
    18741911
     
    18901927                }
    18911928            }
    1892             psTrace(TRACE, 6, "Lower bound for symmetric range: %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binS), binS);
    1893             psTrace(TRACE, 6, "Upper bound for symmetric range: %f (%ld)\n", PS_BIN_MIDPOINT(histogram, binE), binE);
     1929            psTrace(TRACE, 6, "Lower bound for symmetric range: %f (%ld)\n",
     1930                    PS_BIN_MIDPOINT(histogram, binS), binS);
     1931            psTrace(TRACE, 6, "Upper bound for symmetric range: %f (%ld)\n",
     1932                    PS_BIN_MIDPOINT(histogram, binE), binE);
    18941933
    18951934            psVector *y = psVectorAllocEmpty(binE - binS, PS_TYPE_F32); // Vector of coordinates
     
    19001939                    continue;
    19011940                x->data.F32[j] = PS_BIN_MIDPOINT(histogram, i);
    1902                 y->data.F32[j] = log(histogram->nums->data.F32[i]); // note this is the natural log: expected distribution is A exp(-(x-xo)^2/2sigma^2)
     1941                // note this is the natural log: expected distribution is A exp(-(x-xo)^2/2sigma^2)
     1942                y->data.F32[j] = log(histogram->nums->data.F32[i]);
    19031943                j++;
    19041944            }
     
    19211961
    19221962            // calculate upper mean & stdev from parabolic fit -- ignore this value
    1923 #ifndef PS_NO_TRACE
    19241963            guessStdev = sqrt(-0.5/poly->coeff[2]);
    19251964            guessMean = poly->coeff[1]*PS_SQR(guessStdev);
    1926             psTrace(TRACE, 6, "Parabolic Symmetric fit results: %f + %f x + %f x^2\n", poly->coeff[0], poly->coeff[1], poly->coeff[2]);
     1965#ifndef PS_NO_TRACE
     1966            psTrace(TRACE, 6, "Parabolic Symmetric fit results: %f + %f x + %f x^2\n",
     1967                    poly->coeff[0], poly->coeff[1], poly->coeff[2]);
    19271968            psTrace(TRACE, 6, "The symmetric mean  is %f.\n", guessMean);
    19281969            psTrace(TRACE, 6, "The symmetric stdev is %f.\n", guessStdev);
    19291970#endif
    19301971
    1931             // if we converge on a solution outside the range binMin - binMax, use a more conservative range
    1932             float minValueSym = PS_BIN_MIDPOINT(histogram, binS);
    1933             float maxValueSym = PS_BIN_MIDPOINT(histogram, binE - 1);
    1934 
    1935             // saturate on min or max value
    1936             if (guessMean < minValueSym) {
    1937                 guessMean = minValueSym;
    1938                 psTrace(TRACE, 6, "The symmetric mean is out of bounds, saturating to %f.\n", guessMean);
    1939             }
    1940                
    1941             // saturate on min or max value
    1942             if (guessMean > maxValueSym) {
    1943                 guessMean = maxValueSym;
    1944                 psTrace(TRACE, 6, "The symmetric mean is out of bounds, saturating to %f.\n", guessMean);
    1945             }
     1972            // if we converge on a solution outside the range binMin - binMax, use a more conservative range
     1973            float minValueSym = PS_BIN_MIDPOINT(histogram, binS);
     1974            float maxValueSym = PS_BIN_MIDPOINT(histogram, binE - 1);
     1975
     1976            // saturate on min or max value
     1977            if (guessMean < minValueSym) {
     1978                guessMean = minValueSym;
     1979                psTrace(TRACE, 6, "The symmetric mean is out of bounds, saturating to %f.\n", guessMean);
     1980            }
     1981
     1982            // saturate on min or max value
     1983            if (guessMean > maxValueSym) {
     1984                guessMean = maxValueSym;
     1985                psTrace(TRACE, 6, "The symmetric mean is out of bounds, saturating to %f.\n", guessMean);
     1986            }
    19461987
    19471988            psFree (poly);
     
    19962037        // However, it is also not used anywhere, yet.
    19972038        //
    1998         psWarning("WARNING: vectorSmoothHistGaussian() on non-uniform histograms has not been tested or used.\n");
     2039        psWarning("WARNING: vectorSmoothHistGaussian() on non-uniform "
     2040                  "histograms has not been tested or used.\n");
    19992041
    20002042        for (long i = 0; i < numBins; i++) {
     
    26072649        if (!psVectorFitPolynomial1D(myPoly, NULL, 0, y, NULL, x)) {
    26082650            psError(PS_ERR_UNEXPECTED_NULL, false,
    2609                     _("Failed to fit a 1-dimensional polynomial to the three specified data points.  Returning NAN."));
     2651                    _("Failed to fit a 1-dimensional polynomial to the three specified data points.  "
     2652                      "Returning NAN."));
    26102653            psFree(x);
    26112654            psFree(y);
Note: See TracChangeset for help on using the changeset viewer.