[Yt-svn] yt: 3 new changesets
hg at spacepope.org
hg at spacepope.org
Tue Nov 9 16:20:32 PST 2010
hg Repository: yt
details: yt/rev/c7947fef16ac
changeset: 3515:c7947fef16ac
user: Sam Skillman <sam.skillman at gmail.com>
date:
Tue Nov 09 17:00:35 2010 -0700
description:
Adding AMR kd-Tree framework and rendering, as well as a few more
options in the rendering process.
This involves several pieces of yt, detailed below:
yt/utilities/amr_kdtree/amr_kdtree.py - Can be used to partition an
amr hierarchy into kd-Tree whose leaves make up a homogenized,
non-overlapping, volume. Details and usage can be found in the
docstrings for the AMRKDTree class. The main object is the
AMRKDTree.tree, which is a dictionary whose keys are unique binary
tree ids that point to a position in the kd-Tree. The kd-Tree
framework allows for fast location of positions, and
viewpoint-independent ordering of homogenized bricks, ideal for the
volume rendering.
yt/visualization/volume_rendering/camera.py - This camera object has
been updated to, by default, use the kd-Tree for volume decomposition
and ray casting. Basic usage of the camera object has remained the
same, though there are many new options that can be seen in the
docstrings for the Camera class. An email accompanying this commit
will detail the advantages, options, and several examples.
yt/data_objects/grid_patch.py - Modified to accept a keyword in
get_vertex_centered_data(..., no_ghost=False). When no_ghost=True,
instead of gathering the ghost zones from surrounding grids, the ghost
outer vertices are filled by extrapolation from within the grid. This
leads to very large (up to 300x in some cases) speedup, but is not as
accurate. Details of when this option should be used are in the
Camera class docs.
yt/visualization/volume_rendering/transfer_functions.py - Added a
keyword col_bounds=None that specifies the colorbar range in the
transfer function when using sample_colormap or add_layers. This
confines the colormap to a range of values that can differ from the
transfer function bounds themselves. It has made it easier to sample
a full colormap over a reduced set of values.
hg Repository: yt
details: yt/rev/2ff8e32f52d0
changeset: 3516:2ff8e32f52d0
user: Sam Skillman <sam.skillman at gmail.com>
date:
Tue Nov 09 17:01:03 2010 -0700
description:
merging
hg Repository: yt
details: yt/rev/4b7fca0be218
changeset: 3517:4b7fca0be218
user: Sam Skillman <sam.skillman at gmail.com>
date:
Tue Nov 09 17:07:24 2010 -0700
description:
Adding better docstrings to camera.
diffstat:
yt/analysis_modules/halo_finding/halo_objects.py | 2 +-
yt/data_objects/grid_patch.py | 51 +-
yt/frontends/enzo/fields.py | 19 +-
yt/utilities/_amr_utils/RayIntegrators.pyx | 34 +-
yt/utilities/amr_kdtree/__init__.py | 4 +
yt/utilities/amr_kdtree/amr_kdtree.py | 1159 +++
yt/utilities/amr_utils.c | 4789 ++++++--------
yt/utilities/setup.py | 1 +
yt/visualization/volume_rendering/camera.py | 164 +-
yt/visualization/volume_rendering/transfer_functions.py | 30 +-
10 files changed, 3573 insertions(+), 2680 deletions(-)
diffs (truncated from 17077 to 300 lines):
diff -r f1b6415806b3 -r 4b7fca0be218 yt/analysis_modules/halo_finding/halo_objects.py
--- a/yt/analysis_modules/halo_finding/halo_objects.py Mon Nov 08 09:37:17 2010 -0800
+++ b/yt/analysis_modules/halo_finding/halo_objects.py Tue Nov 09 17:07:24 2010 -0700
@@ -1604,7 +1604,7 @@
bin_width = base_padding
num_bins = int(math.ceil(width / bin_width))
bins = na.arange(num_bins+1, dtype='float64') * bin_width + self._data_source.left_edge[dim]
- counts, bins = na.histogram(data, bins, new=True)
+ counts, bins = na.histogram(data, bins)
# left side.
start = 0
count = counts[0]
diff -r f1b6415806b3 -r 4b7fca0be218 yt/data_objects/grid_patch.py
--- a/yt/data_objects/grid_patch.py Mon Nov 08 09:37:17 2010 -0800
+++ b/yt/data_objects/grid_patch.py Tue Nov 09 17:07:24 2010 -0700
@@ -446,18 +446,43 @@
level, new_left_edge, **kwargs)
return cube
- def get_vertex_centered_data(self, field, smoothed=True):
- cg = self.retrieve_ghost_zones(1, field, smoothed=smoothed)
- # We have two extra zones in every direction
- new_field = na.zeros(self.ActiveDimensions + 1, dtype='float64')
- na.add(new_field, cg[field][1: ,1: ,1: ], new_field)
- na.add(new_field, cg[field][:-1,1: ,1: ], new_field)
- na.add(new_field, cg[field][1: ,:-1,1: ], new_field)
- na.add(new_field, cg[field][1: ,1: ,:-1], new_field)
- na.add(new_field, cg[field][:-1,1: ,:-1], new_field)
- na.add(new_field, cg[field][1: ,:-1,:-1], new_field)
- na.add(new_field, cg[field][:-1,:-1,1: ], new_field)
- na.add(new_field, cg[field][:-1,:-1,:-1], new_field)
- na.multiply(new_field, 0.125, new_field)
+ def get_vertex_centered_data(self, field, smoothed=True,
+ no_ghost=False):
+ if not no_ghost:
+ cg = self.retrieve_ghost_zones(1, field, smoothed=smoothed)
+ # We have two extra zones in every direction
+ new_field = na.zeros(self.ActiveDimensions + 1, dtype='float64')
+ na.add(new_field, cg[field][1: ,1: ,1: ], new_field)
+ na.add(new_field, cg[field][:-1,1: ,1: ], new_field)
+ na.add(new_field, cg[field][1: ,:-1,1: ], new_field)
+ na.add(new_field, cg[field][1: ,1: ,:-1], new_field)
+ na.add(new_field, cg[field][:-1,1: ,:-1], new_field)
+ na.add(new_field, cg[field][1: ,:-1,:-1], new_field)
+ na.add(new_field, cg[field][:-1,:-1,1: ], new_field)
+ na.add(new_field, cg[field][:-1,:-1,:-1], new_field)
+ na.multiply(new_field, 0.125, new_field)
+ else:
+ new_field = na.zeros(self.ActiveDimensions + 1, dtype='float64')
+ of = self[field]
+ new_field[:-1,:-1,:-1] += of
+ new_field[:-1,:-1,1:] += of
+ new_field[:-1,1:,:-1] += of
+ new_field[:-1,1:,1:] += of
+ new_field[1:,:-1,:-1] += of
+ new_field[1:,:-1,1:] += of
+ new_field[1:,1:,:-1] += of
+ new_field[1:,1:,1:] += of
+ na.multiply(new_field, 0.125, new_field)
+ new_field = na.log10(new_field)
+
+ new_field[:,:, -1] = 2.0*new_field[:,:,-2] - new_field[:,:,-3]
+ new_field[:,:, 0] = 2.0*new_field[:,:,1] - new_field[:,:,2]
+
+ new_field[:,-1, :] = 2.0*new_field[:,-2,:] - new_field[:,-3,:]
+ new_field[:,0, :] = 2.0*new_field[:,1,:] - new_field[:,2,:]
+
+ new_field[-1,:,:] = 2.0*new_field[-2,:,:] - new_field[-3,:,:]
+ new_field[0,:,:] = 2.0*new_field[1,:,:] - new_field[2,:,:]
+ na.power(10.0,new_field,new_field)
return new_field
diff -r f1b6415806b3 -r 4b7fca0be218 yt/frontends/enzo/fields.py
--- a/yt/frontends/enzo/fields.py Mon Nov 08 09:37:17 2010 -0800
+++ b/yt/frontends/enzo/fields.py Tue Nov 09 17:07:24 2010 -0700
@@ -212,7 +212,8 @@
_default_fields = ["Density","Temperature",
"x-velocity","y-velocity","z-velocity",
- "x-momentum","y-momentum","z-momentum"]
+ "x-momentum","y-momentum","z-momentum",
+ "Bx", "By", "Bz"]
# else:
# _default_fields = ["Density","Temperature","Gas_Energy","Total_Energy",
# "x-velocity","y-velocity","z-velocity"]
@@ -382,6 +383,22 @@
add_field('IsStarParticle', function=_IsStarParticle,
particle_type = True)
+def _convertBfield(data):
+ return na.sqrt(4*na.pi*data.convert("Density")*data.convert("x-velocity")**2)
+for field in ['Bx','By','Bz']:
+ f = EnzoFieldInfo[field]
+ f._convert_function=_convertBfield
+ f._units=r"\mathrm{Gau\ss}"
+ f.take_log=False
+
+def _Bmag(field, data):
+ """ magnitude of bvec
+ """
+ return na.sqrt(data['Bx']**2 + data['By']**2 + data['Bz']**2)
+
+add_field("Bmag", function=_Bmag,display_name=r"|B|",units=r"\mathrm{Gau\ss}")
+
+
#
# Now we do overrides for 2D fields
#
diff -r f1b6415806b3 -r 4b7fca0be218 yt/utilities/_amr_utils/RayIntegrators.pyx
--- a/yt/utilities/_amr_utils/RayIntegrators.pyx Mon Nov 08 09:37:17 2010 -0800
+++ b/yt/utilities/_amr_utils/RayIntegrators.pyx Tue Nov 09 17:07:24 2010 -0700
@@ -138,23 +138,29 @@
# Do left edge then right edge in each dim
cdef int i, x, y
cdef np.float64_t tl, tr, intersect_t, enter_t, exit_t, dt_tolerance
- cdef np.ndarray[np.int64_t, ndim=1] step = np.empty((3,), dtype=np.int64)
- cdef np.ndarray[np.int64_t, ndim=1] cur_ind = np.empty((3,), dtype=np.int64)
- cdef np.ndarray[np.float64_t, ndim=1] tdelta = np.empty((3,), dtype=np.float64)
- cdef np.ndarray[np.float64_t, ndim=1] tmax = np.empty((3,), dtype=np.float64)
- cdef np.ndarray[np.float64_t, ndim=1] intersect = np.empty((3,), dtype=np.float64)
+ cdef np.float64_t iv_dir[3], tdelta[3], tmax[3], intersect[3]
+ cdef np.int64_t step[3], cur_ind[3]
intersect_t = 1
dt_tolerance = 1e-6
# recall p = v * t + u
# where p is position, v is our vector, u is the start point
for i in range(3):
# As long as we're iterating, set some other stuff, too
- if(v[i] < 0): step[i] = -1
- else: step[i] = 1
+ if(v[i] < 0):
+ step[i] = -1
+ elif (v[i] == 0):
+ step[i] = 1
+ tmax[i] = 1e60
+ iv_dir[i] = 1e60
+ tdelta[i] = 1e-60
+ continue
+ else:
+ step[i] = 1
x = (i+1)%3
y = (i+2)%3
- tl = (left_edge[i] - u[i])/v[i]
- tr = (right_edge[i] - u[i])/v[i]
+ iv_dir[i] = 1.0/v[i]
+ tl = (left_edge[i] - u[i])*iv_dir[i]
+ tr = (right_edge[i] - u[i])*iv_dir[i]
if (left_edge[x] <= (u[x] + tl*v[x]) <= right_edge[x]) and \
(left_edge[y] <= (u[y] + tl*v[y]) <= right_edge[y]) and \
(0.0 <= tl < intersect_t):
@@ -170,15 +176,15 @@
intersect_t = 0.0
if not (0 <= intersect_t <= 1): return
# Now get the indices of the intersection
- intersect = u + intersect_t * v
for i in range(3):
+ intersect[i] = u[i] + intersect_t * v[i]
cur_ind[i] = np.floor((intersect[i] + 1e-8*dx[i] - left_edge[i])/dx[i])
- tmax[i] = (((cur_ind[i]+step[i])*dx[i])+left_edge[i]-u[i])/v[i]
+ tmax[i] = (((cur_ind[i]+step[i])*dx[i])+left_edge[i]-u[i])*iv_dir[i]
if cur_ind[i] == grid_mask.shape[i] and step[i] < 0:
cur_ind[i] = grid_mask.shape[i] - 1
- if step[i] > 0: tmax[i] = (((cur_ind[i]+1)*dx[i])+left_edge[i]-u[i])/v[i]
- if step[i] < 0: tmax[i] = (((cur_ind[i]+0)*dx[i])+left_edge[i]-u[i])/v[i]
- tdelta[i] = (dx[i]/v[i])
+ if step[i] > 0: tmax[i] = (((cur_ind[i]+1)*dx[i])+left_edge[i]-u[i])*iv_dir[i]
+ if step[i] < 0: tmax[i] = (((cur_ind[i]+0)*dx[i])+left_edge[i]-u[i])*iv_dir[i]
+ tdelta[i] = (dx[i]*iv_dir[i])
if tdelta[i] < 0: tdelta[i] *= -1
# The variable intersect contains the point we first pierce the grid
enter_t = intersect_t
diff -r f1b6415806b3 -r 4b7fca0be218 yt/utilities/amr_kdtree/__init__.py
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/yt/utilities/amr_kdtree/__init__.py Tue Nov 09 17:07:24 2010 -0700
@@ -0,0 +1,4 @@
+"""
+Initialize amr_kdtree
+"""
+from amr_kdtree import *
diff -r f1b6415806b3 -r 4b7fca0be218 yt/utilities/amr_kdtree/amr_kdtree.py
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/yt/utilities/amr_kdtree/amr_kdtree.py Tue Nov 09 17:07:24 2010 -0700
@@ -0,0 +1,1159 @@
+"""
+AMR kD-Tree Framework
+
+Authors: Samuel Skillman <samskillman at gmail.com>
+Affiliation: University of Colorado at Boulder
+Wil St. Charles <fallen751 at gmail.com>
+Affiliation: University of Colorado at Boulder
+
+Homepage: http://yt.enzotools.org/
+License:
+ Copyright (C) 2010 Samuel Skillman. All Rights Reserved.
+
+ This file is part of yt.
+
+ yt is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 3 of the License, or
+ (at your option) any later version.
+
+ This program is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program. If not, see <http://www.gnu.org/licenses/>.
+"""
+import numpy as na
+from mpi4py import MPI
+from yt.funcs import *
+from yt.visualization.volume_rendering.grid_partitioner import HomogenizedVolume
+from yt.utilities.amr_utils import PartitionedGrid
+from yt.utilities.performance_counters import yt_counters, time_function
+import yt.utilities.parallel_tools.parallel_analysis_interface as PT
+
+from time import time
+import h5py
+comm = MPI.COMM_WORLD
+my_rank = comm.rank
+nprocs = comm.size
+
+def corner_bounds(split_dim, split, current_left = None, current_right = None):
+ r"""
+ Given a kd-Tree split dimension and position and bound to be
+ modified, returns the new bound.
+
+ A simple function that replaces the `split_dim` dimension of the
+ current left or right bound with `split`. Left or Right bound is
+ chosen by specifying the `current_left` or `current_right`.
+ """
+ if(current_left is not None):
+ new_left = na.array([current_left[0],current_left[1],current_left[2]])
+ new_left[split_dim] = split
+ return new_left
+ elif(current_right is not None):
+ new_right = na.array([current_right[0],current_right[1],current_right[2]])
+ new_right[split_dim] = split
+ return new_right
+
+def _lchild_id(id): return (id<<1) + 1
+def _rchild_id(id): return (id<<1) + 2
+def _parent_id(id): return (id-1)>>1
+
+class MasterNode(object):
+ r"""
+ A MasterNode object is the building block of the AMR kd-Tree.
+ Used during the construction to act as both dividing nodes and
+ leaf nodes.
+ """
+ def __init__(self):
+ self.grids = None
+ self.parent = None
+ self.parent_grid = None
+ self.l_corner = None
+ self.r_corner = None
+
+ self.split_ax = None
+ self.split_pos = None
+
+ self.left_children = None
+ self.right_children = None
+ self.cost = 0
+ self.owner = -1
+
+ self.id = None
+
+ self.grid = None
+ self.brick = None
+ self.li = None
+ self.ri = None
+ self.dims = None
+
+ self.done = 0
+ self.cast_done = 0
+
+def set_leaf(thisnode, grid_id, leaf_l_corner, leaf_r_corner):
+ r"""
+ Sets leaf properties.
+
+ Parameters
+ ----------
+ thisnode : `MasterNode`
+ AMR kd-Tree node to be modified.
+ grid_id : `~yt.data_objects.grid_patch`
+ A grid patch that contains the data spanned by this
+ kd-Tree leaf.
+ leaf_l_corner: array_like, dimension 3
+ The left corner of the volume spanned by this leaf.
+ leaf_r_corner: array_like, dimension 3
+ The right corner of the volume spanned by this leaf.
+
+ Returns
+ -------
+ None
+ """
+ thisnode.grid = grid_id.id
+ thisnode.l_corner = leaf_l_corner
More information about the yt-svn
mailing list