[yt-svn] commit/yt: brittonsmith: Merged in MatthewTurk/yt/yt-3.0 (pull request #1126)
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Tue Sep 16 14:09:18 PDT 2014
1 new commit in yt:
https://bitbucket.org/yt_analysis/yt/commits/6f0db65e21bc/
Changeset: 6f0db65e21bc
Branch: yt
User: brittonsmith
Date: 2014-09-16 23:09:04
Summary: Merged in MatthewTurk/yt/yt-3.0 (pull request #1126)
Address some shortcomings of non-Cartesian plots
Affected #: 26 files
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc doc/source/examining/Loading_Spherical_Data.ipynb
--- /dev/null
+++ b/doc/source/examining/Loading_Spherical_Data.ipynb
@@ -0,0 +1,188 @@
+{
+ "metadata": {
+ "name": "",
+ "signature": "sha256:88ed88ce8d8f4a359052f287aea17a7cbed435ff960e195097b440191ce6c2ab"
+ },
+ "nbformat": 3,
+ "nbformat_minor": 0,
+ "worksheets": [
+ {
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# Loading Spherical Data\n",
+ "\n",
+ "With version 3.0 of yt, it has gained the ability to load data from non-Cartesian systems. This support is still being extended, but here is an example of how to load spherical data from a regularly-spaced grid. For irregularly spaced grids, a similar setup can be used, but the `load_hexahedral_mesh` method will have to be used instead.\n",
+ "\n",
+ "Note that in yt, \"spherical\" means that it is ordered $r$, $\\theta$, $\\phi$, where $\\theta$ is the declination from the azimuth (running from $0$ to $\\pi$ and $\\phi$ is the angle around the zenith (running from $0$ to $2\\pi$).\n",
+ "\n",
+ "We first start out by loading yt."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "import numpy as np\n",
+ "import yt"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Now, we create a few derived fields. The first three are just straight translations of the Cartesian coordinates, so that we can see where we are located in the data, and understand what we're seeing. The final one is just a fun field that is some combination of the three coordinates, and will vary in all dimensions."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "@yt.derived_field(name = \"sphx\", units = \"cm\", take_log=False)\n",
+ "def sphx(field, data):\n",
+ " return np.cos(data[\"phi\"]) * np.sin(data[\"theta\"])*data[\"r\"]\n",
+ "@yt.derived_field(name = \"sphy\", units = \"cm\", take_log=False)\n",
+ "def sphy(field, data):\n",
+ " return np.sin(data[\"phi\"]) * np.sin(data[\"theta\"])*data[\"r\"]\n",
+ "@yt.derived_field(name = \"sphz\", units = \"cm\", take_log=False)\n",
+ "def sphz(field, data):\n",
+ " return np.cos(data[\"theta\"])*data[\"r\"]\n",
+ "@yt.derived_field(name = \"funfield\", units=\"cm\", take_log=False)\n",
+ "def funfield(field, data):\n",
+ " return (np.sin(data[\"phi\"])**2 + np.cos(data[\"theta\"])**2) * (1.0*data[\"r\"].uq+data[\"r\"])"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Loading Data\n",
+ "\n",
+ "Now we can actually load our data. We use the `load_uniform_grid` function here. Normally, the first argument would be a dictionary of field data, where the keys were the field names and the values the field data arrays. Here, we're just going to look at derived fields, so we supply an empty one.\n",
+ "\n",
+ "The next few arguments are the number of dimensions, the bounds, and we then specify the geometry as spherical."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "ds = yt.load_uniform_grid({}, [128, 128, 128],\n",
+ " bbox=np.array([[0.0, 1.0], [0.0, np.pi], [0.0, 2*np.pi]]),\n",
+ " geometry=\"spherical\")"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Looking at Data\n",
+ "\n",
+ "Now we can take slices. The first thing we will try is making a slice of data along the \"phi\" axis, here $\\pi/2$, which will be along the y axis in the positive direction. We use the `.slice` attribute, which creates a slice, and then we convert this into a plot window. Note that here 2 is used to indicate the third axis (0-indexed) which for spherical data is $\\phi$.\n",
+ "\n",
+ "This is the manual way of creating a plot -- below, we'll use the standard, automatic ways. Note that the coordinates run from $-r$ to $r$ along the $z$ axis and from $0$ to $r$ along the $R$ axis. We use the capital $R$ to indicate that it's the $R$ along the $x-y$ plane."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "s = ds.slice(2, np.pi/2)\n",
+ "p = s.to_pw(\"funfield\", origin=\"native\")\n",
+ "p.set_zlim(\"all\", 0.0, 4.0)\n",
+ "p.show()"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We can also slice along $r$. For now, this creates a regular grid with *incorrect* units for phi and theta. We are currently exploring two other options -- a simple aitoff projection, and fixing it to use the correct units as-is."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "s = yt.SlicePlot(ds, \"r\", \"funfield\")\n",
+ "s.set_zlim(\"all\", 0.0, 4.0)\n",
+ "s.show()"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We can also slice at constant $\\theta$. But, this is a weird thing! We're slicing at a constant declination from the azimuth. What this means is that when thought of in a Cartesian domain, this slice is actually a cone. The axes have been labeled appropriately, to indicate that these are not exactly the $x$ and $y$ axes, but instead differ by a factor of $\\sin(\\theta))$."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "s = yt.SlicePlot(ds, \"theta\", \"funfield\")\n",
+ "s.set_zlim(\"all\", 0.0, 4.0)\n",
+ "s.show()"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We've seen lots of the `funfield` plots, but we can also look at the Cartesian axes. This next plot plots the Cartesian $x$, $y$ and $z$ values on a $\\theta$ slice. Because we're not supplying an argument to the `center` parameter, yt will place it at the center of the $\\theta$ axis, which will be at $\\pi/2$, where it will be aligned with the $x-y$ plane. The slight change in `sphz` results from the cells themselves migrating, and plotting the center of those cells."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": false,
+ "input": [
+ "s = yt.SlicePlot(ds, \"theta\", [\"sphx\", \"sphy\", \"sphz\"])\n",
+ "s.show()"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We can do the same with the $\\phi$ axis."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "collapsed": true,
+ "input": [
+ "s = yt.SlicePlot(ds, \"phi\", [\"sphx\", \"sphy\", \"sphz\"])\n",
+ "s.show()"
+ ],
+ "language": "python",
+ "metadata": {},
+ "outputs": []
+ }
+ ],
+ "metadata": {}
+ }
+ ]
+}
\ No newline at end of file
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc doc/source/examining/index.rst
--- a/doc/source/examining/index.rst
+++ b/doc/source/examining/index.rst
@@ -9,4 +9,5 @@
loading_data
generic_array_data
generic_particle_data
+ spherical_data
low_level_inspection
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc doc/source/examining/spherical_data.rst
--- /dev/null
+++ b/doc/source/examining/spherical_data.rst
@@ -0,0 +1,6 @@
+.. _loading-spherical-data:
+
+Loading Spherical Data
+======================
+
+.. notebook:: Loading_Spherical_Data.ipynb
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/data_objects/data_containers.py
--- a/yt/data_objects/data_containers.py
+++ b/yt/data_objects/data_containers.py
@@ -804,7 +804,8 @@
skip += list(set(frb._exclude_fields).difference(set(self._key_fields)))
self.fields = ensure_list(fields) + \
[k for k in self.field_data if k not in skip]
- (bounds, center) = get_window_parameters(axis, center, width, self.ds)
+ (bounds, center, display_center) = \
+ get_window_parameters(axis, center, width, self.ds)
pw = PWViewerMPL(self, bounds, fields=self.fields, origin=origin,
frb_generator=frb, plot_type=plot_type)
pw._setup_plots()
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/data_objects/static_output.py
--- a/yt/data_objects/static_output.py
+++ b/yt/data_objects/static_output.py
@@ -45,17 +45,12 @@
YTArray, \
YTQuantity
-from yt.geometry.cartesian_coordinates import \
- CartesianCoordinateHandler
-from yt.geometry.polar_coordinates import \
- PolarCoordinateHandler
-from yt.geometry.cylindrical_coordinates import \
- CylindricalCoordinateHandler
-from yt.geometry.spherical_coordinates import \
- SphericalCoordinateHandler
-from yt.geometry.geographic_coordinates import \
- GeographicCoordinateHandler
-from yt.geometry.spec_cube_coordinates import \
+from yt.geometry.coordinates.api import \
+ CartesianCoordinateHandler, \
+ PolarCoordinateHandler, \
+ CylindricalCoordinateHandler, \
+ SphericalCoordinateHandler, \
+ GeographicCoordinateHandler, \
SpectralCubeCoordinateHandler
# We want to support the movie format in the future.
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/funcs.py
--- a/yt/funcs.py
+++ b/yt/funcs.py
@@ -758,3 +758,9 @@
if os.path.isfile(_fn):
mylog.info("Loading plugins from %s", _fn)
execfile(_fn, yt.__dict__)
+
+def fix_unitary(u):
+ if u == '1':
+ return 'unitary'
+ else:
+ return u
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/cartesian_coordinates.py
--- a/yt/geometry/cartesian_coordinates.py
+++ /dev/null
@@ -1,122 +0,0 @@
-"""
-Cartesian fields
-
-
-
-
-"""
-
-#-----------------------------------------------------------------------------
-# Copyright (c) 2013, yt Development Team.
-#
-# Distributed under the terms of the Modified BSD License.
-#
-# The full license is in the file COPYING.txt, distributed with this software.
-#-----------------------------------------------------------------------------
-
-import numpy as np
-from .coordinate_handler import \
- CoordinateHandler, \
- _unknown_coord, \
- _get_coord_fields
-import yt.visualization._MPL as _MPL
-
-class CartesianCoordinateHandler(CoordinateHandler):
-
- def __init__(self, ds):
- super(CartesianCoordinateHandler, self).__init__(ds)
-
- def setup_fields(self, registry):
- for axi, ax in enumerate('xyz'):
- f1, f2 = _get_coord_fields(axi)
- registry.add_field(("index", "d%s" % ax), function = f1,
- display_field = False,
- units = "code_length")
- registry.add_field(("index", "%s" % ax), function = f2,
- display_field = False,
- units = "code_length")
- def _cell_volume(field, data):
- rv = data["index", "dx"].copy(order='K')
- rv *= data["index", "dy"]
- rv *= data["index", "dz"]
- return rv
- registry.add_field(("index", "cell_volume"), function=_cell_volume,
- display_field=False, units = "code_length**3")
- registry.check_derived_fields(
- [("index", "dx"), ("index", "dy"), ("index", "dz"),
- ("index", "x"), ("index", "y"), ("index", "z"),
- ("index", "cell_volume")])
-
- def pixelize(self, dimension, data_source, field, bounds, size,
- antialias = True, periodic = True):
- if dimension < 3:
- return self._ortho_pixelize(data_source, field, bounds, size,
- antialias, dimension, periodic)
- else:
- return self._oblique_pixelize(data_source, field, bounds, size,
- antialias)
-
- def _ortho_pixelize(self, data_source, field, bounds, size, antialias,
- dim, periodic):
- # We should be using fcoords
- period = self.period[:2].copy() # dummy here
- period[0] = self.period[self.x_axis[dim]]
- period[1] = self.period[self.y_axis[dim]]
- if hasattr(period, 'in_units'):
- period = period.in_units("code_length").d
- buff = _MPL.Pixelize(data_source['px'], data_source['py'],
- data_source['pdx'], data_source['pdy'],
- data_source[field], size[0], size[1],
- bounds, int(antialias),
- period, int(periodic)).transpose()
- return buff
-
- def _oblique_pixelize(self, data_source, field, bounds, size, antialias):
- indices = np.argsort(data_source['dx'])[::-1]
- buff = _MPL.CPixelize(data_source['x'], data_source['y'],
- data_source['z'], data_source['px'],
- data_source['py'], data_source['pdx'],
- data_source['pdy'], data_source['pdz'],
- data_source.center, data_source._inv_mat, indices,
- data_source[field], size[0], size[1], bounds).transpose()
- return buff
-
- def convert_from_cartesian(self, coord):
- return coord
-
- def convert_to_cartesian(self, coord):
- return coord
-
- def convert_to_cylindrical(self, coord):
- center = self.ds.domain_center
- return cartesian_to_cylindrical(coord, center)
-
- def convert_from_cylindrical(self, coord):
- center = self.ds.domain_center
- return cylindrical_to_cartesian(coord, center)
-
- def convert_to_spherical(self, coord):
- raise NotImplementedError
-
- def convert_from_spherical(self, coord):
- raise NotImplementedError
-
- # Despite being mutables, we uses these here to be clear about how these
- # are generated and to ensure that they are not re-generated unnecessarily
- axis_name = { 0 : 'x', 1 : 'y', 2 : 'z',
- 'x' : 'x', 'y' : 'y', 'z' : 'z',
- 'X' : 'x', 'Y' : 'y', 'Z' : 'z'}
-
- axis_id = { 'x' : 0, 'y' : 1, 'z' : 2,
- 0 : 0, 1 : 1, 2 : 2}
-
- x_axis = { 'x' : 1, 'y' : 2, 'z' : 0,
- 0 : 1, 1 : 2, 2 : 0}
-
- y_axis = { 'x' : 2, 'y' : 0, 'z' : 1,
- 0 : 2, 1 : 0, 2 : 1}
-
- @property
- def period(self):
- return self.ds.domain_width
-
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinate_handler.py
--- a/yt/geometry/coordinate_handler.py
+++ /dev/null
@@ -1,116 +0,0 @@
-"""
-Coordinate handler base class.
-
-
-
-
-"""
-
-#-----------------------------------------------------------------------------
-# Copyright (c) 2013, yt Development Team.
-#
-# Distributed under the terms of the Modified BSD License.
-#
-# The full license is in the file COPYING.txt, distributed with this software.
-#-----------------------------------------------------------------------------
-
-import numpy as np
-import abc
-import weakref
-
-from yt.funcs import *
-from yt.fields.field_info_container import \
- NullFunc, FieldInfoContainer
-from yt.utilities.io_handler import io_registry
-from yt.utilities.logger import ytLogger as mylog
-from yt.utilities.parallel_tools.parallel_analysis_interface import \
- ParallelAnalysisInterface
-from yt.utilities.lib.misc_utilities import \
- pixelize_cylinder
-import yt.visualization._MPL as _MPL
-
-def _unknown_coord(field, data):
- raise YTCoordinateNotImplemented
-
-def _get_coord_fields(axi, units = "code_length"):
- def _dds(field, data):
- rv = data.ds.arr(data.fwidth[...,axi].copy(), units)
- return data._reshape_vals(rv)
- def _coords(field, data):
- rv = data.ds.arr(data.fcoords[...,axi].copy(), units)
- return data._reshape_vals(rv)
- return _dds, _coords
-
-class CoordinateHandler(object):
-
- def __init__(self, ds):
- self.ds = weakref.proxy(ds)
-
- def setup_fields(self):
- # This should return field definitions for x, y, z, r, theta, phi
- raise NotImplementedError
-
- def pixelize(self, dimension, data_source, field, bounds, size, antialias = True):
- # This should *actually* be a pixelize call, not just returning the
- # pixelizer
- raise NotImplementedError
-
- def distance(self, start, end):
- p1 = self.convert_to_cartesian(start)
- p2 = self.convert_to_cartesian(end)
- return np.sqrt(((p1-p2)**2.0).sum())
-
- def convert_from_cartesian(self, coord):
- raise NotImplementedError
-
- def convert_to_cartesian(self, coord):
- raise NotImplementedError
-
- def convert_to_cylindrical(self, coord):
- raise NotImplementedError
-
- def convert_from_cylindrical(self, coord):
- raise NotImplementedError
-
- def convert_to_spherical(self, coord):
- raise NotImplementedError
-
- def convert_from_spherical(self, coord):
- raise NotImplementedError
-
- @property
- def axis_name(self):
- raise NotImplementedError
-
- @property
- def axis_id(self):
- raise NotImplementedError
-
- @property
- def x_axis(self):
- raise NotImplementedError
-
- @property
- def y_axis(self):
- raise NotImplementedError
-
- @property
- def period(self):
- raise NotImplementedError
-
-def cartesian_to_cylindrical(coord, center = (0,0,0)):
- c2 = np.zeros_like(coord)
- c2[...,0] = ((coord[...,0] - center[0])**2.0
- + (coord[...,1] - center[1])**2.0)**0.5
- c2[...,1] = coord[...,2] # rzt
- c2[...,2] = np.arctan2(coord[...,1] - center[1],
- coord[...,0] - center[0])
- return c2
-
-def cylindrical_to_cartesian(coord, center = (0,0,0)):
- c2 = np.zeros_like(coord)
- c2[...,0] = np.cos(coord[...,0]) * coord[...,1] + center[0]
- c2[...,1] = np.sin(coord[...,0]) * coord[...,1] + center[1]
- c2[...,2] = coord[...,2]
- return c2
-
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinates/api.py
--- /dev/null
+++ b/yt/geometry/coordinates/api.py
@@ -0,0 +1,29 @@
+"""
+API for coordinate handlers
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2014, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+from .coordinate_handler import \
+ CoordinateHandler
+
+from .cartesian_coordinates import \
+ CartesianCoordinateHandler
+from .polar_coordinates import \
+ PolarCoordinateHandler
+from .cylindrical_coordinates import \
+ CylindricalCoordinateHandler
+from .spherical_coordinates import \
+ SphericalCoordinateHandler
+from .geographic_coordinates import \
+ GeographicCoordinateHandler
+from .spec_cube_coordinates import \
+ SpectralCubeCoordinateHandler
+
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinates/cartesian_coordinates.py
--- /dev/null
+++ b/yt/geometry/coordinates/cartesian_coordinates.py
@@ -0,0 +1,122 @@
+"""
+Cartesian fields
+
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import numpy as np
+from .coordinate_handler import \
+ CoordinateHandler, \
+ _unknown_coord, \
+ _get_coord_fields
+import yt.visualization._MPL as _MPL
+
+class CartesianCoordinateHandler(CoordinateHandler):
+
+ def __init__(self, ds):
+ super(CartesianCoordinateHandler, self).__init__(ds)
+
+ def setup_fields(self, registry):
+ for axi, ax in enumerate('xyz'):
+ f1, f2 = _get_coord_fields(axi)
+ registry.add_field(("index", "d%s" % ax), function = f1,
+ display_field = False,
+ units = "code_length")
+ registry.add_field(("index", "%s" % ax), function = f2,
+ display_field = False,
+ units = "code_length")
+ def _cell_volume(field, data):
+ rv = data["index", "dx"].copy(order='K')
+ rv *= data["index", "dy"]
+ rv *= data["index", "dz"]
+ return rv
+ registry.add_field(("index", "cell_volume"), function=_cell_volume,
+ display_field=False, units = "code_length**3")
+ registry.check_derived_fields(
+ [("index", "dx"), ("index", "dy"), ("index", "dz"),
+ ("index", "x"), ("index", "y"), ("index", "z"),
+ ("index", "cell_volume")])
+
+ def pixelize(self, dimension, data_source, field, bounds, size,
+ antialias = True, periodic = True):
+ if dimension < 3:
+ return self._ortho_pixelize(data_source, field, bounds, size,
+ antialias, dimension, periodic)
+ else:
+ return self._oblique_pixelize(data_source, field, bounds, size,
+ antialias)
+
+ def _ortho_pixelize(self, data_source, field, bounds, size, antialias,
+ dim, periodic):
+ # We should be using fcoords
+ period = self.period[:2].copy() # dummy here
+ period[0] = self.period[self.x_axis[dim]]
+ period[1] = self.period[self.y_axis[dim]]
+ if hasattr(period, 'in_units'):
+ period = period.in_units("code_length").d
+ buff = _MPL.Pixelize(data_source['px'], data_source['py'],
+ data_source['pdx'], data_source['pdy'],
+ data_source[field], size[0], size[1],
+ bounds, int(antialias),
+ period, int(periodic)).transpose()
+ return buff
+
+ def _oblique_pixelize(self, data_source, field, bounds, size, antialias):
+ indices = np.argsort(data_source['dx'])[::-1]
+ buff = _MPL.CPixelize(data_source['x'], data_source['y'],
+ data_source['z'], data_source['px'],
+ data_source['py'], data_source['pdx'],
+ data_source['pdy'], data_source['pdz'],
+ data_source.center, data_source._inv_mat, indices,
+ data_source[field], size[0], size[1], bounds).transpose()
+ return buff
+
+ def convert_from_cartesian(self, coord):
+ return coord
+
+ def convert_to_cartesian(self, coord):
+ return coord
+
+ def convert_to_cylindrical(self, coord):
+ center = self.ds.domain_center
+ return cartesian_to_cylindrical(coord, center)
+
+ def convert_from_cylindrical(self, coord):
+ center = self.ds.domain_center
+ return cylindrical_to_cartesian(coord, center)
+
+ def convert_to_spherical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_spherical(self, coord):
+ raise NotImplementedError
+
+ # Despite being mutables, we uses these here to be clear about how these
+ # are generated and to ensure that they are not re-generated unnecessarily
+ axis_name = { 0 : 'x', 1 : 'y', 2 : 'z',
+ 'x' : 'x', 'y' : 'y', 'z' : 'z',
+ 'X' : 'x', 'Y' : 'y', 'Z' : 'z'}
+
+ axis_id = { 'x' : 0, 'y' : 1, 'z' : 2,
+ 0 : 0, 1 : 1, 2 : 2}
+
+ x_axis = { 'x' : 1, 'y' : 2, 'z' : 0,
+ 0 : 1, 1 : 2, 2 : 0}
+
+ y_axis = { 'x' : 2, 'y' : 0, 'z' : 1,
+ 0 : 2, 1 : 0, 2 : 1}
+
+ @property
+ def period(self):
+ return self.ds.domain_width
+
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinates/coordinate_handler.py
--- /dev/null
+++ b/yt/geometry/coordinates/coordinate_handler.py
@@ -0,0 +1,212 @@
+"""
+Coordinate handler base class.
+
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import numpy as np
+import abc
+import weakref
+from numbers import Number
+
+from yt.funcs import *
+from yt.fields.field_info_container import \
+ NullFunc, FieldInfoContainer
+from yt.utilities.io_handler import io_registry
+from yt.utilities.logger import ytLogger as mylog
+from yt.utilities.parallel_tools.parallel_analysis_interface import \
+ ParallelAnalysisInterface
+from yt.utilities.lib.misc_utilities import \
+ pixelize_cylinder
+import yt.visualization._MPL as _MPL
+from yt.units.yt_array import \
+ YTArray, YTQuantity
+
+def _unknown_coord(field, data):
+ raise YTCoordinateNotImplemented
+
+def _get_coord_fields(axi, units = "code_length"):
+ def _dds(field, data):
+ rv = data.ds.arr(data.fwidth[...,axi].copy(), units)
+ return data._reshape_vals(rv)
+ def _coords(field, data):
+ rv = data.ds.arr(data.fcoords[...,axi].copy(), units)
+ return data._reshape_vals(rv)
+ return _dds, _coords
+
+def validate_iterable_width(width, ds, unit=None):
+ if isinstance(width[0], tuple) and isinstance(width[1], tuple):
+ validate_width_tuple(width[0])
+ validate_width_tuple(width[1])
+ return (ds.quan(width[0][0], fix_unitary(width[0][1])),
+ ds.quan(width[1][0], fix_unitary(width[1][1])))
+ elif isinstance(width[0], Number) and isinstance(width[1], Number):
+ return (ds.quan(width[0], 'code_length'),
+ ds.quan(width[1], 'code_length'))
+ elif isinstance(width[0], YTQuantity) and isinstance(width[1], YTQuantity):
+ return (ds.quan(width[0]), ds.quan(width[1]))
+ else:
+ validate_width_tuple(width)
+ # If width and unit are both valid width tuples, we
+ # assume width controls x and unit controls y
+ try:
+ validate_width_tuple(unit)
+ return (ds.quan(width[0], fix_unitary(width[1])),
+ ds.quan(unit[0], fix_unitary(unit[1])))
+ except YTInvalidWidthError:
+ return (ds.quan(width[0], fix_unitary(width[1])),
+ ds.quan(width[0], fix_unitary(width[1])))
+
+class CoordinateHandler(object):
+
+ def __init__(self, ds):
+ self.ds = weakref.proxy(ds)
+
+ def setup_fields(self):
+ # This should return field definitions for x, y, z, r, theta, phi
+ raise NotImplementedError
+
+ def pixelize(self, dimension, data_source, field, bounds, size, antialias = True):
+ # This should *actually* be a pixelize call, not just returning the
+ # pixelizer
+ raise NotImplementedError
+
+ def distance(self, start, end):
+ p1 = self.convert_to_cartesian(start)
+ p2 = self.convert_to_cartesian(end)
+ return np.sqrt(((p1-p2)**2.0).sum())
+
+ def convert_from_cartesian(self, coord):
+ raise NotImplementedError
+
+ def convert_to_cartesian(self, coord):
+ raise NotImplementedError
+
+ def convert_to_cylindrical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_cylindrical(self, coord):
+ raise NotImplementedError
+
+ def convert_to_spherical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_spherical(self, coord):
+ raise NotImplementedError
+
+ @property
+ def axis_name(self):
+ raise NotImplementedError
+
+ @property
+ def image_axis_name(self):
+ # Default
+ rv = {}
+ for i in range(3):
+ rv[i] = (self.axis_name[self.x_axis[i]],
+ self.axis_name[self.y_axis[i]])
+ rv[self.axis_name[i]] = rv[i]
+ rv[self.axis_name[i].upper()] = rv[i]
+ return rv
+
+ @property
+ def axis_id(self):
+ raise NotImplementedError
+
+ @property
+ def x_axis(self):
+ raise NotImplementedError
+
+ @property
+ def y_axis(self):
+ raise NotImplementedError
+
+ @property
+ def period(self):
+ raise NotImplementedError
+
+ def sanitize_width(self, axis, width, depth):
+ if width is None:
+ # Default to code units
+ if not iterable(axis):
+ xax = self.x_axis[axis]
+ yax = self.y_axis[axis]
+ w = self.ds.domain_width[[xax, yax]]
+ else:
+ # axis is actually the normal vector
+ # for an off-axis data object.
+ mi = np.argmin(self.ds.domain_width)
+ w = self.ds.domain_width[[mi,mi]]
+ width = (w[0], w[1])
+ elif iterable(width):
+ width = validate_iterable_width(width, self.ds)
+ elif isinstance(width, YTQuantity):
+ width = (width, width)
+ elif isinstance(width, Number):
+ width = (self.ds.quan(width, 'code_length'),
+ self.ds.quan(width, 'code_length'))
+ else:
+ raise YTInvalidWidthError(width)
+ if depth is not None:
+ if iterable(depth):
+ validate_width_tuple(depth)
+ depth = (self.ds.quan(depth[0], fix_unitary(depth[1])), )
+ elif isinstance(depth, Number):
+ depth = (self.ds.quan(depth, 'code_length',
+ registry = self.ds.unit_registry), )
+ elif isinstance(depth, YTQuantity):
+ depth = (depth, )
+ else:
+ raise YTInvalidWidthError(depth)
+ return width + depth
+ return width
+
+ def sanitize_center(self, center, axis):
+ if isinstance(center, basestring):
+ if center.lower() == "m" or center.lower() == "max":
+ v, center = self.ds.find_max(("gas", "density"))
+ center = self.ds.arr(center, 'code_length')
+ elif center.lower() == "c" or center.lower() == "center":
+ center = (self.ds.domain_left_edge + self.ds.domain_right_edge) / 2
+ else:
+ raise RuntimeError('center keyword \"%s\" not recognized' % center)
+ elif isinstance(center, YTArray):
+ return self.ds.arr(center), self.convert_to_cartesian(center)
+ elif iterable(center):
+ if iterable(center[0]) and isinstance(center[1], basestring):
+ center = self.ds.arr(center[0], center[1])
+ else:
+ center = self.ds.arr(center, 'code_length')
+ else:
+ raise RuntimeError("center keyword \"%s\" not recognized" % center)
+ # This has to return both a center and a display_center
+ display_center = self.convert_to_cartesian(center)
+ return center, display_center
+
+
+def cartesian_to_cylindrical(coord, center = (0,0,0)):
+ c2 = np.zeros_like(coord)
+ c2[...,0] = ((coord[...,0] - center[0])**2.0
+ + (coord[...,1] - center[1])**2.0)**0.5
+ c2[...,1] = coord[...,2] # rzt
+ c2[...,2] = np.arctan2(coord[...,1] - center[1],
+ coord[...,0] - center[0])
+ return c2
+
+def cylindrical_to_cartesian(coord, center = (0,0,0)):
+ c2 = np.zeros_like(coord)
+ c2[...,0] = np.cos(coord[...,0]) * coord[...,1] + center[0]
+ c2[...,1] = np.sin(coord[...,0]) * coord[...,1] + center[1]
+ c2[...,2] = coord[...,2]
+ return c2
+
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinates/cylindrical_coordinates.py
--- /dev/null
+++ b/yt/geometry/coordinates/cylindrical_coordinates.py
@@ -0,0 +1,162 @@
+"""
+Cylindrical fields
+
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import numpy as np
+from yt.units.yt_array import YTArray
+from .coordinate_handler import \
+ CoordinateHandler, \
+ _unknown_coord, \
+ _get_coord_fields
+import yt.visualization._MPL as _MPL
+from yt.utilities.lib.misc_utilities import \
+ pixelize_cylinder
+#
+# Cylindrical fields
+#
+
+class CylindricalCoordinateHandler(CoordinateHandler):
+
+ def __init__(self, ds, ordering = 'rzt'):
+ if ordering != 'rzt': raise NotImplementedError
+ super(CylindricalCoordinateHandler, self).__init__(ds)
+
+ def setup_fields(self, registry):
+ # return the fields for r, z, theta
+ registry.add_field(("index", "dx"), function=_unknown_coord)
+ registry.add_field(("index", "dy"), function=_unknown_coord)
+ registry.add_field(("index", "x"), function=_unknown_coord)
+ registry.add_field(("index", "y"), function=_unknown_coord)
+ f1, f2 = _get_coord_fields(0)
+ registry.add_field(("index", "dr"), function = f1,
+ display_field = False,
+ units = "code_length")
+ registry.add_field(("index", "r"), function = f2,
+ display_field = False,
+ units = "code_length")
+
+ f1, f2 = _get_coord_fields(1)
+ registry.add_field(("index", "dz"), function = f1,
+ display_field = False,
+ units = "code_length")
+ registry.add_field(("index", "z"), function = f2,
+ display_field = False,
+ units = "code_length")
+
+ f1, f2 = _get_coord_fields(2, "")
+ registry.add_field(("index", "dtheta"), function = f1,
+ display_field = False,
+ units = "")
+ registry.add_field(("index", "theta"), function = f2,
+ display_field = False,
+ units = "")
+
+ def _CylindricalVolume(field, data):
+ return data["index", "dtheta"] \
+ * data["index", "r"] \
+ * data["index", "dr"] \
+ * data["index", "dz"]
+ registry.add_field(("index", "cell_volume"),
+ function=_CylindricalVolume,
+ units = "code_length**3")
+
+
+ def pixelize(self, dimension, data_source, field, bounds, size,
+ antialias = True, periodic = True):
+ ax_name = self.axis_name[dimension]
+ if ax_name in ('r', 'theta'):
+ return self._ortho_pixelize(data_source, field, bounds, size,
+ antialias, dimension, periodic)
+ elif ax_name == "z":
+ return self._cyl_pixelize(data_source, field, bounds, size,
+ antialias)
+ else:
+ # Pixelizing along a cylindrical surface is a bit tricky
+ raise NotImplementedError
+
+ def _ortho_pixelize(self, data_source, field, bounds, size, antialias,
+ dim, periodic):
+ period = self.period[:2].copy() # dummy here
+ period[0] = self.period[self.x_axis[dim]]
+ period[1] = self.period[self.y_axis[dim]]
+ if hasattr(period, 'in_units'):
+ period = period.in_units("code_length").d
+ buff = _MPL.Pixelize(data_source['px'], data_source['py'],
+ data_source['pdx'], data_source['pdy'],
+ data_source[field], size[0], size[1],
+ bounds, int(antialias),
+ period, int(periodic)).transpose()
+ return buff
+
+ def _cyl_pixelize(self, data_source, field, bounds, size, antialias):
+ buff = pixelize_cylinder(data_source['r'],
+ data_source['dr'],
+ data_source['theta'],
+ data_source['dtheta']/2.0, # half-widths
+ size, data_source[field], bounds)
+ return buff
+
+ axis_name = { 0 : 'r', 1 : 'z', 2 : 'theta',
+ 'r' : 'r', 'z' : 'z', 'theta' : 'theta',
+ 'R' : 'r', 'Z' : 'z', 'Theta' : 'theta'}
+
+ axis_id = { 'r' : 0, 'z' : 1, 'theta' : 2,
+ 0 : 0, 1 : 1, 2 : 2}
+
+ x_axis = { 'r' : 1, 'z' : 0, 'theta' : 0,
+ 0 : 1, 1 : 0, 2 : 0}
+
+ y_axis = { 'r' : 2, 'z' : 2, 'theta' : 1,
+ 0 : 2, 1 : 2, 2 : 1}
+
+ _image_axis_name = None
+
+ @property
+ def image_axis_name(self):
+ if self._image_axis_name is not None:
+ return self._image_axis_name
+ # This is the x and y axes labels that get displayed. For
+ # non-Cartesian coordinates, we usually want to override these for
+ # Cartesian coordinates, since we transform them.
+ rv = {0: ('z', 'theta'),
+ 1: ('x', 'y'),
+ 2: ('r', 'z')}
+ for i in rv.keys():
+ rv[self.axis_name[i]] = rv[i]
+ rv[self.axis_name[i].upper()] = rv[i]
+ self._image_axis_name = rv
+ return rv
+
+ def convert_from_cartesian(self, coord):
+ return cartesian_to_cylindrical(coord)
+
+ def convert_to_cartesian(self, coord):
+ return cylindrical_to_cartesian(coord)
+
+ def convert_to_cylindrical(self, coord):
+ return coord
+
+ def convert_from_cylindrical(self, coord):
+ return coord
+
+ def convert_to_spherical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_spherical(self, coord):
+ raise NotImplementedError
+
+ @property
+ def period(self):
+ return np.array([0.0, 0.0, 2.0*np.pi])
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinates/geographic_coordinates.py
--- /dev/null
+++ b/yt/geometry/coordinates/geographic_coordinates.py
@@ -0,0 +1,194 @@
+"""
+Geographic fields
+
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import numpy as np
+from .coordinate_handler import \
+ CoordinateHandler, \
+ _unknown_coord, \
+ _get_coord_fields
+import yt.visualization._MPL as _MPL
+from yt.utilities.lib.misc_utilities import \
+ pixelize_cylinder, pixelize_aitoff
+
+class GeographicCoordinateHandler(CoordinateHandler):
+
+ def __init__(self, ds, ordering = 'latlonalt'):
+ if ordering != 'latlonalt': raise NotImplementedError
+ super(GeographicCoordinateHandler, self).__init__(ds)
+
+ def setup_fields(self, registry):
+ # return the fields for r, z, theta
+ registry.add_field(("index", "dx"), function=_unknown_coord)
+ registry.add_field(("index", "dy"), function=_unknown_coord)
+ registry.add_field(("index", "dz"), function=_unknown_coord)
+ registry.add_field(("index", "x"), function=_unknown_coord)
+ registry.add_field(("index", "y"), function=_unknown_coord)
+ registry.add_field(("index", "z"), function=_unknown_coord)
+ f1, f2 = _get_coord_fields(0, "")
+ registry.add_field(("index", "dlatitude"), function = f1,
+ display_field = False,
+ units = "")
+ registry.add_field(("index", "latitude"), function = f2,
+ display_field = False,
+ units = "")
+
+ f1, f2 = _get_coord_fields(1, "")
+ registry.add_field(("index", "dlongitude"), function = f1,
+ display_field = False,
+ units = "")
+ registry.add_field(("index", "longitude"), function = f2,
+ display_field = False,
+ units = "")
+
+ f1, f2 = _get_coord_fields(2)
+ registry.add_field(("index", "daltitude"), function = f1,
+ display_field = False,
+ units = "code_length")
+ registry.add_field(("index", "altitude"), function = f2,
+ display_field = False,
+ units = "code_length")
+
+ def _SphericalVolume(field, data):
+ # r**2 sin theta dr dtheta dphi
+ # We can use the transformed coordinates here.
+ vol = data["index", "r"]**2.0
+ vol *= data["index", "dr"]
+ vol *= np.sin(data["index", "theta"])
+ vol *= data["index", "dtheta"]
+ vol *= data["index", "dphi"]
+ return vol
+ registry.add_field(("index", "cell_volume"),
+ function=_SphericalVolume,
+ units = "code_length**3")
+
+ # Altitude is the radius from the central zone minus the radius of the
+ # surface.
+ def _altitude_to_radius(field, data):
+ surface_height = data.get_field_parameter("surface_height")
+ if surface_height is None:
+ surface_height = getattr(data.ds, "surface_height", 0.0)
+ return data["altitude"] + surface_height
+ registry.add_field(("index", "r"),
+ function=_altitude_to_radius,
+ units = "code_length")
+ registry.alias(("index", "dr"), ("index", "daltitude"))
+
+ def _longitude_to_theta(field, data):
+ # longitude runs from -180 to 180.
+ return (data["longitude"] + 180) * np.pi/180.0
+ registry.add_field(("index", "theta"),
+ function = _longitude_to_theta,
+ units = "")
+ def _dlongitude_to_dtheta(field, data):
+ return data["dlongitude"] * np.pi/180.0
+ registry.add_field(("index", "dtheta"),
+ function = _dlongitude_to_dtheta,
+ units = "")
+
+ def _latitude_to_phi(field, data):
+ # latitude runs from -90 to 90
+ return (data["latitude"] + 90) * np.pi/180.0
+ registry.add_field(("index", "phi"),
+ function = _latitude_to_phi,
+ units = "")
+ def _dlatitude_to_dphi(field, data):
+ return data["dlatitude"] * np.pi/180.0
+ registry.add_field(("index", "dphi"),
+ function = _dlatitude_to_dphi,
+ units = "")
+
+ def pixelize(self, dimension, data_source, field, bounds, size,
+ antialias = True, periodic = True):
+ if dimension in (0, 1):
+ return self._cyl_pixelize(data_source, field, bounds, size,
+ antialias, dimension)
+ elif dimension == 2:
+ return self._ortho_pixelize(data_source, field, bounds, size,
+ antialias, dimension, periodic)
+ else:
+ raise NotImplementedError
+
+ def _ortho_pixelize(self, data_source, field, bounds, size, antialias,
+ dim, periodic):
+ buff = pixelize_aitoff(data_source["theta"], data_source["dtheta"]/2.0,
+ data_source["phi"], data_source["dphi"]/2.0,
+ size, data_source[field], None,
+ None).transpose()
+ return buff
+
+ def _cyl_pixelize(self, data_source, field, bounds, size, antialias,
+ dimension):
+ if dimension == 0:
+ buff = pixelize_cylinder(data_source['r'],
+ data_source['dr'] / 2.0,
+ data_source['theta'],
+ data_source['dtheta'] / 2.0, # half-widths
+ size, data_source[field], bounds)
+ elif dimension == 1:
+ buff = pixelize_cylinder(data_source['r'],
+ data_source['dr'] / 2.0,
+ data_source['phi'],
+ data_source['dphi'] / 2.0, # half-widths
+ size, data_source[field], bounds)
+ else:
+ raise RuntimeError
+ return buff
+
+
+ def convert_from_cartesian(self, coord):
+ raise NotImplementedError
+
+ def convert_to_cartesian(self, coord):
+ raise NotImplementedError
+
+ def convert_to_cylindrical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_cylindrical(self, coord):
+ raise NotImplementedError
+
+ def convert_to_spherical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_spherical(self, coord):
+ raise NotImplementedError
+
+ # Despite being mutables, we uses these here to be clear about how these
+ # are generated and to ensure that they are not re-generated unnecessarily
+ axis_name = { 0 : 'latitude', 1 : 'longitude', 2 : 'altitude',
+ 'latitude' : 'latitude',
+ 'longitude' : 'longitude',
+ 'altitude' : 'altitude',
+ 'Latitude' : 'latitude',
+ 'Longitude' : 'longitude',
+ 'Altitude' : 'altitude',
+ 'lat' : 'latitude',
+ 'lon' : 'longitude',
+ 'alt' : 'altitude' }
+
+ axis_id = { 'latitude' : 0, 'longitude' : 1, 'altitude' : 2,
+ 0 : 0, 1 : 1, 2 : 2}
+
+ x_axis = { 'latitude' : 1, 'longitude' : 0, 'altitude' : 0,
+ 0 : 1, 1 : 0, 2 : 0}
+
+ y_axis = { 'latitude' : 2, 'longitude' : 2, 'altitude' : 1,
+ 0 : 2, 1 : 2, 2 : 1}
+
+ @property
+ def period(self):
+ return self.ds.domain_width
+
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinates/polar_coordinates.py
--- /dev/null
+++ b/yt/geometry/coordinates/polar_coordinates.py
@@ -0,0 +1,149 @@
+"""
+Polar fields
+
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import numpy as np
+from .coordinate_handler import \
+ CoordinateHandler, \
+ _unknown_coord, \
+ _get_coord_fields
+from yt.utilities.lib.misc_utilities import \
+ pixelize_cylinder
+
+class PolarCoordinateHandler(CoordinateHandler):
+
+ def __init__(self, ds, ordering = 'rtz'):
+ if ordering != 'rtz': raise NotImplementedError
+ super(PolarCoordinateHandler, self).__init__(ds)
+
+ def setup_fields(self, registry):
+ # return the fields for r, z, theta
+ registry.add_field("dx", function=_unknown_coord)
+ registry.add_field("dy", function=_unknown_coord)
+ registry.add_field("x", function=_unknown_coord)
+ registry.add_field("y", function=_unknown_coord)
+
+ f1, f2 = _get_coord_fields(0)
+ registry.add_field(("index", "dr"), function = f1,
+ display_field = False,
+ units = "code_length")
+ registry.add_field(("index", "r"), function = f2,
+ display_field = False,
+ units = "code_length")
+
+ f1, f2 = _get_coord_fields(1, "")
+ registry.add_field(("index", "dtheta"), function = f1,
+ display_field = False,
+ units = "")
+ registry.add_field(("index", "theta"), function = f2,
+ display_field = False,
+ units = "")
+
+ f1, f2 = _get_coord_fields(2)
+ registry.add_field(("index", "dz"), function = f1,
+ display_field = False,
+ units = "code_length")
+ registry.add_field(("index", "z"), function = f2,
+ display_field = False,
+ units = "code_length")
+
+
+ def _CylindricalVolume(field, data):
+ return data["dtheta"] * data["r"] * data["dr"] * data["dz"]
+ registry.add_field("CellVolume", function=_CylindricalVolume)
+
+ def pixelize(self, dimension, data_source, field, bounds, size, antialias = True):
+ ax_name = self.axis_name[dimension]
+ if ax_name in ('r', 'theta'):
+ return self._ortho_pixelize(data_source, field, bounds, size,
+ antialias)
+ elif ax_name == "z":
+ return self._polar_pixelize(data_source, field, bounds, size,
+ antialias)
+ else:
+ # Pixelizing along a cylindrical surface is a bit tricky
+ raise NotImplementedError
+
+
+ def _ortho_pixelize(self, data_source, field, bounds, size, antialias):
+ buff = _MPL.Pixelize(data_source['px'], data_source['py'],
+ data_source['pdx'], data_source['pdy'],
+ data_source[field], size[0], size[1],
+ bounds, int(antialias),
+ True, self.period).transpose()
+ return buff
+
+ def _polar_pixelize(self, data_source, field, bounds, size, antialias):
+ # Out bounds here will *always* be what plot window thinks are x0, x1,
+ # y0, y1, but which will actually be rmin, rmax, thetamin, thetamax.
+ buff = pixelize_cylinder(data_source['r'],
+ data_source['dr'],
+ data_source['theta'],
+ data_source['dtheta'] / 2.0, # half-widths
+ size, data_source[field], bounds)
+ return buff
+
+ axis_name = { 0 : 'r', 1 : 'theta', 2 : 'z',
+ 'r' : 'r', 'theta' : 'theta', 'z' : 'z',
+ 'R' : 'r', 'Theta' : 'theta', 'Z' : 'z'}
+
+ axis_id = { 'r' : 0, 'theta' : 1, 'z' : 2,
+ 0 : 0, 1 : 1, 2 : 2}
+
+ x_axis = { 'r' : 1, 'theta' : 0, 'z' : 0,
+ 0 : 1, 1 : 0, 2 : 0}
+
+ y_axis = { 'r' : 2, 'theta' : 2, 'z' : 1,
+ 0 : 2, 1 : 2, 2 : 1}
+
+ _image_axis_name = None
+ @property
+ def image_axis_name(self):
+ if self._image_axis_name is not None:
+ return self._image_axis_name
+ # This is the x and y axes labels that get displayed. For
+ # non-Cartesian coordinates, we usually want to override these for
+ # Cartesian coordinates, since we transform them.
+ rv = {0: ('theta', 'z'),
+ 1: ('x', 'y'),
+ 2: ('r', 'z')}
+ for i in rv.keys():
+ rv[self.axis_name[i]] = rv[i]
+ rv[self.axis_name[i].upper()] = rv[i]
+ self._image_axis_name = rv
+ return rv
+
+ def convert_from_cartesian(self, coord):
+ return cartesian_to_cylindrical(coord)
+
+ def convert_to_cartesian(self, coord):
+ return cylindrical_to_cartesian(coord)
+
+ def convert_to_cylindrical(self, coord):
+ return coord
+
+ def convert_from_cylindrical(self, coord):
+ return coord
+
+ def convert_to_spherical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_spherical(self, coord):
+ raise NotImplementedError
+
+ @property
+ def period(self):
+ return np.array([0.0, 0.0, 2.0*np.pi])
+
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinates/spec_cube_coordinates.py
--- /dev/null
+++ b/yt/geometry/coordinates/spec_cube_coordinates.py
@@ -0,0 +1,65 @@
+"""
+Cartesian fields
+
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import numpy as np
+from .cartesian_coordinates import \
+ CartesianCoordinateHandler
+
+class SpectralCubeCoordinateHandler(CartesianCoordinateHandler):
+
+ def __init__(self, ds):
+ super(SpectralCubeCoordinateHandler, self).__init__(ds)
+
+ self.axis_name = {}
+ self.axis_id = {}
+
+ for axis, axis_name in zip([ds.lon_axis, ds.lat_axis, ds.spec_axis],
+ ["Image\ x", "Image\ y", ds.spec_name]):
+ lower_ax = "xyz"[axis]
+ upper_ax = lower_ax.upper()
+
+ self.axis_name[axis] = axis_name
+ self.axis_name[lower_ax] = axis_name
+ self.axis_name[upper_ax] = axis_name
+ self.axis_name[axis_name] = axis_name
+
+ self.axis_id[lower_ax] = axis
+ self.axis_id[axis] = axis
+ self.axis_id[axis_name] = axis
+
+ self.default_unit_label = {}
+ self.default_unit_label[ds.lon_axis] = "pixel"
+ self.default_unit_label[ds.lat_axis] = "pixel"
+ self.default_unit_label[ds.spec_axis] = ds.spec_unit
+
+ def _spec_axis(ax, x, y):
+ p = (x,y)[ax]
+ return [self.ds.pixel2spec(pp).v for pp in p]
+
+ self.axis_field = {}
+ self.axis_field[self.ds.spec_axis] = _spec_axis
+
+ def convert_to_cylindrical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_cylindrical(self, coord):
+ raise NotImplementedError
+
+ x_axis = { 'x' : 1, 'y' : 0, 'z' : 0,
+ 0 : 1, 1 : 0, 2 : 0}
+
+ y_axis = { 'x' : 2, 'y' : 2, 'z' : 1,
+ 0 : 2, 1 : 2, 2 : 1}
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/coordinates/spherical_coordinates.py
--- /dev/null
+++ b/yt/geometry/coordinates/spherical_coordinates.py
@@ -0,0 +1,221 @@
+"""
+Spherical fields
+
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import numpy as np
+from .coordinate_handler import \
+ CoordinateHandler, \
+ _unknown_coord, \
+ _get_coord_fields
+import yt.visualization._MPL as _MPL
+from yt.utilities.lib.misc_utilities import \
+ pixelize_cylinder, pixelize_aitoff
+
+class SphericalCoordinateHandler(CoordinateHandler):
+
+ def __init__(self, ds, ordering = 'rtp'):
+ if ordering != 'rtp': raise NotImplementedError
+ super(SphericalCoordinateHandler, self).__init__(ds)
+
+ def setup_fields(self, registry):
+ # return the fields for r, z, theta
+ registry.add_field(("index", "dx"), function=_unknown_coord)
+ registry.add_field(("index", "dy"), function=_unknown_coord)
+ registry.add_field(("index", "dz"), function=_unknown_coord)
+ registry.add_field(("index", "x"), function=_unknown_coord)
+ registry.add_field(("index", "y"), function=_unknown_coord)
+ registry.add_field(("index", "z"), function=_unknown_coord)
+ f1, f2 = _get_coord_fields(0)
+ registry.add_field(("index", "dr"), function = f1,
+ display_field = False,
+ units = "code_length")
+ registry.add_field(("index", "r"), function = f2,
+ display_field = False,
+ units = "code_length")
+
+ f1, f2 = _get_coord_fields(1, "")
+ registry.add_field(("index", "dtheta"), function = f1,
+ display_field = False,
+ units = "")
+ registry.add_field(("index", "theta"), function = f2,
+ display_field = False,
+ units = "")
+
+ f1, f2 = _get_coord_fields(2, "")
+ registry.add_field(("index", "dphi"), function = f1,
+ display_field = False,
+ units = "")
+ registry.add_field(("index", "phi"), function = f2,
+ display_field = False,
+ units = "")
+
+ def _SphericalVolume(field, data):
+ # r**2 sin theta dr dtheta dphi
+ vol = data["index", "r"]**2.0
+ vol *= data["index", "dr"]
+ vol *= np.sin(data["index", "theta"])
+ vol *= data["index", "dtheta"]
+ vol *= data["index", "dphi"]
+ return vol
+ registry.add_field(("index", "cell_volume"),
+ function=_SphericalVolume,
+ units = "code_length**3")
+
+ def pixelize(self, dimension, data_source, field, bounds, size,
+ antialias = True, periodic = True):
+ self.period
+ if dimension == 0:
+ return self._ortho_pixelize(data_source, field, bounds, size,
+ antialias, dimension, periodic)
+ elif dimension in (1, 2):
+ return self._cyl_pixelize(data_source, field, bounds, size,
+ antialias, dimension)
+ else:
+ raise NotImplementedError
+
+ def _ortho_pixelize(self, data_source, field, bounds, size, antialias,
+ dim, periodic):
+ # We should be using fcoords
+ period = self.period[:2].copy() # dummy here
+ period[0] = self.period[self.x_axis[dim]]
+ period[1] = self.period[self.y_axis[dim]]
+ period = period.in_units("code_length").d
+ buff = _MPL.Pixelize(data_source['px'], data_source['py'],
+ data_source['pdx'], data_source['pdy'],
+ data_source[field], size[0], size[1],
+ bounds, int(antialias),
+ period, int(periodic)).transpose()
+ return buff
+
+ def _cyl_pixelize(self, data_source, field, bounds, size, antialias,
+ dimension):
+ if dimension == 1:
+ buff = pixelize_cylinder(data_source['r'],
+ data_source['dr'] / 2.0,
+ data_source['phi'],
+ data_source['dphi'] / 2.0, # half-widths
+ size, data_source[field], bounds)
+ elif dimension == 2:
+ buff = pixelize_cylinder(data_source['r'],
+ data_source['dr'] / 2.0,
+ data_source['theta'],
+ data_source['dtheta'] / 2.0, # half-widths
+ size, data_source[field], bounds)
+ buff = buff.transpose()
+ else:
+ raise RuntimeError
+ return buff
+
+
+ def convert_from_cartesian(self, coord):
+ raise NotImplementedError
+
+ def convert_to_cartesian(self, coord):
+ if isinstance(coord, np.ndarray) and len(coord.shape) > 1:
+ r = coord[:,0]
+ theta = coord[:,1]
+ phi = coord[:,2]
+ nc = np.zeros_like(coord)
+ # r, theta, phi
+ nc[:,0] = np.cos(phi) * np.sin(theta)*r
+ nc[:,1] = np.sin(phi) * np.sin(theta)*r
+ nc[:,2] = np.cos(theta) * r
+ else:
+ r, theta, phi = coord
+ nc = (np.cos(phi) * np.sin(theta)*r,
+ np.sin(phi) * np.sin(theta)*r,
+ np.cos(theta) * r)
+ return nc
+
+ def convert_to_cylindrical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_cylindrical(self, coord):
+ raise NotImplementedError
+
+ def convert_to_spherical(self, coord):
+ raise NotImplementedError
+
+ def convert_from_spherical(self, coord):
+ raise NotImplementedError
+
+ # Despite being mutables, we uses these here to be clear about how these
+ # are generated and to ensure that they are not re-generated unnecessarily
+ axis_name = { 0 : 'r', 1 : 'theta', 2 : 'phi',
+ 'r' : 'r', 'theta' : 'theta', 'phi' : 'phi',
+ 'R' : 'r', 'Theta' : 'theta', 'Phi' : 'phi'}
+
+ _image_axis_name = None
+ @property
+ def image_axis_name(self):
+ if self._image_axis_name is not None:
+ return self._image_axis_name
+ # This is the x and y axes labels that get displayed. For
+ # non-Cartesian coordinates, we usually want to override these for
+ # Cartesian coordinates, since we transform them.
+ rv = {0: ('theta', 'phi'),
+ 1: ('x / \\sin(\\theta)', 'y / \\sin(\\theta)'),
+ 2: ('R', 'z')}
+ for i in rv.keys():
+ rv[self.axis_name[i]] = rv[i]
+ rv[self.axis_name[i].upper()] = rv[i]
+ self._image_axis_name = rv
+ return rv
+
+
+ axis_id = { 'r' : 0, 'theta' : 1, 'phi' : 2,
+ 0 : 0, 1 : 1, 2 : 2}
+
+ x_axis = { 'r' : 1, 'theta' : 0, 'phi' : 0,
+ 0 : 1, 1 : 0, 2 : 0}
+
+ y_axis = { 'r' : 2, 'theta' : 2, 'phi' : 1,
+ 0 : 2, 1 : 2, 2 : 1}
+
+ @property
+ def period(self):
+ return self.ds.domain_width
+
+ def sanitize_center(self, center, axis):
+ center, display_center = super(
+ SphericalCoordinateHandler, self).sanitize_center(center, axis)
+ if axis == 0:
+ display_center = center
+ elif axis == 1:
+ display_center = (0.0 * display_center[0],
+ 0.0 * display_center[1],
+ 0.0 * display_center[2])
+ elif axis ==2:
+ display_center = (self.ds.domain_width[0]/2.0,
+ 0.0 * display_center[1],
+ 0.0 * display_center[2])
+ return center, display_center
+
+ def sanitize_width(self, axis, width, depth):
+ if width is not None:
+ width = super(SphericalCoordinateHandler, self).sanitize_width(
+ axis, width, depth)
+ elif axis == 0:
+ width = [self.ds.domain_width[self.x_axis[0]],
+ self.ds.domain_width[self.y_axis[0]]]
+ elif axis == 1:
+ # Remember, in spherical coordinates when we cut in theta,
+ # we create a conic section
+ width = [2.0*self.ds.domain_width[0],
+ 2.0*self.ds.domain_width[0]]
+ elif axis == 2:
+ width = [self.ds.domain_width[0],
+ 2.0*self.ds.domain_width[0]]
+ return width
diff -r 66bbf028d950fca7b523c06712c4f155bfbc15a4 -r 6f0db65e21bcabd8283cadb52e30b6c1ed9a4afc yt/geometry/cylindrical_coordinates.py
--- a/yt/geometry/cylindrical_coordinates.py
+++ /dev/null
@@ -1,145 +0,0 @@
-"""
-Cylindrical fields
-
-
-
-
-"""
-
-#-----------------------------------------------------------------------------
-# Copyright (c) 2013, yt Development Team.
-#
-# Distributed under the terms of the Modified BSD License.
-#
-# The full license is in the file COPYING.txt, distributed with this software.
-#-----------------------------------------------------------------------------
-
-import numpy as np
-from yt.units.yt_array import YTArray
-from .coordinate_handler import \
- CoordinateHandler, \
- _unknown_coord, \
- _get_coord_fields
-import yt.visualization._MPL as _MPL
-from yt.utilities.lib.misc_utilities import \
- pixelize_cylinder
-#
-# Cylindrical fields
-#
-
-class CylindricalCoordinateHandler(CoordinateHandler):
-
- def __init__(self, ds, ordering = 'rzt'):
- if ordering != 'rzt': raise NotImplementedError
- super(CylindricalCoordinateHandler, self).__init__(ds)
-
- def setup_fields(self, registry):
- # return the fields for r, z, theta
- registry.add_field(("index", "dx"), function=_unknown_coord)
- registry.add_field(("index", "dy"), function=_unknown_coord)
- registry.add_field(("index", "x"), function=_unknown_coord)
- registry.add_field(("index", "y"), function=_unknown_coord)
- f1, f2 = _get_coord_fields(0)
- registry.add_field(("index", "dr"), function = f1,
- display_field = False,
- units = "code_length")
- registry.add_field(("index", "r"), function = f2,
- display_field = False,
- units = "code_length")
-
- f1, f2 = _get_coord_fields(1)
- registry.add_field(("index", "dz"), function = f1,
- display_field = False,
- units = "code_length")
- registry.add_field(("index", "z"), function = f2,
- display_field = False,
- units = "code_length")
-
- f1, f2 = _get_coord_fields(2, "")
- registry.add_field(("index", "dtheta"), function = f1,
- display_field = False,
- units = "")
- registry.add_field(("index", "theta"), function = f2,
- display_field = False,
- units = "")
-
- def _CylindricalVolume(field, data):
- return data["index", "dtheta"] \
- * data["index", "r"] \
- * data["index", "dr"] \
- * data["index", "dz"]
- registry.add_field(("index", "cell_volume"),
- function=_CylindricalVolume,
- units = "code_length**3")
-
-
- def pixelize(self, dimension, data_source, field, bounds, size,
- antialias = True, periodic = True):
- ax_name = self.axis_name[dimension]
- if ax_name in ('r', 'theta'):
- return self._ortho_pixelize(data_source, field, bounds, size,
- antialias, dimension, periodic)
- elif ax_name == "z":
- return self._cyl_pixelize(data_source, field, bounds, size,
- antialias)
- else:
- # Pixelizing along a cylindrical surface is a bit tricky
- raise NotImplementedError
-
- def _ortho_pixelize(self, data_source, field, bounds, size, antialias,
- dim, periodic):
- period = self.period[:2].copy() # dummy here
- period[0] = self.period[self.x_axis[dim]]
- period[1] = self.period[self.y_axis[dim]]
- if hasattr(period, 'in_units'):
- period = period.in_units("code_length").d
- buff = _MPL.Pixelize(data_source['px'], data_source['py'],
- data_source['pdx'], data_source['pdy'],
- data_source[field], size[0], size[1],
- bounds, int(antialias),
- period, int(periodic)).transpose()
- return buff
-
- def _cyl_pixelize(self, data_source, field, bounds, size, antialias):
- buff = pixelize_cylinder(data_source['r'],
- data_source['dr'],
- data_source['theta'],
- data_source['dtheta']/2.0, # half-widths
- size, data_source[field], bounds)
- return buff
-
- axis_name = { 0 : 'r', 1 : 'z', 2 : 'theta',
- 'r' : 'r', 'z' : 'z', 'theta' : 'theta',
- 'R' : 'r', 'Z' : 'z', 'Theta' : 'theta'}
-
- axis_id = { 'r' : 0, 'z' : 1, 'theta' : 2,
- 0 : 0, 1 : 1, 2 : 2}
-
- x_axis = { 'r' : 1, 'z' : 0, 'theta' : 0,
- 0 : 1, 1 : 0, 2 : 0}
-
- y_axis = { 'r' : 2, 'z' : 2, 'theta' : 1,
- 0 : 2, 1 : 2, 2 : 1}
-
- def convert_from_cartesian(self, coord):
- return cartesian_to_cylindrical(coord)
-
- def convert_to_cartesian(self, coord):
- return cylindrical_to_cartesian(coord)
-
- def convert_to_cylindrical(self, coord):
- return coord
-
- def convert_from_cylindrical(self, coord):
- return coord
-
- def convert_to_spherical(self, coord):
- raise NotImplementedError
-
- def convert_from_spherical(self, coord):
- raise NotImplementedError
-
- @property
- def period(self):
- return np.array([0.0, 0.0, 2.0*np.pi])
-
This diff is so big that we needed to truncate the remainder.
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