IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 14863


Ignore:
Timestamp:
Sep 16, 2007, 3:15:55 PM (19 years ago)
Author:
magnier
Message:

fixed up tests for psImageMapFit

Location:
branches/eam_branch_20070830/psLib/test/imageops
Files:
5 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branch_20070830/psLib/test/imageops

    • Property svn:ignore

      •  

        old new  
        3434tap_psImageStructManip
        3535tap_psImageMap
         36tap_psImageMapFit
         37tap_psImageMapFit2

  • branches/eam_branch_20070830/psLib/test/imageops/.cvsignore

    r14777 r14863  
    3434tap_psImageStructManip
    3535tap_psImageMap
     36
     37tap_psImageMapFit
     38tap_psImageMapFit2

  • branches/eam_branch_20070830/psLib/test/imageops/Makefile.am

    r14859 r14863  
    2424        tap_psImageMap \
    2525        tap_psImageMapFit \
     26        tap_psImageMapFit2 \
    2627        tap_psImageMaskOps
    2728

  • branches/eam_branch_20070830/psLib/test/imageops/tap_psImageMapFit.c

    r14860 r14863  
    66#include "pstap.h"
    77
     8// save function used to dump out test images while debugging algorithm
    89int SaveImage (psMetadata *header, psImage *image, char *filename) {
    910
    … …  
    1415}
    1516
    16 # define TEST_4PT_0 0
    17 # define TEST_4PT_1 0
    18 # define TEST_4PT_2 0
    19 # define TEST_4PT_3 0
    20 
    21 # define TEST_6PT_0 0
    22 
    23 # define TEST_9PT_0 0
    24 # define TEST_9PT_1 0
     17# define TEST_4PT_0 1
     18# define TEST_4PT_1 1
     19# define TEST_4PT_2 1
     20# define TEST_4PT_3 1
     21# define TEST_6PT_0 1
     22# define TEST_9PT_0 1
     23# define TEST_9PT_1 1
    2524# define TEST_9PT_2 1
    26 
    27 # define TEST1 0
    28 # define TEST2 0
    29 # define TEST3 0
    30 
    31 # define C00 +0.0
    32 # define C01 +0.0
    33 # define C10 +1.0
     25# define TEST_9PT_3 1
     26# define TEST_9PT_4 1
    3427
    3528int main (void)
    … …  
    4134
    4235    // test for more points: 3x3 grid fitted to 9 points with simple slope and scale difference
     36    # if (TEST_9PT_4)
     37    {
     38        psMemId id = psMemGetId();
     39
     40        psImageBinning *binning = psImageBinningAlloc();
     41        binning->nXfine = 6;
     42        binning->nYfine = 6;
     43        binning->nXruff = 3;
     44        binning->nYruff = 3;
     45
     46        // generate a grid of test data points
     47        psVector *x = psVectorAlloc (50, PS_TYPE_F32);
     48        psVector *y = psVectorAlloc (50, PS_TYPE_F32);
     49        psVector *f = psVectorAlloc (50, PS_TYPE_F32);
     50
     51        // the underlying field is f = ix + iy, where ix,iy are fine pixel coordinates
     52        // place the measurement points exactly on the ruff reference pixel centers
     53        int n = 0;
     54        for (float ix = 0.5; ix < 6.0; ix += 1.0) {
     55            for (float iy = 0.5; iy < 6.0; iy += 1.0) {
     56                x->data.F32[n] = ix;
     57                y->data.F32[n] = iy;
     58                f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
     59                n++;
     60            }
     61        }
     62        x->n = n;
     63        y->n = n;
     64        f->n = n;
     65
     66        psImage *field = psImageAlloc(6, 6, PS_TYPE_F32);
     67        for (int ix = 0; ix < 6; ix++) {
     68            for (int iy = 0; iy < 6; iy++) {
     69                field->data.F32[iy][ix] = (ix + 0.5) + (iy + 0.5);
     70            }
     71        }
     72
     73        psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
     74
     75        // scale defines both field and map image sizes (nXfine, nXruff)
     76        psImageMap *map = psImageMapAlloc (NULL, binning, stats);
     77
     78        // fit the data to the map
     79        psImageMapFit (map, x, y, f, NULL);
     80
     81        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
     82        for (int ix = 0; ix < model->numCols; ix++) {
     83            for (int iy = 0; iy < model->numRows; iy++) {
     84                model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
     85                is_float_tol (model->data.F32[iy][ix], field->data.F32[iy][ix], 1e-5, "model matches inputs");
     86            }
     87        }
     88
     89        // SaveImage (NULL, map->map, "map.fits");
     90        // SaveImage (NULL, field, "field.fits");
     91        // SaveImage (NULL, model, "model.fits");
     92
     93        psFree (model);
     94        psFree (binning);
     95        psFree (map);
     96        psFree (stats);
     97        psFree (field);
     98        psFree (x);
     99        psFree (y);
     100        psFree (f);
     101        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     102    }
     103    # endif
     104
     105    // test for more points: 3x3 grid fitted to 9 points with simple slope and scale difference
     106    # if (TEST_9PT_3)
     107    {
     108        psMemId id = psMemGetId();
     109
     110        psImageBinning *binning = psImageBinningAlloc();
     111        binning->nXfine = 6;
     112        binning->nYfine = 6;
     113        binning->nXruff = 3;
     114        binning->nYruff = 3;
     115
     116        // generate a grid of test data points
     117        psVector *x = psVectorAlloc (50, PS_TYPE_F32);
     118        psVector *y = psVectorAlloc (50, PS_TYPE_F32);
     119        psVector *f = psVectorAlloc (50, PS_TYPE_F32);
     120
     121        // the underlying field is f = ix + iy, where ix,iy are fine pixel coordinates
     122        // place the measurement points exactly on the ruff reference pixel centers
     123        int n = 0;
     124        for (float ix = 1.0; ix < 6; ix += 2.0) {
     125            for (float iy = 1.0; iy < 6; iy += 2.0) {
     126                y->data.F32[n] = iy;
     127                x->data.F32[n] = ix;
     128                if ((ix == 1.0) && (iy == 1.0)) {
     129                    x->data.F32[n] = ix - 0.1;
     130                }
     131                if ((ix == 3.0) && (iy == 1.0)) {
     132                    y->data.F32[n] = iy - 0.1;
     133                }
     134                if ((ix == 5.0) && (iy == 3.0)) {
     135                    x->data.F32[n] = ix + 0.1;
     136                }
     137                if ((ix == 3.0) && (iy == 5.0)) {
     138                    y->data.F32[n] = iy + 0.1;
     139                }
     140                f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
     141                n++;
     142            }
     143        }
     144        x->n = n;
     145        y->n = n;
     146        f->n = n;
     147
     148        psImage *field = psImageAlloc(6, 6, PS_TYPE_F32);
     149        for (int ix = 0; ix < 6; ix++) {
     150            for (int iy = 0; iy < 6; iy++) {
     151                field->data.F32[iy][ix] = (ix + 0.5) + (iy + 0.5);
     152            }
     153        }
     154
     155        psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
     156
     157        // scale defines both field and map image sizes (nXfine, nXruff)
     158        psImageMap *map = psImageMapAlloc (NULL, binning, stats);
     159
     160        // fit the data to the map
     161        psImageMapFit (map, x, y, f, NULL);
     162
     163        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
     164        for (int ix = 0; ix < model->numCols; ix++) {
     165            for (int iy = 0; iy < model->numRows; iy++) {
     166                model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
     167                is_float_tol (model->data.F32[iy][ix], field->data.F32[iy][ix], 1e-5, "model matches inputs");
     168            }
     169        }
     170
     171        // SaveImage (NULL, map->map, "map.fits");
     172        // SaveImage (NULL, field, "field.fits");
     173        // SaveImage (NULL, model, "model.fits");
     174       
     175        psFree (model);
     176        psFree (binning);
     177        psFree (map);
     178        psFree (stats);
     179        psFree (field);
     180        psFree (x);
     181        psFree (y);
     182        psFree (f);
     183        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     184    }
     185    # endif
     186
     187    // test for more points: 3x3 grid fitted to 9 points with simple slope and scale difference
     188# if (TEST_9PT_2)
     189    {
     190        psMemId id = psMemGetId();
     191
     192        psImageBinning *binning = psImageBinningAlloc();
     193        binning->nXfine = 3;
     194        binning->nYfine = 3;
     195        binning->nXruff = 3;
     196        binning->nYruff = 3;
     197
     198        // generate a grid of test data points
     199        psVector *x = psVectorAlloc (50, PS_TYPE_F32);
     200        psVector *y = psVectorAlloc (50, PS_TYPE_F32);
     201        psVector *f = psVectorAlloc (50, PS_TYPE_F32);
     202
     203        // the underlying field is f = ix + iy, where ix,iy are fine pixel coordinates
     204        // place the measurement points exactly on the ruff reference pixel centers
     205        int n = 0;
     206        for (float ix = 0.0; ix < 3; ix += 1.0) {
     207            for (float iy = 0.0; iy < 3; iy += 1.0) {
     208                x->data.F32[n] = ix + 0.5;
     209                y->data.F32[n] = iy + 0.5;
     210
     211                if (ix == 0.0) {
     212                    x->data.F32[n] -= 0.1;
     213                }
     214                if (iy == 0.0) {
     215                    y->data.F32[n] -= 0.1;
     216                }
     217                if (ix == 2.0) {
     218                    x->data.F32[n] += 0.1;
     219                }
     220                if (iy == 2.0) {
     221                    y->data.F32[n] += 0.1;
     222                }
     223                f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
     224                n++;
     225            }
     226        }
     227        x->n = n;
     228        y->n = n;
     229        f->n = n;
     230
     231        psImage *field = psImageAlloc(3, 3, PS_TYPE_F32);
     232        for (int ix = 0; ix < 3; ix++) {
     233            for (int iy = 0; iy < 3; iy++) {
     234                field->data.F32[iy][ix] = (ix + 0.5) + (iy + 0.5);
     235            }
     236        }
     237
     238        psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
     239
     240        // scale defines both field and map image sizes (nXfine, nXruff)
     241        psImageMap *map = psImageMapAlloc (NULL, binning, stats);
     242
     243        // fit the data to the map
     244        psImageMapFit (map, x, y, f, NULL);
     245
     246        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
     247        for (int ix = 0; ix < model->numCols; ix++) {
     248            for (int iy = 0; iy < model->numRows; iy++) {
     249                model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
     250                is_float_tol (model->data.F32[iy][ix], field->data.F32[iy][ix], 1e-5, "model matches inputs");
     251            }
     252        }
     253        // SaveImage (NULL, map->map, "map.fits");
     254        // SaveImage (NULL, field, "field.fits");
     255        // SaveImage (NULL, model, "model.fits");
     256
     257        psFree (model);
     258        psFree (binning);
     259        psFree (map);
     260        psFree (stats);
     261        psFree (field);
     262        psFree (x);
     263        psFree (y);
     264        psFree (f);
     265        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     266    }
     267# endif
     268
     269    // still a simple test: 3x3 grid fitted to 9 points with simple slope and scale difference
     270    # if (TEST_9PT_1)
     271    {
     272        psMemId id = psMemGetId();
     273
     274        psImageBinning *binning = psImageBinningAlloc();
     275        binning->nXfine = 6;
     276        binning->nYfine = 6;
     277        binning->nXruff = 3;
     278        binning->nYruff = 3;
     279
     280        // generate a grid of test data points
     281        psVector *x = psVectorAlloc (9, PS_TYPE_F32);
     282        psVector *y = psVectorAlloc (9, PS_TYPE_F32);
     283        psVector *f = psVectorAlloc (9, PS_TYPE_F32);
     284
     285        // the underlying field is f = ix + iy, where ix,iy are fine pixel coordinates
     286        // place the measurement points exactly on the ruff reference pixel centers
     287        int n = 0;
     288        for (int ix = 1; ix < 6; ix += 2) {
     289            for (int iy = 1; iy < 6; iy += 2) {
     290                x->data.F32[n] = ix;
     291                y->data.F32[n] = iy;
     292                f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
     293                n++;
     294            }
     295        }
     296
     297        psImage *field = psImageAlloc(6, 6, PS_TYPE_F32);
     298        for (int ix = 0; ix < 6; ix++) {
     299            for (int iy = 0; iy < 6; iy++) {
     300                field->data.F32[iy][ix] = (ix + 0.5) + (iy + 0.5);
     301            }
     302        }
     303
     304        psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
     305
     306        // scale defines both field and map image sizes (nXfine, nXruff)
     307        psImageMap *map = psImageMapAlloc (NULL, binning, stats);
     308
     309        // fit the data to the map
     310        psImageMapFit (map, x, y, f, NULL);
     311
     312        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
     313        for (int ix = 0; ix < model->numCols; ix++) {
     314            for (int iy = 0; iy < model->numRows; iy++) {
     315                model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
     316                is_float_tol (model->data.F32[iy][ix], field->data.F32[iy][ix], FLT_EPSILON, "model matches inputs");
     317            }
     318        }
     319
     320        // SaveImage (NULL, map->map, "map.fits");
     321        // SaveImage (NULL, field, "field.fits");
     322        // SaveImage (NULL, model, "model.fits");
     323
     324        psFree (model);
     325        psFree (binning);
     326        psFree (map);
     327        psFree (stats);
     328        psFree (field);
     329        psFree (x);
     330        psFree (y);
     331        psFree (f);
     332        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     333    }
     334    # endif
     335
     336    // very simple test: 3x3 grid fitted to 9 points with simple slope
     337    # if (TEST_9PT_0)
     338    {
     339        // function is defined over the range 0-1000, 0-1000
     340        psMemId id = psMemGetId();
     341
     342        psImageBinning *binning = psImageBinningAlloc();
     343        binning->nXfine = 3;
     344        binning->nYfine = 3;
     345        binning->nXruff = 3;
     346        binning->nYruff = 3;
     347
     348        // generate a grid of test data points
     349        psVector *x = psVectorAlloc (9, PS_TYPE_F32);
     350        psVector *y = psVectorAlloc (9, PS_TYPE_F32);
     351        psVector *f = psVectorAlloc (9, PS_TYPE_F32);
     352
     353        int n = 0;
     354        psImage *field = psImageAlloc(3, 3, PS_TYPE_F32);
     355        for (int ix = 0; ix < 3; ix++) {
     356            for (int iy = 0; iy < 3; iy++) {
     357                x->data.F32[n] = ix + 0.5;
     358                y->data.F32[n] = iy + 0.5;
     359                f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
     360                field->data.F32[iy][ix] = f->data.F32[n];
     361                n++;
     362            }
     363        }
     364
     365        psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
     366
     367        // scale defines both field and map image sizes (nXfine, nXruff)
     368        psImageMap *map = psImageMapAlloc (NULL, binning, stats);
     369
     370        // fit the data to the map
     371        psImageMapFit (map, x, y, f, NULL);
     372
     373        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
     374        for (int ix = 0; ix < model->numCols; ix++) {
     375            for (int iy = 0; iy < model->numRows; iy++) {
     376                model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
     377                is_float_tol (model->data.F32[iy][ix], field->data.F32[iy][ix], FLT_EPSILON, "model matches inputs");
     378            }
     379        }
     380
     381        // SaveImage (NULL, map->map, "map.fits");
     382        // SaveImage (NULL, field, "field.fits");
     383        // SaveImage (NULL, model, "model.fits");
     384
     385        psFree (model);
     386        psFree (binning);
     387        psFree (map);
     388        psFree (stats);
     389        psFree (field);
     390        psFree (x);
     391        psFree (y);
     392        psFree (f);
     393        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     394    }
     395    # endif
     396
     397    // test for more points: 3x3 grid fitted to 9 points with simple slope and scale difference
    43398    # if (TEST_6PT_0)
    44399    {
     400        psMemId id = psMemGetId();
     401
    45402        psImageBinning *binning = psImageBinningAlloc();
    46403        binning->nXfine = 3;
    … …  
    94451        // fit the data to the map
    95452        psImageMapFit (map, x, y, f, NULL);
    96         psFree (binning);
    97 
    98         SaveImage (NULL, map->map, "map.fits");
    99        
    100         SaveImage (NULL, field, "field.fits");
    101453
    102454        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
    … …  
    107459            }
    108460        }
    109         SaveImage (NULL, model, "model.fits");
    110         // ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    111     }
    112     # endif
    113 
    114     // test for more points: 3x3 grid fitted to 9 points with simple slope and scale difference
    115     # if (TEST_9PT_2)
    116     {
    117         psImageBinning *binning = psImageBinningAlloc();
    118         binning->nXfine = 6;
    119         binning->nYfine = 6;
    120         binning->nXruff = 3;
    121         binning->nYruff = 3;
    122 
    123         // generate a grid of test data points
    124         psVector *x = psVectorAlloc (50, PS_TYPE_F32);
    125         psVector *y = psVectorAlloc (50, PS_TYPE_F32);
    126         psVector *f = psVectorAlloc (50, PS_TYPE_F32);
    127 
    128         // the underlying field is f = ix + iy, where ix,iy are fine pixel coordinates
    129         // place the measurement points exactly on the ruff reference pixel centers
    130         int n = 0;
    131         for (float ix = 1.0; ix < 6; ix += 2.0) {
    132             for (float iy = 1.0; iy < 6; iy += 2.0) {
    133                 if ((ix == 3.0) && (iy == 1.0)) {
    134                     x->data.F32[n] = ix + 0.0;
    135                     y->data.F32[n] = iy + 0.1; // add in both points. 
    136                     f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
    137                     n++;
    138                     x->data.F32[n] = ix + 0.0;
    139                     y->data.F32[n] = iy - 0.1;
    140                 } else {
    141                     x->data.F32[n] = ix;
    142                     y->data.F32[n] = iy;
    143                 }
    144                 f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
    145                 n++;
    146             }
    147         }
    148         x->n = n;
    149         y->n = n;
    150         f->n = n;
    151 
    152         psImage *field = psImageAlloc(6, 6, PS_TYPE_F32);
    153         for (int ix = 0; ix < 6; ix++) {
    154             for (int iy = 0; iy < 6; iy++) {
    155                 field->data.F32[iy][ix] = (ix + 0.5) + (iy + 0.5);
    156             }
    157         }
    158 
    159         psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
    160 
    161         // scale defines both field and map image sizes (nXfine, nXruff)
    162         psImageMap *map = psImageMapAlloc (NULL, binning, stats);
    163 
    164         // fit the data to the map
    165         psImageMapFit (map, x, y, f, NULL);
    166         psFree (binning);
    167 
    168         SaveImage (NULL, map->map, "map.fits");
    169        
    170         SaveImage (NULL, field, "field.fits");
    171 
    172         psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
    173         for (int ix = 0; ix < model->numCols; ix++) {
    174             for (int iy = 0; iy < model->numRows; iy++) {
    175                 model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
    176                 is_float_tol (model->data.F32[iy][ix], field->data.F32[iy][ix], FLT_EPSILON, "model matches inputs");
    177             }
    178         }
    179         SaveImage (NULL, model, "model.fits");
    180         // ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    181     }
    182     # endif
    183 
    184     // still a simple test: 3x3 grid fitted to 9 points with simple slope and scale difference
    185     # if (TEST_9PT_1)
    186     {
    187         psImageBinning *binning = psImageBinningAlloc();
    188         binning->nXfine = 6;
    189         binning->nYfine = 6;
    190         binning->nXruff = 3;
    191         binning->nYruff = 3;
    192 
    193         // generate a grid of test data points
    194         psVector *x = psVectorAlloc (9, PS_TYPE_F32);
    195         psVector *y = psVectorAlloc (9, PS_TYPE_F32);
    196         psVector *f = psVectorAlloc (9, PS_TYPE_F32);
    197 
    198         // the underlying field is f = ix + iy, where ix,iy are fine pixel coordinates
    199         // place the measurement points exactly on the ruff reference pixel centers
    200         int n = 0;
    201         for (int ix = 1; ix < 6; ix += 2) {
    202             for (int iy = 1; iy < 6; iy += 2) {
    203                 x->data.F32[n] = ix;
    204                 y->data.F32[n] = iy;
    205                 f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
    206                 n++;
    207             }
    208         }
    209 
    210         psImage *field = psImageAlloc(6, 6, PS_TYPE_F32);
    211         for (int ix = 0; ix < 6; ix++) {
    212             for (int iy = 0; iy < 6; iy++) {
    213                 field->data.F32[iy][ix] = (ix + 0.5) + (iy + 0.5);
    214             }
    215         }
    216 
    217         psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
    218 
    219         // scale defines both field and map image sizes (nXfine, nXruff)
    220         psImageMap *map = psImageMapAlloc (NULL, binning, stats);
    221 
    222         // fit the data to the map
    223         psImageMapFit (map, x, y, f, NULL);
    224         psFree (binning);
    225 
    226         SaveImage (NULL, map->map, "map.fits");
    227        
    228         SaveImage (NULL, field, "field.fits");
    229 
    230         psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
    231         for (int ix = 0; ix < model->numCols; ix++) {
    232             for (int iy = 0; iy < model->numRows; iy++) {
    233                 model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
    234                 is_float_tol (model->data.F32[iy][ix], field->data.F32[iy][ix], FLT_EPSILON, "model matches inputs");
    235             }
    236         }
    237         SaveImage (NULL, model, "model.fits");
    238         // ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    239     }
    240     # endif
    241 
    242     // very simple test: 3x3 grid fitted to 9 points with simple slope
    243     # if (TEST_9PT_0)
    244     {
    245         // function is defined over the range 0-1000, 0-1000
    246         psImageBinning *binning = psImageBinningAlloc();
    247         binning->nXfine = 3;
    248         binning->nYfine = 3;
    249         binning->nXruff = 3;
    250         binning->nYruff = 3;
    251 
    252         // generate a grid of test data points
    253         psVector *x = psVectorAlloc (9, PS_TYPE_F32);
    254         psVector *y = psVectorAlloc (9, PS_TYPE_F32);
    255         psVector *f = psVectorAlloc (9, PS_TYPE_F32);
    256 
    257         int n = 0;
    258         psImage *field = psImageAlloc(3, 3, PS_TYPE_F32);
    259         for (int ix = 0; ix < 3; ix++) {
    260             for (int iy = 0; iy < 3; iy++) {
    261                 x->data.F32[n] = ix + 0.5;
    262                 y->data.F32[n] = iy + 0.5;
    263                 f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
    264                 field->data.F32[iy][ix] = f->data.F32[n];
    265                 n++;
    266             }
    267         }
    268 
    269         psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
    270 
    271         // scale defines both field and map image sizes (nXfine, nXruff)
    272         psImageMap *map = psImageMapAlloc (NULL, binning, stats);
    273 
    274         // fit the data to the map
    275         psImageMapFit (map, x, y, f, NULL);
    276         psFree (binning);
    277 
    278         SaveImage (NULL, map->map, "map.fits");
    279        
    280         SaveImage (NULL, field, "field.fits");
    281 
    282         psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
    283         for (int ix = 0; ix < model->numCols; ix++) {
    284             for (int iy = 0; iy < model->numRows; iy++) {
    285                 model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
    286                 is_float_tol (model->data.F32[iy][ix], field->data.F32[iy][ix], FLT_EPSILON, "model matches inputs");
    287             }
    288         }
    289         SaveImage (NULL, model, "model.fits");
    290         // ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     461
     462        // SaveImage (NULL, map->map, "map.fits");
     463        // SaveImage (NULL, field, "field.fits");
     464        // SaveImage (NULL, model, "model.fits");
     465
     466        psFree (model);
     467        psFree (binning);
     468        psFree (map);
     469        psFree (stats);
     470        psFree (field);
     471        psFree (x);
     472        psFree (y);
     473        psFree (f);
     474        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    291475    }
    292476    # endif
    … …  
    298482    {
    299483        // function is defined over the range 0-1000, 0-1000
     484        psMemId id = psMemGetId();
     485
    300486        psImageBinning *binning = psImageBinningAlloc();
    301487        binning->nXfine = 4;
    … …  
    334520        // fit the data to the map
    335521        psImageMapFit (map, x, y, f, NULL);
    336         psFree (binning);
    337 
    338         SaveImage (NULL, map->map, "map.fits");
    339        
    340         SaveImage (NULL, field, "field.fits");
    341522
    342523        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
    … …  
    347528            }
    348529        }
    349         SaveImage (NULL, model, "model.fits");
    350         // ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     530
     531        // SaveImage (NULL, map->map, "map.fits");
     532        // SaveImage (NULL, field, "field.fits");
     533        // SaveImage (NULL, model, "model.fits");
     534
     535        psFree (model);
     536        psFree (binning);
     537        psFree (map);
     538        psFree (stats);
     539        psFree (field);
     540        psFree (x);
     541        psFree (y);
     542        psFree (f);
     543        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    351544    }
    352545    # endif
    … …  
    358551    {
    359552        // function is defined over the range 0-1000, 0-1000
     553        psMemId id = psMemGetId();
     554
    360555        psImageBinning *binning = psImageBinningAlloc();
    361556        binning->nXfine = 2;
    … …  
    365560
    366561        // generate a grid of test data points
    367         psVector *x = psVectorAlloc (4, PS_TYPE_F32);
    368         psVector *y = psVectorAlloc (4, PS_TYPE_F32);
    369         psVector *f = psVectorAlloc (4, PS_TYPE_F32);
     562        psVector *x = psVectorAlloc (5, PS_TYPE_F32);
     563        psVector *y = psVectorAlloc (5, PS_TYPE_F32);
     564        psVector *f = psVectorAlloc (5, PS_TYPE_F32);
    370565
    371566        int n = 0;
    … …  
    374569        for (int ix = 0; ix < 2; ix++) {
    375570            for (int iy = 0; iy < 2; iy++) {
    376                 // add in one extra point on the surface
     571                # if (1)
    377572                if (!ix && !iy) {
    378                     x->data.F32[n] = ix + 0.5;
    379                     y->data.F32[n] = iy + 0.5 + 0.1;
     573                    x->data.F32[n] = ix + 0.5 - 0.1;
     574                    y->data.F32[n] = iy + 0.5 - 0.0;
     575//                  f->data.F32[n] = x->data.F32[n] + y->data.F32[n];
     576//                  n++;
     577//                  x->data.F32[n] = ix + 0.5;
     578//                  y->data.F32[n] = iy + 0.5 - 0.1;
    380579                } else {
    381580                    x->data.F32[n] = ix + 0.5;
    382581                    y->data.F32[n] = iy + 0.5;
    383582                }
     583                # endif
    384584                # if (0)
     585                // offset points from centers
    385586                if (ix) {
    386587                    x->data.F32[n] = ix + 0.5 - 0.1;
    … …  
    399600            }
    400601        }
    401 
    402         psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
    403 
    404         // scale defines both field and map image sizes (nXfine, nXruff)
    405         psImageMap *map = psImageMapAlloc (NULL, binning, stats);
    406 
    407         // fit the data to the map
    408         psImageMapFit (map, x, y, f, NULL);
    409         psFree (binning);
    410 
    411         SaveImage (NULL, map->map, "map.fits");
    412        
    413         SaveImage (NULL, field, "field.fits");
     602        x->n = n;
     603        y->n = n;
     604        f->n = n;
     605
     606        psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
     607
     608        // scale defines both field and map image sizes (nXfine, nXruff)
     609        psImageMap *map = psImageMapAlloc (NULL, binning, stats);
     610
     611        // fit the data to the map
     612        psImageMapFit (map, x, y, f, NULL);
    414613
    415614        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
    … …  
    420619            }
    421620        }
    422         SaveImage (NULL, model, "model.fits");
    423         // ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     621
     622        // SaveImage (NULL, map->map, "map.fits");
     623        // SaveImage (NULL, field, "field.fits");
     624        // SaveImage (NULL, model, "model.fits");
     625
     626        psFree (model);
     627        psFree (binning);
     628        psFree (map);
     629        psFree (stats);
     630        psFree (field);
     631        psFree (x);
     632        psFree (y);
     633        psFree (f);
     634        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    424635    }
    425636    # endif
    … …  
    429640    {
    430641        // function is defined over the range 0-1000, 0-1000
     642        psMemId id = psMemGetId();
     643
    431644        psImageBinning *binning = psImageBinningAlloc();
    432645        binning->nXfine = 4;
    … …  
    465678        // fit the data to the map
    466679        psImageMapFit (map, x, y, f, NULL);
    467         psFree (binning);
    468 
    469         SaveImage (NULL, map->map, "map.fits");
    470        
    471         SaveImage (NULL, field, "field.fits");
    472680
    473681        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
    … …  
    478686            }
    479687        }
    480         SaveImage (NULL, model, "model.fits");
    481         // ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     688
     689        // SaveImage (NULL, map->map, "map.fits");
     690        // SaveImage (NULL, field, "field.fits");
     691        // SaveImage (NULL, model, "model.fits");
     692
     693        psFree (model);
     694        psFree (binning);
     695        psFree (map);
     696        psFree (stats);
     697        psFree (field);
     698        psFree (x);
     699        psFree (y);
     700        psFree (f);
     701        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    482702    }
    483703    # endif
    … …  
    487707    {
    488708        // function is defined over the range 0-1000, 0-1000
     709        psMemId id = psMemGetId();
     710
    489711        psImageBinning *binning = psImageBinningAlloc();
    490712        binning->nXfine = 2;
    … …  
    517739        // fit the data to the map
    518740        psImageMapFit (map, x, y, f, NULL);
    519         psFree (binning);
    520 
    521         SaveImage (NULL, map->map, "map.fits");
    522        
    523         SaveImage (NULL, field, "field.fits");
    524741
    525742        psImage *model = psImageAlloc(field->numCols, field->numRows, PS_TYPE_F32);
    … …  
    530747            }
    531748        }
    532         SaveImage (NULL, model, "model.fits");
    533         // ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     749
     750        // SaveImage (NULL, map->map, "map.fits");
     751        // SaveImage (NULL, field, "field.fits");
     752        // SaveImage (NULL, model, "model.fits");
     753
     754        psFree (model);
     755        psFree (binning);
     756        psFree (map);
     757        psFree (stats);
     758        psFree (field);
     759        psFree (x);
     760        psFree (y);
     761        psFree (f);
     762        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    534763    }
    535764    # endif

  • branches/eam_branch_20070830/psLib/test/imageops/tap_psImageMapFit2.c

    r14859 r14863  
     1#include <stdio.h>
     2#include <string.h>
     3#include <pslib.h>
     4
     5#include "tap.h"
     6#include "pstap.h"
     7
     8int SaveImage (psMetadata *header, psImage *image, char *filename) {
     9
     10    psFits *fits = psFitsOpen (filename, "w");
     11    psFitsWriteImage (fits, NULL, image, 0, NULL);
     12    psFitsClose (fits);
     13    return (TRUE);
     14}
     15
     16# define TEST1 0
     17# define TEST2 1
     18# define TEST3 0
     19
     20# define C00 +0.0
     21# define C01 -2.0
     22# define C10 +1.0
     23
     24int main (void)
     25{
     26
     27    plan_tests(1);
    128
    229    // make a model for a well-sampled field of a simple function (f = ax + by + c)
    330    # if (TEST1)
    431    {
     32        psMemId id = psMemGetId();
     33
    534        // function is defined over the range 0-1000, 0-1000
    635        psImageBinning *binning = psImageBinningAlloc();
    … …  
    3362        // fit the data to the map
    3463        psImageMapFit (map, x, y, f, NULL);
    35         psFree (binning);
    36 
    37         // fprintf (stderr, "nGood: %d\n", map->nGood);
    38         // fprintf (stderr, "nPoor: %d\n", map->nPoor);
    39         // fprintf (stderr, "nBad:  %d\n", map->nBad);
    40        
    41         SaveImage (NULL, map->map, "map.fits");
    42        
     64
    4365        psImage *field = psImageAlloc(1000, 1000, PS_TYPE_F32);
    4466        for (int ix = 0; ix < 1000; ix++) {
    … …  
    4769            }
    4870        }
    49         SaveImage (NULL, field, "field.fits");
    5071
    5172        // measure difference between model (map) and data (field)
    … …  
    5475                int xo = psImageBinningGetFineX(map->binning, ix + 0.5);
    5576                int yo = psImageBinningGetFineY(map->binning, iy + 0.5);
    56                 float df = field->data.F32[yo][xo] - map->map->data.F32[iy][ix];
    57                 fprintf (stderr, "%d %d -> %d %d  :  %f = %f - %f\n", ix, iy, xo, yo, df, field->data.F32[yo][xo], map->map->data.F32[iy][ix]);
    58             }
    59         }
    60 
     77                is_float_tol (field->data.F32[yo][xo], map->map->data.F32[iy][ix], 1e-3, "model matches data");
     78            }
     79        }
     80
     81        // XXX test on the model or don't bother
    6182        psImage *model = psImageAlloc(1000, 1000, PS_TYPE_F32);
    6283        for (int ix = 0; ix < 1000; ix++) {
    … …  
    6586            }
    6687        }
    67         SaveImage (NULL, model, "model.fits");
     88
     89        // SaveImage (NULL, model, "model.fits");
     90        // SaveImage (NULL, field, "field.fits");
     91        // SaveImage (NULL, map->map, "map.fits");
     92       
     93        psFree (model);
     94        psFree (binning);
     95        psFree (map);
     96        psFree (stats);
     97        psFree (field);
     98        psFree (x);
     99        psFree (y);
     100        psFree (f);
     101        ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    68102    }
    69103    # endif
    … …  
    90124        for (int ix = 0; ix < 1000; ix += 20) {
    91125            for (int iy = 0; iy < 1000; iy += 20) {
    92                 x->data.F32[x->n] = ix;
    93                 y->data.F32[y->n] = iy;
    94                 f->data.F32[f->n] = C10*PS_SQR(ix-500) + C01*PS_SQR(iy-500);
     126                x->data.F32[x->n] = ix + 0,5;
     127                y->data.F32[y->n] = iy + 0.5;
     128                f->data.F32[f->n] = PS_SQR(x->data.F32[x->n]) + PS_SQR(y->data.F32[x->n]);
    95129                psVectorExtend (x, 100, 1);
    96130                psVectorExtend (y, 100, 1);
    … …  
    107141        // XXX this function needs to correct for the mean superpixel position
    108142        // which is sampled...
    109         psImageMapGenerate (map, x, y, f, 0.1);
     143        psImageMapFit (map, x, y, f, NULL);
    110144        psFree (binning);
    111145
    112         fprintf (stderr, "nGood: %d\n", map->nGood);
    113         fprintf (stderr, "nPoor: %d\n", map->nPoor);
    114         fprintf (stderr, "nBad:  %d\n", map->nBad);
    115        
    116146        SaveImage (NULL, map->map, "map.fits");
    117147       
    … …  
    119149        for (int ix = 0; ix < 1000; ix++) {
    120150            for (int iy = 0; iy < 1000; iy++) {
    121                 field->data.F32[iy][ix] = C10*PS_SQR(ix-500) + C01*PS_SQR(iy-500);
     151                field->data.F32[iy][ix] = PS_SQR(ix+0.5) + PS_SQR(iy+0.5);
    122152            }
    123153        }
    … …  
    133163            }
    134164        }
     165
     166        // XXX test on the model or don't bother
     167        psImage *model = psImageAlloc(1000, 1000, PS_TYPE_F32);
     168        for (int ix = 0; ix < 1000; ix++) {
     169            for (int iy = 0; iy < 1000; iy++) {
     170                model->data.F32[iy][ix] = psImageUnbinPixel_V2 (ix + 0.5, iy + 0.5, map->map, map->binning);
     171            }
     172        }
     173        SaveImage (NULL, model, "model.fits");
    135174    }
    136175    # endif
    … …  
    215254    // choose a model scale for a poorly-sampled field of a non-polynomial function (f = (a(x-xo)^2 + b(y-yo)^2)^-1)
    216255   
     256    return exit_status();
     257}
Note: See TracChangeset for help on using the changeset viewer.