Changeset 14726 for branches/eam_branch_20070830/psLib/test/imageops

Timestamp:

Sep 3, 2007, 10:29:07 AM (19 years ago)

Author:

magnier

Message:

updates to image map code, tests

File:

• branches/eam_branch_20070830/psLib/test/imageops/tap_psImageMap.c (modified) (1 diff)

branches/eam_branch_20070830/psLib/test/imageops/tap_psImageMap.c

@@ -92 +92 @@
 
     // make a model for a well-sampled field of a non-polynomial function (f = (a(x-xo)^2 + b(y-yo)^2)^-1)
+    {
+        // function is defined over the range 0-1000, 0-1000
+        psImageBinning *binning = psImageBinningAlloc();
+        binning->nXfine = 1000;
+        binning->nYfine = 1000;
+        binning->nXruff = 5;
+        binning->nYruff = 5;
+
+        // generate a grid of test data points
+        psVector *x = psVectorAllocEmpty (100, PS_TYPE_F32);
+        psVector *y = psVectorAllocEmpty (100, PS_TYPE_F32);
+        psVector *f = psVectorAllocEmpty (100, PS_TYPE_F32);
+
+        for (int ix = 0; ix < 1000; ix += 10) {
+            for (int iy = 0; iy < 1000; iy += 10) {
+                x->data.F32[x->n] = ix;
+                y->data.F32[y->n] = iy;
+                f->data.F32[f->n] = C00 / (C10*PS_SQR(ix-500) + C01*PS_SQR(iy-500) + 0.1);
+                psVectorExtend (x, 100, 1);
+                psVectorExtend (y, 100, 1);
+                psVectorExtend (f, 100, 1);
+            }
+        }
+
+        psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEDIAN);
+
+        // scale defines both field and map image sizes (nXfine, nXruff)
+        psImageMap *map = psImageMapAlloc (NULL, binning, stats);
+
+        // allocate a map, but we have no field image to supply
+        // XXX this function needs to correct for the mean superpixel position 
+        // which is sampled...
+        psImageMapGenerate (map, x, y, f, 0.1);
+        psFree (binning);
+
+        fprintf (stderr, "nGood: %d\n", map->nGood);
+        fprintf (stderr, "nPoor: %d\n", map->nPoor);
+        fprintf (stderr, "nBad:  %d\n", map->nBad);
+        
+        SaveImage (NULL, map->map, "map.fits");
+        
+        psImage *field = psImageAlloc(1000, 1000, PS_TYPE_F32);
+        for (int ix = 0; ix < 1000; ix++) {
+            for (int iy = 0; iy < 1000; iy++) {
+                field->data.F32[iy][ix] = C00 / (C10*PS_SQR(ix-500) + C01*PS_SQR(iy-500) + 0.1);
+            }
+        }
+        SaveImage (NULL, field, "field.fits");
+
+        // measure difference between model (map) and data (field)
+        for (int ix = 0; ix < map->map->numCols; ix++) {
+            for (int iy = 0; iy < map->map->numRows; iy++) {
+              
+                int xo = (ix + 0.5)*map->binning->nXbin;
+                int yo = (iy + 0.5)*map->binning->nYbin;
+                float df = field->data.F32[yo][xo] - map->map->data.F32[iy][ix];
+                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]);
+            }
+        }
+    }
 
     // make a model for a poorly-sampled field of a non-polynomial function (f = (a(x-xo)^2 + b(y-yo)^2)^-1)