IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 7105


Ignore:
Timestamp:
May 10, 2006, 3:19:57 AM (20 years ago)
Author:
Paul Price
Message:

Fixing tests: standardised inputs, distributed inputs around [-1,1] so as to avoid numerical effects. tst_psPolyFit1D required special treatment for chebyshev polynomial fitting when x==NULL (indices used as the ordinate)

Location:
trunk/psLib/test/math
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/psLib/test/math/tst_psPolyFit1D.c

    r6484 r7105  
    1313#include "pslib.h"
    1414#include "psTest.h"
    15 #define NUM_DATA 35
     15#define NUM_DATA 100
    1616#define POLY_ORDER 5
    1717#define A 2.0
    … …  
    5252}
    5353
    54 psS32 genericTest(
     54bool genericTest(
    5555    psU32 flags,
    5656    psS32 polyOrder,
    … …  
    5959{
    6060    psS32 currentId = psMemGetId();
    61     psS32 testStatus = true;
     61    bool testStatus = true;
    6262    psS32 memLeaks = 0;
    6363    psPolynomial1D *myPoly = NULL;
    … …  
    9494    if (flags & TS00_X_NULL) {
    9595        printf(" using a NULL x vector\n");
    96     }
     96        numData = 30;
     97    }
     98
     99    psVector *xTruth = psVectorAlloc(numData, PS_TYPE_F64);
     100    psVector *fTruth = psVectorAlloc(numData, PS_TYPE_F64);
     101    xTruth->n = numData;
     102    fTruth->n = numData;
     103    psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS, 1); // Using a known seed
     104    for (int i = 0; i < numData; i++) {
     105        xTruth->data.F64[i] = (flags & TS00_X_NULL) ? i : 2.0*psRandomUniform(rng) - 1.0;
     106        fTruth->data.F64[i] = setData(xTruth->data.F64[i]);
     107    }
     108    if (flags & TS00_X_NULL && flags & TS00_POLY_CHEB) {
     109        // Renormalise the indices
     110        p_psNormalizeVectorRange(xTruth, -1.0, 1.0);
     111    }
     112    psFree(rng);
     113    #if VERBOSE
     114
     115    for (int i = 0; i < numData; i++) {
     116        printf("Original %d: %f\t%f\n", i, xTruth->data.F64[i], fTruth->data.F64[i]);
     117    }
     118    #endif
    97119
    98120    if (flags & TS00_X_F32) {
    99121        printf(" using a psF32 x vector\n");
    100         x = psVectorAlloc(numData, PS_TYPE_F32);
    101         x->n = numData;
    102         for (psS32 i=0;i<numData;i++) {
    103             x->data.F32[i] = (psF32) i;
    104         }
     122        x = psVectorCopy(NULL, xTruth, PS_TYPE_F32);
     123
     124        #if 0
    105125
    106126        if (flags & TS00_POLY_CHEB) {
    107127            p_psNormalizeVectorRange(x, -1.0, 1.0);
    108128        }
     129        #endif
     130
    109131    }
    110132
    111133    if (flags & TS00_X_S32) {
    112134        printf(" using a psS32 x vector\n");
    113         x = psVectorAlloc(numData, PS_TYPE_S32);
    114         x->n = numData;
    115         for (psS32 i=0;i<numData;i++) {
    116             x->data.S32[i] = (psS32) i;
    117         }
     135        x = psVectorCopy(NULL, xTruth, PS_TYPE_S32);
     136
     137        #if 0
    118138
    119139        if (flags & TS00_POLY_CHEB) {
    120140            p_psNormalizeVectorRange(x, -1, 1);
    121141        }
     142        #endif
     143
    122144    }
    123145
    124146    if (flags & TS00_X_F64) {
    125147        printf(" using a psF64 x vector\n");
    126         x = psVectorAlloc(numData, PS_TYPE_F64);
    127         x->n = numData;
    128         for (psS32 i=0;i<numData;i++) {
    129             x->data.F64[i] = (psF64) i;
    130         }
     148        x = psVectorCopy(NULL, xTruth, PS_TYPE_F64);
     149
     150        #if 0
    131151
    132152        if (flags & TS00_POLY_CHEB) {
    133153            p_psNormalizeVectorRange(x, -1.0, 1.0);
    134154        }
     155        #endif
     156
    135157    }
    136158
    … …  
    141163    if (flags & TS00_F_F32) {
    142164        printf(" using a psF32 f vector\n");
    143         f = psVectorAlloc(numData, PS_TYPE_F32);
    144         f->n = numData;
    145         for (psS32 i=0;i<numData;i++) {
    146             f->data.F32[i] = setData((psF32) i);
    147         }
     165        f = psVectorCopy(NULL, fTruth, PS_TYPE_F32);
    148166        // Set a few outliers in the data.
    149167        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    152170            f->data.F32[3*numData/4]*= 2.0;
    153171        }
    154 
    155         if (VERBOSE) {
    156             for (psS32 i=0;i<numData;i++) {
    157                 printf("Original data %d: (%.1f %.1f)\n", i, (psF32) i, f->data.F32[i]);
    158             }
    159         }
    160172    }
    161173
    162174    if (flags & TS00_F_S32) {
    163175        printf(" using a psS32 f vector\n");
    164         f = psVectorAlloc(numData, PS_TYPE_S32);
    165         f->n = numData;
    166         for (psS32 i=0;i<numData;i++) {
    167             f->data.S32[i] = (psS32) setData((psF32) i);
    168         }
     176        f = psVectorCopy(NULL, fTruth, PS_TYPE_S32);
    169177        // Set a few outliers in the data.
    170178        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    174182        }
    175183
    176         if (VERBOSE) {
    177             for (psS32 i=0;i<numData;i++) {
    178                 printf("Original data %d: (%.1f %d)\n", i, (psF32) i, f->data.S32[i]);
    179             }
    180         }
    181184    }
    182185
    183186    if (flags & TS00_F_F64) {
    184187        printf(" using a psF64 f vector\n");
    185         f = psVectorAlloc(numData, PS_TYPE_F64);
    186         f->n = numData;
    187         for (psS32 i=0;i<numData;i++) {
    188             f->data.F64[i] = (psF64) setData((psF32) i);
    189         }
     188        f = psVectorCopy(NULL, fTruth, PS_TYPE_F64);
    190189        // Set a few outliers in the data.
    191190        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    193192            f->data.F64[numData/2]*= 2.0;
    194193            f->data.F64[3*numData/4]*= 2.0;
    195         }
    196 
    197         if (VERBOSE) {
    198             for (psS32 i=0;i<numData;i++) {
    199                 printf("Original data %d: (%.1f %.1f)\n", i, (psF32) i, f->data.F64[i]);
    200             }
    201194        }
    202195    }
    … …  
    279272        }
    280273
    281         for (psS32 i=0 ;i<numData; i++) {
     274        psVector *result = psPolynomial1DEvalVector(myPoly, xTruth);
     275        for (psS32 i=0; i<numData; i++) {
    282276            // Skip the outliers.
    283277            if ((i == numData/4) || (i == numData/2) || (i == 3*numData/4)) {
    284278                continue;
    285279            }
    286             psF32 expectData;
    287             psF32 xData;
    288             if (flags & TS00_F_F32) {
    289                 expectData = f->data.F32[i];
    290             } else if (flags & TS00_F_F64) {
    291                 expectData = (psF32) f->data.F64[i];
    292             } else if (flags & TS00_F_S32) {
    293                 expectData = (psF32) f->data.S32[i];
    294             }
    295 
    296             if (flags & TS00_X_F32) {
    297                 xData = x->data.F32[i];
    298             } else if (flags & TS00_X_F64) {
    299                 xData = (psF32) x->data.F64[i];
    300             } else if (flags & TS00_X_S32) {
    301                 xData = (psF32) x->data.S32[i];
    302             } else if (flags & TS00_X_NULL) {
    303                 if (flags & TS00_POLY_ORD) {
    304                     xData = (psF32) i;
    305                 } else if (flags & TS00_POLY_CHEB) {
    306                     xData = ((2.0 / ((psF32) (numData - 1))) * ((psF32) i)) - 1.0;
    307                 }
    308             }
    309 
    310             psF32 actualData = psPolynomial1DEval(myPoly, xData);
    311 
     280            psF32 expectData = fTruth->data.F64[i];
     281            psF32 actualData = result->data.F64[i];
    312282            if (fabs(actualData-expectData) > fabs(ERROR_TOLERANCE * expectData)) {
    313                 printf("TEST ERROR: Fitted data %d: (%.1f %.1f), expected was (%.1f)\n",
    314                        i, xData, actualData, expectData);
     283                printf("TEST ERROR: Fitted data %d: %.1f --> %.1f vs %.1f\n",
     284                       i, xTruth->data.F64[i], actualData, expectData);
    315285                testStatus = false;
    316286            } else {
    317287                if (VERBOSE) {
    318                     printf("GOOD: Fitted data %d: (%.1f %.1f), expected was (%.1f)\n",
    319                            i, xData, actualData, expectData);
     288                    printf("GOOD: Fitted data %d: %1.f --> %.1f vs %.1f\n",
     289                           i, xTruth->data.F64[i], actualData, expectData);
    320290                }
    321291            }
    322292        }
     293        psFree(result);
    323294    }
    324295
    … …  
    328299    psFree(x);
    329300    psFree(f);
     301    psFree(xTruth);
     302    psFree(fTruth);
    330303    psFree(fErr);
    331304    psFree(stats);
    … …  
    353326    F64 tests: Chebyshev polys, clip fit
    354327 *****************************************************************************/
    355 psS32 main()
     328int main()
    356329{
    357330    psBool testStatus = true;

  • trunk/psLib/test/math/tst_psPolyFit2D.c

    r6484 r7105  
    1111#include "pslib.h"
    1212#include "psTest.h"
    13 #define NUM_DATA 35
     13#define NUM_DATA 100
    1414#define POLY_ORDER_X 2
    1515#define POLY_ORDER_Y 3
    … …  
    9393    }
    9494
     95
     96    psVector *xTruth = psVectorAlloc(numData, PS_TYPE_F64);
     97    psVector *yTruth = psVectorAlloc(numData, PS_TYPE_F64);
     98    psVector *fTruth = psVectorAlloc(numData, PS_TYPE_F64);
     99    xTruth->n = numData;
     100    yTruth->n = numData;
     101    fTruth->n = numData;
     102    psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS, 1); // Using an RNG with a known seed
     103    for (int i = 0; i < numData; i++) {
     104        xTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     105        yTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     106        fTruth->data.F64[i] = setData(xTruth->data.F64[i], yTruth->data.F64[i]);
     107    }
     108    psFree(rng);
     109
    95110    if (flags & TS00_X_NULL) {
    96111        printf(" using a NULL x vector\n");
    … …  
    99114    if (flags & TS00_X_F32) {
    100115        printf(" using a psF32 x vector\n");
    101         x = psVectorAlloc(numData, PS_TYPE_F32);
    102         x->n = numData;
    103         for (psS32 i=0;i<numData;i++) {
    104             x->data.F32[i] = (psF32) i;
    105         }
     116        x = psVectorCopy(NULL, xTruth, PS_TYPE_F32);
    106117    }
    107118
    108119    if (flags & TS00_X_S32) {
    109120        printf(" using a psS32 x vector\n");
    110         x = psVectorAlloc(numData, PS_TYPE_S32);
    111         x->n = numData;
    112         for (psS32 i=0;i<numData;i++) {
    113             x->data.S32[i] = (psS32) i;
    114         }
     121        x = psVectorCopy(NULL, xTruth, PS_TYPE_S32);
    115122    }
    116123
    117124    if (flags & TS00_X_F64) {
    118125        printf(" using a psF64 x vector\n");
    119         x = psVectorAlloc(numData, PS_TYPE_F64);
    120         x->n = numData;
    121         for (psS32 i=0;i<numData;i++) {
    122             x->data.F64[i] = (psF64) i;
    123         }
    124     }
    125 
     126        x = psVectorCopy(NULL, xTruth, PS_TYPE_F64);
     127    }
    126128
    127129    if (flags & TS00_Y_NULL) {
    … …  
    131133    if (flags & TS00_Y_F32) {
    132134        printf(" using a psF32 y vector\n");
    133         y = psVectorAlloc(numData, PS_TYPE_F32);
    134         y->n = numData;
    135         for (psS32 i=0;i<numData;i++) {
    136             y->data.F32[i] = (psF32) i;
    137         }
     135        y = psVectorCopy(NULL, yTruth, PS_TYPE_F32);
    138136    }
    139137
    140138    if (flags & TS00_Y_S32) {
    141139        printf(" using a psS32 y vector\n");
    142         y = psVectorAlloc(numData, PS_TYPE_S32);
    143         y->n = numData;
    144         for (psS32 i=0;i<numData;i++) {
    145             y->data.S32[i] = (psS32) i;
    146         }
     140        y = psVectorCopy(NULL, yTruth, PS_TYPE_S32);
    147141    }
    148142
    149143    if (flags & TS00_Y_F64) {
    150144        printf(" using a psF64 y vector\n");
    151         y = psVectorAlloc(numData, PS_TYPE_F64);
    152         y->n = numData;
    153         for (psS32 i=0;i<numData;i++) {
    154             y->data.F64[i] = (psF64) i;
    155         }
     145        y = psVectorCopy(NULL, yTruth, PS_TYPE_F64);
    156146    }
    157147
    … …  
    162152    if (flags & TS00_F_F32) {
    163153        printf(" using a psF32 f vector\n");
    164         f = psVectorAlloc(numData, PS_TYPE_F32);
    165         f->n = numData;
    166         for (psS32 i=0;i<numData;i++) {
    167             f->data.F32[i] = setData((psF32) i, (psF32) i);
    168         }
     154        f = psVectorCopy(NULL, fTruth, PS_TYPE_F32);
    169155        // Set a few outliers in the data.
    170156        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    183169    if (flags & TS00_F_S32) {
    184170        printf(" using a psS32 f vector\n");
    185         f = psVectorAlloc(numData, PS_TYPE_S32);
    186         f->n = numData;
    187         for (psS32 i=0;i<numData;i++) {
    188             f->data.S32[i] = (psS32) setData((psF32) i, (psF32) i);
    189         }
     171        f = psVectorCopy(NULL, fTruth, PS_TYPE_S32);
    190172        // Set a few outliers in the data.
    191173        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    204186    if (flags & TS00_F_F64) {
    205187        printf(" using a psF64 f vector\n");
    206         f = psVectorAlloc(numData, PS_TYPE_F64);
    207         f->n = numData;
    208         for (psS32 i=0;i<numData;i++) {
    209             f->data.F64[i] = (psF64) setData((psF32) i, (psF32) i);
    210         }
     188        f = psVectorCopy(NULL, fTruth, PS_TYPE_F64);
    211189        // Set a few outliers in the data.
    212190        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    301279        }
    302280
     281        psVector *result = psPolynomial2DEvalVector(myPoly, xTruth, yTruth);
    303282        for (psS32 i=0 ;i<numData; i++) {
    304283            // Skip the outliers.
    … …  
    306285                continue;
    307286            }
    308             psF32 expectData;
    309             psF32 xData;
    310             psF32 yData;
    311             if (flags & TS00_F_F32) {
    312                 expectData = f->data.F32[i];
    313             } else if (flags & TS00_F_F64) {
    314                 expectData = (psF32) f->data.F64[i];
    315             } else if (flags & TS00_F_S32) {
    316                 expectData = (psF32) f->data.S32[i];
    317             }
    318 
    319             if (flags & TS00_X_F32) {
    320                 xData = x->data.F32[i];
    321             } else if (flags & TS00_X_F64) {
    322                 xData = (psF32) x->data.F64[i];
    323             } else if (flags & TS00_X_S32) {
    324                 xData = (psF32) x->data.S32[i];
    325             } else if (flags & TS00_X_NULL) {
    326                 xData = (psF32) i;
    327             }
    328 
    329             if (flags & TS00_Y_F32) {
    330                 yData = y->data.F32[i];
    331             } else if (flags & TS00_Y_F64) {
    332                 yData = (psF32) y->data.F64[i];
    333             } else if (flags & TS00_Y_S32) {
    334                 yData = (psF32) y->data.S32[i];
    335             } else if (flags & TS00_Y_NULL) {
    336                 yData = (psF32) i;
    337             }
    338 
    339             psF32 actualData = psPolynomial2DEval(myPoly, xData, yData);
     287            psF64 actualData = result->data.F64[i];
     288            psF64 expectData = fTruth->data.F64[i];
    340289
    341290            if (fabs(actualData-expectData) > fabs(ERROR_TOLERANCE * expectData)) {
    342                 printf("TEST ERROR: Fitted data %d: (%.1f %.1f), expected was (%.1f)\n",
    343                        i, xData, actualData, expectData);
     291                printf("TEST ERROR: Fitted data %d: (%.1f), expected was (%.1f)\n",
     292                       i, actualData, expectData);
    344293                testStatus = false;
    345294            } else {
    346295                if (VERBOSE) {
    347                     printf("GOOD: Fitted data %d: (%.1f %.1f), expected was (%.1f)\n",
    348                            i, xData, actualData, expectData);
     296                    printf("GOOD: Fitted data %d: (%.1f), expected was (%.1f)\n",
     297                           i, actualData, expectData);
    349298                }
    350299            }
    351300        }
     301        psFree(result);
    352302    }
    353303
    … …  
    355305    psFree(myPoly);
    356306    psFree(mask);
     307    psFree(xTruth);
     308    psFree(yTruth);
     309    psFree(fTruth);
    357310    psFree(x);
    358311    psFree(y);

  • trunk/psLib/test/math/tst_psPolyFit3D.c

    r6484 r7105  
    1111#include "pslib.h"
    1212#include "psTest.h"
    13 #define NUM_DATA 35
     13#define NUM_DATA 30
    1414#define POLY_ORDER_X 2
    1515#define POLY_ORDER_Y 3
    … …  
    2727#define ERROR_TOLERANCE 0.10
    2828#define YERR 10.0
    29 #define VERBOSE 0
     29#define VERBOSE 1
    3030#define NUM_ITERATIONS 5
    3131#define CLIP_SIGMA 4.0
    … …  
    9595    printPositiveTestHeader(stdout, "psMinimize functions", "3D Polynomial Fitting Functions");
    9696
     97    psVector *xTruth = psVectorAlloc(numData, PS_TYPE_F64);
     98    psVector *yTruth = psVectorAlloc(numData, PS_TYPE_F64);
     99    psVector *zTruth = psVectorAlloc(numData, PS_TYPE_F64);
     100    psVector *fTruth = psVectorAlloc(numData, PS_TYPE_F64);
     101    xTruth->n = numData;
     102    yTruth->n = numData;
     103    zTruth->n = numData;
     104    fTruth->n = numData;
     105    psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS, 1); // Using known seed
     106    for (int i = 0; i < numData; i++) {
     107        xTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     108        yTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     109        zTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     110        fTruth->data.F64[i] = setData(xTruth->data.F64[i], yTruth->data.F64[i], zTruth->data.F64[i]);
     111    }
     112    psFree(rng);
     113
    97114    if (expectedRC == false) {
    98115        printf("This test should generate an error message, and return NULL.\n");
    … …  
    116133    if (flags & TS00_X_F32) {
    117134        printf(" using a psF32 x vector\n");
    118         x = psVectorAlloc(numData, PS_TYPE_F32);
    119         x->n = numData;
    120         for (psS32 i=0;i<numData;i++) {
    121             x->data.F32[i] = (psF32) i;
    122         }
     135        x = psVectorCopy(NULL, xTruth, PS_TYPE_F32);
    123136    }
    124137
    125138    if (flags & TS00_X_S32) {
    126139        printf(" using a psS32 x vector\n");
    127         x = psVectorAlloc(numData, PS_TYPE_S32);
    128         x->n = numData;
    129         for (psS32 i=0;i<numData;i++) {
    130             x->data.S32[i] = (psS32) i;
    131         }
     140        x = psVectorCopy(NULL, xTruth, PS_TYPE_S32);
    132141    }
    133142
    134143    if (flags & TS00_X_F64) {
    135144        printf(" using a psF64 x vector\n");
    136         x = psVectorAlloc(numData, PS_TYPE_F64);
    137         x->n = numData;
    138         for (psS32 i=0;i<numData;i++) {
    139             x->data.F64[i] = (psF64) i;
    140         }
     145        x = psVectorCopy(NULL, xTruth, PS_TYPE_F64);
    141146    }
    142147
    … …  
    148153    if (flags & TS00_Y_F32) {
    149154        printf(" using a psF32 y vector\n");
    150         y = psVectorAlloc(numData, PS_TYPE_F32);
    151         y->n = numData;
    152         for (psS32 i=0;i<numData;i++) {
    153             y->data.F32[i] = (psF32) i;
    154         }
     155        y = psVectorCopy(NULL, yTruth, PS_TYPE_F32);
    155156    }
    156157
    157158    if (flags & TS00_Y_S32) {
    158159        printf(" using a psS32 y vector\n");
    159         y = psVectorAlloc(numData, PS_TYPE_S32);
    160         y->n = numData;
    161         for (psS32 i=0;i<numData;i++) {
    162             y->data.S32[i] = (psS32) i;
    163         }
     160        y = psVectorCopy(NULL, yTruth, PS_TYPE_S32);
    164161    }
    165162
    166163    if (flags & TS00_Y_F64) {
    167164        printf(" using a psF64 y vector\n");
    168         y = psVectorAlloc(numData, PS_TYPE_F64);
    169         y->n = numData;
    170         for (psS32 i=0;i<numData;i++) {
    171             y->data.F64[i] = (psF64) i;
    172         }
     165        y = psVectorCopy(NULL, yTruth, PS_TYPE_F64);
    173166    }
    174167
    … …  
    179172    if (flags & TS00_Z_F32) {
    180173        printf(" using a psF32 z vector\n");
    181         z = psVectorAlloc(numData, PS_TYPE_F32);
    182         z->n = numData;
    183         for (psS32 i=0;i<numData;i++) {
    184             z->data.F32[i] = (psF32) i;
    185         }
     174        z = psVectorCopy(NULL, zTruth, PS_TYPE_F32);
    186175    }
    187176
    188177    if (flags & TS00_Z_S32) {
    189178        printf(" using a psS32 z vector\n");
    190         z = psVectorAlloc(numData, PS_TYPE_S32);
    191         z->n = numData;
    192         for (psS32 i=0;i<numData;i++) {
    193             z->data.S32[i] = (psS32) i;
    194         }
     179        z = psVectorCopy(NULL, zTruth, PS_TYPE_S32);
    195180    }
    196181
    197182    if (flags & TS00_Z_F64) {
    198183        printf(" using a psF64 z vector\n");
    199         z = psVectorAlloc(numData, PS_TYPE_F64);
    200         z->n = numData;
    201         for (psS32 i=0;i<numData;i++) {
    202             z->data.F64[i] = (psF64) i;
    203         }
     184        z = psVectorCopy(NULL, zTruth, PS_TYPE_F64);
    204185    }
    205186
    … …  
    211192    if (flags & TS00_F_F32) {
    212193        printf(" using a psF32 f vector\n");
    213         f = psVectorAlloc(numData, PS_TYPE_F32);
    214         f->n = numData;
    215         for (psS32 i=0;i<numData;i++) {
    216             f->data.F32[i] = setData((psF32) i, (psF32) i, (psF32) i);
    217         }
     194        f = psVectorCopy(NULL, fTruth, PS_TYPE_F32);
    218195        // Set a few outliers in the data.
    219196        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    225202        if (VERBOSE) {
    226203            for (psS32 i=0;i<numData;i++) {
    227                 printf("Original data %d: (%.1f %.1f)\n", i, (psF32) i, f->data.F32[i]);
     204                printf("Original data %d: (%.1f)\n", i, f->data.F32[i]);
    228205            }
    229206        }
    … …  
    232209    if (flags & TS00_F_S32) {
    233210        printf(" using a psS32 f vector\n");
    234         f = psVectorAlloc(numData, PS_TYPE_S32);
    235         f->n = numData;
    236         for (psS32 i=0;i<numData;i++) {
    237             f->data.S32[i] = (psS32) setData((psF32) i, (psF32) i, (psF32) i);
    238         }
     211        f = psVectorCopy(NULL, fTruth, PS_TYPE_S32);
    239212        // Set a few outliers in the data.
    240213        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    246219        if (VERBOSE) {
    247220            for (psS32 i=0;i<numData;i++) {
    248                 printf("Original data %d: (%.1f %d)\n", i, (psF32) i, f->data.S32[i]);
     221                printf("Original data %d: (%d)\n", i, f->data.S32[i]);
    249222            }
    250223        }
    … …  
    253226    if (flags & TS00_F_F64) {
    254227        printf(" using a psF64 f vector\n");
    255         f = psVectorAlloc(numData, PS_TYPE_F64);
    256         f->n = numData;
    257         for (psS32 i=0;i<numData;i++) {
    258             f->data.F64[i] = (psF64) setData((psF32) i, (psF32) i, (psF32) i);
    259         }
     228        f = psVectorCopy(NULL, fTruth, PS_TYPE_F64);
    260229        // Set a few outliers in the data.
    261230        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    267236        if (VERBOSE) {
    268237            for (psS32 i=0;i<numData;i++) {
    269                 printf("Original data %d: (%.1f %.1f)\n", i, (psF32) i, f->data.F64[i]);
     238                printf("Original data %d: (%.1f)\n", i, f->data.F64[i]);
    270239            }
    271240        }
    … …  
    353322        }
    354323
     324        psVector *result = psPolynomial3DEvalVector(myPoly, xTruth, yTruth, zTruth);
    355325        for (psS32 i=0 ;i<numData; i++) {
    356326            // Skip the outliers.
    … …  
    358328                continue;
    359329            }
    360             psF32 expectData;
    361             psF32 xData;
    362             psF32 yData;
    363             psF32 zData;
    364             if (flags & TS00_F_F32) {
    365                 expectData = f->data.F32[i];
    366             } else if (flags & TS00_F_F64) {
    367                 expectData = (psF32) f->data.F64[i];
    368             } else if (flags & TS00_F_S32) {
    369                 expectData = (psF32) f->data.S32[i];
    370             }
    371 
    372             if (flags & TS00_X_F32) {
    373                 xData = x->data.F32[i];
    374             } else if (flags & TS00_X_F64) {
    375                 xData = (psF32) x->data.F64[i];
    376             } else if (flags & TS00_X_S32) {
    377                 xData = (psF32) x->data.S32[i];
    378             } else if (flags & TS00_X_NULL) {
    379                 xData = (psF32) i;
    380             }
    381 
    382             if (flags & TS00_Y_F32) {
    383                 yData = y->data.F32[i];
    384             } else if (flags & TS00_Y_F64) {
    385                 yData = (psF32) y->data.F64[i];
    386             } else if (flags & TS00_Y_S32) {
    387                 yData = (psF32) y->data.S32[i];
    388             } else if (flags & TS00_Y_NULL) {
    389                 yData = (psF32) i;
    390             }
    391 
    392             if (flags & TS00_Z_F32) {
    393                 zData = z->data.F32[i];
    394             } else if (flags & TS00_Z_F64) {
    395                 zData = (psF32) z->data.F64[i];
    396             } else if (flags & TS00_Z_S32) {
    397                 zData = (psF32) z->data.S32[i];
    398             } else if (flags & TS00_Z_NULL) {
    399                 zData = (psF32) i;
    400             }
    401 
    402             psF32 actualData = psPolynomial3DEval(myPoly, xData, yData, zData);
     330            psF32 expectData = fTruth->data.F64[i];
     331            psF32 actualData = result->data.F64[i];
    403332
    404333            if (fabs(actualData-expectData) > fabs(ERROR_TOLERANCE * expectData)) {
    405                 printf("TEST ERROR: Fitted data %d: (%.1f %.1f), expected was (%.1f)\n",
    406                        i, xData, actualData, expectData);
     334                printf("TEST ERROR: Fitted data %d: (%.1f), expected was (%.1f)\n",
     335                       i, actualData, expectData);
    407336                testStatus = false;
    408337            } else {
    409338                if (VERBOSE) {
    410                     printf("GOOD: Fitted data %d: (%.1f %.1f), expected was (%.1f)\n",
    411                            i, xData, actualData, expectData);
     339                    printf("GOOD: Fitted data %d: (%.1f), expected was (%.1f)\n",
     340                           i, actualData, expectData);
    412341                }
    413342            }
    414343        }
     344        psFree(result);
    415345    }
    416346
    … …  
    422352    psFree(z);
    423353    psFree(f);
     354    psFree(xTruth);
     355    psFree(yTruth);
     356    psFree(zTruth);
     357    psFree(fTruth);
    424358    psFree(fErr);
    425359    psFree(stats);

  • trunk/psLib/test/math/tst_psPolyFit4D.c

    r6527 r7105  
    1111#include "pslib.h"
    1212#include "psTest.h"
    13 #define NUM_DATA 35
     13#define NUM_DATA 30
    1414#define POLY_ORDER_X 2
    15 #define POLY_ORDER_Y 3
     15#define POLY_ORDER_Y 2
    1616#define POLY_ORDER_Z 2
    1717#define POLY_ORDER_T 2
    … …  
    3333#define ERROR_TOLERANCE 0.10
    3434#define YERR 10.0
    35 #define VERBOSE 0
     35#define VERBOSE 1
    3636#define NUM_ITERATIONS 5
    3737#define CLIP_SIGMA 4.0
    … …  
    8282}
    8383
    84 psS32 genericTest(
     84bool genericTest(
    8585    psU32 flags,
    8686    psS32 polyOrderX,
    … …  
    9292{
    9393    psS32 currentId = psMemGetId();
    94     psS32 testStatus = true;
     94    bool testStatus = true;
    9595    psS32 memLeaks = 0;
    9696    psPolynomial4D *myPoly = NULL;
    … …  
    108108    printPositiveTestHeader(stdout, "psMinimize functions", "4D Polynomial Fitting Functions");
    109109
     110    psVector *xTruth = psVectorAlloc(numData, PS_TYPE_F64);
     111    psVector *yTruth = psVectorAlloc(numData, PS_TYPE_F64);
     112    psVector *zTruth = psVectorAlloc(numData, PS_TYPE_F64);
     113    psVector *tTruth = psVectorAlloc(numData, PS_TYPE_F64);
     114    psVector *fTruth = psVectorAlloc(numData, PS_TYPE_F64);
     115    xTruth->n = numData;
     116    yTruth->n = numData;
     117    zTruth->n = numData;
     118    tTruth->n = numData;
     119    fTruth->n = numData;
     120    psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS, 1); // Using known seed
     121    for (int i = 0; i < numData; i++) {
     122        xTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     123        yTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     124        zTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     125        tTruth->data.F64[i] = 2.0*psRandomUniform(rng) - 1.0;
     126        fTruth->data.F64[i] = setData(xTruth->data.F64[i], yTruth->data.F64[i],
     127                                      zTruth->data.F64[i], tTruth->data.F64[i]);
     128    }
     129    psFree(rng);
     130
    110131    if (expectedRC == false) {
    111132        printf("This test should generate an error message, and return NULL.\n");
    … …  
    129150    if (flags & TS00_X_F32) {
    130151        printf(" using a psF32 x vector\n");
    131         x = psVectorAlloc(numData, PS_TYPE_F32);
    132         x->n = numData;
    133         for (psS32 i=0;i<numData;i++) {
    134             x->data.F32[i] = (psF32) i;
    135         }
     152        x = psVectorCopy(NULL, xTruth, PS_TYPE_F32);
    136153    }
    137154
    138155    if (flags & TS00_X_S32) {
    139156        printf(" using a psS32 x vector\n");
    140         x = psVectorAlloc(numData, PS_TYPE_S32);
    141         x->n = numData;
    142         for (psS32 i=0;i<numData;i++) {
    143             x->data.S32[i] = (psS32) i;
    144         }
     157        x = psVectorCopy(NULL, xTruth, PS_TYPE_S32);
    145158    }
    146159
    147160    if (flags & TS00_X_F64) {
    148161        printf(" using a psF64 x vector\n");
    149         x = psVectorAlloc(numData, PS_TYPE_F64);
    150         x->n = numData;
    151         for (psS32 i=0;i<numData;i++) {
    152             x->data.F64[i] = (psF64) i;
    153         }
     162        x = psVectorCopy(NULL, xTruth, PS_TYPE_F64);
    154163    }
    155164
    … …  
    161170    if (flags & TS00_Y_F32) {
    162171        printf(" using a psF32 y vector\n");
    163         y = psVectorAlloc(numData, PS_TYPE_F32);
    164         y->n = numData;
    165         for (psS32 i=0;i<numData;i++) {
    166             y->data.F32[i] = (psF32) i;
    167         }
     172        y = psVectorCopy(NULL, yTruth, PS_TYPE_F32);
    168173    }
    169174
    170175    if (flags & TS00_Y_S32) {
    171176        printf(" using a psS32 y vector\n");
    172         y = psVectorAlloc(numData, PS_TYPE_S32);
    173         y->n = numData;
    174         for (psS32 i=0;i<numData;i++) {
    175             y->data.S32[i] = (psS32) i;
    176         }
     177        y = psVectorCopy(NULL, yTruth, PS_TYPE_S32);
    177178    }
    178179
    179180    if (flags & TS00_Y_F64) {
    180181        printf(" using a psF64 y vector\n");
    181         y = psVectorAlloc(numData, PS_TYPE_F64);
    182         y->n = numData;
    183         for (psS32 i=0;i<numData;i++) {
    184             y->data.F64[i] = (psF64) i;
    185         }
     182        y = psVectorCopy(NULL, yTruth, PS_TYPE_F64);
    186183    }
    187184
    … …  
    192189    if (flags & TS00_Z_F32) {
    193190        printf(" using a psF32 z vector\n");
    194         z = psVectorAlloc(numData, PS_TYPE_F32);
    195         z->n = numData;
    196         for (psS32 i=0;i<numData;i++) {
    197             z->data.F32[i] = (psF32) i;
    198         }
     191        z = psVectorCopy(NULL, zTruth, PS_TYPE_F32);
    199192    }
    200193
    201194    if (flags & TS00_Z_S32) {
    202195        printf(" using a psS32 z vector\n");
    203         z = psVectorAlloc(numData, PS_TYPE_S32);
    204         z->n = numData;
    205         for (psS32 i=0;i<numData;i++) {
    206             z->data.S32[i] = (psS32) i;
    207         }
     196        z = psVectorCopy(NULL, zTruth, PS_TYPE_S32);
    208197    }
    209198
    210199    if (flags & TS00_Z_F64) {
    211200        printf(" using a psF64 z vector\n");
    212         z = psVectorAlloc(numData, PS_TYPE_F64);
    213         z->n = numData;
    214         for (psS32 i=0;i<numData;i++) {
    215             z->data.F64[i] = (psF64) i;
    216         }
     201        z = psVectorCopy(NULL, zTruth, PS_TYPE_F64);
    217202    }
    218203
    … …  
    223208    if (flags & TS00_T_F32) {
    224209        printf(" using a psF32 t vector\n");
    225         t = psVectorAlloc(numData, PS_TYPE_F32);
    226         t->n = numData;
    227         for (psS32 i=0;i<numData;i++) {
    228             t->data.F32[i] = (psF32) i;
    229         }
     210        t = psVectorCopy(NULL, tTruth, PS_TYPE_F32);
    230211    }
    231212
    232213    if (flags & TS00_T_S32) {
    233214        printf(" using a psS32 t vector\n");
    234         t = psVectorAlloc(numData, PS_TYPE_S32);
    235         t->n = numData;
    236         for (psS32 i=0;i<numData;i++) {
    237             t->data.S32[i] = (psS32) i;
    238         }
     215        t = psVectorCopy(NULL, tTruth, PS_TYPE_S32);
    239216    }
    240217
    241218    if (flags & TS00_T_F64) {
    242219        printf(" using a psF64 t vector\n");
    243         t = psVectorAlloc(numData, PS_TYPE_F64);
    244         t->n = numData;
    245         for (psS32 i=0;i<numData;i++) {
    246             t->data.F64[i] = (psF64) i;
    247         }
     220        t = psVectorCopy(NULL, tTruth, PS_TYPE_F64);
    248221    }
    249222
    … …  
    255228    if (flags & TS00_F_F32) {
    256229        printf(" using a psF32 f vector\n");
    257         f = psVectorAlloc(numData, PS_TYPE_F32);
    258         f->n = numData;
    259         for (psS32 i=0;i<numData;i++) {
    260             f->data.F32[i] = setData((psF32) i, (psF32) i, (psF32) i, (psF32) i);
    261         }
     230        f = psVectorCopy(NULL, fTruth, PS_TYPE_F32);
    262231        // Set a few outliers in the data.
    263232        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    266235            f->data.F32[3*numData/4]*= 2.0;
    267236        }
    268 
    269         if (VERBOSE) {
    270             for (psS32 i=0;i<numData;i++) {
    271                 printf("Original data %d: (%.1f %.1f)\n", i, (psF32) i, f->data.F32[i]);
    272             }
    273         }
    274237    }
    275238
    276239    if (flags & TS00_F_S32) {
    277240        printf(" using a psS32 f vector\n");
    278         f = psVectorAlloc(numData, PS_TYPE_S32);
    279         f->n = numData;
    280         for (psS32 i=0;i<numData;i++) {
    281             f->data.S32[i] = (psS32) setData((psF32) i, (psF32) i, (psF32) i, (psF32) i);
    282         }
     241        f = psVectorCopy(NULL, fTruth, PS_TYPE_S32);
    283242        // Set a few outliers in the data.
    284243        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    287246            f->data.S32[3*numData/4]*= 2.0;
    288247        }
    289 
    290         if (VERBOSE) {
    291             for (psS32 i=0;i<numData;i++) {
    292                 printf("Original data %d: (%.1f %d)\n", i, (psF32) i, f->data.S32[i]);
    293             }
    294         }
    295248    }
    296249
    297250    if (flags & TS00_F_F64) {
    298251        printf(" using a psF64 f vector\n");
    299         f = psVectorAlloc(numData, PS_TYPE_F64);
    300         f->n = numData;
    301         for (psS32 i=0;i<numData;i++) {
    302             f->data.F64[i] = (psF64) setData((psF32) i, (psF32) i, (psF32) i, (psF32) i);
    303         }
     252        f = psVectorCopy(NULL, fTruth, PS_TYPE_F64);
    304253        // Set a few outliers in the data.
    305254        if (OUTLIERS && (flags & TS00_CLIP_FIT)) {
    … …  
    307256            f->data.F64[numData/2]*= 2.0;
    308257            f->data.F64[3*numData/4]*= 2.0;
    309         }
    310 
    311         if (VERBOSE) {
    312             for (psS32 i=0;i<numData;i++) {
    313                 printf("Original data %d: (%.1f %.1f)\n", i, (psF32) i, f->data.F64[i]);
    314             }
    315258        }
    316259    }
    … …  
    399342        }
    400343
    401         for (psS32 i=0 ;i<numData; i++) {
     344        psVector *result = psPolynomial4DEvalVector(myPoly, xTruth, yTruth, zTruth, tTruth);
     345        for (int i=0;i<numData; i++) {
    402346            // Skip the outliers.
    403347            if ((i == numData/4) || (i == numData/2) || (i == 3*numData/4)) {
    404348                continue;
    405349            }
    406             psF32 expectData;
    407             psF32 xData;
    408             psF32 yData;
    409             psF32 zData;
    410             psF32 tData;
    411             if (flags & TS00_F_F32) {
    412                 expectData = f->data.F32[i];
    413             } else if (flags & TS00_F_F64) {
    414                 expectData = (psF32) f->data.F64[i];
    415             } else if (flags & TS00_F_S32) {
    416                 expectData = (psF32) f->data.S32[i];
    417             }
    418 
    419             if (flags & TS00_X_F32) {
    420                 xData = x->data.F32[i];
    421             } else if (flags & TS00_X_F64) {
    422                 xData = (psF32) x->data.F64[i];
    423             } else if (flags & TS00_X_S32) {
    424                 xData = (psF32) x->data.S32[i];
    425             } else if (flags & TS00_X_NULL) {
    426                 xData = (psF32) i;
    427             }
    428 
    429             if (flags & TS00_Y_F32) {
    430                 yData = y->data.F32[i];
    431             } else if (flags & TS00_Y_F64) {
    432                 yData = (psF32) y->data.F64[i];
    433             } else if (flags & TS00_Y_S32) {
    434                 yData = (psF32) y->data.S32[i];
    435             } else if (flags & TS00_Y_NULL) {
    436                 yData = (psF32) i;
    437             }
    438 
    439             if (flags & TS00_Z_F32) {
    440                 zData = z->data.F32[i];
    441             } else if (flags & TS00_Z_F64) {
    442                 zData = (psF32) z->data.F64[i];
    443             } else if (flags & TS00_Z_S32) {
    444                 zData = (psF32) z->data.S32[i];
    445             } else if (flags & TS00_Z_NULL) {
    446                 zData = (psF32) i;
    447             }
    448 
    449             if (flags & TS00_T_F32) {
    450                 tData = t->data.F32[i];
    451             } else if (flags & TS00_T_F64) {
    452                 tData = (psF32) t->data.F64[i];
    453             } else if (flags & TS00_T_S32) {
    454                 tData = (psF32) t->data.S32[i];
    455             } else if (flags & TS00_T_NULL) {
    456                 tData = (psF32) i;
    457             }
    458 
    459             psF32 actualData = psPolynomial4DEval(myPoly, xData, yData, zData, tData);
     350            psF32 expectData = fTruth->data.F64[i];
     351            psF32 actualData = result->data.F64[i];
    460352
    461353            if (fabs(actualData-expectData) > fabs(ERROR_TOLERANCE * expectData)) {
    462                 printf("TEST ERROR: Fitted data %d: (%.1f %.1f), expected was (%.1f)\n",
    463                        i, xData, actualData, expectData);
     354                printf("TEST ERROR: Fitted data %d: (%.1f), expected was (%.1f)\n",
     355                       i, actualData, expectData);
    464356                testStatus = false;
    465357            } else {
    466358                if (VERBOSE) {
    467                     printf("GOOD: Fitted data %d: (%.1f %.1f), expected was (%.1f)\n",
    468                            i, xData, actualData, expectData);
     359                    printf("GOOD: Fitted data %d: (%.1f), expected was (%.1f)\n",
     360                           i, actualData, expectData);
    469361                }
    470362            }
    471363        }
     364        psFree(result);
    472365    }
    473366
    … …  
    480373    psFree(t);
    481374    psFree(f);
     375    psFree(xTruth);
     376    psFree(yTruth);
     377    psFree(zTruth);
     378    psFree(tTruth);
     379    psFree(fTruth);
    482380    psFree(fErr);
    483381    psFree(stats);
    … …  
    505403    F64 tests: Chebyshev polys, clip fit
    506404 *****************************************************************************/
    507 psS32 main()
     405int main()
    508406{
    509407    psBool testStatus = true;
    … …  
    611509    testStatus &= genericTest(TS00_MASK_S32 | TS00_FERR_F64 | TS00_X_F64 | TS00_Y_F64 | TS00_Z_F64 | TS00_T_F64 | TS00_F_F64 | TS00_POLY_ORD | TS00_CLIP_FIT, POLY_ORDER_X, POLY_ORDER_Y, POLY_ORDER_Z, POLY_ORDER_T, NUM_DATA, false);
    612510
     511    printf("Status: %x\n", testStatus);
     512
    613513    printFooter(stdout, "psMinimize functions: 4D Polynomial Fitting Functions", "", testStatus);
    614514}
Note: See TracChangeset for help on using the changeset viewer.