IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 6484 for trunk/psLib/src/math


Ignore:
Timestamp:
Feb 24, 2006, 1:43:16 PM (20 years ago)
Author:
drobbin
Message:

Changed psArray & psVector to zero 'n' upon allocation. Used to set vector->n = nalloc.

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

Legend:

Unmodified
Added
Removed

  • trunk/psLib/src/math/psConstants.h

    r6437 r6484  
    66 *  @author GLG, MHPCC
    77 *
    8  *  @version $Revision: 1.85 $ $Name: not supported by cvs2svn $
    9  *  @date $Date: 2006-02-17 00:56:48 $
     8 *  @version $Revision: 1.86 $ $Name: not supported by cvs2svn $
     9 *  @date $Date: 2006-02-24 23:43:15 $
    1010 *
    1111 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    312312if (VEC->n != SIZE) { \
    313313    psError(PS_ERR_BAD_PARAMETER_SIZE, true, \
    314             "psVector %s has size %d, should be %d." \
     314            "psVector %s has size %d, should be %d.", \
    315315            #VEC, VEC->n, SIZE); \
    316316    return(RVAL); \

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

    r6346 r6484  
    1010 *  @author EAM, IfA
    1111 *
    12  *  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
    13  *  @date $Date: 2006-02-07 23:14:21 $
     12 *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
     13 *  @date $Date: 2006-02-24 23:43:15 $
    1414 *
    1515 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    239239    psF64 ymodel;
    240240    psVector *deriv = psVectorAlloc(params->n, PS_TYPE_F32);
     241    deriv->n = deriv->nalloc;
    241242
    242243    // zero alpha and beta for summing below
    … …  
    370371    psVector *Beta   = psVectorAlloc(params->n, PS_TYPE_F64);
    371372    psVector *Params = psVectorAlloc(params->n, PS_TYPE_F32);
     373    beta->n = beta->nalloc;
     374    Beta->n = Beta->nalloc;
     375    Params->n = Params->nalloc;
    372376    psVector *dy     = NULL;
    373377    psF64 Chisq = 0.0;
    … …  
    385389        param_min = psVectorAlloc(params->n, PS_TYPE_F32);
    386390        param_max = psVectorAlloc(params->n, PS_TYPE_F32);
     391        beta_lim->n = beta_lim->nalloc;
     392        param_min->n = param_min->nalloc;
     393        param_max->n = param_max->nalloc;
    387394        for (int i = 0; i < params->n; i++) {
    388395            beta_lim->data.F32[i] = covar->data.F64[0][i];
    … …  
    402409    } else {
    403410        dy = psVectorAlloc(y->n, PS_TYPE_F32);
     411        dy->n = dy->nalloc;
    404412        psVectorInit(dy, 1.0);
    405413    }

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

    r6305 r6484  
    1010 *  @author EAM, IfA
    1111 *
    12  *  @version $Revision: 1.7 $ $Name: not supported by cvs2svn $
    13  *  @date $Date: 2006-02-02 21:09:07 $
     12 *  @version $Revision: 1.8 $ $Name: not supported by cvs2svn $
     13 *  @date $Date: 2006-02-24 23:43:15 $
    1414 *
    1515 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    4444    for (psS32 i = 0 ; i < SIZE ; i++) { \
    4545        VEC->data.F64[i] = ((2.0 / ((psF64) (SIZE - 1))) * ((psF64) i)) - 1.0; \
     46        VEC->n++; \
    4647    }\
    4748} else if (TYPE == PS_TYPE_F32){ \
    4849    for (psS32 i = 0 ; i < SIZE ; i++) { \
    4950        VEC->data.F32[i] = ((2.0 / ((psF32) (SIZE - 1))) * ((psF32) i)) - 1.0; \
     51        VEC->n++; \
    5052    }\
    5153}\
    … …  
    8789    if (sums == NULL) {
    8890        sums = psVectorAlloc(nSum, PS_TYPE_F64);
     91        sums->n = sums->nalloc;
    8992    } else if (nSum > sums->n) {
    9093        sums = psVectorRealloc(sums, nSum);
     94        sums->n = sums->nalloc;
    9195    }
    9296
    … …  
    305309    linear equations which can be easily solved.  The resulting vector is the
    306310    coefficients of the Chebyshev polys.
    307     
     311 
    308312    This method is significantly slower than the standard NR algorithm.  It
    309313    was explicitly requested that we not use the NR algorithm.
    … …  
    371375        }
    372376        B->data.F64[i] = ordPoly->coeff[i];
     377        B->n++;
    373378    }
    374379
    … …  
    534539    // Compute the B vector
    535540    psVector *B = psVectorAlloc(NUM_POLY, PS_TYPE_F64);
     541    B->n = B->nalloc;
    536542    for (psS32 i = 0 ; i < NUM_POLY ; i++) {
    537543        B->data.F64[i] = 0.0;
    … …  
    780786    A = psImageAlloc(nTerm, nTerm, PS_TYPE_F64);
    781787    B = psVectorAlloc(nTerm, PS_TYPE_F64);
     788    B->n = B->nalloc;
    782789    // Initialize data structures.
    783790    if (!psImageInit(A, 0.0) || !psVectorInit(B, 0.0)) {
    … …  
    10451052    psVector *fit   = NULL;
    10461053    psVector *resid = psVectorAlloc(f->n, PS_TYPE_F64);
     1054    resid->n = resid->nalloc;
    10471055
    10481056    // eventual expansion: user supplies one of various stats option pairs,
    … …  
    12051213    A = psImageAlloc(nTerm, nTerm, PS_TYPE_F64);
    12061214    B = psVectorAlloc(nTerm, PS_TYPE_F64);
     1215    B->n = B->nalloc;
    12071216    // Initialize data structures.
    12081217    if (!psImageInit(A, 0.0) || !psVectorInit(B, 0.0)) {
    … …  
    14621471    }
    14631472    psVector *resid = psVectorAlloc(f->n, PS_TYPE_F64);
     1473    resid->n = resid->nalloc;
    14641474
    14651475    // eventual expansion: user supplies one of various stats option pairs,
    … …  
    16301640    A = psImageAlloc(nTerm, nTerm, PS_TYPE_F64);
    16311641    B = psVectorAlloc(nTerm, PS_TYPE_F64);
     1642    B->n = B->nalloc;
    16321643    // Initialize data structures.
    16331644    if (!psImageInit(A, 0.0) || !psVectorInit(B, 0.0)) {
    … …  
    20022013    psVector *fit   = NULL;
    20032014    psVector *resid = psVectorAlloc(f->n, PS_TYPE_F64);
     2015    resid->n = resid->nalloc;
    20042016
    20052017    // eventual expansion: user supplies one of various stats option pairs,
    … …  
    21742186    A = psImageAlloc(nTerm, nTerm, PS_TYPE_F64);
    21752187    B = psVectorAlloc(nTerm, PS_TYPE_F64);
     2188    B->n = B->nalloc;
    21762189    // Initialize data structures.
    21772190    if (!psImageInit(A, 0.0) || !psVectorInit(B, 0.0)) {
    … …  
    25962609    psVector *fit   = NULL;
    25972610    psVector *resid = psVectorAlloc(f->n, PS_TYPE_F64);
     2611    resid->n = resid->nalloc;
    25982612
    25992613    // eventual expansion: user supplies one of various stats option pairs,

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

    r6449 r6484  
    1111 *  NOTE: XXX: The SDR is silent about data types.  F32 is implemented here.
    1212 *
    13  *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
    14  *  @date $Date: 2006-02-18 00:56:44 $
     13 *  @version $Revision: 1.7 $ $Name: not supported by cvs2svn $
     14 *  @date $Date: 2006-02-24 23:43:15 $
    1515 *
    1616 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    525525        p_psMemSetPersistent(myParamMask, true);
    526526        p_psMemSetPersistent(myParamMask->data.U8, true);
    527         for (i=0;i<myParamMask->n;i++) {
     527        for (i=0;i<myParamMask->nalloc;i++) {
    528528            myParamMask->data.U8[i] = 0;
     529            myParamMask->n++;
    529530        }
    530531    } else {
    … …  
    543544                ((psVector *) (v->data[i]))->data.F32[j] = 0.0;
    544545            }
     546            ((psVector *)(v->data[i]))->n++;
    545547        }
    546548        v->n++;
    … …  
    550552    for (i=0;i<numDims;i++) {
    551553        Q->data.F32[i] = params->data.F32[i];
     554        Q->n++;
    552555    }
    553556
    … …  
    595598            if (myParamMask->data.U8[i] == 0) {
    596599                u->data.F32[i] = Q->data.F32[i] - params->data.F32[i];
     600                u->n++;
    597601
    598602                psTrace(__func__, 6, "u[i]=Q[i]-P[i] (%f = %f - %f)\n", u->data.F32[i],
    … …  
    602606            } else {
    603607                u->data.F32[i] = 0.0;
     608                u->n++;
    604609            }
    605610        }

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

    r6437 r6484  
    77*  polynomials.  It also contains a Gaussian functions.
    88*
    9 *  @version $Revision: 1.143 $ $Name: not supported by cvs2svn $
    10 *  @date $Date: 2006-02-17 00:56:48 $
     9*  @version $Revision: 1.144 $ $Name: not supported by cvs2svn $
     10*  @date $Date: 2006-02-24 23:43:15 $
    1111*
    1212*  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    589589    for (unsigned int i = 0; i < Npts; i++) {
    590590        gauss->data.F32[i] = mean + p_psRandomGaussian(r, sigma);
     591        gauss->n++;
    591592    }
    592593    psFree(r);
    … …  
    802803        for (unsigned int i=0;i<x->n;i++) {
    803804            tmp->data.F64[i] = psPolynomial1DEval(poly, x->data.F64[i]);
     805            tmp->n++;
    804806        }
    805807        break;
    … …  
    808810        for (unsigned int i=0;i<x->n;i++) {
    809811            tmp->data.F32[i] = psPolynomial1DEval(poly, x->data.F32[i]);
     812            tmp->n++;
    810813        }
    811814        break;
    … …  
    869872        for (unsigned int i=0; i<vecLen; i++) {
    870873            tmp->data.F32[i] = psPolynomial2DEval(poly,x->data.F32[i],y->data.F32[i]);
     874            tmp->n++;
    871875        }
    872876        break;
    … …  
    883887        for (unsigned int i=0; i<vecLen; i++) {
    884888            tmp->data.F64[i] = psPolynomial2DEval(poly,x->data.F64[i],y->data.F64[i]);
     889            tmp->n++;
    885890        }
    886891        break;
    … …  
    944949    // Evaluate polynomial
    945950    // XXX: Consult with IfA: is this how they want to handle multiple data types?
    946     if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F32)) {
    947         for (unsigned int i = 0; i < vecLen; i++) {
    948             tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F32[i], y->data.F32[i], z->data.F32[i]);
    949         }
    950     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F64)) {
    951         for (unsigned int i = 0; i < vecLen; i++) {
    952             tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F32[i], y->data.F32[i], z->data.F64[i]);
    953         }
    954     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F32)) {
    955         for (unsigned int i = 0; i < vecLen; i++) {
    956             tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F32[i], y->data.F64[i], z->data.F32[i]);
    957         }
    958     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F64)) {
    959         for (unsigned int i = 0; i < vecLen; i++) {
    960             tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F32[i], y->data.F64[i], z->data.F64[i]);
    961         }
    962     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F32)) {
    963         for (unsigned int i = 0; i < vecLen; i++) {
    964             tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F64[i], y->data.F32[i], z->data.F32[i]);
    965         }
    966     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F64)) {
    967         for (unsigned int i = 0; i < vecLen; i++) {
    968             tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F64[i], y->data.F32[i], z->data.F64[i]);
    969         }
    970     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F32)) {
    971         for (unsigned int i = 0; i < vecLen; i++) {
    972             tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F64[i], y->data.F64[i], z->data.F32[i]);
    973         }
    974     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F64)) {
    975         for (unsigned int i = 0; i < vecLen; i++) {
    976             tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F64[i], y->data.F64[i], z->data.F64[i]);
     951    if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32)
     952            && (z->type.type == PS_TYPE_F32)) {
     953        for (unsigned int i = 0; i < vecLen; i++) {
     954            tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F32[i],
     955                                                  y->data.F32[i], z->data.F32[i]);
     956            tmp->n++;
     957        }
     958    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32)
     959               && (z->type.type == PS_TYPE_F64)) {
     960        for (unsigned int i = 0; i < vecLen; i++) {
     961            tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F32[i],
     962                                                  y->data.F32[i], z->data.F64[i]);
     963            tmp->n++;
     964        }
     965    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64)
     966               && (z->type.type == PS_TYPE_F32)) {
     967        for (unsigned int i = 0; i < vecLen; i++) {
     968            tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F32[i],
     969                                                  y->data.F64[i], z->data.F32[i]);
     970            tmp->n++;
     971        }
     972    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64)
     973               && (z->type.type == PS_TYPE_F64)) {
     974        for (unsigned int i = 0; i < vecLen; i++) {
     975            tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F32[i],
     976                                                  y->data.F64[i], z->data.F64[i]);
     977            tmp->n++;
     978        }
     979    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32)
     980               && (z->type.type == PS_TYPE_F32)) {
     981        for (unsigned int i = 0; i < vecLen; i++) {
     982            tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F64[i],
     983                                                  y->data.F32[i], z->data.F32[i]);
     984            tmp->n++;
     985        }
     986    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32)
     987               && (z->type.type == PS_TYPE_F64)) {
     988        for (unsigned int i = 0; i < vecLen; i++) {
     989            tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F64[i],
     990                                                  y->data.F32[i], z->data.F64[i]);
     991            tmp->n++;
     992        }
     993    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64)
     994               && (z->type.type == PS_TYPE_F32)) {
     995        for (unsigned int i = 0; i < vecLen; i++) {
     996            tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F64[i],
     997                                                  y->data.F64[i], z->data.F32[i]);
     998            tmp->n++;
     999        }
     1000    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64)
     1001               && (z->type.type == PS_TYPE_F64)) {
     1002        for (unsigned int i = 0; i < vecLen; i++) {
     1003            tmp->data.F64[i] = psPolynomial3DEval(poly, x->data.F64[i],
     1004                                                  y->data.F64[i], z->data.F64[i]);
     1005            tmp->n++;
    9771006        }
    9781007    }
    … …  
    10401069    // Evaluate polynomial
    10411070    // XXX: Consult with IfA: is this how they want to handle multiple data types?
    1042     if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F32)) {
    1043         for (unsigned int i = 0; i < vecLen; i++) {
    1044             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i], y->data.F32[i], z->data.F32[i], t->data.F32[i]);
    1045         }
    1046     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F64)) {
    1047         for (unsigned int i = 0; i < vecLen; i++) {
    1048             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i], y->data.F32[i], z->data.F32[i], t->data.F64[i]);
    1049         }
    1050     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F32)) {
    1051         for (unsigned int i = 0; i < vecLen; i++) {
    1052             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i], y->data.F32[i], z->data.F64[i], t->data.F32[i]);
    1053         }
    1054     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F64)) {
    1055         for (unsigned int i = 0; i < vecLen; i++) {
    1056             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i], y->data.F32[i], z->data.F64[i], t->data.F64[i]);
    1057         }
    1058     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F32)) {
    1059         for (unsigned int i = 0; i < vecLen; i++) {
    1060             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i], y->data.F64[i], z->data.F32[i], t->data.F32[i]);
    1061         }
    1062     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F64)) {
    1063         for (unsigned int i = 0; i < vecLen; i++) {
    1064             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i], y->data.F64[i], z->data.F32[i], t->data.F64[i]);
    1065         }
    1066     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F32)) {
    1067         for (unsigned int i = 0; i < vecLen; i++) {
    1068             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i], y->data.F64[i], z->data.F64[i], t->data.F32[i]);
    1069         }
    1070     } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F64)) {
    1071         for (unsigned int i = 0; i < vecLen; i++) {
    1072             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i], y->data.F64[i], z->data.F64[i], t->data.F64[i]);
    1073         }
    1074     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F32)) {
    1075         for (unsigned int i = 0; i < vecLen; i++) {
    1076             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i], y->data.F32[i], z->data.F32[i], t->data.F32[i]);
    1077         }
    1078     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F64)) {
    1079         for (unsigned int i = 0; i < vecLen; i++) {
    1080             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i], y->data.F32[i], z->data.F32[i], t->data.F64[i]);
    1081         }
    1082     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F32)) {
    1083         for (unsigned int i = 0; i < vecLen; i++) {
    1084             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i], y->data.F32[i], z->data.F64[i], t->data.F32[i]);
    1085         }
    1086     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32) && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F64)) {
    1087         for (unsigned int i = 0; i < vecLen; i++) {
    1088             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i], y->data.F32[i], z->data.F64[i], t->data.F64[i]);
    1089         }
    1090     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F32)) {
    1091         for (unsigned int i = 0; i < vecLen; i++) {
    1092             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i], y->data.F64[i], z->data.F32[i], t->data.F32[i]);
    1093         }
    1094     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F64)) {
    1095         for (unsigned int i = 0; i < vecLen; i++) {
    1096             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i], y->data.F64[i], z->data.F32[i], t->data.F64[i]);
    1097         }
    1098     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F32)) {
    1099         for (unsigned int i = 0; i < vecLen; i++) {
    1100             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i], y->data.F64[i], z->data.F64[i], t->data.F32[i]);
    1101         }
    1102     } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64) && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F64)) {
    1103         for (unsigned int i = 0; i < vecLen; i++) {
    1104             tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i], y->data.F64[i], z->data.F64[i],
     1071    if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32)
     1072            && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F32)) {
     1073        for (unsigned int i = 0; i < vecLen; i++) {
     1074            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i],
     1075                                                  y->data.F32[i], z->data.F32[i], t->data.F32[i]);
     1076            tmp->n++;
     1077        }
     1078    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32)
     1079               && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F64)) {
     1080        for (unsigned int i = 0; i < vecLen; i++) {
     1081            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i],
     1082                                                  y->data.F32[i], z->data.F32[i], t->data.F64[i]);
     1083            tmp->n++;
     1084        }
     1085    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32)
     1086               && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F32)) {
     1087        for (unsigned int i = 0; i < vecLen; i++) {
     1088            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i],
     1089                                                  y->data.F32[i], z->data.F64[i], t->data.F32[i]);
     1090            tmp->n++;
     1091        }
     1092    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F32)
     1093               && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F64)) {
     1094        for (unsigned int i = 0; i < vecLen; i++) {
     1095            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i],
     1096                                                  y->data.F32[i], z->data.F64[i], t->data.F64[i]);
     1097            tmp->n++;
     1098        }
     1099    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64)
     1100               && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F32)) {
     1101        for (unsigned int i = 0; i < vecLen; i++) {
     1102            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i],
     1103                                                  y->data.F64[i], z->data.F32[i], t->data.F32[i]);
     1104            tmp->n++;
     1105        }
     1106    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64)
     1107               && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F64)) {
     1108        for (unsigned int i = 0; i < vecLen; i++) {
     1109            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i],
     1110                                                  y->data.F64[i], z->data.F32[i], t->data.F64[i]);
     1111            tmp->n++;
     1112        }
     1113    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64)
     1114               && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F32)) {
     1115        for (unsigned int i = 0; i < vecLen; i++) {
     1116            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i],
     1117                                                  y->data.F64[i], z->data.F64[i], t->data.F32[i]);
     1118            tmp->n++;
     1119        }
     1120    } else if ((x->type.type == PS_TYPE_F32) && (y->type.type == PS_TYPE_F64)
     1121               && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F64)) {
     1122        for (unsigned int i = 0; i < vecLen; i++) {
     1123            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F32[i],
     1124                                                  y->data.F64[i], z->data.F64[i], t->data.F64[i]);
     1125            tmp->n++;
     1126        }
     1127    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32)
     1128               && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F32)) {
     1129        for (unsigned int i = 0; i < vecLen; i++) {
     1130            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i],
     1131                                                  y->data.F32[i], z->data.F32[i], t->data.F32[i]);
     1132            tmp->n++;
     1133        }
     1134    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32)
     1135               && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F64)) {
     1136        for (unsigned int i = 0; i < vecLen; i++) {
     1137            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i],
     1138                                                  y->data.F32[i], z->data.F32[i], t->data.F64[i]);
     1139            tmp->n++;
     1140        }
     1141    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32)
     1142               && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F32)) {
     1143        for (unsigned int i = 0; i < vecLen; i++) {
     1144            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i],
     1145                                                  y->data.F32[i], z->data.F64[i], t->data.F32[i]);
     1146            tmp->n++;
     1147        }
     1148    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F32)
     1149               && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F64)) {
     1150        for (unsigned int i = 0; i < vecLen; i++) {
     1151            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i],
     1152                                                  y->data.F32[i], z->data.F64[i], t->data.F64[i]);
     1153            tmp->n++;
     1154        }
     1155    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64)
     1156               && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F32)) {
     1157        for (unsigned int i = 0; i < vecLen; i++) {
     1158            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i],
     1159                                                  y->data.F64[i], z->data.F32[i], t->data.F32[i]);
     1160            tmp->n++;
     1161        }
     1162    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64)
     1163               && (z->type.type == PS_TYPE_F32) && (t->type.type == PS_TYPE_F64)) {
     1164        for (unsigned int i = 0; i < vecLen; i++) {
     1165            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i],
     1166                                                  y->data.F64[i], z->data.F32[i], t->data.F64[i]);
     1167            tmp->n++;
     1168        }
     1169    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64)
     1170               && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F32)) {
     1171        for (unsigned int i = 0; i < vecLen; i++) {
     1172            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i],
     1173                                                  y->data.F64[i], z->data.F64[i], t->data.F32[i]);
     1174            tmp->n++;
     1175        }
     1176    } else if ((x->type.type == PS_TYPE_F64) && (y->type.type == PS_TYPE_F64)
     1177               && (z->type.type == PS_TYPE_F64) && (t->type.type == PS_TYPE_F64)) {
     1178        for (unsigned int i = 0; i < vecLen; i++) {
     1179            tmp->data.F64[i] = psPolynomial4DEval(poly, x->data.F64[i],
     1180                                                  y->data.F64[i], z->data.F64[i],
    11051181                                                  t->data.F64[i]);
     1182            tmp->n++;
    11061183        }
    11071184    }

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

    r6437 r6484  
    11/** @file psSpline.c
    22*
    3 *  @brief Contains basic spline allocation, deallocation, fitting, 
     3*  @brief Contains basic spline allocation, deallocation, fitting,
    44*         and evaluation routines.
    55*
    66*  This file contains the routines that allocate, free, and evaluate splines.
    77*
    8 *  @version $Revision: 1.136 $ $Name: not supported by cvs2svn $
    9 *  @date $Date: 2006-02-17 00:56:48 $
     8*  @version $Revision: 1.137 $ $Name: not supported by cvs2svn $
     9*  @date $Date: 2006-02-24 23:43:15 $
    1010*
    1111*  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    239239            for (psS32 i = 0 ; i < x->n ; i++) {
    240240                spline->knots->data.F32[i] = x->data.F32[i];
     241                spline->knots->n++;
    241242            }
    242243        } else if (x->type.type == PS_TYPE_F64) {
    243244            for (psS32 i = 0 ; i < x->n ; i++) {
    244245                spline->knots->data.F32[i] = (psF32) x->data.F64[i];
     246                spline->knots->n++;
    245247            }
    246248        }
    … …  
    248250        for (psS32 i = 0 ; i < y->n ; i++) {
    249251            spline->knots->data.F32[i] = (psF32) i;
     252            spline->knots->n++;
    250253        }
    251254    }
    … …  
    412415        for (psS32 i=0;i<x->n;i++) {
    413416            tmpVector->data.F32[i] = psSpline1DEval(spline, x->data.F32[i]);
     417            tmpVector->n++;
    414418        }
    415419    } else if (x->type.type == PS_TYPE_F64) {
    416420        for (psS32 i=0;i<x->n;i++) {
    417421            tmpVector->data.F32[i] = psSpline1DEval(spline, (psF32) x->data.F64[i]);
     422            tmpVector->n++;
    418423        }
    419424    }

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

    r6449 r6484  
    1616 * use ->min and ->max (PS_STAT_USE_RANGE)
    1717 *
    18  *  @version $Revision: 1.169 $ $Name: not supported by cvs2svn $
    19  *  @date $Date: 2006-02-18 00:56:44 $
     18 *  @version $Revision: 1.170 $ $Name: not supported by cvs2svn $
     19 *  @date $Date: 2006-02-24 23:43:15 $
    2020 *
    2121 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    569569    // Allocate temporary vectors for the data.
    570570    unsortedVector = psVectorAlloc(nValues, PS_TYPE_F32);
     571    unsortedVector->n = nValues;
    571572
    572573    // Determine if we must only use data points within a min/max range.
    … …  
    657658    psS32 numBins = histogram->nums->n;
    658659    psVector *smooth = psVectorAlloc(numBins, PS_TYPE_F32);
     660    smooth->n = smooth->nalloc;
    659661    psS32 jMin = 0;
    660662    psS32 jMax = 0;
    … …  
    783785    // Allocate temporary vectors for the data.
    784786    unsortedVector = psVectorAlloc(nValues, PS_TYPE_F32);
     787    unsortedVector->n = unsortedVector->nalloc;
    785788
    786789    if (stats->options & PS_STAT_USE_RANGE) {
    … …  
    10371040    // However, we do no want to modify the original mask vector.
    10381041    tmpMask = psVectorAlloc(myVector->n, PS_TYPE_U8);
     1042    tmpMask->n = tmpMask->nalloc;
    10391043
    10401044    // If we were called with a mask vector, then initialize the temporary
    … …  
    11981202    psVector *x = psVectorAlloc(3, PS_TYPE_F64);
    11991203    psVector *y = psVectorAlloc(3, PS_TYPE_F64);
     1204    x->n = 3;
     1205    y->n = 3;
    12001206    psF32 tmpFloat = 0.0f;
    12011207
    … …  
    13131319    if (deriv == NULL) {
    13141320        deriv = psVectorAlloc(2, PS_TYPE_F32);
     1321        deriv->n = 2;
    13151322    } else {
    13161323        PS_ASSERT_VECTOR_SIZE(deriv, 2, NAN);
    … …  
    13651372    psS32 rcBool = false;
    13661373    psVector *tmpMaskVec = psVectorAlloc(myVector->n, PS_TYPE_U8);
     1374    tmpMaskVec->n = tmpMaskVec->nalloc;
    13671375    if (maskVector != NULL) {
    13681376        for (psS32 i = 0 ; i < myVector->n ; i++) {
    … …  
    18501858        psArray *x = psArrayAlloc((1 + (binMax - binMin)));
    18511859        psS32 j = 0;
     1860        y->n = y->nalloc;
    18521861
    18531862        for (psS32 i = binMin ; i <= binMax ; i++) {
    … …  
    18561865            x->n++;
    18571866            ((psVector *) x->data[j])->data.F32[0] = PS_BIN_MIDPOINT(newHistogram, i);
     1867            ((psVector *) x->data[j])->n++;
    18581868            j++;
    18591869        }
    … …  
    18861896        psMinimization *min = psMinimizationAlloc(100, 0.01);
    18871897        psVector *params = psVectorAlloc(2, PS_TYPE_F32);
     1898        params->n = params->nalloc;
    18881899        // Initial guess for the mean ([0]) and standard dev.
    18891900        params->data.F32[0] = stats->robustMedian;
    … …  
    20492060    // Allocate the bins, and initialize them to zero.
    20502061    newHist->nums = psVectorAlloc(n, PS_TYPE_F32);
    2051     for (i = 0; i < newHist->nums->n; i++) {
     2062    for (i = 0; i < newHist->nums->nalloc; i++) {
    20522063        newHist->nums->data.F32[i] = 0.0;
     2064        newHist->nums->n++;
    20532065    }
    20542066
    … …  
    20942106    // then there are N-1 bins.
    20952107    newHist->nums = psVectorAlloc((bounds->n) - 1, PS_TYPE_F32);
    2096     for (i = 0; i < newHist->nums->n; i++) {
     2108    for (i = 0; i < newHist->nums->nalloc; i++) {
    20972109        newHist->nums->data.F32[i] = 0.0;
     2110        newHist->nums->n++;
    20982111    }
    20992112
    … …  
    23662379        for (i = 0; i < in->n; i++) {
    23672380            tmp->data.F32[i] = (psF32)in->data.S8[i];
     2381            tmp->n++;
    23682382        }
    23692383    } else if (in->type.type == PS_TYPE_S16) {
    … …  
    23712385        for (i = 0; i < in->n; i++) {
    23722386            tmp->data.F32[i] = (psF32) in->data.S16[i];
     2387            tmp->n++;
    23732388        }
    23742389    } else if (in->type.type == PS_TYPE_S32) {
    … …  
    23762391        for (i = 0; i < in->n; i++) {
    23772392            tmp->data.F32[i] = (psF32)in->data.S32[i];
     2393            tmp->n++;
    23782394        }
    23792395    } else if (in->type.type == PS_TYPE_S64) {
    … …  
    23812397        for (i = 0; i < in->n; i++) {
    23822398            tmp->data.F32[i] = (psF32)in->data.S64[i];
     2399            tmp->n++;
    23832400        }
    23842401    } else if (in->type.type == PS_TYPE_U8) {
    … …  
    23862403        for (i = 0; i < in->n; i++) {
    23872404            tmp->data.F32[i] = (psF32)in->data.U8[i];
     2405            tmp->n++;
    23882406        }
    23892407    } else if (in->type.type == PS_TYPE_U16) {
    … …  
    23912409        for (i = 0; i < in->n; i++) {
    23922410            tmp->data.F32[i] = (psF32)in->data.U16[i];
     2411            tmp->n++;
    23932412        }
    23942413    } else if (in->type.type == PS_TYPE_U32) {
    … …  
    23962415        for (i = 0; i < in->n; i++) {
    23972416            tmp->data.F32[i] = (psF32)in->data.U32[i];
     2417            tmp->n++;
    23982418        }
    23992419    } else if (in->type.type == PS_TYPE_U64) {
    … …  
    24012421        for (i = 0; i < in->n; i++) {
    24022422            tmp->data.F32[i] = (psF32)in->data.U64[i];
     2423            tmp->n++;
    24032424        }
    24042425    } else if (in->type.type == PS_TYPE_F64) {
    … …  
    24062427        for (i = 0; i < in->n; i++) {
    24072428            tmp->data.F32[i] = (psF32)in->data.F64[i];
     2429            tmp->n++;
    24082430        }
    24092431    } else if (in->type.type == PS_TYPE_F32) {
    … …  
    24632485    inF32 = p_psConvertToF32((psVector *) in);
    24642486    if (inF32 == NULL) {
    2465         inF32 = (psVector *) in;
    2466         mustFreeVectorIn = 0;
    2467     }
     2487        //        printf("\n\ninF32 is NULL\n\n");
     2488        //        inF32 = (psVector *) in;
     2489        //        mustFreeVectorIn = 0;
     2490        inF32 = psVectorCopy(inF32, in, PS_TYPE_F32);
     2491        inF32->n = inF32->nalloc;
     2492        mustFreeVectorIn = 1;
     2493    }
     2494    //    else printf("\ninF32 has n=%ld, nalloc=%ld\n", inF32->n, inF32->nalloc);
    24682495    errorsF32 = p_psConvertToF32((psVector *) errors);
    24692496    if (errorsF32 == NULL) {
    … …  
    24712498        mustFreeVectorErrors = 0;
    24722499    }
    2473 
     2500    //    inF32->n = inF32->nalloc;
     2501    //    errorsF32->n = errorsF32->nalloc;
    24742502    if ((stats->options & PS_STAT_USE_RANGE) && (stats->min >= stats->max)) {
    24752503        PS_ASSERT_FLOAT_LARGER_THAN_OR_EQUAL(stats->max, stats->min, stats);
Note: See TracChangeset for help on using the changeset viewer.