IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 6322


Ignore:
Timestamp:
Feb 3, 2006, 12:05:22 PM (20 years ago)
Author:
gusciora
Message:

Misc code cleaning

Location:
trunk/psLib/src/math
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/psLib/src/math/psMinimizeLMM.c

    r6226 r6322  
    1010 *  @author EAM, IfA
    1111 *
    12  *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
    13  *  @date $Date: 2006-01-27 20:08:58 $
     12 *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
     13 *  @date $Date: 2006-02-03 22:05:22 $
    1414 *
    1515 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    4848
    4949// XXX EAM : can we use static copies of LUv, LUm, A?
    50 // XXX: Add trace messages, check return codes.
    5150psBool p_psMinLM_GuessABP(
    5251    psImage  *Alpha,
    … …  
    7574    A = psImageCopy(NULL, alpha, PS_TYPE_F64);
    7675    for (int j = 0; j < params->n; j++) {
    77         if ((paramMask != NULL) && (paramMask->data.U8[j]))
     76        if ((paramMask != NULL) && (paramMask->data.U8[j])) {
    7877            continue;
     78        }
    7979        A->data.F64[j][j] = alpha->data.F64[j][j] * (1.0 + lambda);
    8080    }
    … …  
    8383    // these operations do not modify the input values (creates LUm, LUv)
    8484    LUm   = psMatrixLUD(NULL, &LUv, A);
     85    if (LUm == NULL) {
     86        psError(PS_ERR_UNKNOWN, false, "psMatrixLUD() returned NULL\n");
     87    }
    8588    Beta  = psMatrixLUSolve(Beta, LUm, beta, LUv);
     89    if (Beta == NULL) {
     90        psError(PS_ERR_UNKNOWN, false, "psMatrixLUSolve() returned NULL\n");
     91    }
    8692    Alpha = psMatrixInvert(Alpha, A, &det);
    87 
     93    psFree(LUm);
     94    psFree(LUv);
     95    psFree(A);
     96
     97    if (Alpha == NULL) {
     98        psError(PS_ERR_UNKNOWN, false, "psMatrixInvert() returned NULL\n");
     99        return(false);
     100    }
    88101    # else
    89102        // gauss-jordan version
    … …  
    94107    Alpha = psImageCopy(Alpha, alpha, PS_TYPE_F64);
    95108    for (int j = 0; j < params->n; j++) {
    96         if ((paramMask != NULL) && (paramMask->data.U8[j]))
     109        if ((paramMask != NULL) && (paramMask->data.U8[j])) {
    97110            continue;
     111        }
    98112        Alpha->data.F64[j][j] = alpha->data.F64[j][j] * (1.0 + lambda);
    99113    }
    100114
    101     // XXX: Check error codes!
    102     psGaussJordan(Alpha, Beta);
    103     psFree(A);
    104     psFree(LUm);
    105     psFree(LUv);
     115    if (false == psGaussJordan(Alpha, Beta)) {
     116        psError(PS_ERR_UNKNOWN, false, "psMatrixInvert() returned NULL\n");
     117        return(false);
     118    }
    106119    # endif
    107120
    108121    // apply Beta to get new Params values
    109122    for (int j = 0; j < params->n; j++) {
    110         if ((paramMask != NULL) && (paramMask->data.U8[j]))
     123        if ((paramMask != NULL) && (paramMask->data.U8[j])) {
    111124            continue;
     125        }
    112126        // Params->data.F32[j] = params->data.F32[j] - Beta->data.F64[j];
    113127        // compare Beta to beta limits
    … …  
    126140        }
    127141    }
    128     # if (USE_LU_DECOMP)
    129         psFree(A);
    130     psFree(LUm);
    131     psFree(LUv);
    132     # endif
    133142
    134143    return(true);
    … …  
    136145
    137146
    138 // XXX: Add trace messages, check return codes.
    139147bool psMinimizeGaussNewtonDelta(
    140148    psVector *delta,
    … …  
    146154    psMinimizeLMChi2Func func)
    147155{
     156    psTrace(__func__, 3, "---- %s() begin ----\n", __func__);
    148157    // allocate internal arrays (current vs Guess)
    149158    psImage  *alpha  = psImageAlloc (params->n, params->n, PS_TYPE_F64);
    … …  
    152161    psVector *Params = psVectorAlloc(params->n, PS_TYPE_F64);
    153162    psVector *dy     = NULL;
     163    psBool rc = true;
    154164
    155165    // the user provides the error or NULL.  we need to convert
    … …  
    162172    }
    163173
    164     p_psMinLM_SetABX(alpha, beta, params, paramMask, x, y, dy, func);
    165     p_psMinLM_GuessABP(Alpha, delta, Params, alpha, beta, params, paramMask, NULL, NULL, NULL, 0.0);
     174    psF64 rcF64 = p_psMinLM_SetABX(alpha, beta, params, paramMask, x, y, dy, func);
     175    if (isnan(rcF64)) {
     176        psError(PS_ERR_UNKNOWN, false, "p_psMinLM_SetABX() retruned a NAN.\n");
     177        rc = false;
     178    }
     179    psTrace(__func__, 5, "p_psMinLM_SetABX() was succesful\n", __func__);
     180
     181    psBool rcBool = p_psMinLM_GuessABP(Alpha, delta, Params, alpha, beta, params, paramMask, NULL, NULL, NULL, 0.0);
     182    if (rcBool == false) {
     183        psError(PS_ERR_UNKNOWN, false, "p_psMinLM_GuessABP() retruned FALSE.\n");
     184        rc = false;
     185    }
     186    psTrace(__func__, 5, "p_psMinLM_GuessABP() was succesful\n", __func__);
    166187
    167188    psFree(alpha);
    … …  
    172193        psFree(dy);
    173194    }
    174     return (true);
     195    psTrace(__func__, 3, "---- %s() end ----\n", __func__);
     196    return(rc);
    175197}
    176198
    … …  
    235257
    236258        for (psS32 j = 0; j < params->n; j++) {
    237             if ((paramMask != NULL) && (paramMask->data.U8[j]))
     259            if ((paramMask != NULL) && (paramMask->data.U8[j])) {
    238260                continue;
     261            }
    239262            weight = deriv->data.F32[j] * dy->data.F32[i];
    240263            for (psS32 k = 0; k <= j; k++) {
    241                 if ((paramMask != NULL) && (paramMask->data.U8[k]))
     264                if ((paramMask != NULL) && (paramMask->data.U8[k])) {
    242265                    continue;
     266                }
    243267                alpha->data.F64[j][k] += weight * deriv->data.F32[k];
    244268            }
    … …  
    394418    }
    395419    if (psTraceGetLevel (__func__) >= 6) {
    396         //XXX:  p_psVectorPrintRow(psTraceGetDestination(), Params, "params guess");
     420        psTrace(__func__, 6, "The current Param vector: \n");
     421        for (psS32 i = 0 ; i < Params->n ; i++) {
     422            psTrace(__func__, 6, "Params[%d] is %f\n", Params->data.F32[i]);
     423        }
    397424    }
    398425
    … …  
    416443        }
    417444        if (psTraceGetLevel(__func__) >= 6) {
    418             //XXX: p_psVectorPrintRow(psTraceGetDestination(), Params, "params guess");
     445            if (psTraceGetLevel (__func__) >= 6) {
     446                psTrace(__func__, 6, "The current Param vector: \n");
     447                for (psS32 i = 0 ; i < Params->n ; i++) {
     448                    psTrace(__func__, 6, "Params[%d] is %f\n", Params->data.F32[i]);
     449                }
     450            }
    419451        }
    420452
    … …  
    628660bool psMemCheckConstrain(psPtr tmp)
    629661{
    630     return( psMemGetDeallocator(tmp) == (psFreeFunc)constrainFree );
    631 }
     662    return(psMemGetDeallocator(tmp) == (psFreeFunc) constrainFree);
     663}

  • trunk/psLib/src/math/psMinimizePowell.c

    r6186 r6322  
    99 *  @author GLG, MHPCC
    1010 *
    11  *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
    12  *  @date $Date: 2006-01-23 22:25:31 $
     11 *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
     12 *  @date $Date: 2006-02-03 22:05:22 $
    1313 *
    1414 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    15  *
    16  *  XXX: must follow coding name standards on local functions.
    17  *  XXX: put local functions in front.
    1815 *
    1916 */
    … …  
    6966/*****************************************************************************/
    7067/* GLOBAL VARIABLES                                                          */
    71 /* XXX: Do these conform to code standard?         */
    7268/*****************************************************************************/
    7369static psMinimizeChi2PowellFunc Chi2PowellFunc = NULL;
    74 static psVector *myValue;
    75 static psVector *myError;
     70static psVector *PowellValue;
     71static psVector *PowellError;
    7672/*****************************************************************************/
    7773/* FILE STATIC VARIABLES                                                     */
    … …  
    9793Repeat this process until a local minimum is found.
    9894 
    99 XXX:
    100 new algorithm: 
     95XXX: new algorithm: 
    10196start at x=0, expand in one direction until the function
    10297decreases.  Then you have two points in the bracket.  Keep going until it
    … …  
    10499direction.
    105100 
    106 XXX:
    107 This is F32 only.  Must add
    108 F64 support (actually, make the defaults F64,
     101XXX: This is F32 only.  Must add F64 support (actually, make the defaults F64,
    109102and convert F32 vectors to F64).
    110103 
    111 XXX:
    112 output bracket vector should be an input as well.
     104XXX: output bracket vector should be an input as well.
    113105*****************************************************************************/
    114106psVector *p_psDetermineBracket(
    … …  
    355347 *****************************************************************************/
    356348#define PS_LINEMIN_MAX_ITERATIONS 30
    357 psF32 p_psLineMin(psMinimization *min,
    358                   psVector *params,
    359                   psVector *line,
    360                   const psVector *paramMask,
    361                   const psArray *coords,
    362                   psMinimizePowellFunc func)
     349psF32 p_psLineMin(
     350    psMinimization *min,
     351    psVector *params,
     352    psVector *line,
     353    const psVector *paramMask,
     354    const psArray *coords,
     355    psMinimizePowellFunc func)
    363356{
    364357    PS_ASSERT_PTR_NON_NULL(min, NAN);
    … …  
    502495#define PS_MINIMIZE_POWELL_LINEMIN_ERROR_TOLERANCE 0.01
    503496
    504 bool psMinimizePowell(psMinimization *min,
    505                       psVector *params,
    506                       const psVector *paramMask,
    507                       const psArray *coords,
    508                       psMinimizePowellFunc func)
     497bool psMinimizePowell(
     498    psMinimization *min,
     499    psVector *params,
     500    const psVector *paramMask,
     501    const psArray *coords,
     502    psMinimizePowellFunc func)
    509503{
    510504    PS_ASSERT_PTR_NON_NULL(min, NULL);
    … …  
    696690This functions uses global variables to receive the function pointer, the
    697691data values, and the data errors.
     692 
    698693XXX: This is F32 only.  Must implement F64.
    699694 *****************************************************************************/
    700 static psF32 myPowellChi2Func(const psVector *params,
    701                               const psArray *coords)
     695static psF32 myPowellChi2Func(
     696    const psVector *params,
     697    const psArray *coords)
    702698{
    703699    psTrace(__func__, 4, "---- myPowellChi2Func() begin ----\n");
    704700    PS_ASSERT_VECTOR_NON_NULL(params, NAN);
    705701    PS_ASSERT_VECTOR_NON_EMPTY(params, NAN);
    706     PS_ASSERT_VECTOR_NON_NULL(myValue, NAN);
    707     PS_ASSERT_VECTOR_NON_EMPTY(myValue, NAN);
     702    PS_ASSERT_VECTOR_NON_NULL(PowellValue, NAN);
     703    PS_ASSERT_VECTOR_NON_EMPTY(PowellValue, NAN);
    708704    PS_ASSERT_PTR_NON_NULL(coords, NAN);
    709705
    … …  
    714710
    715711    tmp = Chi2PowellFunc(params, coords);
    716     if (myError == NULL) {
     712    if (PowellError == NULL) {
    717713        for (i=0;i<coords->n;i++) {
    718             d = (tmp->data.F32[i] - myValue->data.F32[i]);
     714            d = (tmp->data.F32[i] - PowellValue->data.F32[i]);
    719715            chi2+= d * d;
    720716        }
    721717    } else {
    722718        for (i=0;i<coords->n;i++) {
    723             d = (tmp->data.F32[i] - myValue->data.F32[i]) / myError->data.F32[i];
     719            d = (tmp->data.F32[i] - PowellValue->data.F32[i]) / PowellError->data.F32[i];
    724720            chi2+= d * d;
    725721        }
    … …  
    741737psMinimizePowell().
    742738 *****************************************************************************/
    743 bool psMinimizeChi2Powell(psMinimization *min,
    744                           psVector *params,
    745                           const psVector *paramMask,
    746                           const psArray *coords,
    747                           const psVector *value,
    748                           const psVector *error,
    749                           psMinimizeChi2PowellFunc model)
     739bool psMinimizeChi2Powell(
     740    psMinimization *min,
     741    psVector *params,
     742    const psVector *paramMask,
     743    const psArray *coords,
     744    const psVector *value,
     745    const psVector *error,
     746    psMinimizeChi2PowellFunc model)
    750747{
    751     myValue = (psVector *) value;
    752     myError = (psVector *) error;
     748    PowellValue = (psVector *) value;
     749    PowellError = (psVector *) error;
    753750
    754751    Chi2PowellFunc = model;

  • trunk/psLib/src/math/psStats.c

    r6315 r6322  
    1616 * use ->min and ->max (PS_STAT_USE_RANGE)
    1717 *
    18  *  @version $Revision: 1.166 $ $Name: not supported by cvs2svn $
    19  *  @date $Date: 2006-02-03 01:30:14 $
     18 *  @version $Revision: 1.167 $ $Name: not supported by cvs2svn $
     19 *  @date $Date: 2006-02-03 22:05:22 $
    2020 *
    2121 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    303303max of the input vector.  If there was a problem with the max calculation,
    304304this routine sets stats->max to NAN.
    305  
    306 XXX: Do we need to factor errors into it?
    307305 *****************************************************************************/
    308306psS32 p_psVectorMax(const psVector* myVector,
    … …  
    549547median of the input vector.  Returns true on success (including if there were
    550548no valid input vector elements).
    551  
    552 XXX: Use static vectors for sort arrays.
    553549 *****************************************************************************/
    554550bool p_psVectorSampleMedian(const psVector* myVector,
    … …  
    783779{
    784780    psTrace(__func__, 4, "---- %s() begin ----\n", __func__);
    785     psVector* unsortedVector = NULL;    // Temporary vector
    786     psVector* sortedVector = NULL;      // Temporary vector
    787     psS32 i = 0;                  // Loop index variable
    788     psS32 count = 0;              // # of points in this mean.
    789     psS32 nValues = 0;            // # data points
     781    psVector* unsortedVector = NULL;
     782    psVector* sortedVector = NULL;
     783    psS32 i = 0;
     784    psS32 count = 0;
     785    psS32 nValues = 0;
    790786
    791787    // Determine how many data points fit inside this min/max range
    792     // and are not maxed, IF the maskVector is not NULL.
     788    // and are not masked, if the maskVector is not NULL.
    793789    nValues = p_psVectorNValues(myVector, maskVector, maskVal, stats);
    794790
    … …  
    796792    unsortedVector = psVectorAlloc(nValues, PS_TYPE_F32);
    797793
    798     // Determine if we must only use data points within a min/max range.
    799794    if (stats->options & PS_STAT_USE_RANGE) {
    800795        // Store all non-masked data points within the min/max range
    … …  
    1001996    psTrace(__func__, 4, "---- %s() begin ----\n", __func__);
    1002997    psTrace(__func__, 4, "Trace level is %d\n", psTraceGetLevel(__func__));
    1003     psF32 clippedMean = 0.0;    // self-explanatory
    1004     psF32 clippedStdev = 0.0;   // self-explanatory
    1005     psVector* tmpMask = NULL;   // Temporary vector for masks during iterations.
    1006     static psStats *statsTmp = NULL;   // Temporary psStats struct.
    1007     psS32 rc = 0;               // Return code.
     998    psF32 clippedMean = 0.0;
     999    psF32 clippedStdev = 0.0;
     1000    psVector* tmpMask = NULL;
     1001    static psStats *statsTmp = NULL;
     1002    psS32 rc = 0;
    10081003
    10091004    // Ensure that stats->clipIter is within the proper range.
    … …  
    11891184parameter).
    11901185 
    1191 XXX: After you fit the polynomial, solve for X analytically.
    1192  
    1193 XXX: the vectors do not have to be the same length.  Must insert the proper
    1194 tests to ensure that binNum is within acceptable ranges for both vectors.
    1195  
    1196 XXX: This currently assumes that the three points are monotonically increasing
    1197 or decreasing: so, it works for the cumulative histogram vectors, but not for
    1198 arbitrary vectors.  We should probably test that condition.
    11991186*****************************************************************************/
    12001187psF32 fitQuadraticSearchForYThenReturnX(
    … …  
    12211208    psVector *x = psVectorAlloc(3, PS_TYPE_F64);
    12221209    psVector *y = psVectorAlloc(3, PS_TYPE_F64);
    1223     psVector *yErr = psVectorAlloc(3, PS_TYPE_F64);
    1224 
    12251210    psF32 tmpFloat = 0.0f;
    12261211
    1227     if ((binNum > 0) && (binNum < (yVec->n - 2))) {
     1212    if ((binNum >= 1) && (binNum < (yVec->n - 2)) && (binNum < (xVec->n - 2))) {
    12281213        // The general case.  We have all three points.
    12291214        x->data.F64[0] = (psF64) (0.5 * (xVec->data.F32[binNum - 1] + xVec->data.F32[binNum]));
    … …  
    12341219        y->data.F64[2] = yVec->data.F32[binNum + 1];
    12351220        psTrace(__func__, 6, "x vec (orig) is (%f %f %f %f)\n", xVec->data.F32[binNum - 1],
    1236                 xVec->data.F32[binNum],
    1237                 xVec->data.F32[binNum+1],
    1238                 xVec->data.F32[binNum+2]);
    1239         psTrace(__func__, 6, "x vec is (%f %f %f)\n", x->data.F64[0], x->data.F64[1], x->data.F64[2]);
    1240         psTrace(__func__, 6, "y vec is (%f %f %f)\n", y->data.F64[0], y->data.F64[1], y->data.F64[2]);
     1221                xVec->data.F32[binNum], xVec->data.F32[binNum+1], xVec->data.F32[binNum+2]);
     1222        psTrace(__func__, 6, "x data is (%f %f %f)\n", x->data.F64[0], x->data.F64[1], x->data.F64[2]);
     1223        psTrace(__func__, 6, "y data is (%f %f %f)\n", y->data.F64[0], y->data.F64[1], y->data.F64[2]);
    12411224
    12421225        //
    … …  
    12461229        // so that the following checks are not necessary.
    12471230        //
    1248         if (y->data.F64[0] < y->data.F64[1]) {
    1249             if (!(y->data.F64[1] <= y->data.F64[2])) {
    1250                 psError(PS_ERR_UNKNOWN, true, "This routine must be called with montically increasing or decreasing data points.\n");
    1251                 psFree(x);
    1252                 psFree(y);
    1253                 psFree(yErr);
    1254                 psTrace(__func__, 5, "---- %s(NAN) end ----\n", __func__);
    1255                 return(NAN);
    1256             }
    1257             // Ensure that yVal is within the range of the bins we are using.
    1258             if (!((y->data.F64[0] <= yVal) && (yVal <= y->data.F64[2]))) {
    1259                 psError(PS_ERR_BAD_PARAMETER_VALUE, true,
    1260                         PS_ERRORTEXT_psStats_YVAL_OUT_OF_RANGE,
    1261                         (psF64)yVal, y->data.F64[2], y->data.F64[0]);
    1262             }
    1263         } else if (y->data.F64[0] > y->data.F64[1]) {
    1264             if (!(y->data.F64[1] >= y->data.F64[2])) {
    1265                 psError(PS_ERR_UNKNOWN, true, "This routine must be called with montically increasing or decreasing data points.\n");
    1266                 psFree(x);
    1267                 psFree(y);
    1268                 psFree(yErr);
    1269                 psTrace(__func__, 5, "---- %s(NAN) end ----\n", __func__);
    1270                 return(NAN);
    1271             }
    1272             // Ensure that yVal is within the range of the bins we are using.
    1273             if (!((y->data.F64[2] <= yVal) && (yVal <= y->data.F64[0]))) {
    1274                 psError(PS_ERR_BAD_PARAMETER_VALUE, true,
    1275                         PS_ERRORTEXT_psStats_YVAL_OUT_OF_RANGE,
    1276                         (psF64)yVal, y->data.F64[2], y->data.F64[0]);
    1277             }
    1278         }
    1279 
    1280         yErr->data.F64[0] = 1.0;
    1281         yErr->data.F64[1] = 1.0;
    1282         yErr->data.F64[2] = 1.0;
     1231        if (! (((y->data.F64[0] <= yVal) && (yVal <= y->data.F64[2])) ||
     1232                ((y->data.F64[2] <= yVal) && (yVal <= y->data.F64[0]))) ) {
     1233            psError(PS_ERR_BAD_PARAMETER_VALUE, true,
     1234                    PS_ERRORTEXT_psStats_YVAL_OUT_OF_RANGE,
     1235                    (psF64)yVal, y->data.F64[0], y->data.F64[2]);
     1236        }
     1237
     1238        if (((y->data.F64[0] < y->data.F64[1]) && !(y->data.F64[1] <= y->data.F64[2])) ||
     1239                ((y->data.F64[0] > y->data.F64[1]) && !(y->data.F64[1] >= y->data.F64[2]))) {
     1240            psError(PS_ERR_UNKNOWN, true, "This routine must be called with monotonically increasing or decreasing data points.\n");
     1241            psFree(x);
     1242            psFree(y);
     1243            psTrace(__func__, 5, "---- %s() end ----\n", __func__);
     1244            return(NAN);
     1245        }
    12831246
    12841247        // Determine the coefficients of the polynomial.
    12851248        psPolynomial1D *myPoly = psPolynomial1DAlloc(PS_POLYNOMIAL_ORD, 2);
    1286         myPoly = psVectorFitPolynomial1D(myPoly, NULL, 0, y, yErr, x);
     1249        myPoly = psVectorFitPolynomial1D(myPoly, NULL, 0, y, NULL, x);
    12871250        if (myPoly == NULL) {
    1288             psError(PS_ERR_UNEXPECTED_NULL,
    1289                     false,
     1251            psError(PS_ERR_UNEXPECTED_NULL, false,
    12901252                    PS_ERRORTEXT_psStats_STATS_FIT_QUADRATIC_POLYNOMIAL_1D_FIT);
    12911253            psFree(x);
    12921254            psFree(y);
    1293             psFree(yErr);
    12941255            psTrace(__func__, 5, "---- %s(NAN) end ----\n", __func__);
    12951256            return(NAN);
    … …  
    13091270        if (isnan(tmpFloat)) {
    13101271            psError(PS_ERR_UNEXPECTED_NULL,
    1311                     false,
    1312                     PS_ERRORTEXT_psStats_STATS_FIT_QUADRATIC_POLY_MEDIAN);
     1272                    false, PS_ERRORTEXT_psStats_STATS_FIT_QUADRATIC_POLY_MEDIAN);
    13131273            psFree(x);
    13141274            psFree(y);
    1315             psFree(yErr);
    13161275            psTrace(__func__, 5, "---- %s(NAN) end ----\n", __func__);
    13171276            return(NAN);
    13181277        }
    1319 
    13201278    } else {
    1321         // The special case where we have two points only at the beginning of
    1322         // the vectors x and y.
     1279        // These are special cases where the bin is at the beginning or end of the vector.
    13231280        if (binNum == 0) {
     1281            // We have two points only at the beginning of the vectors x and y.
    13241282            tmpFloat = 0.5 * (xVec->data.F32[binNum] +
    13251283                              xVec->data.F32[binNum + 1]);
    … …  
    13391297    psFree(x);
    13401298    psFree(y);
    1341     psFree(yErr);
    13421299
    13431300    psTrace(__func__, 5, "---- %s(%f) end ----\n", __func__, tmpFloat);
    … …  
    15421499        // ADD: Step 3.
    15431500        // Interpolate to the exact 50% position: this is the robust histogram median.
    1544         // XXX: Check return codes.
    15451501        //
    15461502        stats->robustMedian = fitQuadraticSearchForYThenReturnX(
    … …  
    15491505                                  binMedian,
    15501506                                  totalDataPoints/2.0);
     1507        if (isnan(stats->robustMedian)) {
     1508            psError(PS_ERR_UNKNOWN, false, "Failed to fit a quadratic and calculate the 50-percent position.\n");
     1509            psFree(tmpStatsMinMax);
     1510            psFree(robustHistogram);
     1511            psFree(cumulativeRobustHistogram);
     1512            psFree(tmpMaskVec);
     1513            psTrace(__func__, 4, "---- %s(1) end  ----\n", __func__);
     1514            return(1);
     1515        }
    15511516        psTrace(__func__, 6, "Current robust median is %f\n", stats->robustMedian);
    15521517
    … …  
    23242289                    binNum = (psS32)((inF32->data.F32[i] - out->bounds->data.F32[0]) / binSize);
    23252290                    if (errorsF32 != NULL) {
    2326                         // XXX: Check return codes.
    2327                         UpdateHistogramBins(binNum, out,
    2328                                             inF32->data.F32[i],
    2329                                             errorsF32->data.F32[i]);
     2291                        psS32 rc = UpdateHistogramBins(binNum, out, inF32->data.F32[i],
     2292                                                       errorsF32->data.F32[i]);
     2293                        if (rc < 0) {
     2294                            psLogMsg(__func__, PS_LOG_WARN, "WARNING: Failed to update the histogram bins with the errors vector.\n");
     2295                        }
    23302296                    } else {
    23312297                        // XXX: This if-statement really shouldn't be necessary.
    … …  
    23482314                    } else {
    23492315                        if (errorsF32 != NULL) {
    2350                             // XXX: Check return codes.
    2351                             UpdateHistogramBins(binNum, out,
    2352                                                 inF32->data.F32[i],
    2353                                                 errors->data.F32[i]);
     2316                            psS32 rc = UpdateHistogramBins(binNum, out,
     2317                                                           inF32->data.F32[i],
     2318                                                           errors->data.F32[i]);
     2319                            if (rc < 0) {
     2320                                psLogMsg(__func__, PS_LOG_WARN, "WARNING: Failed to update the histogram bins with the errors vector.\n");
     2321                            }
    23542322                        } else {
    23552323                            (out->nums->data.F32[binNum])+= 1.0;
Note: See TracChangeset for help on using the changeset viewer.