IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 5873


Ignore:
Timestamp:
Jan 3, 2006, 11:12:54 AM (21 years ago)
Author:
bswift
Message:

Added more polynomial tests.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/pstest/src/tst_psPolynomial.c

    r5613 r5873  
    11#include <stdio.h>
     2#include <math.h>
    23#include "pslib_strict.h"
    34#include "psTest.h"
    45
    56static psS32 testPolynomial1DEval_01(void);
     7static psS32 testPolynomial1DEvalVector_01(void);
     8static psS32 testPolynomial1DEvalVector_02(void);
     9static psS32 testVectorFitPolynomial1D_01(void);
     10static psS32 testVectorFitPolynomial1D_02(void);
     11static psS32 testVectorClipFitPolynomial1D_01(void);
    612
    713/* list of test functions in test suite */
    814testDescription tests[] = {
    915                              {testPolynomial1DEval_01,1,"psPolynomial1DEval_01",0,false},
     16                              {testPolynomial1DEvalVector_01,2,"psPolynomial1DEvalVector_01",0,false},
     17                              {testPolynomial1DEvalVector_02,2,"psPolynomial1DEvalVector_02",0,false},
     18                              {testVectorFitPolynomial1D_01,2,"psVectorFitPolynomial1D_01",0,false},
     19                              {testVectorFitPolynomial1D_02,2,"psVectorFitPolynomial1D_02",0,false},
     20                              {testVectorClipFitPolynomial1D_01,2,"psVectorFitClipPolynomial1D_01",0,false},
    1021                              {NULL}
    1122                          };
    … …  
    1425int main(int argc, char* argv[])
    1526{
     27    //psTraceSetLevel("VectorFitPolynomial1DOrd", 6);
    1628    psLogSetFormat("HLNM");
    1729    psLogSetLevel(PS_LOG_INFO);
    … …  
    2436# define X3 -2.0
    2537
    26 static psPolynomial1D *makePolynomial1D(int Nord) {
     38
     39static psPolynomial1D *makePolynomial1D(int Nord, psPolynomialType type, int ordMask) {
    2740
    2841    // define the polynomial
    29     psPolynomial1D *poly = psPolynomial1DAlloc(Nord, PS_POLYNOMIAL_ORD);
     42    psPolynomial1D *poly = psPolynomial1DAlloc(Nord, type);
    3043    switch (Nord) {
    3144      case 3:
    32         poly->coeff[3] = X3;
     45        poly->coeff[3] = (ordMask==3) ?  0:X3;
    3346      case 2:
    34         poly->coeff[2] = X2;
     47        poly->coeff[2] = (ordMask==2) ?  0:X2;
    3548      case 1:
    36         poly->coeff[1] = X1;
     49        poly->coeff[1] = (ordMask==1) ?  0:X1;
    3750      case 0:
    38         poly->coeff[0] = X0;
     51        poly->coeff[0] = (ordMask==0) ?  0:X0;
    3952        break;
    4053      default:
    … …  
    4457}
    4558
    46 static void fillVectors (psVector *x, psVector *y, int Nord, int Npts) {
    47 
    48     psF32 X;
     59static void fillVectors(psVector *x, psVector *y, int Nord, int Npts, psPolynomialType type) {
     60
     61    psF32 X;//, P1, P2, P3;
     62    double c[] = {(double)X0,(double)X1,(double)X2,(double)X3};
     63    //double sv,d,dd;
     64    double sum, Tn, Tn_1, Tn_2;
    4965
    5066    for (int i = 0; i < Npts; i++) {
    51         X = i - 0.5*Npts;
     67        X = (float)i/Npts - 0.5;
    5268        x->data.F32[i] = X;
    53         switch (Nord) {
    54           case 0:
    55             y->data.F32[i] = X0;
     69        switch (type) {
     70          case PS_POLYNOMIAL_ORD:
     71            switch (Nord) {
     72              case 0:
     73                y->data.F32[i] = c[0];
     74                break;
     75              case 1:
     76                y->data.F32[i] = c[0] + c[1]*X;
     77                break;
     78              case 2:
     79                y->data.F32[i] = c[0] + c[1]*X + c[2]*X*X;
     80                break;
     81              case 3:
     82                y->data.F32[i] = c[0] + c[1]*X + c[2]*X*X + c[3]*X*X*X;
     83                break;
     84              default:
     85                y->data.F32[i] = 0;
     86            }
    5687            break;
    57           case 1:
    58             y->data.F32[i] = X0 + X*X1;
     88          case PS_POLYNOMIAL_CHEB:
     89            switch (Nord) {
     90              case 0:
     91                y->data.F32[i] = c[0];
     92                break;
     93              case 1:
     94                y->data.F32[i] =c[0] + c[1]*X;
     95                break;
     96              default:
     97                //Clenshaw's formula approach (NR 5.5)
     98                //d=0.0;
     99                //dd=0.0;
     100                //for (int j=Nord; j>=1; j--) {
     101                //    sv = d;
     102                //    d = (double)2.0*(double)X*d - dd + c[j];
     103                //    dd = sv;
     104                //}
     105                //y->data.F32[i] = (double)X*d-dd+(double)0.5*c[0];
     106
     107                //Straightforward approach
     108                Tn_1=X; //T1
     109                Tn_2=1; //T0
     110                sum=0.;
     111                for (int j=2; j<=Nord; j++) {
     112                    Tn = 2.0*X*Tn_1 - Tn_2;
     113                    sum = sum + c[j]*Tn;
     114                    Tn_2 = Tn_1;
     115                    Tn_1 = Tn;
     116                }
     117                y->data.F32[i] = sum + c[1]*X + c[0] - 0.5*c[0];
     118            }
    59119            break;
    60           case 2:
    61             y->data.F32[i] = X0 + X*X1 + X*X*X2;
    62             break;
    63           case 3:
    64             y->data.F32[i] = X0 + X*X1 + X*X*X2 + X*X*X*X3;
    65             break;
    66           default:
    67             y->data.F32[i] = 0;
    68120        }
    69121    }
    … …  
    71123}
    72124
     125# define NPTS 10
     126
    73127int testPolynomial1DEval_01(void) {
    74 
    75 # define NPTS 100
    76 # define NORD 3
    77 
    78     psPolynomial1D *poly = makePolynomial1D(NORD);
    79 
    80     // create x vector, known y vector:
    81     psVector *x = psVectorAlloc (NPTS, PS_TYPE_F32);
    82     psVector *y = psVectorAlloc (NPTS, PS_TYPE_F32);
    83     fillVectors (x, y, NORD, NPTS);
    84 
    85     psF32 xV, yV;
    86 
    87     int Ngood = 0;
    88     int Nbad = 0;
    89 
    90     for (int i = 0; i < NPTS; i++) {
    91         xV = x->data.F32[i];
    92         yV = psPolynomial1DEval (poly, xV);
    93         if (yV == y->data.F32[i]) {
    94             Ngood ++;
    95         } else {
    96             Nbad ++;
    97         }
    98     }
    99     psFree (x);
    100     psFree (y);
    101     psFree (poly);
    102 
    103     if (Nbad > 0) {
    104         psError(PS_ERR_UNKNOWN,true,"%d incorrect polynomial evaluations\n", Nbad);
    105         return 1;
    106     }
    107     return 0;
    108 }
     128   
     129    //# define NORD 3
     130    int NORD;
     131    int maxNORD = 3;
     132    int retval = 0;
     133    psPolynomialType type[] = {PS_POLYNOMIAL_ORD, PS_POLYNOMIAL_CHEB};
     134
     135    for (int k = 0; k<1; k++) {  //loop over the 2 polynomial types
     136        for (int j = 1; j <= maxNORD; j++) {  //loop over the polynomial orders
     137            int nGood = 0, nBad = 0;
     138            NORD = j;
     139            //printf("Testing NORD=%d\n", NORD);
     140            psPolynomial1D *poly = makePolynomial1D(NORD, type[k], -1);
     141           
     142            //create x vector, known y vector:
     143            psVector *x = psVectorAlloc (NPTS, PS_TYPE_F32);
     144            psVector *y = psVectorAlloc (NPTS, PS_TYPE_F32);
     145            fillVectors (x, y, NORD, NPTS, type[k]);
     146           
     147            psF32 xV, yV;
     148           
     149            for (int i = 0; i < NPTS; i++) {
     150                xV = x->data.F32[i];
     151                yV = psPolynomial1DEval (poly, xV);
     152                //Change some values around to make it fail!
     153                //if (i>90 && NORD==2) {yV = yV + 0.1;}
     154                if (yV == y->data.F32[i]) {
     155                    nGood++;
     156                } else {
     157                    nBad++;
     158                    printf("%.3f   -->  %.10f  & %.10f\tdel=%e\n",xV,yV,y->data.F32[i],y->data.F32[i]-yV);
     159                }
     160            }
     161            psFree (x);
     162            psFree (y);
     163            psFree (poly);
     164           
     165            if (nBad > 0) {
     166                psError(PS_ERR_UNKNOWN,true,"%d incorrect polynomial evaluations for NORD=%d\n", nBad, NORD);
     167                retval = retval + pow(2,NORD);
     168            }
     169        }
     170    }
     171    return retval;
     172}
     173
     174
     175int testPolynomial1DEvalVector_01(void) {
     176   
     177    //# define NORD 3
     178    int NORD;
     179    int maxNORD = 3;
     180    int retval = 0;
     181    psPolynomialType type[] = {PS_POLYNOMIAL_ORD, PS_POLYNOMIAL_CHEB};
     182   
     183    for (int k = 0; k<1; k++) {  //loop over the 2 polynomial types
     184        for (int j = 1; j <= maxNORD; j++) {  //loop over the polynomial orders
     185            int nGood = 0, nBad = 0;
     186            NORD = j;
     187            //printf("Testing NORD=%d\n", NORD);
     188            psPolynomial1D *poly = makePolynomial1D(NORD, type[k], -1);
     189           
     190            // create x vector, known y vector:
     191            psVector *x = psVectorAlloc (NPTS, PS_TYPE_F32);
     192            psVector *y = psVectorAlloc (NPTS, PS_TYPE_F32);
     193            fillVectors (x, y, NORD, NPTS, type[k]);
     194           
     195            psVector *xV = x;
     196            psVector *yV = psPolynomial1DEvalVector(poly, xV);
     197           
     198            for (int i = 0; i < NPTS; i++) {
     199                //Change some values around to make it fail!
     200                //if (i>90 && NORD==2) {yV->data.F32[i] = yV->data.F32[i] + 0.1;}
     201                if (yV->data.F32[i] == y->data.F32[i]) {
     202                    nGood++;
     203                } else {
     204                    nBad++;
     205                    printf("%.3f   -->  %.3f  & %.3f\tdel=%e\n",xV->data.F32[i],yV->data.F32[i],y->data.F32[i],y->data.F32[i]-yV->data.F32[i]);
     206                }
     207            }
     208            psFree (x);
     209            psFree (y);
     210            psFree (yV);
     211            psFree (poly);
     212           
     213            if (nBad > 0) {
     214                psError(PS_ERR_UNKNOWN,true,"%d incorrect polynomial evaluations for NORD=%d\n", nBad, NORD);
     215                retval = retval + pow(2,NORD);
     216            }
     217        }
     218    }
     219    return retval;
     220}
     221int testPolynomial1DEvalVector_02(void) {
     222    //  check that psPolynomial1DEvalVector obeys the psPolynomial1D.mask
     223    int NORD;
     224    int maxNORD = 3;
     225    int retval = 0;
     226    psPolynomialType type[] = {PS_POLYNOMIAL_ORD, PS_POLYNOMIAL_CHEB};
     227   
     228    for (int k = 0; k<1; k++) {  //loop over the 2 polynomial types
     229        for (int imask=0; imask<=maxNORD; imask++) {  // mask one order at a time
     230            int j= maxNORD;
     231            NORD = j;
     232            int nGood = 0, nBad = 0;
     233            psPolynomial1D *poly       = makePolynomial1D(NORD, type[k], -1);
     234            psPolynomial1D *polynomask = makePolynomial1D(NORD, type[k], imask);
     235            poly->mask[imask]=1;
     236
     237            // create x vector, known y vector:
     238            psVector *x = psVectorAlloc (NPTS, PS_TYPE_F32);
     239            psVector *y = psVectorAlloc (NPTS, PS_TYPE_F32);
     240            fillVectors (x, y, NORD, NPTS, type[k]);
     241           
     242            psVector *xV = x;
     243            psVector *yV = psPolynomial1DEvalVector(poly, xV);
     244            psVector *yVnm = psPolynomial1DEvalVector(polynomask, xV);
     245           
     246            for (int i = 0; i < NPTS; i++) {
     247                if (yV->data.F32[i] == yVnm->data.F32[i]) {
     248                    nGood++;
     249                } else {
     250                    nBad++;
     251                    printf("x:%.3f  y:%.3f , %.3f , %.3f\tdel=%e\n",xV->data.F32[i],y->data.F32[i],yV->data.F32[i],yVnm->data.F32[i],yVnm->data.F32[i]-yV->data.F32[i]);
     252                }
     253            }
     254            psFree (x);
     255            psFree (y);
     256            psFree (yV);
     257            psFree (yVnm);
     258            psFree (poly);
     259            psFree (polynomask);
     260           
     261            if (nBad > 0) {
     262                psError(PS_ERR_UNKNOWN,true,"%d incorrect polynomial evaluations for NORD=%d\n", nBad, NORD);
     263                retval = retval + pow(2,NORD);
     264            }
     265        }
     266    }
     267    return retval;
     268}
     269
     270
     271
     272int testVectorFitPolynomial1D_01(void) {
     273    int retval=0;
     274    int NORD, maxNORD=3;
     275    float tolerance = 1.0E-5;
     276    psPolynomialType type[] = {PS_POLYNOMIAL_ORD, PS_POLYNOMIAL_CHEB};
     277   
     278    for (int k = 0; k<1; k++) {  //loop over the 2 polynomial types   XXX: not using Chebyshevs right now!!!
     279        for (int j = 1; j <= maxNORD; j++) {  //loop over the polynomial orders
     280            int nGood = 0, nBad = 0;
     281            NORD = j;
     282            printf("Testing NORD=%d\n", NORD);
     283            psPolynomial1D *poly = makePolynomial1D(NORD, type[k], -1);
     284            psPolynomial1D *polycopy = makePolynomial1D(NORD, type[k], -1);
     285            psVector *f = psVectorAlloc(NPTS, PS_TYPE_F32);
     286            psVector *fErr = psVectorAlloc(NPTS, PS_TYPE_F32);
     287            psVector *x = psVectorAlloc(NPTS, PS_TYPE_F32);
     288           
     289            fillVectors(x, f, NORD, NPTS, type[k]);
     290            psTime *now = psTimeGetNow(PS_TIME_UTC);
     291            const psRandom *r = psRandomAlloc(PS_RANDOM_TAUS, now->sec);
     292            for (int i = 0; i < NPTS; i++) {
     293                double errgen = psRandomGaussian(r);
     294                double amp = 1.0E-3;
     295                fErr->data.F32[i] = amp*errgen;
     296                //fErr->data.F32[i] = 0.;
     297            }
     298            psFree(now);
     299            psFree(r);
     300
     301            //f->data.F32[NPTS-1] = X0 + 10.*X1*x->data.F32[NPTS-1];
     302            //fErr->data.F32[NPTS-1] = 0.0000001;
     303
     304
     305            psVectorFitPolynomial1D(poly, NULL, 0, f, fErr, x);
     306            for (int i = 0; i <= NORD; i++) {
     307                if (poly->coeff[i] <= polycopy->coeff[i]+tolerance
     308                      &&
     309                    poly->coeff[i] >= polycopy->coeff[i]-tolerance) {
     310                    nGood++;
     311                } else {
     312                    nBad++;
     313                    printf("polycopy->coeff[%d]=%.20f    poly->coeff[%d]=%.20f\n",i,polycopy->coeff[i],i,poly->coeff[i]);
     314                }
     315                if (polycopy->coeffErr[i] != poly->coeffErr[i]) {
     316                    printf("polycopy->coeffErr[%d]=%.20f    poly->coeffErr[%d]=%.20f\n",i,polycopy->coeffErr[i],i,poly->coeffErr[i]);
     317                }               
     318            }
     319
     320            if (nBad > 0) {
     321                psError(PS_ERR_UNKNOWN,true,"%d incorrect polynomial evaluations for NORD=%d\n", nBad, NORD);
     322                retval = retval + pow(2,NORD);
     323            }
     324
     325            printf("NORD: %d   retval: %d\n",NORD,retval);
     326            psFree(poly);
     327            psFree(polycopy);
     328            psFree(f);
     329            psFree(fErr);
     330            psFree(x);     
     331        }
     332    }
     333    return retval;
     334}
     335
     336
     337int testVectorFitPolynomial1D_02(void) {
     338    // this time test that the masking in psPolynomial1D works correctly
     339    int retval=0;
     340    int NORD, maxNORD=3;
     341    float tolerance = 1.0E-5;
     342    psPolynomialType type[] = {PS_POLYNOMIAL_ORD, PS_POLYNOMIAL_CHEB};
     343   
     344    for (int k = 0; k<1; k++) {  //loop over the 2 polynomial types   XXX: not using Chebyshevs right now!!!
     345        for (int imask=0; imask<=maxNORD; imask++) {  // mask one order at a time
     346            int j= maxNORD;
     347            NORD = j;
     348            int nGood = 0, nBad = 0;
     349            psPolynomial1D *poly       = makePolynomial1D(NORD, type[k], -1);
     350            psPolynomial1D *polynomask = makePolynomial1D(NORD, type[k], imask);
     351            poly->mask[imask]=1;
     352            //poly and polynomask look the same now.
     353
     354            psVector *mask = psVectorAlloc(NPTS, PS_TYPE_U8);
     355            psVector *f = psVectorAlloc(NPTS, PS_TYPE_F32);
     356            psVector *fErr = psVectorAlloc(NPTS, PS_TYPE_F32);
     357            psVector *x = psVectorAlloc(NPTS, PS_TYPE_F32);
     358           
     359            fillVectors(x, f, NORD, NPTS, type[k]);
     360           
     361            // add random errors to the vector
     362            psTime *now = psTimeGetNow(PS_TIME_UTC);
     363            const psRandom *r = psRandomAlloc(PS_RANDOM_TAUS, now->sec);
     364            for (int i = 0; i < NPTS; i++) {
     365                double errgen = psRandomGaussian(r);
     366                double amp = 1.0E-3;
     367                fErr->data.F32[i] = amp*errgen;
     368                //fErr->data.F32[i] = 0.;
     369            }
     370            psFree(now);
     371            psFree(r);
     372           
     373           
     374            //mask->data.U8[imask] = 0xff;
     375            //for (int ii=0; ii<=mask->n; ii++) {
     376            //printf("  %d",mask->data.U8[ii]);
     377            //}
     378            //printf("\n");
     379
     380            poly->mask[imask] = 1;
     381            for (int ii=0; ii<=poly->nX; ii++) {
     382                printf("%x  ", poly->mask[ii]);
     383            }
     384            printf("\n");
     385
     386            psVectorFitPolynomial1D(poly, NULL, 0xff, f, fErr, x);
     387            psVectorFitPolynomial1D(polynomask, NULL, 0xff, f, fErr, x);
     388
     389            psVector *xx = psVectorAlloc(NPTS, PS_TYPE_F32);
     390            psVector *yy = psPolynomial1DEvalVector(poly, xx);
     391            psVector *del;
     392            psBinaryOp(del, f, "-", yy);
     393            psPolynomial1D *polydel = makePolynomial1D(NORD, type[k], -1);
     394            psVectorFitPolynomial1D(polydel, NULL,0xff, f, fErr, x);
     395           
     396
     397
     398            //mask->data.U8[imask] = 0;
     399            poly->mask[imask] = 0;
     400            for (int i = 0; i <= NORD; i++) {
     401                printf("polynomask->coeff[%d]=%.8f    poly->coeff[%d]=%.8f    polydel->coeff[%d]=%.8f\n",i,polynomask->coeff[i],i,poly->coeff[i],i,polydel->coeff[i]);
     402                if (poly->coeff[i] <= polynomask->coeff[i]+tolerance
     403                    &&
     404                    poly->coeff[i] >= polynomask->coeff[i]-tolerance) {
     405                    nGood++;
     406                } else {
     407                    nBad++;
     408                    //printf("polynomask->coeff[%d]=%.8f    poly->coeff[%d]=%.8f\n",i,polynomask->coeff[i],i,poly->coeff[i]);
     409                }
     410                if (polynomask->coeffErr[i] != poly->coeffErr[i]) {
     411                    //printf("polynomask->coeffErr[%d]=%.8f    poly->coeffErr[%d]=%.8f\n",i,polynomask->coeffErr[i],i,poly->coeffErr[i]);
     412                }               
     413            }
     414           
     415            if (nBad > 0) {
     416                psError(PS_ERR_UNKNOWN,true,"%d incorrect polynomial evaluations for NORD=%d\n", nBad, NORD);
     417                retval = retval + pow(2,NORD);
     418            }
     419           
     420            printf("Mask Ord: %d   retval: %d\n", imask, retval);
     421            psFree(poly);
     422            psFree(polynomask);
     423            psFree(mask);
     424            psFree(f);
     425            psFree(fErr);
     426            psFree(x);     
     427            psFree(xx);
     428            psFree(yy);
     429            psFree(del);
     430            psFree(polydel);
     431
     432        }
     433    }
     434    return retval;
     435}
     436
     437
     438
     439
     440
     441int testVectorClipFitPolynomial1D_01(void) {
     442    int retval=0;
     443    int NORD, maxNORD=3;
     444    double Nsigma = 3.;
     445    int Niter= 3.;
     446    float tolerance = 1.0E-5;
     447    psPolynomialType type[] = {PS_POLYNOMIAL_ORD, PS_POLYNOMIAL_CHEB};
     448   
     449    for (int k = 0; k<1; k++) {  //loop over the 2 polynomial types  XXX: psVectorClipFirPolynomial1D does not
     450                                  //                                  work with Chebyshev now
     451        for (int j = 1; j <= maxNORD; j++) {  //loop over the polynomial orders
     452            int nGood = 0, nBad = 0;
     453            NORD = j;
     454            printf("Testing NORD=%d\n", NORD);
     455            psPolynomial1D *poly = makePolynomial1D(NORD, type[k], -1);
     456            psPolynomial1D *polycopy = makePolynomial1D(NORD, type[k], -1);
     457            psStats *stats = psStatsAlloc(PS_STAT_CLIPPED_MEAN | PS_STAT_CLIPPED_STDEV);
     458            stats->clipSigma = Nsigma;
     459            stats->clipIter  = Niter;
     460            psVector *f = psVectorAlloc(NPTS, PS_TYPE_F32);
     461            psVector *fErr = psVectorAlloc(NPTS, PS_TYPE_F32);
     462            psVector *x = psVectorAlloc(NPTS, PS_TYPE_F32);
     463           
     464            fillVectors(x, f, NORD, NPTS, type[k]);
     465            psTime *now = psTimeGetNow(PS_TIME_UTC);
     466            const psRandom *r = psRandomAlloc(PS_RANDOM_TAUS, now->sec);
     467            for (int i = 0; i < NPTS; i++) {
     468                double errgen = psRandomGaussian(r);
     469                double amp = 1.0E-3;
     470                fErr->data.F32[i] = amp*errgen;
     471                //fErr->data.F32[i] = 0.;
     472            }
     473            psFree(now);
     474            psFree(r);
     475
     476            //f->data.F32[NPTS-1] = X0 + 10.*X1*x->data.F32[NPTS-1];
     477            //fErr->data.F32[NPTS-1] = 0.0000001;
     478
     479
     480            psVectorClipFitPolynomial1D(poly, stats, NULL, 2, f, fErr, x);
     481            for (int i = 0; i <= NORD; i++) {
     482                if (poly->coeff[i] <= polycopy->coeff[i]+tolerance
     483                      &&
     484                    poly->coeff[i] >= polycopy->coeff[i]-tolerance) {
     485                    nGood++;
     486                } else {
     487                    nBad++;
     488                    printf("polycopy->coeff[%d]=%.20f    poly->coeff[%d]=%.20f\n",i,polycopy->coeff[i],i,poly->coeff[i]);
     489                }
     490                if (polycopy->coeffErr[i] != poly->coeffErr[i]) {
     491                    printf("polycopy->coeffErr[%d]=%.20f    poly->coeffErr[%d]=%.20f\n",i,polycopy->coeffErr[i],i,poly->coeffErr[i]);
     492                }               
     493            }
     494
     495            if (nBad > 0) {
     496                psError(PS_ERR_UNKNOWN,true,"%d incorrect polynomial evaluations for NORD=%d\n", nBad, NORD);
     497                retval = retval + pow(2,NORD);
     498            }
     499
     500            printf("NORD: %d   retval: %d\n",NORD,retval);
     501            psFree(poly);
     502            psFree(polycopy);
     503            psFree(stats);
     504            psFree(f);
     505            psFree(fErr);
     506            psFree(x);     
     507        }
     508    }
     509    return retval;
     510}
     511
     512
     513
     514
     515
     516
     517
     518//test psStats
     519
     520//test psGaussian
     521
     522//test psMinimization (and psMinimizeLMChi2 and psMinimizeChi2Powell?)
     523
     524//test psRandomUniform / psRandomGausian / psRandomPoisson
     525
     526//test psSpline?
Note: See TracChangeset for help on using the changeset viewer.