[yt-svn] commit/yt: MatthewTurk: Merged in drudd/yt/yt-3.0 (pull request #1015)
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Mon Jul 14 11:02:48 PDT 2014
1 new commit in yt:
https://bitbucket.org/yt_analysis/yt/commits/561d2e56371c/
Changeset: 561d2e56371c
Branch: yt-3.0
User: MatthewTurk
Date: 2014-07-14 20:02:36
Summary: Merged in drudd/yt/yt-3.0 (pull request #1015)
Point search refactoring (take 2)
Affected #: 17 files
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/analysis_modules/particle_trajectories/particle_trajectories.py
--- a/yt/analysis_modules/particle_trajectories/particle_trajectories.py
+++ b/yt/analysis_modules/particle_trajectories/particle_trajectories.py
@@ -225,7 +225,8 @@
x = self["particle_position_x"][:,step].ndarray_view()
y = self["particle_position_y"][:,step].ndarray_view()
z = self["particle_position_z"][:,step].ndarray_view()
- particle_grids, particle_grid_inds = ds.index.find_points(x,y,z)
+ # This will fail for non-grid index objects
+ particle_grids, particle_grid_inds = ds.index._find_points(x,y,z)
for grid in particle_grids:
cube = grid.retrieve_ghost_zones(1, [fd])
CICSample_3(x,y,z,pfield,
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/data_objects/data_containers.py
--- a/yt/data_objects/data_containers.py
+++ b/yt/data_objects/data_containers.py
@@ -725,6 +725,12 @@
self.index._identify_base_chunk(self)
return self._current_chunk.fwidth
+class YTSelectionContainer0D(YTSelectionContainer):
+ _spatial = False
+ def __init__(self, pf, field_parameters):
+ super(YTSelectionContainer0D, self).__init__(
+ pf, field_parameters)
+
class YTSelectionContainer1D(YTSelectionContainer):
_spatial = False
def __init__(self, pf, field_parameters):
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/data_objects/selection_data_containers.py
--- a/yt/data_objects/selection_data_containers.py
+++ b/yt/data_objects/selection_data_containers.py
@@ -21,7 +21,8 @@
from yt.utilities.lib.alt_ray_tracers import cylindrical_ray_trace
from yt.utilities.orientation import Orientation
from .data_containers import \
- YTSelectionContainer1D, YTSelectionContainer2D, YTSelectionContainer3D
+ YTSelectionContainer0D, YTSelectionContainer1D, \
+ YTSelectionContainer2D, YTSelectionContainer3D
from yt.data_objects.derived_quantities import \
DerivedQuantityCollection
from yt.utilities.exceptions import \
@@ -34,6 +35,31 @@
from yt.utilities.math_utils import get_rotation_matrix
from yt.units.yt_array import YTQuantity
+
+class YTPointBase(YTSelectionContainer0D):
+ """
+ A 0-dimensional object defined by a single point
+
+ Parameters
+ ----------
+ p: array_like
+ A points defined within the domain. If the domain is
+ periodic its position will be corrected to lie inside
+ the range [DLE,DRE) to ensure one and only one cell may
+ match that point
+
+ Examples
+ --------
+ >>> pf = load("DD0010/moving7_0010")
+ >>> c = [0.5,0.5,0.5]
+ >>> point = pf.point(c)
+ """
+ _type_name = "point"
+ _con_args = ('p',)
+ def __init__(self, p, pf = None, field_parameters = None):
+ super(YTPointBase, self).__init__(pf, field_parameters)
+ self.p = p
+
class YTOrthoRayBase(YTSelectionContainer1D):
"""
This is an orthogonal ray cast through the entire domain, at a specific
@@ -529,7 +555,7 @@
class YTSphereBase(YTSelectionContainer3D):
"""
- A sphere f points defined by a *center* and a *radius*.
+ A sphere of points defined by a *center* and a *radius*.
Parameters
----------
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/data_objects/static_output.py
--- a/yt/data_objects/static_output.py
+++ b/yt/data_objects/static_output.py
@@ -541,6 +541,34 @@
min_val, mx, my, mz)
return min_val, self.arr([mx, my, mz], 'code_length', dtype="float64")
+ def find_field_values_at_point(self, fields, coord):
+ """
+ Returns the values [field1, field2,...] of the fields at the given
+ coordinates. Returns a list of field values in the same order as
+ the input *fields*.
+ """
+ return self.point(coords)[fields]
+
+ def find_field_values_at_points(self, fields, coords):
+ """
+ Returns the values [field1, field2,...] of the fields at the given
+ [(x1, y1, z2), (x2, y2, z2),...] points. Returns a list of field
+ values in the same order as the input *fields*.
+
+ This is quite slow right now as it creates a new data object for each
+ point. If an optimized version exists on the Index object we'll use
+ that instead.
+ """
+ if hasattr(self,"index") and \
+ hasattr(self.index,"_find_field_values_at_points"):
+ return self.index._find_field_values_at_points(fields,coords)
+
+ fields = ensure_list(fields)
+ out = np.zeros((len(fields),len(coords)), dtype=np.float64)
+ for i,coord in enumerate(coords):
+ out[:][i] = self.point(coord)[fields]
+ return out
+
# Now all the object related stuff
def all_data(self, find_max=False):
if find_max: c = self.find_max("density")[1]
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/data_objects/tests/test_points.py
--- /dev/null
+++ b/yt/data_objects/tests/test_points.py
@@ -0,0 +1,10 @@
+from yt.testing import *
+import numpy as np
+
+def setup():
+ from yt.config import ytcfg
+ ytcfg["yt","__withintesting"] = "True"
+
+def test_domain_point():
+ pf = fake_random_pf(16, fields = ("density"))
+ p = pf.point(pf.domain_center)
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/frontends/artio/_artio_caller.pyx
--- a/yt/frontends/artio/_artio_caller.pyx
+++ b/yt/frontends/artio/_artio_caller.pyx
@@ -300,6 +300,10 @@
cdef int num_fields = len(fields)
cdef np.float64_t pos[3]
+ # since RuntimeErrors are not fatal, ensure no artio_particles* functions
+ # called if fileset lacks particles
+ if not self.has_particles: return
+
data = {}
accessed_species = np.zeros( self.num_species, dtype="int")
selected_mass = [ None for i in range(self.num_species)]
@@ -587,15 +591,25 @@
self.handle = artio_handle.handle
self.oct_count = None
self.root_mesh_data = NULL
- self.pcount = <np.int64_t **> malloc(sizeof(np.int64_t*)
- * artio_handle.num_species)
- self.nvars[0] = artio_handle.num_species
- self.nvars[1] = artio_handle.num_grid_variables
- for i in range(artio_handle.num_species):
- self.pcount[i] = <np.int64_t*> malloc(sizeof(np.int64_t)
- * (self.sfc_end - self.sfc_start + 1))
- for sfc in range(self.sfc_end - self.sfc_start + 1):
- self.pcount[i][sfc] = 0
+ self.pcount = NULL
+
+ if artio_handle.has_particles:
+ self.pcount = <np.int64_t **> malloc(sizeof(np.int64_t*)
+ * artio_handle.num_species)
+ self.nvars[0] = artio_handle.num_species
+ for i in range(artio_handle.num_species):
+ self.pcount[i] = <np.int64_t*> malloc(sizeof(np.int64_t)
+ * (self.sfc_end - self.sfc_start + 1))
+ for sfc in range(self.sfc_end - self.sfc_start + 1):
+ self.pcount[i][sfc] = 0
+ else:
+ self.nvars[0] = 0
+
+ if artio_handle.has_grid:
+ self.nvars[1] = artio_handle.num_grid_variables
+ else:
+ self.nvars[1] = 0
+
for i in range(3):
self.dims[i] = domain_dimensions[i]
self.DLE[i] = domain_left_edge[i]
@@ -604,12 +618,15 @@
def __dealloc__(self):
cdef int i
- for i in range(self.nvars[0]):
- free(self.pcount[i])
- for i in range(self.nvars[1]):
- free(self.root_mesh_data[i])
- free(self.pcount)
- free(self.root_mesh_data)
+ if self.artio_handle.has_particles:
+ for i in range(self.nvars[0]):
+ free(self.pcount[i])
+ free(self.pcount)
+ if self.artio_handle.has_grid:
+ if self.root_mesh_data != NULL:
+ for i in range(self.nvars[1]):
+ free(self.root_mesh_data[i])
+ free(self.root_mesh_data)
@cython.boundscheck(False)
@cython.wraparound(False)
@@ -623,8 +640,7 @@
cdef int *num_octs_per_level = <int *>malloc(
(max_level + 1)*sizeof(int))
cdef int num_species = self.artio_handle.num_species
- cdef int *num_particles_per_species = <int *>malloc(
- sizeof(int)*num_species)
+ cdef int *num_particles_per_species
cdef ARTIOOctreeContainer octree
ngv = self.nvars[1]
cdef float *grid_variables = <float *>malloc(
@@ -648,10 +664,10 @@
status = artio_grid_read_root_cell_begin( self.handle,
sfc, dpos, grid_variables, &num_oct_levels,
num_octs_per_level)
+ check_artio_status(status)
for i in range(ngv):
self.root_mesh_data[i][sfc - self.sfc_start] = \
grid_variables[i]
- check_artio_status(status)
if num_oct_levels > 0:
oc = 0
for level in range(num_oct_levels):
@@ -660,33 +676,39 @@
oct_count[sfc - self.sfc_start] = oc
octree.initialize_local_mesh(oc, num_oct_levels,
num_octs_per_level, sfc)
- status = artio_grid_read_root_cell_end( self.handle )
+ status = artio_grid_read_root_cell_end(self.handle)
check_artio_status(status)
- status = artio_grid_clear_sfc_cache( self.handle)
+ status = artio_grid_clear_sfc_cache(self.handle)
check_artio_status(status)
- # Now particles
- status = artio_particle_cache_sfc_range(self.handle, self.sfc_start,
+
+ if self.artio_handle.has_particles:
+ num_particles_per_species = <int *>malloc(
+ sizeof(int)*num_species)
+
+ # Now particles
+ status = artio_particle_cache_sfc_range(self.handle, self.sfc_start,
self.sfc_end)
- check_artio_status(status)
- for sfc in range(self.sfc_start, self.sfc_end + 1):
- # Now particles
- status = artio_particle_read_root_cell_begin( self.handle,
- sfc, num_particles_per_species )
+ check_artio_status(status)
+ for sfc in range(self.sfc_start, self.sfc_end + 1):
+ # Now particles
+ status = artio_particle_read_root_cell_begin(self.handle,
+ sfc, num_particles_per_species)
+ check_artio_status(status)
+
+ for i in range(num_species):
+ self.pcount[i][sfc - self.sfc_start] = \
+ num_particles_per_species[i]
+
+ status = artio_particle_read_root_cell_end(self.handle)
+ check_artio_status(status)
+
+ status = artio_particle_clear_sfc_cache(self.handle)
check_artio_status(status)
- for i in range(num_species):
- self.pcount[i][sfc - self.sfc_start] = \
- num_particles_per_species[i]
-
- status = artio_particle_read_root_cell_end( self.handle)
- check_artio_status(status)
-
- status = artio_particle_clear_sfc_cache( self.handle)
- check_artio_status(status)
+ free(num_particles_per_species)
free(grid_variables)
free(num_octs_per_level)
- free(num_particles_per_species)
self.oct_count = oct_count
self.doct_count = <np.int64_t *> oct_count.data
self.root_mesh_handler = ARTIORootMeshContainer(self)
@@ -861,7 +883,7 @@
cdef np.int64_t sfc_start, sfc_end
sfc_start = self.domains[0].con_id
sfc_end = self.domains[self.num_domains - 1].con_id
- status = artio_grid_cache_sfc_range(handle, sfc_start, sfc_end )
+ status = artio_grid_cache_sfc_range(handle, sfc_start, sfc_end)
check_artio_status(status)
cdef np.int64_t offset = 0
for si in range(self.num_domains):
@@ -884,7 +906,7 @@
status = artio_grid_read_level_end(handle)
check_artio_status(status)
lp += num_octs_per_level[level]
- status = artio_grid_read_root_cell_end( handle )
+ status = artio_grid_read_root_cell_end(handle)
check_artio_status(status)
# Now we have all our sources.
for j in range(nf):
@@ -943,6 +965,10 @@
cdef np.ndarray[np.int64_t, ndim=1] npi64arr
cdef np.ndarray[np.float64_t, ndim=1] npf64arr
+ if not artio_handle.has_particles:
+ raise RuntimeError("Attempted to read non-existent particles in ARTIO")
+ return
+
# Now we set up our field pointers
params = artio_handle.parameters
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/frontends/artio/data_structures.py
--- a/yt/frontends/artio/data_structures.py
+++ b/yt/frontends/artio/data_structures.py
@@ -182,29 +182,6 @@
def convert(self, unit):
return self.parameter_file.conversion_factors[unit]
- def find_max(self, field, finest_levels=3):
- """
- Returns (value, center) of location of maximum for a given field.
- """
- if (field, finest_levels) in self._max_locations:
- return self._max_locations[(field, finest_levels)]
- mv, pos = self.find_max_cell_location(field, finest_levels)
- self._max_locations[(field, finest_levels)] = (mv, pos)
- return mv, pos
-
- def find_max_cell_location(self, field, finest_levels=3):
- source = self.all_data()
- if finest_levels is not False:
- source.min_level = self.max_level - finest_levels
- mylog.debug("Searching for maximum value of %s", field)
- max_val, maxi, mx, my, mz = \
- source.quantities["MaxLocation"](field)
- mylog.info("Max Value is %0.5e at %0.16f %0.16f %0.16f",
- max_val, mx, my, mz)
- self.pf.parameters["Max%sValue" % (field)] = max_val
- self.pf.parameters["Max%sPos" % (field)] = "%s" % ((mx, my, mz),)
- return max_val, np.array((mx, my, mz), dtype='float64')
-
def _detect_output_fields(self):
self.fluid_field_list = self._detect_fluid_fields()
self.particle_field_list = self._detect_particle_fields()
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/frontends/boxlib/data_structures.py
--- a/yt/frontends/boxlib/data_structures.py
+++ b/yt/frontends/boxlib/data_structures.py
@@ -343,7 +343,6 @@
self._parallel_locking = False
self._data_file = None
self._data_mode = None
- self._max_locations = {}
def _detect_output_fields(self):
# This is all done in _parse_header_file
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/frontends/stream/data_structures.py
--- a/yt/frontends/stream/data_structures.py
+++ b/yt/frontends/stream/data_structures.py
@@ -402,7 +402,7 @@
def assign_particle_data(pf, pdata) :
"""
- Assign particle data to the grids using find_points. This
+ Assign particle data to the grids using MatchPointsToGrids. This
will overwrite any existing particle data, so be careful!
"""
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/geometry/geometry_handler.py
--- a/yt/geometry/geometry_handler.py
+++ b/yt/geometry/geometry_handler.py
@@ -72,7 +72,6 @@
self._parallel_locking = False
self._data_file = None
self._data_mode = None
- self._max_locations = {}
self.num_grids = None
def _initialize_data_storage(self):
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/geometry/grid_geometry_handler.py
--- a/yt/geometry/grid_geometry_handler.py
+++ b/yt/geometry/grid_geometry_handler.py
@@ -201,32 +201,36 @@
for item in ("Mpc", "pc", "AU", "cm"):
print "\tWidth: %0.3e %s" % (dx.in_units(item), item)
- def find_max(self, field, finest_levels = 3):
+ def _find_field_values_at_points(self, fields, coords):
+ r"""Find the value of fields at a set of coordinates.
+
+ Returns the values [field1, field2,...] of the fields at the given
+ (x, y, z) points. Returns a numpy array of field values cross coords
"""
- Returns (value, center) of location of maximum for a given field.
- """
- if (field, finest_levels) in self._max_locations:
- return self._max_locations[(field, finest_levels)]
- mv, pos = self.find_max_cell_location(field, finest_levels)
- self._max_locations[(field, finest_levels)] = (mv, pos)
- return mv, pos
+ coords = YTArray(ensure_numpy_array(coords),'code_length', registry=self.pf.unit_registry)
+ grids = self._find_points(coords[:,0], coords[:,1], coords[:,2])[0]
+ fields = ensure_list(fields)
+ mark = np.zeros(3, dtype=np.int)
+ out = []
- def find_max_cell_location(self, field, finest_levels = 3):
- if finest_levels is not False:
- gi = (self.grid_levels >= self.max_level - finest_levels).ravel()
- source = self.data_collection([0.0]*3, self.grids[gi])
- else:
- source = self.all_data()
- mylog.debug("Searching for maximum value of %s", field)
- max_val, maxi, mx, my, mz = \
- source.quantities["MaxLocation"](field)
- mylog.info("Max Value is %0.5e at %0.16f %0.16f %0.16f",
- max_val, mx, my, mz)
- self.parameters["Max%sValue" % (field)] = max_val
- self.parameters["Max%sPos" % (field)] = "%s" % ((mx,my,mz),)
- return max_val, np.array((mx,my,mz), dtype='float64')
+ # create point -> grid mapping
+ grid_index = {}
+ for coord_index, grid in enumerate(grids):
+ if not grid_index.has_key(grid):
+ grid_index[grid] = []
+ grid_index[grid].append(coord_index)
- def find_points(self, x, y, z) :
+ out = np.zeros((len(fields),len(coords)), dtype=np.float64)
+ for grid in grid_index:
+ cellwidth = (grid.RightEdge - grid.LeftEdge) / grid.ActiveDimensions
+ for field in fields:
+ for coord_index in grid_index[grid]:
+ mark = ((coords[coord_index,:] - grid.LeftEdge) / cellwidth).astype('int')
+ out[:,coord_index] = grid[field][mark[0],mark[1],mark[2]]
+ return out
+
+
+ def _find_points(self, x, y, z) :
"""
Returns the (objects, indices) of leaf grids containing a number of (x,y,z) points
"""
@@ -236,46 +240,12 @@
if not len(x) == len(y) == len(z):
raise AssertionError("Arrays of indices must be of the same size")
- grid_tree = self.get_grid_tree()
+ grid_tree = self._get_grid_tree()
pts = MatchPointsToGrids(grid_tree, len(x), x, y, z)
ind = pts.find_points_in_tree()
return self.grids[ind], ind
- def find_field_value_at_point(self, fields, coord):
- r"""Find the value of fields at a coordinate.
-
- Returns the values [field1, field2,...] of the fields at the given
- (x, y, z) points. Returns a list of field values in the same order as
- the input *fields*.
-
- Parameters
- ----------
- fields : string or list of strings
- The field(s) that will be returned.
-
- coord : list or array of coordinates
- The location for which field values will be returned.
-
- Examples
- --------
- >>> pf.h.find_field_value_at_point(['Density', 'Temperature'],
- [0.4, 0.3, 0.8])
- [2.1489e-24, 1.23843e4]
- """
- this = self.find_points(*coord)[0][-1]
- cellwidth = (this.RightEdge - this.LeftEdge) / this.ActiveDimensions
- mark = np.zeros(3).astype('int')
- # Find the index for the cell containing this point.
- for dim in xrange(len(coord)):
- mark[dim] = int((coord[dim] - this.LeftEdge[dim]) / cellwidth[dim])
- out = []
- fields = ensure_list(fields)
- # Pull out the values and add it to the out list.
- for field in fields:
- out.append(this[field][mark[0], mark[1], mark[2]])
- return out
-
- def get_grid_tree(self) :
+ def _get_grid_tree(self) :
left_edge = self.pf.arr(np.zeros((self.num_grids, 3)),
'code_length')
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/geometry/oct_geometry_handler.py
--- a/yt/geometry/oct_geometry_handler.py
+++ b/yt/geometry/oct_geometry_handler.py
@@ -49,27 +49,3 @@
def convert(self, unit):
return self.parameter_file.conversion_factors[unit]
-
- def find_max(self, field, finest_levels = 3):
- """
- Returns (value, center) of location of maximum for a given field.
- """
- if (field, finest_levels) in self._max_locations:
- return self._max_locations[(field, finest_levels)]
- mv, pos = self.find_max_cell_location(field, finest_levels)
- self._max_locations[(field, finest_levels)] = (mv, pos)
- return mv, pos
-
- def find_max_cell_location(self, field, finest_levels = 3):
- source = self.all_data()
- if finest_levels is not False:
- source.min_level = self.max_level - finest_levels
- mylog.debug("Searching for maximum value of %s", field)
- max_val, maxi, mx, my, mz = \
- source.quantities["MaxLocation"](field)
- mylog.info("Max Value is %0.5e at %0.16f %0.16f %0.16f",
- max_val, mx, my, mz)
- self.pf.parameters["Max%sValue" % (field,)] = max_val
- self.pf.parameters["Max%sPos" % (field,)] = "%s" % ((mx,my,mz),)
- return max_val, np.array((mx,my,mz), dtype='float64')
-
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/geometry/selection_routines.pyx
--- a/yt/geometry/selection_routines.pyx
+++ b/yt/geometry/selection_routines.pyx
@@ -534,6 +534,63 @@
self.domain_width[1],
self.domain_width[2])
+
+cdef class PointSelector(SelectorObject):
+ cdef np.float64_t p[3]
+
+ def __init__(self, dobj):
+ for i in range(3):
+ self.p[i] = _ensure_code(dobj.p[i])
+
+ # ensure the point lies in the domain
+ if self.periodicity[i]:
+ self.p[i] = np.fmod(self.p[i], self.domain_width[i])
+ if self.p[i] < 0.0:
+ self.p[i] += self.domain_width[i]
+
+ @cython.boundscheck(False)
+ @cython.wraparound(False)
+ @cython.cdivision(True)
+ cdef int select_cell(self, np.float64_t pos[3], np.float64_t dds[3]) nogil:
+ if (pos[0] - 0.5*dds[0] <= self.p[0] < pos[0]+0.5*dds[0] and
+ pos[1] - 0.5*dds[1] <= self.p[1] < pos[1]+0.5*dds[1] and
+ pos[2] - 0.5*dds[2] <= self.p[2] < pos[2]+0.5*dds[2]):
+ return 1
+ else:
+ return 0
+
+ @cython.boundscheck(False)
+ @cython.wraparound(False)
+ @cython.cdivision(True)
+ cdef int select_sphere(self, np.float64_t pos[3], np.float64_t radius) nogil:
+ cdef int i
+ cdef np.float64_t dist, dist2 = 0
+ for i in range(3):
+ dist = self.difference(pos[i], self.p[i], i)
+ dist2 += dist*dist
+ if dist2 <= radius*radius: return 1
+ return 0
+
+ @cython.boundscheck(False)
+ @cython.wraparound(False)
+ @cython.cdivision(True)
+ cdef int select_bbox(self, np.float64_t left_edge[3],
+ np.float64_t right_edge[3]) nogil:
+ cdef int i
+ # point definitely can only be in one cell
+ if (left_edge[0] <= self.p[0] < right_edge[0] and
+ left_edge[1] <= self.p[1] < right_edge[1] and
+ left_edge[2] <= self.p[2] < right_edge[2]):
+ return 1
+ else:
+ return 0
+
+ def _hash_vals(self):
+ return (self.p[0], self.p[1], self.p[2])
+
+point_selector = PointSelector
+
+
cdef class SphereSelector(SelectorObject):
cdef np.float64_t radius
cdef np.float64_t radius2
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/utilities/lib/tests/test_grid_tree.py
--- a/yt/utilities/lib/tests/test_grid_tree.py
+++ b/yt/utilities/lib/tests/test_grid_tree.py
@@ -50,7 +50,7 @@
def test_grid_tree():
"""Main test suite for GridTree"""
- grid_tree = test_pf.index.get_grid_tree()
+ grid_tree = test_pf.index._get_grid_tree()
indices, levels, nchild, children = grid_tree.return_tree_info()
grid_levels = [grid.Level for grid in test_pf.index.grids]
@@ -67,53 +67,3 @@
grid_children = np.array([child.id - child._id_offset
for child in grid.Children])
yield assert_equal, grid_children, children[i]
-
-
-def test_find_points():
- """Main test suite for MatchPoints"""
- num_points = 100
- randx = np.random.uniform(low=test_pf.domain_left_edge[0],
- high=test_pf.domain_right_edge[0],
- size=num_points)
- randy = np.random.uniform(low=test_pf.domain_left_edge[1],
- high=test_pf.domain_right_edge[1],
- size=num_points)
- randz = np.random.uniform(low=test_pf.domain_left_edge[2],
- high=test_pf.domain_right_edge[2],
- size=num_points)
-
- point_grids, point_grid_inds = test_pf.index.find_points(randx, randy, randz)
-
- grid_inds = np.zeros((num_points), dtype='int64')
-
- for ind, ixx, iyy, izz in zip(range(num_points), randx, randy, randz):
-
- pos = np.array([ixx, iyy, izz])
- pt_level = -1
-
- for grid in test_pf.index.grids:
-
- if np.all(pos >= grid.LeftEdge) and \
- np.all(pos <= grid.RightEdge) and \
- grid.Level > pt_level:
- pt_level = grid.Level
- grid_inds[ind] = grid.id - grid._id_offset
-
- yield assert_equal, point_grid_inds, grid_inds
-
- # Test wheter find_points works for lists
- point_grids, point_grid_inds = test_pf.index.find_points(randx.tolist(),
- randy.tolist(),
- randz.tolist())
- yield assert_equal, point_grid_inds, grid_inds
-
- # Test if find_points works for scalar
- ind = random.randint(0, num_points - 1)
- point_grids, point_grid_inds = test_pf.index.find_points(randx[ind],
- randy[ind],
- randz[ind])
- yield assert_equal, point_grid_inds, grid_inds[ind]
-
- # Test if find_points fails properly for non equal indices' array sizes
- yield assert_raises, AssertionError, test_pf.index.find_points, \
- [0], 1.0, [2, 3]
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/utilities/particle_generator.py
--- a/yt/utilities/particle_generator.py
+++ b/yt/utilities/particle_generator.py
@@ -99,7 +99,7 @@
Assigns grids to particles and sets up particle positions. *setup_fields* is
a dict of fields other than the particle positions to set up.
"""
- particle_grids, particle_grid_inds = self.pf.index.find_points(x,y,z)
+ particle_grids, particle_grid_inds = self.pf.index._find_points(x,y,z)
idxs = np.argsort(particle_grid_inds)
self.particles[:,self.posx_index] = x[idxs]
self.particles[:,self.posy_index] = y[idxs]
diff -r ef7c988b774a159c854d7a2fefc000e767007eda -r 561d2e56371ca0053ee093df7968f59ae978481a yt/utilities/tests/test_selectors.py
--- a/yt/utilities/tests/test_selectors.py
+++ b/yt/utilities/tests/test_selectors.py
@@ -1,7 +1,7 @@
import numpy as np
from yt.testing import \
fake_random_pf, assert_equal, assert_array_less, \
- YTArray
+ YTArray, fake_amr_pf
from yt.utilities.math_utils import periodic_dist
@@ -9,6 +9,26 @@
from yt.config import ytcfg
ytcfg["yt", "__withintesting"] = "True"
+def test_point_selector():
+ # generate fake amr data
+ pf = fake_random_pf(64, nprocs=51)
+ assert(all(pf.periodicity))
+
+ dd = pf.h.all_data()
+ positions = np.array([dd[ax] for ax in 'xyz'])
+ delta = 0.5*np.array([dd['d'+ax] for ax in 'xyz'])
+ # ensure cell centers and corners always return one and
+ # only one point object
+ for p in positions:
+ data = pf.point(p)
+ assert_equal(data["ones"].shape[0], 1)
+ for p in positions - delta:
+ data = pf.point(p)
+ assert_equal(data["ones"].shape[0], 1)
+ for p in positions + delta:
+ data = pf.point(p)
+ assert_equal(data["ones"].shape[0], 1)
+
def test_sphere_selector():
# generate fake data with a number of non-cubical grids
pf = fake_random_pf(64, nprocs=51)
Repository URL: https://bitbucket.org/yt_analysis/yt/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the yt-svn
mailing list