IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Jun 28, 2005, 2:43:46 PM (21 years ago)
Author:
drobbin
Message:

made requested revisions to psPolynomials per apidelta-report-cycle6

File:
1 edited

Legend:

Unmodified
Added
Removed

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

    r4405 r4422  
    77 *  polynomials.  It also contains a Gaussian functions.
    88 *
    9  *  @version $Revision: 1.112 $ $Name: not supported by cvs2svn $
    10  *  @date $Date: 2005-06-28 00:53:53 $
     9 *  @version $Revision: 1.113 $ $Name: not supported by cvs2svn $
     10 *  @date $Date: 2005-06-29 00:43:46 $
    1111 *
    1212 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    4646/* TYPE DEFINITIONS                                                          */
    4747/*****************************************************************************/
    48 static void polynomial1DFree(psPolynomial1D* myPoly);
    49 static void polynomial2DFree(psPolynomial2D* myPoly);
    50 static void polynomial3DFree(psPolynomial3D* myPoly);
    51 static void polynomial4DFree(psPolynomial4D* myPoly);
    52 static void dPolynomial1DFree(psDPolynomial1D* myPoly);
    53 static void dPolynomial2DFree(psDPolynomial2D* myPoly);
    54 static void dPolynomial3DFree(psDPolynomial3D* myPoly);
    55 static void dPolynomial4DFree(psDPolynomial4D* myPoly);
     48static void polynomial1DFree(psPolynomial1D* poly);
     49static void polynomial2DFree(psPolynomial2D* poly);
     50static void polynomial3DFree(psPolynomial3D* poly);
     51static void polynomial4DFree(psPolynomial4D* poly);
     52static void dPolynomial1DFree(psDPolynomial1D* poly);
     53static void dPolynomial2DFree(psDPolynomial2D* poly);
     54static void dPolynomial3DFree(psDPolynomial3D* poly);
     55static void dPolynomial4DFree(psDPolynomial4D* poly);
    5656static void spline1DFree(psSpline1D *tmpSpline);
    5757static psS32 vectorBinDisectF32(psF32 *bins,psS32 numBins,psF32 x);
    … …  
    9797}
    9898
    99 static void polynomial1DFree(psPolynomial1D* myPoly)
    100 {
    101     psFree(myPoly->coeff);
    102     psFree(myPoly->coeffErr);
    103     psFree(myPoly->mask);
    104 }
    105 
    106 static void polynomial2DFree(psPolynomial2D* myPoly)
    107 {
    108     psS32 x = 0;
    109 
    110     for (x = 0; x < myPoly->nX; x++) {
    111         psFree(myPoly->coeff[x]);
    112         psFree(myPoly->coeffErr[x]);
    113         psFree(myPoly->mask[x]);
    114     }
    115     psFree(myPoly->coeff);
    116     psFree(myPoly->coeffErr);
    117     psFree(myPoly->mask);
    118 }
    119 
    120 static void polynomial3DFree(psPolynomial3D* myPoly)
    121 {
    122     psS32 x = 0;
    123     psS32 y = 0;
    124 
    125     for (x = 0; x < myPoly->nX; x++) {
    126         for (y = 0; y < myPoly->nY; y++) {
    127             psFree(myPoly->coeff[x][y]);
    128             psFree(myPoly->coeffErr[x][y]);
    129             psFree(myPoly->mask[x][y]);
    130         }
    131         psFree(myPoly->coeff[x]);
    132         psFree(myPoly->coeffErr[x]);
    133         psFree(myPoly->mask[x]);
    134     }
    135 
    136     psFree(myPoly->coeff);
    137     psFree(myPoly->coeffErr);
    138     psFree(myPoly->mask);
    139 }
    140 
    141 static void polynomial4DFree(psPolynomial4D* myPoly)
    142 {
    143     psS32 w = 0;
    144     psS32 x = 0;
    145     psS32 y = 0;
    146 
    147     for (w = 0; w < myPoly->nW; w++) {
    148         for (x = 0; x < myPoly->nX; x++) {
    149             for (y = 0; y < myPoly->nY; y++) {
    150                 psFree(myPoly->coeff[w][x][y]);
    151                 psFree(myPoly->coeffErr[w][x][y]);
    152                 psFree(myPoly->mask[w][x][y]);
     99static void polynomial1DFree(psPolynomial1D* poly)
     100{
     101    psFree(poly->coeff);
     102    psFree(poly->coeffErr);
     103    psFree(poly->mask);
     104}
     105
     106static void polynomial2DFree(psPolynomial2D* poly)
     107{
     108    unsigned int x = 0;
     109
     110    for (x = 0; x < poly->nX; x++) {
     111        psFree(poly->coeff[x]);
     112        psFree(poly->coeffErr[x]);
     113        psFree(poly->mask[x]);
     114    }
     115    psFree(poly->coeff);
     116    psFree(poly->coeffErr);
     117    psFree(poly->mask);
     118}
     119
     120static void polynomial3DFree(psPolynomial3D* poly)
     121{
     122    unsigned int x = 0;
     123    unsigned int y = 0;
     124
     125    for (x = 0; x < poly->nX; x++) {
     126        for (y = 0; y < poly->nY; y++) {
     127            psFree(poly->coeff[x][y]);
     128            psFree(poly->coeffErr[x][y]);
     129            psFree(poly->mask[x][y]);
     130        }
     131        psFree(poly->coeff[x]);
     132        psFree(poly->coeffErr[x]);
     133        psFree(poly->mask[x]);
     134    }
     135
     136    psFree(poly->coeff);
     137    psFree(poly->coeffErr);
     138    psFree(poly->mask);
     139}
     140
     141static void polynomial4DFree(psPolynomial4D* poly)
     142{
     143    unsigned int x = 0;
     144    unsigned int y = 0;
     145    unsigned int z = 0;
     146
     147    for (x = 0; x < poly->nX; x++) {
     148        for (y = 0; y < poly->nY; y++) {
     149            for (z = 0; z < poly->nZ; z++) {
     150                psFree(poly->coeff[x][y][z]);
     151                psFree(poly->coeffErr[x][y][z]);
     152                psFree(poly->mask[x][y][z]);
    153153            }
    154             psFree(myPoly->coeff[w][x]);
    155             psFree(myPoly->coeffErr[w][x]);
    156             psFree(myPoly->mask[w][x]);
    157         }
    158         psFree(myPoly->coeff[w]);
    159         psFree(myPoly->coeffErr[w]);
    160         psFree(myPoly->mask[w]);
    161     }
    162 
    163     psFree(myPoly->coeff);
    164     psFree(myPoly->coeffErr);
    165     psFree(myPoly->mask);
    166 }
    167 
    168 static void dPolynomial1DFree(psDPolynomial1D* myPoly)
    169 {
    170     psFree(myPoly->coeff);
    171     psFree(myPoly->coeffErr);
    172     psFree(myPoly->mask);
    173 }
    174 
    175 static void dPolynomial2DFree(psDPolynomial2D* myPoly)
    176 {
    177     for (psS32 x = 0; x < myPoly->nX; x++) {
    178         psFree(myPoly->coeff[x]);
    179         psFree(myPoly->coeffErr[x]);
    180         psFree(myPoly->mask[x]);
    181     }
    182     psFree(myPoly->coeff);
    183     psFree(myPoly->coeffErr);
    184     psFree(myPoly->mask);
    185 }
    186 
    187 static void dPolynomial3DFree(psDPolynomial3D* myPoly)
    188 {
    189     psS32 x = 0;
    190     psS32 y = 0;
    191 
    192     for (x = 0; x < myPoly->nX; x++) {
    193         for (y = 0; y < myPoly->nY; y++) {
    194             psFree(myPoly->coeff[x][y]);
    195             psFree(myPoly->coeffErr[x][y]);
    196             psFree(myPoly->mask[x][y]);
    197         }
    198         psFree(myPoly->coeff[x]);
    199         psFree(myPoly->coeffErr[x]);
    200         psFree(myPoly->mask[x]);
    201     }
    202 
    203     psFree(myPoly->coeff);
    204     psFree(myPoly->coeffErr);
    205     psFree(myPoly->mask);
    206 }
    207 
    208 static void dPolynomial4DFree(psDPolynomial4D* myPoly)
    209 {
    210     psS32 w = 0;
    211     psS32 x = 0;
    212     psS32 y = 0;
    213 
    214     for (w = 0; w < myPoly->nW; w++) {
    215         for (x = 0; x < myPoly->nX; x++) {
    216             for (y = 0; y < myPoly->nY; y++) {
    217                 psFree(myPoly->coeff[w][x][y]);
    218                 psFree(myPoly->coeffErr[w][x][y]);
    219                 psFree(myPoly->mask[w][x][y]);
     154            psFree(poly->coeff[x][y]);
     155            psFree(poly->coeffErr[x][y]);
     156            psFree(poly->mask[x][y]);
     157        }
     158        psFree(poly->coeff[x]);
     159        psFree(poly->coeffErr[x]);
     160        psFree(poly->mask[x]);
     161    }
     162
     163    psFree(poly->coeff);
     164    psFree(poly->coeffErr);
     165    psFree(poly->mask);
     166}
     167
     168static void dPolynomial1DFree(psDPolynomial1D* poly)
     169{
     170    psFree(poly->coeff);
     171    psFree(poly->coeffErr);
     172    psFree(poly->mask);
     173}
     174
     175static void dPolynomial2DFree(psDPolynomial2D* poly)
     176{
     177    for (unsigned int x = 0; x < poly->nX; x++) {
     178        psFree(poly->coeff[x]);
     179        psFree(poly->coeffErr[x]);
     180        psFree(poly->mask[x]);
     181    }
     182    psFree(poly->coeff);
     183    psFree(poly->coeffErr);
     184    psFree(poly->mask);
     185}
     186
     187static void dPolynomial3DFree(psDPolynomial3D* poly)
     188{
     189    unsigned int x = 0;
     190    unsigned int y = 0;
     191
     192    for (x = 0; x < poly->nX; x++) {
     193        for (y = 0; y < poly->nY; y++) {
     194            psFree(poly->coeff[x][y]);
     195            psFree(poly->coeffErr[x][y]);
     196            psFree(poly->mask[x][y]);
     197        }
     198        psFree(poly->coeff[x]);
     199        psFree(poly->coeffErr[x]);
     200        psFree(poly->mask[x]);
     201    }
     202
     203    psFree(poly->coeff);
     204    psFree(poly->coeffErr);
     205    psFree(poly->mask);
     206}
     207
     208static void dPolynomial4DFree(psDPolynomial4D* poly)
     209{
     210    unsigned int x = 0;
     211    unsigned int y = 0;
     212    unsigned int z = 0;
     213
     214    for (x = 0; x < poly->nX; x++) {
     215        for (y = 0; y < poly->nY; y++) {
     216            for (z = 0; z < poly->nZ; z++) {
     217                psFree(poly->coeff[x][y][z]);
     218                psFree(poly->coeffErr[x][y][z]);
     219                psFree(poly->mask[x][y][z]);
    220220            }
    221             psFree(myPoly->coeff[w][x]);
    222             psFree(myPoly->coeffErr[w][x]);
    223             psFree(myPoly->mask[w][x]);
    224         }
    225         psFree(myPoly->coeff[w]);
    226         psFree(myPoly->coeffErr[w]);
    227         psFree(myPoly->mask[w]);
    228     }
    229 
    230     psFree(myPoly->coeff);
    231     psFree(myPoly->coeffErr);
    232     psFree(myPoly->mask);
     221            psFree(poly->coeff[x][y]);
     222            psFree(poly->coeffErr[x][y]);
     223            psFree(poly->mask[x][y]);
     224        }
     225        psFree(poly->coeff[x]);
     226        psFree(poly->coeffErr[x]);
     227        psFree(poly->mask[x]);
     228    }
     229
     230    psFree(poly->coeff);
     231    psFree(poly->coeffErr);
     232    psFree(poly->mask);
    233233}
    234234
    … …  
    280280    Polynomial coefficients will be accessed in [w][x][y][z] fashion.
    281281 *****************************************************************************/
    282 static psF32 ordPolynomial1DEval(psF32 x, const psPolynomial1D* myPoly)
     282static psF32 ordPolynomial1DEval(psF32 x, const psPolynomial1D* poly)
    283283{
    284284    psS32 loop_x = 0;
    … …  
    289289            "---- Calling ordPolynomial1DEval(%f)\n", x);
    290290    psTrace(".psLib.dataManip.psFunctions.ordPolynomial1DEval", 4,
    291             "Polynomial order is %d\n", myPoly->n);
    292     for (loop_x = 0; loop_x < myPoly->n; loop_x++) {
     291            "Polynomial order is %d\n", poly->n);
     292    for (loop_x = 0; loop_x < poly->n; loop_x++) {
    293293        psTrace(".psLib.dataManip.psFunctions.ordPolynomial1DEval", 4,
    294                 "Polynomial coeff[%d] is %f\n", loop_x, myPoly->coeff[loop_x]);
    295     }
    296 
    297     for (loop_x = 0; loop_x < myPoly->n; loop_x++) {
    298         if (myPoly->mask[loop_x] == 0) {
     294                "Polynomial coeff[%d] is %f\n", loop_x, poly->coeff[loop_x]);
     295    }
     296
     297    for (loop_x = 0; loop_x < poly->n; loop_x++) {
     298        if (poly->mask[loop_x] == 0) {
    299299            psTrace(".psLib.dataManip.psFunctions.ordPolynomial1DEval", 10,
    300                     "polysum+= sum*coeff [%f+= (%f * %f)\n", polySum, xSum, myPoly->coeff[loop_x]);
    301             polySum += xSum * myPoly->coeff[loop_x];
     300                    "polysum+= sum*coeff [%f+= (%f * %f)\n", polySum, xSum, poly->coeff[loop_x]);
     301            polySum += xSum * poly->coeff[loop_x];
    302302        }
    303303        xSum *= x;
    … …  
    310310// XXX: How does the mask vector effect Crenshaw's formula?
    311311// XXX: We assume that x is scaled between -1.0 and 1.0;
    312 static psF32 chebPolynomial1DEval(psF32 x, const psPolynomial1D* myPoly)
     312static psF32 chebPolynomial1DEval(psF32 x, const psPolynomial1D* poly)
    313313{
    314314    PS_ASSERT_FLOAT_WITHIN_RANGE(x, -1.0, 1.0, 0.0);
    315315    // XXX: Create a macro for this in psConstants.h
    316     if (myPoly->n < 1) {
    317         psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Error: Chebyshev polynomial as order %d.", myPoly->n);
     316    if (poly->n < 1) {
     317        psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Error: Chebyshev polynomial as order %d.", poly->n);
    318318        return(NAN);
    319319    }
    320320    psVector *d;
    321     psS32 n = myPoly->n;
     321    psS32 n = poly->n;
    322322    psS32 i;
    323323    psF32 tmp = 0.0;
    … …  
    325325    // Special case where the Chebyshev poly is constant.
    326326    if (n == 1) {
    327         if (myPoly->mask[0] == 0) {
    328             tmp += myPoly->coeff[0];
     327        if (poly->mask[0] == 0) {
     328            tmp += poly->coeff[0];
    329329        }
    330330        return(tmp);
    … …  
    333333    // Special case where the Chebyshev poly is linear.
    334334    if (n == 2) {
    335         if (myPoly->mask[0] == 0) {
    336             tmp+= myPoly->coeff[0];
    337         }
    338         if (myPoly->mask[1] == 0) {
    339             tmp+= myPoly->coeff[1] * x;
     335        if (poly->mask[0] == 0) {
     336            tmp+= poly->coeff[0];
     337        }
     338        if (poly->mask[1] == 0) {
     339            tmp+= poly->coeff[1] * x;
    340340        }
    341341        return(tmp);
    … …  
    344344    // General case where the Chebyshev poly has 2 or more terms.
    345345    d = psVectorAlloc(n, PS_TYPE_F32);
    346     if(myPoly->mask[n-1] == 0) {
    347         d->data.F32[n-1] = myPoly->coeff[n-1];
     346    if(poly->mask[n-1] == 0) {
     347        d->data.F32[n-1] = poly->coeff[n-1];
    348348    } else {
    349349        d->data.F32[n-1] = 0.0;
    … …  
    351351
    352352    d->data.F32[n-2] = (2.0 * x * d->data.F32[n-1]);
    353     if(myPoly->mask[n-2] == 0) {
    354         d->data.F32[n-2] += myPoly->coeff[n-2];
     353    if(poly->mask[n-2] == 0) {
     354        d->data.F32[n-2] += poly->coeff[n-2];
    355355    }
    356356
    … …  
    358358        d->data.F32[i] = (2.0 * x * d->data.F32[i+1]) -
    359359                         (d->data.F32[i+2]);
    360         if(myPoly->mask[i] == 0) {
    361             d->data.F32[i] += myPoly->coeff[i];
     360        if(poly->mask[i] == 0) {
     361            d->data.F32[i] += poly->coeff[i];
    362362        }
    363363    }
    … …  
    365365    tmp = (x * d->data.F32[1]) -
    366366          (d->data.F32[2]);
    367     if(myPoly->mask[0] == 0) {
    368         tmp += (0.5 * myPoly->coeff[0]);
     367    if(poly->mask[0] == 0) {
     368        tmp += (0.5 * poly->coeff[0]);
    369369    }
    370370    psFree(d);
    … …  
    378378    psPolynomial1D **chebPolys = NULL;
    379379
    380     n = myPoly->n;
     380    n = poly->n;
    381381    chebPolys = createChebyshevPolys(n);
    382382
    383383    tmp = 0.0;
    384     for (i=0;i<myPoly->n;i++) {
    385         tmp+= (myPoly->coeff[i] * psPolynomial1DEval(x, chebPolys[i]));
    386     }
    387     tmp-= (myPoly->coeff[0]/2.0);
     384    for (i=0;i<poly->n;i++) {
     385        tmp+= (poly->coeff[i] * psPolynomial1DEval(x, chebPolys[i]));
     386    }
     387    tmp-= (poly->coeff[0]/2.0);
    388388
    389389
    … …  
    394394static psF32 ordPolynomial2DEval(psF32 x,
    395395                                 psF32 y,
    396                                  const psPolynomial2D* myPoly)
    397 {
    398     PS_ASSERT_POLY_NON_NULL(myPoly, NAN);
     396                                 const psPolynomial2D* poly)
     397{
     398    PS_ASSERT_POLY_NON_NULL(poly, NAN);
    399399
    400400    psS32 loop_x = 0;
    … …  
    404404    psF32 ySum = 1.0;
    405405
    406     for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
     406    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
    407407        ySum = xSum;
    408         for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    409             if (myPoly->mask[loop_x][loop_y] == 0) {
    410                 polySum += ySum * myPoly->coeff[loop_x][loop_y];
     408        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     409            if (poly->mask[loop_x][loop_y] == 0) {
     410                polySum += ySum * poly->coeff[loop_x][loop_y];
    411411            }
    412412            ySum *= y;
    … …  
    418418}
    419419
    420 static psF32 chebPolynomial2DEval(psF32 x, psF32 y, const psPolynomial2D* myPoly)
     420static psF32 chebPolynomial2DEval(psF32 x, psF32 y, const psPolynomial2D* poly)
    421421{
    422422    PS_ASSERT_FLOAT_WITHIN_RANGE(x, -1.0, 1.0, 0.0);
    423423    PS_ASSERT_FLOAT_WITHIN_RANGE(y, -1.0, 1.0, 0.0);
    424     PS_ASSERT_POLY_NON_NULL(myPoly, NAN);
     424    PS_ASSERT_POLY_NON_NULL(poly, NAN);
    425425
    426426    psS32 loop_x = 0;
    … …  
    433433    // Determine how many Chebyshev polynomials
    434434    // are needed, then create them.
    435     maxChebyPoly = myPoly->nX;
    436     if (myPoly->nY > maxChebyPoly) {
    437         maxChebyPoly = myPoly->nY;
     435    maxChebyPoly = poly->nX;
     436    if (poly->nY > maxChebyPoly) {
     437        maxChebyPoly = poly->nY;
    438438    }
    439439    chebPolys = createChebyshevPolys(maxChebyPoly);
    440440
    441     for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
    442         for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    443             if (myPoly->mask[loop_x][loop_y] == 0) {
    444                 polySum += myPoly->coeff[loop_x][loop_y] *
     441    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
     442        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     443            if (poly->mask[loop_x][loop_y] == 0) {
     444                polySum += poly->coeff[loop_x][loop_y] *
    445445                           psPolynomial1DEval(chebPolys[loop_x], x) *
    446446                           psPolynomial1DEval(chebPolys[loop_y], y);
    … …  
    455455}
    456456
    457 static psF32 ordPolynomial3DEval(psF32 x, psF32 y, psF32 z, const psPolynomial3D* myPoly)
     457static psF32 ordPolynomial3DEval(psF32 x, psF32 y, psF32 z, const psPolynomial3D* poly)
    458458{
    459459    psS32 loop_x = 0;
    … …  
    465465    psF32 zSum = 1.0;
    466466
    467     for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
     467    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
    468468        ySum = xSum;
    469         for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
     469        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
    470470            zSum = ySum;
    471             for (loop_z = 0; loop_z < myPoly->nZ; loop_z++) {
    472                 if (myPoly->mask[loop_x][loop_y][loop_z] == 0) {
    473                     polySum += zSum * myPoly->coeff[loop_x][loop_y][loop_z];
     471            for (loop_z = 0; loop_z < poly->nZ; loop_z++) {
     472                if (poly->mask[loop_x][loop_y][loop_z] == 0) {
     473                    polySum += zSum * poly->coeff[loop_x][loop_y][loop_z];
    474474                }
    475475                zSum *= z;
    … …  
    483483}
    484484
    485 static psF32 chebPolynomial3DEval(psF32 x, psF32 y, psF32 z, const psPolynomial3D* myPoly)
     485static psF32 chebPolynomial3DEval(psF32 x, psF32 y, psF32 z, const psPolynomial3D* poly)
    486486{
    487487    PS_ASSERT_FLOAT_WITHIN_RANGE(x, -1.0, 1.0, 0.0);
    … …  
    498498    // Determine how many Chebyshev polynomials
    499499    // are needed, then create them.
    500     maxChebyPoly = myPoly->nX;
    501     if (myPoly->nY > maxChebyPoly) {
    502         maxChebyPoly = myPoly->nY;
    503     }
    504     if (myPoly->nZ > maxChebyPoly) {
    505         maxChebyPoly = myPoly->nZ;
     500    maxChebyPoly = poly->nX;
     501    if (poly->nY > maxChebyPoly) {
     502        maxChebyPoly = poly->nY;
     503    }
     504    if (poly->nZ > maxChebyPoly) {
     505        maxChebyPoly = poly->nZ;
    506506    }
    507507    chebPolys = createChebyshevPolys(maxChebyPoly);
    508508
    509     for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
    510         for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    511             for (loop_z = 0; loop_z < myPoly->nZ; loop_z++) {
    512                 if (myPoly->mask[loop_x][loop_y][loop_z] == 0) {
    513                     polySum += myPoly->coeff[loop_x][loop_y][loop_z] *
     509    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
     510        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     511            for (loop_z = 0; loop_z < poly->nZ; loop_z++) {
     512                if (poly->mask[loop_x][loop_y][loop_z] == 0) {
     513                    polySum += poly->coeff[loop_x][loop_y][loop_z] *
    514514                               psPolynomial1DEval(chebPolys[loop_x], x) *
    515515                               psPolynomial1DEval(chebPolys[loop_y], y) *
    … …  
    527527}
    528528
    529 static psF32 ordPolynomial4DEval(psF32 w, psF32 x, psF32 y, psF32 z, const psPolynomial4D* myPoly)
    530 {
    531     psS32 loop_w = 0;
     529static psF32 ordPolynomial4DEval(psF32 x, psF32 y, psF32 z, psF32 t, const psPolynomial4D* poly)
     530{
    532531    psS32 loop_x = 0;
    533532    psS32 loop_y = 0;
    534533    psS32 loop_z = 0;
     534    psS32 loop_t = 0;
    535535    psF32 polySum = 0.0;
    536     psF32 wSum = 1.0;
    537536    psF32 xSum = 1.0;
    538537    psF32 ySum = 1.0;
    539538    psF32 zSum = 1.0;
    540 
    541     for (loop_w = 0; loop_w < myPoly->nW; loop_w++) {
    542         xSum = wSum;
    543         for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
    544             ySum = xSum;
    545             for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    546                 zSum = ySum;
    547                 for (loop_z = 0; loop_z < myPoly->nZ; loop_z++) {
    548                     if (myPoly->mask[loop_w][loop_x][loop_y][loop_z] == 0) {
    549                         polySum += zSum * myPoly->coeff[loop_w][loop_x][loop_y][loop_z];
     539    psF32 tSum = 1.0;
     540
     541    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
     542        ySum = xSum;
     543        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     544            zSum = ySum;
     545            for (loop_z = 0; loop_z < poly->nZ; loop_z++) {
     546                tSum = zSum;
     547                for (loop_t = 0; loop_t < poly->nT; loop_t++) {
     548                    if (poly->mask[loop_x][loop_y][loop_z][loop_t] == 0) {
     549                        polySum += tSum * poly->coeff[loop_x][loop_y][loop_z][loop_t];
    550550                    }
    551                     zSum *= z;
     551                    tSum *= t;
    552552                }
    553                 ySum *= y;
     553                zSum *= z;
    554554            }
    555             xSum *= x;
    556         }
    557         wSum *= w;
     555            ySum *= y;
     556        }
     557        xSum *= x;
    558558    }
    559559
    … …  
    561561}
    562562
    563 static psF32 chebPolynomial4DEval(psF32 w, psF32 x, psF32 y, psF32 z, const psPolynomial4D* myPoly)
    564 {
    565     PS_ASSERT_FLOAT_WITHIN_RANGE(w, -1.0, 1.0, 0.0);
     563static psF32 chebPolynomial4DEval(psF32 x, psF32 y, psF32 z, psF32 t, const psPolynomial4D* poly)
     564{
    566565    PS_ASSERT_FLOAT_WITHIN_RANGE(x, -1.0, 1.0, 0.0);
    567566    PS_ASSERT_FLOAT_WITHIN_RANGE(y, -1.0, 1.0, 0.0);
    568567    PS_ASSERT_FLOAT_WITHIN_RANGE(z, -1.0, 1.0, 0.0);
    569     psS32 loop_w = 0;
     568    PS_ASSERT_FLOAT_WITHIN_RANGE(t, -1.0, 1.0, 0.0);
    570569    psS32 loop_x = 0;
    571570    psS32 loop_y = 0;
    572571    psS32 loop_z = 0;
     572    psS32 loop_t = 0;
    573573    psS32 i = 0;
    574574    psF32 polySum = 0.0;
    … …  
    578578    // Determine how many Chebyshev polynomials
    579579    // are needed, then create them.
    580     maxChebyPoly = myPoly->nW;
    581     if (myPoly->nX > maxChebyPoly) {
    582         maxChebyPoly = myPoly->nX;
    583     }
    584     if (myPoly->nY > maxChebyPoly) {
    585         maxChebyPoly = myPoly->nY;
    586     }
    587     if (myPoly->nZ > maxChebyPoly) {
    588         maxChebyPoly = myPoly->nZ;
     580    maxChebyPoly = poly->nX;
     581    if (poly->nY > maxChebyPoly) {
     582        maxChebyPoly = poly->nY;
     583    }
     584    if (poly->nZ > maxChebyPoly) {
     585        maxChebyPoly = poly->nZ;
     586    }
     587    if (poly->nT > maxChebyPoly) {
     588        maxChebyPoly = poly->nT;
    589589    }
    590590    chebPolys = createChebyshevPolys(maxChebyPoly);
    591591
    592     for (loop_w = 0; loop_w < myPoly->nW; loop_w++) {
    593         for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
    594             for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    595                 for (loop_z = 0; loop_z < myPoly->nZ; loop_z++) {
    596                     if (myPoly->mask[loop_w][loop_x][loop_y][loop_z] == 0) {
    597                         polySum += myPoly->coeff[loop_w][loop_x][loop_y][loop_z] *
    598                                    psPolynomial1DEval(chebPolys[loop_w], w) *
     592    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
     593        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     594            for (loop_z = 0; loop_z < poly->nZ; loop_z++) {
     595                for (loop_t = 0; loop_t < poly->nT; loop_t++) {
     596                    if (poly->mask[loop_x][loop_y][loop_z][loop_t] == 0) {
     597                        polySum += poly->coeff[loop_x][loop_y][loop_z][loop_t] *
    599598                                   psPolynomial1DEval(chebPolys[loop_x], x) *
    600599                                   psPolynomial1DEval(chebPolys[loop_y], y) *
    601                                    psPolynomial1DEval(chebPolys[loop_z], z);
     600                                   psPolynomial1DEval(chebPolys[loop_z], z) *
     601                                   psPolynomial1DEval(chebPolys[loop_t], t);
    602602                    }
    603603                }
    … …  
    616616    Polynomial coefficients will be accessed in [w][x][y][z] fashion.
    617617 *****************************************************************************/
    618 static psF64 dOrdPolynomial1DEval(psF64 x, const psDPolynomial1D* myPoly)
     618static psF64 dOrdPolynomial1DEval(psF64 x, const psDPolynomial1D* poly)
    619619{
    620620    psS32 loop_x = 0;
    … …  
    622622    psF64 xSum = 1.0;
    623623
    624     for (loop_x = 0; loop_x < myPoly->n; loop_x++) {
    625         if (myPoly->mask[loop_x] == 0) {
    626             polySum += xSum * myPoly->coeff[loop_x];
     624    for (loop_x = 0; loop_x < poly->n; loop_x++) {
     625        if (poly->mask[loop_x] == 0) {
     626            polySum += xSum * poly->coeff[loop_x];
    627627        }
    628628        xSum *= x;
    … …  
    634634// XXX: You can do this without having to psAlloc() vector d.
    635635// XXX: How does the mask vector effect Crenshaw's formula?
    636 static psF64 dChebPolynomial1DEval(psF64 x, const psDPolynomial1D* myPoly)
     636static psF64 dChebPolynomial1DEval(psF64 x, const psDPolynomial1D* poly)
    637637{
    638638    PS_ASSERT_FLOAT_WITHIN_RANGE(x, -1.0, 1.0, 0.0);
    … …  
    642642    psF64 tmp;
    643643
    644     n = myPoly->n;
     644    n = poly->n;
    645645    d = psVectorAlloc(n, PS_TYPE_F64);
    646     if(myPoly->mask[n-1] == 0) {
    647         d->data.F64[n-1] = myPoly->coeff[n-1];
     646    if(poly->mask[n-1] == 0) {
     647        d->data.F64[n-1] = poly->coeff[n-1];
    648648    } else {
    649649        d->data.F64[n-1] = 0.0;
    650650    }
    651651    d->data.F64[n-2] = (2.0 * x * d->data.F64[n-1]);
    652     if(myPoly->mask[n-2] == 0) {
    653         d->data.F64[n-2] += myPoly->coeff[n-2];
     652    if(poly->mask[n-2] == 0) {
     653        d->data.F64[n-2] += poly->coeff[n-2];
    654654    }
    655655    for (i=n-3;i>=1;i--) {
    656656        d->data.F64[i] = (2.0 * x * d->data.F64[i+1]) -
    657657                         (d->data.F64[i+2]);
    658         if(myPoly->mask[i] == 0) {
    659             d->data.F64[i] += myPoly->coeff[i];
     658        if(poly->mask[i] == 0) {
     659            d->data.F64[i] += poly->coeff[i];
    660660        }
    661661    }
    … …  
    663663    tmp = (x * d->data.F64[1]) -
    664664          (d->data.F64[2]);
    665     if(myPoly->mask[0] == 0) {
    666         tmp += (0.5 * myPoly->coeff[0]);
     665    if(poly->mask[0] == 0) {
     666        tmp += (0.5 * poly->coeff[0]);
    667667    }
    668668
    … …  
    673673static psF64 dOrdPolynomial2DEval(psF64 x,
    674674                                  psF64 y,
    675                                   const psDPolynomial2D* myPoly)
     675                                  const psDPolynomial2D* poly)
    676676{
    677677    psS32 loop_x = 0;
    … …  
    681681    psF64 ySum = 1.0;
    682682
    683     for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
     683    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
    684684        ySum = xSum;
    685         for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    686             if (myPoly->mask[loop_x][loop_y] == 0) {
    687                 polySum += ySum * myPoly->coeff[loop_x][loop_y];
     685        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     686            if (poly->mask[loop_x][loop_y] == 0) {
     687                polySum += ySum * poly->coeff[loop_x][loop_y];
    688688            }
    689689            ySum *= y;
    … …  
    695695}
    696696
    697 static psF64 dChebPolynomial2DEval(psF64 x, psF64 y, const psDPolynomial2D* myPoly)
     697static psF64 dChebPolynomial2DEval(psF64 x, psF64 y, const psDPolynomial2D* poly)
    698698{
    699699    PS_ASSERT_FLOAT_WITHIN_RANGE(x, -1.0, 1.0, 0.0);
    … …  
    708708    // Determine how many Chebyshev polynomials
    709709    // are needed, then create them.
    710     maxChebyPoly = myPoly->nX;
    711     if (myPoly->nY > maxChebyPoly) {
    712         maxChebyPoly = myPoly->nY;
     710    maxChebyPoly = poly->nX;
     711    if (poly->nY > maxChebyPoly) {
     712        maxChebyPoly = poly->nY;
    713713    }
    714714    chebPolys = createChebyshevPolys(maxChebyPoly);
    715715
    716     for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
    717         for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    718             if (myPoly->mask[loop_x][loop_y] == 0) {
    719                 polySum += myPoly->coeff[loop_x][loop_y] *
     716    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
     717        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     718            if (poly->mask[loop_x][loop_y] == 0) {
     719                polySum += poly->coeff[loop_x][loop_y] *
    720720                           psPolynomial1DEval(chebPolys[loop_x], x) *
    721721                           psPolynomial1DEval(chebPolys[loop_y], y);
    … …  
    731731}
    732732
    733 static psF64 dOrdPolynomial3DEval(psF64 x, psF64 y, psF64 z, const psDPolynomial3D* myPoly)
     733static psF64 dOrdPolynomial3DEval(psF64 x, psF64 y, psF64 z, const psDPolynomial3D* poly)
    734734{
    735735    psS32 loop_x = 0;
    … …  
    741741    psF64 zSum = 1.0;
    742742
    743     for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
     743    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
    744744        ySum = xSum;
    745         for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
     745        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
    746746            zSum = ySum;
    747             for (loop_z = 0; loop_z < myPoly->nZ; loop_z++) {
    748                 if (myPoly->mask[loop_x][loop_y][loop_z] == 0) {
    749                     polySum += zSum * myPoly->coeff[loop_x][loop_y][loop_z];
     747            for (loop_z = 0; loop_z < poly->nZ; loop_z++) {
     748                if (poly->mask[loop_x][loop_y][loop_z] == 0) {
     749                    polySum += zSum * poly->coeff[loop_x][loop_y][loop_z];
    750750                }
    751751                zSum *= z;
    … …  
    759759}
    760760
    761 static psF64 dChebPolynomial3DEval(psF64 x, psF64 y, psF64 z, const psDPolynomial3D* myPoly)
     761static psF64 dChebPolynomial3DEval(psF64 x, psF64 y, psF64 z, const psDPolynomial3D* poly)
    762762{
    763763    PS_ASSERT_FLOAT_WITHIN_RANGE(x, -1.0, 1.0, 0.0);
    … …  
    774774    // Determine how many Chebyshev polynomials
    775775    // are needed, then create them.
    776     maxChebyPoly = myPoly->nX;
    777     if (myPoly->nY > maxChebyPoly) {
    778         maxChebyPoly = myPoly->nY;
    779     }
    780     if (myPoly->nZ > maxChebyPoly) {
    781         maxChebyPoly = myPoly->nZ;
     776    maxChebyPoly = poly->nX;
     777    if (poly->nY > maxChebyPoly) {
     778        maxChebyPoly = poly->nY;
     779    }
     780    if (poly->nZ > maxChebyPoly) {
     781        maxChebyPoly = poly->nZ;
    782782    }
    783783    chebPolys = createChebyshevPolys(maxChebyPoly);
    784784
    785     for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
    786         for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    787             for (loop_z = 0; loop_z < myPoly->nZ; loop_z++) {
    788                 if (myPoly->mask[loop_x][loop_y][loop_z] == 0) {
    789                     polySum += myPoly->coeff[loop_x][loop_y][loop_z] *
     785    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
     786        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     787            for (loop_z = 0; loop_z < poly->nZ; loop_z++) {
     788                if (poly->mask[loop_x][loop_y][loop_z] == 0) {
     789                    polySum += poly->coeff[loop_x][loop_y][loop_z] *
    790790                               psPolynomial1DEval(chebPolys[loop_x], x) *
    791791                               psPolynomial1DEval(chebPolys[loop_y], y) *
    … …  
    803803}
    804804
    805 static psF64 dOrdPolynomial4DEval(psF64 w, psF64 x, psF64 y, psF64 z, const psDPolynomial4D* myPoly)
    806 {
    807     psS32 loop_w = 0;
     805static psF64 dOrdPolynomial4DEval(psF64 x, psF64 y, psF64 z, psF64 t, const psDPolynomial4D* poly)
     806{
    808807    psS32 loop_x = 0;
    809808    psS32 loop_y = 0;
    810809    psS32 loop_z = 0;
     810    psS32 loop_t = 0;
    811811    psF64 polySum = 0.0;
    812     psF64 wSum = 1.0;
    813812    psF64 xSum = 1.0;
    814813    psF64 ySum = 1.0;
    815814    psF64 zSum = 1.0;
    816 
    817     for (loop_w = 0; loop_w < myPoly->nW; loop_w++) {
    818         xSum = wSum;
    819         for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
    820             ySum = xSum;
    821             for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    822                 zSum = ySum;
    823                 for (loop_z = 0; loop_z < myPoly->nZ; loop_z++) {
    824                     if (myPoly->mask[loop_w][loop_x][loop_y][loop_z] == 0) {
    825                         polySum += zSum * myPoly->coeff[loop_w][loop_x][loop_y][loop_z];
     815    psF64 tSum = 1.0;
     816
     817    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
     818        ySum = xSum;
     819        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     820            zSum = ySum;
     821            for (loop_z = 0; loop_z < poly->nZ; loop_z++) {
     822                tSum = zSum;
     823                for (loop_t = 0; loop_t < poly->nT; loop_t++) {
     824                    if (poly->mask[loop_x][loop_y][loop_z][loop_t] == 0) {
     825                        polySum += tSum * poly->coeff[loop_x][loop_y][loop_z][loop_t];
    826826                    }
    827                     zSum *= z;
     827                    tSum *= t;
    828828                }
    829                 ySum *= y;
     829                zSum *= z;
    830830            }
    831             xSum *= x;
    832         }
    833         wSum *= w;
     831            ySum *= y;
     832        }
     833        xSum *= x;
    834834    }
    835835
    … …  
    837837}
    838838
    839 static psF64 dChebPolynomial4DEval(psF64 w, psF64 x, psF64 y, psF64 z, const psDPolynomial4D* myPoly)
    840 {
    841     PS_ASSERT_FLOAT_WITHIN_RANGE(w, -1.0, 1.0, 0.0);
     839static psF64 dChebPolynomial4DEval(psF64 x, psF64 y, psF64 z, psF64 t, const psDPolynomial4D* poly)
     840{
    842841    PS_ASSERT_FLOAT_WITHIN_RANGE(x, -1.0, 1.0, 0.0);
    843842    PS_ASSERT_FLOAT_WITHIN_RANGE(y, -1.0, 1.0, 0.0);
    844843    PS_ASSERT_FLOAT_WITHIN_RANGE(z, -1.0, 1.0, 0.0);
    845     psS32 loop_w = 0;
     844    PS_ASSERT_FLOAT_WITHIN_RANGE(t, -1.0, 1.0, 0.0);
    846845    psS32 loop_x = 0;
    847846    psS32 loop_y = 0;
    848847    psS32 loop_z = 0;
     848    psS32 loop_t = 0;
    849849    psS32 i = 0;
    850850    psF64 polySum = 0.0;
    … …  
    854854    // Determine how many Chebyshev polynomials
    855855    // are needed, then create them.
    856     maxChebyPoly = myPoly->nW;
    857     if (myPoly->nX > maxChebyPoly) {
    858         maxChebyPoly = myPoly->nX;
    859     }
    860     if (myPoly->nY > maxChebyPoly) {
    861         maxChebyPoly = myPoly->nY;
    862     }
    863     if (myPoly->nZ > maxChebyPoly) {
    864         maxChebyPoly = myPoly->nZ;
     856    maxChebyPoly = poly->nX;
     857    if (poly->nY > maxChebyPoly) {
     858        maxChebyPoly = poly->nY;
     859    }
     860    if (poly->nZ > maxChebyPoly) {
     861        maxChebyPoly = poly->nZ;
     862    }
     863    if (poly->nT > maxChebyPoly) {
     864        maxChebyPoly = poly->nT;
    865865    }
    866866    chebPolys = createChebyshevPolys(maxChebyPoly);
    867867
    868     for (loop_w = 0; loop_w < myPoly->nW; loop_w++) {
    869         for (loop_x = 0; loop_x < myPoly->nX; loop_x++) {
    870             for (loop_y = 0; loop_y < myPoly->nY; loop_y++) {
    871                 for (loop_z = 0; loop_z < myPoly->nZ; loop_z++) {
    872                     if (myPoly->mask[loop_w][loop_x][loop_y][loop_z] == 0) {
    873                         polySum += myPoly->coeff[loop_w][loop_x][loop_y][loop_z] *
    874                                    psPolynomial1DEval(chebPolys[loop_w], w) *
     868    for (loop_x = 0; loop_x < poly->nX; loop_x++) {
     869        for (loop_y = 0; loop_y < poly->nY; loop_y++) {
     870            for (loop_z = 0; loop_z < poly->nZ; loop_z++) {
     871                for (loop_t = 0; loop_t < poly->nT; loop_t++) {
     872                    if (poly->mask[loop_x][loop_y][loop_z][loop_t] == 0) {
     873                        polySum += poly->coeff[loop_x][loop_y][loop_z][loop_t] *
    875874                                   psPolynomial1DEval(chebPolys[loop_x], x) *
    876875                                   psPolynomial1DEval(chebPolys[loop_y], y) *
    877                                    psPolynomial1DEval(chebPolys[loop_z], z);
     876                                   psPolynomial1DEval(chebPolys[loop_z], z) *
     877                                   psPolynomial1DEval(chebPolys[loop_t], t);
    878878                    }
    879879                }
    … …  
    10691069    This routine must allocate memory for the polynomial structures.
    10701070 *****************************************************************************/
    1071 psPolynomial1D* psPolynomial1DAlloc(psS32 n,
     1071psPolynomial1D* psPolynomial1DAlloc(int n,
    10721072                                    psPolynomialType type)
    10731073{
    10741074    PS_ASSERT_INT_POSITIVE(n, NULL);
    10751075
    1076     psS32 i = 0;
     1076    int i = 0;
    10771077    psPolynomial1D* newPoly = NULL;
    10781078
    … …  
    10941094}
    10951095
    1096 psPolynomial2D* psPolynomial2DAlloc(psS32 nX, psS32 nY,
     1096psPolynomial2D* psPolynomial2DAlloc(int nX, int nY,
    10971097                                    psPolynomialType type)
    10981098{
    … …  
    11001100    PS_ASSERT_INT_POSITIVE(nY, NULL);
    11011101
    1102     psS32 x = 0;
    1103     psS32 y = 0;
     1102    int x = 0;
     1103    int y = 0;
    11041104    psPolynomial2D* newPoly = NULL;
    11051105
    … …  
    11301130}
    11311131
    1132 psPolynomial3D* psPolynomial3DAlloc(psS32 nX, psS32 nY, psS32 nZ,
     1132psPolynomial3D* psPolynomial3DAlloc(int nX, int nY, int nZ,
    11331133                                    psPolynomialType type)
    11341134{
    … …  
    11761176}
    11771177
    1178 psPolynomial4D* psPolynomial4DAlloc(psS32 nW, psS32 nX, psS32 nY, psS32 nZ,
     1178psPolynomial4D* psPolynomial4DAlloc(int nX, int nY, int nZ, int nT,
    11791179                                    psPolynomialType type)
    11801180{
    1181     PS_ASSERT_INT_POSITIVE(nW, NULL);
    11821181    PS_ASSERT_INT_POSITIVE(nX, NULL);
    11831182    PS_ASSERT_INT_POSITIVE(nY, NULL);
    11841183    PS_ASSERT_INT_POSITIVE(nZ, NULL);
    1185 
    1186     psS32 w = 0;
     1184    PS_ASSERT_INT_POSITIVE(nT, NULL);
     1185
    11871186    psS32 x = 0;
    11881187    psS32 y = 0;
    11891188    psS32 z = 0;
     1189    psS32 t = 0;
    11901190    psPolynomial4D* newPoly = NULL;
    11911191
    … …  
    11941194
    11951195    newPoly->type = type;
    1196     newPoly->nW = nW;
    11971196    newPoly->nX = nX;
    11981197    newPoly->nY = nY;
    11991198    newPoly->nZ = nZ;
    1200 
    1201     newPoly->coeff = (psF32 ****)psAlloc(nW * sizeof(psF32 ***));
    1202     newPoly->coeffErr = (psF32 ****)psAlloc(nW * sizeof(psF32 ***));
    1203     newPoly->mask = (psU8 ****)psAlloc(nW * sizeof(psU8 ***));
    1204     for (w = 0; w < nW; w++) {
    1205         newPoly->coeff[w] = (psF32 ***)psAlloc(nX * sizeof(psF32 **));
    1206         newPoly->coeffErr[w] = (psF32 ***)psAlloc(nX * sizeof(psF32 **));
    1207         newPoly->mask[w] = (psU8 ***)psAlloc(nX * sizeof(psU8 **));
    1208         for (x = 0; x < nX; x++) {
    1209             newPoly->coeff[w][x] = (psF32 **)psAlloc(nY * sizeof(psF32 *));
    1210             newPoly->coeffErr[w][x] = (psF32 **)psAlloc(nY * sizeof(psF32 *));
    1211             newPoly->mask[w][x] = (psU8 **)psAlloc(nY * sizeof(psU8 *));
    1212             for (y = 0; y < nY; y++) {
    1213                 newPoly->coeff[w][x][y] = (psF32 *)psAlloc(nZ * sizeof(psF32));
    1214                 newPoly->coeffErr[w][x][y] = (psF32 *)psAlloc(nZ * sizeof(psF32));
    1215                 newPoly->mask[w][x][y] = (psU8 *)psAlloc(nZ * sizeof(psU8));
     1199    newPoly->nT = nT;
     1200
     1201    newPoly->coeff = (psF32 ****)psAlloc(nX * sizeof(psF32 ***));
     1202    newPoly->coeffErr = (psF32 ****)psAlloc(nX * sizeof(psF32 ***));
     1203    newPoly->mask = (psU8 ****)psAlloc(nX * sizeof(psU8 ***));
     1204    for (x = 0; x < nX; x++) {
     1205        newPoly->coeff[x] = (psF32 ***)psAlloc(nY * sizeof(psF32 **));
     1206        newPoly->coeffErr[x] = (psF32 ***)psAlloc(nY * sizeof(psF32 **));
     1207        newPoly->mask[x] = (psU8 ***)psAlloc(nY * sizeof(psU8 **));
     1208        for (y = 0; y < nY; y++) {
     1209            newPoly->coeff[x][y] = (psF32 **)psAlloc(nZ * sizeof(psF32 *));
     1210            newPoly->coeffErr[x][y] = (psF32 **)psAlloc(nZ * sizeof(psF32 *));
     1211            newPoly->mask[x][y] = (psU8 **)psAlloc(nZ * sizeof(psU8 *));
     1212            for (z = 0; z < nZ; z++) {
     1213                newPoly->coeff[x][y][z] = (psF32 *)psAlloc(nT * sizeof(psF32));
     1214                newPoly->coeffErr[x][y][z] = (psF32 *)psAlloc(nT * sizeof(psF32));
     1215                newPoly->mask[x][y][z] = (psU8 *)psAlloc(nT * sizeof(psU8));
    12161216            }
    12171217        }
    12181218    }
    1219     for (w = 0; w < nW; w++) {
    1220         for (x = 0; x < nX; x++) {
    1221             for (y = 0; y < nY; y++) {
    1222                 for (z = 0; z < nZ; z++) {
    1223                     newPoly->coeff[w][x][y][z] = 0.0;
    1224                     newPoly->coeffErr[w][x][y][z] = 0.0;
    1225                     newPoly->mask[w][x][y][z] = 0;
     1219    for (x = 0; x < nX; x++) {
     1220        for (y = 0; y < nY; y++) {
     1221            for (z = 0; z < nZ; z++) {
     1222                for (t = 0; t < nT; t++) {
     1223                    newPoly->coeff[x][y][z][t] = 0.0;
     1224                    newPoly->coeffErr[x][y][z][t] = 0.0;
     1225                    newPoly->mask[x][y][z][t] = 0;
    12261226                }
    12271227            }
    … …  
    14311431
    14321432
    1433 psDPolynomial1D* psDPolynomial1DAlloc(psS32 n,
     1433psDPolynomial1D* psDPolynomial1DAlloc(int n,
    14341434                                      psPolynomialType type)
    14351435{
    14361436    PS_ASSERT_INT_POSITIVE(n, NULL);
    14371437
    1438     psS32 i = 0;
     1438    unsigned int i = 0;
    14391439    psDPolynomial1D* newPoly = NULL;
    14401440
    … …  
    14561456}
    14571457
    1458 psDPolynomial2D* psDPolynomial2DAlloc(psS32 nX, psS32 nY,
     1458psDPolynomial2D* psDPolynomial2DAlloc(int nX, int nY,
    14591459                                      psPolynomialType type)
    14601460{
    … …  
    14621462    PS_ASSERT_INT_POSITIVE(nY, NULL);
    14631463
    1464     psS32 x = 0;
    1465     psS32 y = 0;
     1464    unsigned int x = 0;
     1465    unsigned int y = 0;
    14661466    psDPolynomial2D* newPoly = NULL;
    14671467
    … …  
    14921492}
    14931493
    1494 psDPolynomial3D* psDPolynomial3DAlloc(psS32 nX, psS32 nY, psS32 nZ,
     1494psDPolynomial3D* psDPolynomial3DAlloc(int nX, int nY, int nZ,
    14951495                                      psPolynomialType type)
    14961496{
    … …  
    14991499    PS_ASSERT_INT_POSITIVE(nZ, NULL);
    15001500
    1501     psS32 x = 0;
    1502     psS32 y = 0;
    1503     psS32 z = 0;
     1501    unsigned int x = 0;
     1502    unsigned int y = 0;
     1503    unsigned int z = 0;
    15041504    psDPolynomial3D* newPoly = NULL;
    15051505
    … …  
    15381538}
    15391539
    1540 psDPolynomial4D* psDPolynomial4DAlloc(psS32 nW, psS32 nX, psS32 nY, psS32 nZ,
     1540psDPolynomial4D* psDPolynomial4DAlloc(int nX, int nY, int nZ, int nT,
    15411541                                      psPolynomialType type)
    15421542{
    1543     PS_ASSERT_INT_POSITIVE(nW, NULL);
    15441543    PS_ASSERT_INT_POSITIVE(nX, NULL);
    15451544    PS_ASSERT_INT_POSITIVE(nY, NULL);
    15461545    PS_ASSERT_INT_POSITIVE(nZ, NULL);
    1547 
    1548     psS32 w = 0;
    1549     psS32 x = 0;
    1550     psS32 y = 0;
    1551     psS32 z = 0;
     1546    PS_ASSERT_INT_POSITIVE(nT, NULL);
     1547
     1548    unsigned int x = 0;
     1549    unsigned int y = 0;
     1550    unsigned int z = 0;
     1551    unsigned int t = 0;
    15521552    psDPolynomial4D* newPoly = NULL;
    15531553
    … …  
    15561556
    15571557    newPoly->type = type;
    1558     newPoly->nW = nW;
    15591558    newPoly->nX = nX;
    15601559    newPoly->nY = nY;
    15611560    newPoly->nZ = nZ;
    1562 
    1563     newPoly->coeff = (psF64 ****)psAlloc(nW * sizeof(psF64 ***));
    1564     newPoly->coeffErr = (psF64 ****)psAlloc(nW * sizeof(psF64 ***));
    1565     newPoly->mask = (psU8 ****)psAlloc(nW * sizeof(psU8 ***));
    1566     for (w = 0; w < nW; w++) {
    1567         newPoly->coeff[w] = (psF64 ***)psAlloc(nX * sizeof(psF64 **));
    1568         newPoly->coeffErr[w] = (psF64 ***)psAlloc(nX * sizeof(psF64 **));
    1569         newPoly->mask[w] = (psU8 ***)psAlloc(nX * sizeof(psU8 **));
    1570         for (x = 0; x < nX; x++) {
    1571             newPoly->coeff[w][x] = (psF64 **)psAlloc(nY * sizeof(psF64 *));
    1572             newPoly->coeffErr[w][x] = (psF64 **)psAlloc(nY * sizeof(psF64 *));
    1573             newPoly->mask[w][x] = (psU8 **)psAlloc(nY * sizeof(psU8 *));
    1574             for (y = 0; y < nY; y++) {
    1575                 newPoly->coeff[w][x][y] = (psF64 *)psAlloc(nZ * sizeof(psF64));
    1576                 newPoly->coeffErr[w][x][y] = (psF64 *)psAlloc(nZ * sizeof(psF64));
    1577                 newPoly->mask[w][x][y] = (psU8 *)psAlloc(nZ * sizeof(psU8));
     1561    newPoly->nT = nT;
     1562
     1563    newPoly->coeff = (psF64 ****)psAlloc(nX * sizeof(psF64 ***));
     1564    newPoly->coeffErr = (psF64 ****)psAlloc(nX * sizeof(psF64 ***));
     1565    newPoly->mask = (psU8 ****)psAlloc(nX * sizeof(psU8 ***));
     1566    for (x = 0; x < nX; x++) {
     1567        newPoly->coeff[x] = (psF64 ***)psAlloc(nY * sizeof(psF64 **));
     1568        newPoly->coeffErr[x] = (psF64 ***)psAlloc(nY * sizeof(psF64 **));
     1569        newPoly->mask[x] = (psU8 ***)psAlloc(nY * sizeof(psU8 **));
     1570        for (y = 0; y < nY; y++) {
     1571            newPoly->coeff[x][y] = (psF64 **)psAlloc(nZ * sizeof(psF64 *));
     1572            newPoly->coeffErr[x][y] = (psF64 **)psAlloc(nZ * sizeof(psF64 *));
     1573            newPoly->mask[x][y] = (psU8 **)psAlloc(nZ * sizeof(psU8 *));
     1574            for (z = 0; z < nZ; z++) {
     1575                newPoly->coeff[x][y][z] = (psF64 *)psAlloc(nT * sizeof(psF64));
     1576                newPoly->coeffErr[x][y][z] = (psF64 *)psAlloc(nT * sizeof(psF64));
     1577                newPoly->mask[x][y][z] = (psU8 *)psAlloc(nT * sizeof(psU8));
    15781578            }
    15791579        }
    15801580    }
    1581     for (w = 0; w < nW; w++) {
    1582         for (x = 0; x < nX; x++) {
    1583             for (y = 0; y < nY; y++) {
    1584                 for (z = 0; z < nZ; z++) {
    1585                     newPoly->coeff[w][x][y][z] = 0.0;
    1586                     newPoly->coeffErr[w][x][y][z] = 0.0;
    1587                     newPoly->mask[w][x][y][z] = 0;
     1581    for (x = 0; x < nX; x++) {
     1582        for (y = 0; y < nY; y++) {
     1583            for (z = 0; z < nZ; z++) {
     1584                for (t = 0; t < nT; t++) {
     1585                    newPoly->coeff[x][y][z][t] = 0.0;
     1586                    newPoly->coeffErr[x][y][z][t] = 0.0;
     1587                    newPoly->mask[x][y][z][t] = 0;
    15881588                }
    15891589            }
    … …  
    15951595
    15961596
    1597 psF64 psDPolynomial1DEval(const psDPolynomial1D* myPoly, psF64 x)
    1598 {
    1599     PS_ASSERT_POLY_NON_NULL(myPoly, NAN);
    1600 
    1601     if (myPoly->type == PS_POLYNOMIAL_ORD) {
    1602         return(dOrdPolynomial1DEval(x, myPoly));
    1603     } else if (myPoly->type == PS_POLYNOMIAL_CHEB) {
    1604         return(dChebPolynomial1DEval(x, myPoly));
     1597psF64 psDPolynomial1DEval(const psDPolynomial1D* poly, psF64 x)
     1598{
     1599    PS_ASSERT_POLY_NON_NULL(poly, NAN);
     1600
     1601    if (poly->type == PS_POLYNOMIAL_ORD) {
     1602        return(dOrdPolynomial1DEval(x, poly));
     1603    } else if (poly->type == PS_POLYNOMIAL_CHEB) {
     1604        return(dChebPolynomial1DEval(x, poly));
    16051605    } else {
    16061606        psError(PS_ERR_BAD_PARAMETER_TYPE, true,
    16071607                PS_ERRORTEXT_psFunctions_INVALID_POLYNOMIAL_TYPE,
    1608                 myPoly->type);
     1608                poly->type);
    16091609    }
    16101610    return(NAN);
    16111611}
    16121612
    1613 psVector *psDPolynomial1DEvalVector(const psDPolynomial1D *myPoly,
     1613psVector *psDPolynomial1DEvalVector(const psDPolynomial1D *poly,
    16141614                                    const psVector *x)
    16151615
    16161616{
    1617     PS_ASSERT_POLY_NON_NULL(myPoly, NULL);
     1617    PS_ASSERT_POLY_NON_NULL(poly, NULL);
    16181618    PS_ASSERT_VECTOR_NON_NULL(x, NULL);
    16191619    PS_ASSERT_VECTOR_TYPE(x, PS_TYPE_F64, NULL);
    … …  
    16231623    tmp = psVectorAlloc(x->n, PS_TYPE_F64);
    16241624    for (psS32 i=0;i<x->n;i++) {
    1625         tmp->data.F64[i] = psDPolynomial1DEval(myPoly,
     1625        tmp->data.F64[i] = psDPolynomial1DEval(poly,
    16261626                                               x->data.F64[i]);
    16271627    }
    … …  
    16311631
    16321632
    1633 psF64 psDPolynomial2DEval(const psDPolynomial2D* myPoly,
     1633psF64 psDPolynomial2DEval(const psDPolynomial2D* poly,
    16341634                          psF64 x,
    16351635                          psF64 y)
    16361636{
    1637     PS_ASSERT_POLY_NON_NULL(myPoly, NAN);
    1638     if (myPoly->type == PS_POLYNOMIAL_ORD) {
    1639         return(dOrdPolynomial2DEval(x, y, myPoly));
    1640     } else if (myPoly->type == PS_POLYNOMIAL_CHEB) {
    1641         return(dChebPolynomial2DEval(x, y, myPoly));
     1637    PS_ASSERT_POLY_NON_NULL(poly, NAN);
     1638    if (poly->type == PS_POLYNOMIAL_ORD) {
     1639        return(dOrdPolynomial2DEval(x, y, poly));
     1640    } else if (poly->type == PS_POLYNOMIAL_CHEB) {
     1641        return(dChebPolynomial2DEval(x, y, poly));
    16421642    } else {
    16431643        psError(PS_ERR_BAD_PARAMETER_TYPE, true,
    16441644                PS_ERRORTEXT_psFunctions_INVALID_POLYNOMIAL_TYPE,
    1645                 myPoly->type);
     1645                poly->type);
    16461646    }
    16471647    return(NAN);
    16481648}
    16491649
    1650 psVector *psDPolynomial2DEvalVector(const psDPolynomial2D *myPoly,
     1650psVector *psDPolynomial2DEvalVector(const psDPolynomial2D *poly,
    16511651                                    const psVector *x,
    16521652                                    const psVector *y)
    16531653{
    1654     PS_ASSERT_POLY_NON_NULL(myPoly, NULL);
     1654    PS_ASSERT_POLY_NON_NULL(poly, NULL);
    16551655    PS_ASSERT_VECTOR_NON_NULL(x, NULL);
    16561656    PS_ASSERT_VECTOR_TYPE(x, PS_TYPE_F64, NULL);
    … …  
    16711671    // Evaluate the polynomial
    16721672    for (psS32 i = 0; i < vecLen; i++) {
    1673         tmp->data.F64[i] = psDPolynomial2DEval(myPoly,x->data.F64[i],y->data.F64[i]);
     1673        tmp->data.F64[i] = psDPolynomial2DEval(poly,x->data.F64[i],y->data.F64[i]);
    16741674    }
    16751675
    … …  
    16791679
    16801680
    1681 psF64 psDPolynomial3DEval(const psDPolynomial3D* myPoly,
     1681psF64 psDPolynomial3DEval(const psDPolynomial3D* poly,
    16821682                          psF64 x,
    16831683                          psF64 y,
    16841684                          psF64 z)
    16851685{
    1686     PS_ASSERT_POLY_NON_NULL(myPoly, NAN);
    1687 
    1688     if (myPoly->type == PS_POLYNOMIAL_ORD) {
    1689         return(dOrdPolynomial3DEval(x, y, z, myPoly));
    1690     } else if (myPoly->type == PS_POLYNOMIAL_CHEB) {
    1691         return(dChebPolynomial3DEval(x, y, z, myPoly));
     1686    PS_ASSERT_POLY_NON_NULL(poly, NAN);
     1687
     1688    if (poly->type == PS_POLYNOMIAL_ORD) {
     1689        return(dOrdPolynomial3DEval(x, y, z, poly));
     1690    } else if (poly->type == PS_POLYNOMIAL_CHEB) {
     1691        return(dChebPolynomial3DEval(x, y, z, poly));
    16921692    } else {
    16931693        psError(PS_ERR_BAD_PARAMETER_TYPE, true,
    16941694                PS_ERRORTEXT_psFunctions_INVALID_POLYNOMIAL_TYPE,
    1695                 myPoly->type);
     1695                poly->type);
    16961696    }
    16971697    return(NAN);
    16981698}
    16991699
    1700 psVector *psDPolynomial3DEvalVector(const psDPolynomial3D *myPoly,
     1700psVector *psDPolynomial3DEvalVector(const psDPolynomial3D *poly,
    17011701                                    const psVector *x,
    17021702                                    const psVector *y,
    … …  
    17041704
    17051705{
    1706     PS_ASSERT_POLY_NON_NULL(myPoly, NULL);
     1706    PS_ASSERT_POLY_NON_NULL(poly, NULL);
    17071707    PS_ASSERT_VECTOR_NON_NULL(x, NULL);
    17081708    PS_ASSERT_VECTOR_TYPE(x, PS_TYPE_F64, NULL);
    … …  
    17281728    // Evaluate polynomial
    17291729    for (psS32 i = 0; i < vecLen; i++) {
    1730         tmp->data.F64[i] = psDPolynomial3DEval(myPoly,
     1730        tmp->data.F64[i] = psDPolynomial3DEval(poly,
    17311731                                               x->data.F64[i],
    17321732                                               y->data.F64[i],
    … …  
    17381738}
    17391739
    1740 psF64 psDPolynomial4DEval(const psDPolynomial4D* myPoly,
    1741                           psF64 w,
     1740psF64 psDPolynomial4DEval(const psDPolynomial4D* poly,
    17421741                          psF64 x,
    17431742                          psF64 y,
    1744                           psF64 z)
    1745 {
    1746     PS_ASSERT_POLY_NON_NULL(myPoly, NAN);
    1747 
    1748     if (myPoly->type == PS_POLYNOMIAL_ORD) {
    1749         return(dOrdPolynomial4DEval(w,x,y,z, myPoly));
    1750     } else if (myPoly->type == PS_POLYNOMIAL_CHEB) {
    1751         return(dChebPolynomial4DEval(w,x,y,z, myPoly));
     1743                          psF64 z,
     1744                          psF64 t)
     1745{
     1746    PS_ASSERT_POLY_NON_NULL(poly, NAN);
     1747
     1748    if (poly->type == PS_POLYNOMIAL_ORD) {
     1749        return(dOrdPolynomial4DEval(x,y,z,t, poly));
     1750    } else if (poly->type == PS_POLYNOMIAL_CHEB) {
     1751        return(dChebPolynomial4DEval(x,y,z,t, poly));
    17521752    } else {
    17531753        psError(PS_ERR_BAD_PARAMETER_TYPE, true,
    17541754                PS_ERRORTEXT_psFunctions_INVALID_POLYNOMIAL_TYPE,
    1755                 myPoly->type);
     1755                poly->type);
    17561756    }
    17571757    return(NAN);
    17581758}
    17591759
    1760 psVector *psDPolynomial4DEvalVector(const psDPolynomial4D *myPoly,
    1761                                     const psVector *w,
     1760psVector *psDPolynomial4DEvalVector(const psDPolynomial4D *poly,
    17621761                                    const psVector *x,
    17631762                                    const psVector *y,
    1764                                     const psVector *z)
    1765 {
    1766     PS_ASSERT_POLY_NON_NULL(myPoly, NULL);
    1767     PS_ASSERT_VECTOR_NON_NULL(w, NULL);
    1768     PS_ASSERT_VECTOR_TYPE(w, PS_TYPE_F64, NULL);
     1763                                    const psVector *z,
     1764                                    const psVector *t)
     1765{
     1766    PS_ASSERT_POLY_NON_NULL(poly, NULL);
    17691767    PS_ASSERT_VECTOR_NON_NULL(x, NULL);
    17701768    PS_ASSERT_VECTOR_TYPE(x, PS_TYPE_F64, NULL);
    … …  
    17731771    PS_ASSERT_VECTOR_NON_NULL(z, NULL);
    17741772    PS_ASSERT_VECTOR_TYPE(z, PS_TYPE_F64, NULL);
     1773    PS_ASSERT_VECTOR_NON_NULL(t, NULL);
     1774    PS_ASSERT_VECTOR_TYPE(t, PS_TYPE_F64, NULL);
    17751775
    17761776    psVector *tmp;
    1777     psS32 vecLen=w->n;
     1777    psS32 vecLen=x->n;
    17781778
    17791779    // Determine the output vector size from min of input vectors
     1780    if (z->n < vecLen) {
     1781        vecLen = z->n;
     1782    }
    17801783    if (y->n < vecLen) {
    17811784        vecLen = y->n;
    17821785    }
    1783     if (x->n < vecLen) {
    1784         vecLen = x->n;
    1785     }
    1786     if (z->n < vecLen) {
    1787         vecLen = z->n;
     1786    if (t->n < vecLen) {
     1787        vecLen = t->n;
    17881788    }
    17891789
    … …  
    17931793    // Evaluate the polynomial
    17941794    for (psS32 i = 0; i < vecLen; i++) {
    1795         tmp->data.F64[i] = psDPolynomial4DEval(myPoly,
    1796                                                w->data.F64[i],
     1795        tmp->data.F64[i] = psDPolynomial4DEval(poly,
    17971796                                               x->data.F64[i],
    17981797                                               y->data.F64[i],
    1799                                                z->data.F64[i]);
     1798                                               z->data.F64[i],
     1799                                               t->data.F64[i]);
    18001800    }
    18011801
Note: See TracChangeset for help on using the changeset viewer.