[Yt-svn] yt-commit r398 - in trunk: tests tests/regression_scripts yt/lagos

mturk at wrangler.dreamhost.com mturk at wrangler.dreamhost.com
Fri Apr 11 14:23:44 PDT 2008


Author: mturk
Date: Fri Apr 11 14:23:43 2008
New Revision: 398
URL: http://yt.spacepope.org/changeset/398

Log:
Fixed the regressed behavior in the Contour Finder.  The issue came about in
that I was directly accessing the ->data[] attributes of the PyArrayObjects; I
am not sure why this was not working, but it wasn't.  I switched to
PyArray_GETPTR1.  Now I get correct behavior for my regression tests, and all
the units tests passed.  Additionally, I cleaned up some 'type' stuff in the C
code and in the ContourFinder.py file.

Added more unit tests for datacubes and contour finding.  Also made it so that
you could flush back a data field from a datacube even if that field is new.



Modified:
   trunk/tests/regression_scripts/gal.py
   trunk/tests/test_hdf5_reader.py
   trunk/tests/test_lagos.py
   trunk/yt/lagos/BaseDataTypes.py
   trunk/yt/lagos/ContourFinder.py
   trunk/yt/lagos/PointCombine.c

Modified: trunk/tests/regression_scripts/gal.py
==============================================================================
--- trunk/tests/regression_scripts/gal.py	(original)
+++ trunk/tests/regression_scripts/gal.py	Fri Apr 11 14:23:43 2008
@@ -27,6 +27,20 @@
         Contour id 4.0 has: 4.56339e-14 7.18032e-14  (1669) (3 grids, 11039.0 11051.0)
 """
 
+"""With C extension and no ->data[] accesses
+yt.lagos   INFO       2008-04-11 14:19:24,079 Getting field Density from 197
+yt.lagos   INFO       2008-04-11 14:19:25,356 Contouring over 1061337 cells with 30115 candidates
+yt.lagos   INFO       2008-04-11 14:19:58,070 Getting field tempContours from 197
+yt.lagos   INFO       2008-04-11 14:19:58,455 Identified 5 contours between 4.56276e-14 and 2.15510e-13
+yt.lagos   INFO       2008-04-11 14:19:58,462 Getting field Contours from 197
+yt.lagos   INFO       2008-04-11 14:19:59,418 Getting field GridIndices from 197
+        Contour id 0.0 has: 4.56297e-14 9.88512e-14  (2351) (4 grids, 11081.0 11087.0)
+        Contour id 1.0 has: 4.58381e-14 5.41707e-14  (13) (1 grids, 11079.0 11079.0)
+        Contour id 2.0 has: 4.60523e-14 4.79728e-14  (19) (2 grids, 11079.0 11082.0)
+        Contour id 3.0 has: 4.56287e-14 2.15439e-13  (26063) (32 grids, 11040.0 11103.0)
+        Contour id 4.0 has: 4.56339e-14 7.18032e-14  (1669) (3 grids, 11039.0 11051.0)
+"""
+
 
 import sys
 sys.path.insert(0,"/Users/matthewturk/Development/yt/trunk")

Modified: trunk/tests/test_hdf5_reader.py
==============================================================================
--- trunk/tests/test_hdf5_reader.py	(original)
+++ trunk/tests/test_hdf5_reader.py	Fri Apr 11 14:23:43 2008
@@ -17,6 +17,7 @@
         my_table.close()
         recv_array = ReadData("testing_h5lt_io.h5", "/%s" % (self.dtype))
         self.assert_(numpy.all(recv_array == self.rand_array))
+        self.assert_(recv_array.shape == self.rand_array.shape)
     def tearDown(self):
         os.unlink("testing_h5lt_io.h5")
 

Modified: trunk/tests/test_lagos.py
==============================================================================
--- trunk/tests/test_lagos.py	(original)
+++ trunk/tests/test_lagos.py	Fri Apr 11 14:23:43 2008
@@ -228,6 +228,15 @@
         cg2 = self.hierarchy.covering_grid(3, [0.0]*3, [1.0]*3, [64,64,64])
         self.assertTrue(na.all(cg["Ones"] == cg2["Ones"]))
 
+    def testRawFlushBack(self):
+        cg = self.hierarchy.covering_grid(3, [0.0]*3, [1.0]*3, [64,64,64])
+        cg["DensityNew"] = cg["Density"] * 2.111
+        cg.flush_data(field="DensityNew")
+        for g in self.hierarchy.grids:
+            ni = g["DensityNew"] > 0
+            self.assertTrue(na.all(g["DensityNew"][ni]/2.111 
+                                == g["Density"][ni]))
+
     def testAllCover(self):
         cg = self.hierarchy.covering_grid(0, [0.0]*3, [1.0]*3, [32,32,32])
         self.assertTrue(cg["Density"].max() \

Modified: trunk/yt/lagos/BaseDataTypes.py
==============================================================================
--- trunk/yt/lagos/BaseDataTypes.py	(original)
+++ trunk/yt/lagos/BaseDataTypes.py	Fri Apr 11 14:23:43 2008
@@ -1504,6 +1504,8 @@
         g_dx = na.array([grid.dx, grid.dy, grid.dz])
         c_dx = na.array([self.dx, self.dy, self.dz])
         for field in ensure_list(fields):
+            if not grid.has_key(field): grid[field] = \
+               na.zeros(grid.ActiveDimensions, dtype=self[field].dtype)
             PointCombine.DataCubeReplace(
                 grid.LeftEdge, g_dx, grid[field], grid.child_mask,
                 self.left_edge, self.right_edge, c_dx, self[field],

Modified: trunk/yt/lagos/ContourFinder.py
==============================================================================
--- trunk/yt/lagos/ContourFinder.py	(original)
+++ trunk/yt/lagos/ContourFinder.py	Fri Apr 11 14:23:43 2008
@@ -102,7 +102,7 @@
             kk = na.arange(cur_max_id, cur_max_id-local_ind[0].size, -1)
             cg["tempContours"][local_ind] = kk[:]
             cur_max_id -= local_ind[0].size
-        fd = cg["tempContours"]
+        fd = cg["tempContours"].astype('int64')
         fd_original = fd.copy()
         xi_u,yi_u,zi_u = na.where(fd > -1)
         cor_order = na.argsort(-1*fd[(xi_u,yi_u,zi_u)])
@@ -110,7 +110,7 @@
         yi = yi_u[cor_order]
         zi = zi_u[cor_order]
         PointCombine.FindContours(fd, xi, yi, zi)
-        cg["tempContours"] = fd.copy()
+        cg["tempContours"] = fd.copy().astype('float64')
         cg.flush_data("tempContours")
         my_queue.add(cg._grids)
         force_ind = na.unique(cg["GridIndices"][na.where(

Modified: trunk/yt/lagos/PointCombine.c
==============================================================================
--- trunk/yt/lagos/PointCombine.c	(original)
+++ trunk/yt/lagos/PointCombine.c	Fri Apr 11 14:23:43 2008
@@ -585,12 +585,13 @@
 
 static PyObject *_findContoursError;
 
+npy_int64 process_neighbors(PyArrayObject*, npy_int32, npy_int32, npy_int32);
 static PyObject *
 Py_FindContours(PyObject *obj, PyObject *args)
 {
     PyObject *ocon_ids, *oxi, *oyi, *ozi;
     PyArrayObject *con_ids, *xi, *yi, *zi;
-    int i, j, k, n;
+    npy_int32 i, j, k, n;
 
     i = 0;
     if (!PyArg_ParseTuple(args, "OOOO",
@@ -599,7 +600,7 @@
                     "FindContours: Invalid parameters.");
     
     con_ids   = (PyArrayObject *) PyArray_FromAny(ocon_ids,
-                    PyArray_DescrFromType(NPY_FLOAT64), 3, 3,
+                    PyArray_DescrFromType(NPY_INT64), 3, 3,
                     NPY_INOUT_ARRAY | NPY_UPDATEIFCOPY, NULL);
     if((con_ids==NULL) || (con_ids->nd != 3)) {
     PyErr_Format(_findContoursError,
@@ -634,8 +635,10 @@
     goto _fail;
     }
     
-    for(n=0;n<PyArray_SIZE(xi);n++) {
-      i=xi->data[n];j=yi->data[n];k=zi->data[n];
+    for(n=0;n<xi->dimensions[0];n++) {
+      i=*(npy_int32 *)PyArray_GETPTR1(xi,n);
+      j=*(npy_int32 *)PyArray_GETPTR1(yi,n);
+      k=*(npy_int32 *)PyArray_GETPTR1(zi,n);
       process_neighbors(con_ids, i, j, k);
     }
 
@@ -655,38 +658,39 @@
         return NULL;
 }
 
-int process_neighbors(PyArrayObject *con_ids, int i, int j, int k)
+npy_int64 process_neighbors(PyArrayObject *con_ids, npy_int32 i, npy_int32 j, npy_int32 k)
 {
-  int off_i, off_j, off_k, spawn_check;
+  npy_int32 off_i, off_j, off_k;
+  int spawn_check;
   int mi, mj, mk;
-  npy_float64 val1, val2, new_cid;
+  npy_int64 *fd_off, *fd_ijk;
+  npy_int64 new_cid, original_cid;
   mi = con_ids->dimensions[0];
   mj = con_ids->dimensions[1];
   mk = con_ids->dimensions[2];
+  fd_ijk = ((npy_int64*)PyArray_GETPTR3(con_ids, i, j, k));
+  //if(*fd_ijk == -1){return ((npy_int64)*fd_ijk);}
   do {
     spawn_check = 0;
     for (off_i=max(i-1,0);off_i<=min(i+1,mi-1);off_i++)
       for (off_j=max(j-1,0);off_j<=min(j+1,mj-1);off_j++)
         for (off_k=max(k-1,0);off_k<=min(k+1,mk-1);off_k++) {
           if((off_i==i)&&(off_j==j)&&(off_k==k)) continue;
-          val1 = *((npy_float64*)PyArray_GETPTR3(con_ids, off_i, off_j, off_k));
-          val2 = *((npy_float64*)PyArray_GETPTR3(con_ids, i, j, k));
-          if(val1!=-1) {
-            if(val1 > val2){
-              val2 = *((npy_float64*)PyArray_GETPTR3(con_ids, i,j,k)) = val1;
-              spawn_check += 1;
-            }
-            if(val1 < val2){
-              *((npy_float64*)PyArray_GETPTR3(con_ids, off_i,off_j,off_k)) = val2;
-              new_cid = process_neighbors(con_ids, off_i, off_j, off_k);
-              if (new_cid != val2) spawn_check += 1;
-              *((npy_float64*)PyArray_GETPTR3(con_ids, i,j,k)) = new_cid;
-            }
+          fd_off = ((npy_int64*)PyArray_GETPTR3(con_ids, off_i, off_j, off_k));
+          if(*fd_off == -1) continue;
+          if(*fd_off > *fd_ijk){
+            *fd_ijk = *fd_off;
+            spawn_check += 1;
+          }
+          if(*fd_off < *fd_ijk){
+            *fd_off = *fd_ijk;
+            new_cid = process_neighbors(con_ids, off_i, off_j, off_k);
+            if (new_cid != *fd_ijk) spawn_check += 1;
+            *fd_ijk = new_cid;
           }
         }
   } while (spawn_check > 0);
-  val2 = *((npy_float64*)PyArray_GETPTR3(con_ids, i, j, k));
-  return val2;
+  return (npy_int64) *fd_ijk;
 }
 
 static PyObject *_interpolateError;



More information about the yt-svn mailing list