[yt-svn] commit/yt: MatthewTurk: Merged in jisuoqing/yt (pull request #1702)

commits-noreply at bitbucket.org commits-noreply at bitbucket.org
Thu Aug 20 09:17:24 PDT 2015


1 new commit in yt:

https://bitbucket.org/yt_analysis/yt/commits/98b9f2c2d791/
Changeset:   98b9f2c2d791
Branch:      yt
User:        MatthewTurk
Date:        2015-08-20 16:17:12+00:00
Summary:     Merged in jisuoqing/yt (pull request #1702)

[NEW] [DOCS] Add line integral convolution for plot modification callback
Affected #:  6 files

diff -r 1e76952e015dfb14db5d97333c182edc658d4aaa -r 98b9f2c2d791d921eccc8e2edc419d863eb14650 doc/source/visualizing/callbacks.rst
--- a/doc/source/visualizing/callbacks.rst
+++ b/doc/source/visualizing/callbacks.rst
@@ -463,6 +463,35 @@
    s.annotate_streamlines('velocity_x', 'velocity_y')
    s.save()
 
+.. _annotate-line-integral-convolution:
+
+Overplot Line Integral Convolution
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+.. function:: annotate_line_integral_convolution(self, field_x, field_y, \
+                                                 texture=None, kernellen=50., \
+                                                 lim=(0.5,0.6), cmap='binary', \
+                                                 alpha=0.8, const_alpha=False)
+
+   (This is a proxy for
+   :class:`~yt.visualization.plot_modifications.LineIntegralConvolutionCallback`.)
+
+   Add line integral convolution to any plot, using the ``field_x`` and ``field_y`` 
+   from the associated data. A white noise background will be used for ``texture`` 
+   as default. Adjust the bounds of ``lim`` in the range of ``[0, 1]`` which applies 
+   upper and lower bounds to the values of line integral convolution and enhance 
+   the visibility of plots. When ``const_alpha=False``, alpha will be weighted 
+   spatially by the values of line integral convolution; otherwise a constant value 
+   of the given alpha is used.
+
+.. python-script::
+
+   import yt
+   ds = yt.load("IsolatedGalaxy/galaxy0030/galaxy0030")
+   s = yt.SlicePlot(ds, 'z', 'density', center='c', width=(20, 'kpc'))
+   s.annotate_line_integral_convolution('velocity_x', 'velocity_y', lim=(0.5,0.65))
+   s.save()
+
 .. _annotate-text:
 
 Overplot Text

diff -r 1e76952e015dfb14db5d97333c182edc658d4aaa -r 98b9f2c2d791d921eccc8e2edc419d863eb14650 yt/utilities/lib/api.py
--- a/yt/utilities/lib/api.py
+++ b/yt/utilities/lib/api.py
@@ -29,3 +29,4 @@
 from .write_array import *
 from .mesh_utilities import *
 from .ContourFinding import *
+from .line_integral_convolution import *

diff -r 1e76952e015dfb14db5d97333c182edc658d4aaa -r 98b9f2c2d791d921eccc8e2edc419d863eb14650 yt/utilities/lib/line_integral_convolution.pyx
--- /dev/null
+++ b/yt/utilities/lib/line_integral_convolution.pyx
@@ -0,0 +1,109 @@
+"""
+Utilities for line integral convolution annotation
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2015, yt Development Team.
+#
+# Code originally from Scipy Cookbook (http://wiki.scipy.org/Cookbook/LineIntegralConvolution),
+# with bug fixed which leads to crash when non equal-size vector field in two
+# dimensions is provided.
+#
+# 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
+cimport numpy as np
+cimport cython
+
+ at cython.cdivision(True)
+cdef void _advance_2d(double vx, double vy,
+                     int* x, int* y,
+                     double* fx, double* fy,
+                     int w, int h):
+    cdef double tx, ty
+    if vx>=0:
+        tx = (1-fx[0])/vx
+    else:
+        tx = -fx[0]/vx
+    if vy>=0:
+        ty = (1-fy[0])/vy
+    else:
+        ty = -fy[0]/vy
+    if tx<ty:
+        if vx>=0:
+            x[0]+=1
+            fx[0]=0
+        else:
+            x[0]-=1
+            fx[0]=1
+        fy[0]+=tx*vy
+    else:
+        if vy>=0:
+            y[0]+=1
+            fy[0]=0
+        else:
+            y[0]-=1
+            fy[0]=1
+        fx[0]+=ty*vx
+    if x[0]>=w:
+        x[0]=w-1
+    if x[0]<0:
+        x[0]=0
+    if y[0]<0:
+        y[0]=0
+    if y[0]>=h:
+        y[0]=h-1
+
+def line_integral_convolution_2d(
+        np.ndarray[double, ndim=3] vectors,
+        np.ndarray[double, ndim=2] texture,
+        np.ndarray[double, ndim=1] kernel):
+    cdef int i,j,l,x,y
+    cdef int h,w,kernellen
+    cdef int t
+    cdef double fx, fy, tx, ty
+    cdef np.ndarray[double, ndim=2] result
+
+    w = vectors.shape[0]
+    h = vectors.shape[1]
+    t = vectors.shape[2]
+
+    kernellen = kernel.shape[0]
+    result = np.zeros((w,h),dtype=np.double)
+
+    vectors = vectors[...,::-1].copy()
+
+    for i in range(w):
+        for j in range(h):
+            x = i
+            y = j
+            fx = 0.5
+            fy = 0.5
+            
+            l = kernellen//2
+            result[i,j] += kernel[l]*texture[x,y]
+            while l<kernellen-1:
+                _advance_2d(vectors[x,y,0],vectors[x,y,1],
+                        &x, &y, &fx, &fy, w, h)
+                l+=1
+                result[i,j] += kernel[l]*texture[x,y]
+            
+            x = i
+            y = j
+            fx = 0.5
+            fy = 0.5
+            
+            while l>0:
+                _advance_2d(-vectors[x,y,0],-vectors[x,y,1],
+                        &x, &y, &fx, &fy, w, h)
+                l-=1
+                result[i,j] += kernel[l]*texture[x,y]
+                    
+    return result

diff -r 1e76952e015dfb14db5d97333c182edc658d4aaa -r 98b9f2c2d791d921eccc8e2edc419d863eb14650 yt/utilities/lib/setup.py
--- a/yt/utilities/lib/setup.py
+++ b/yt/utilities/lib/setup.py
@@ -155,6 +155,8 @@
     config.add_extension("amr_kdtools", 
                          ["yt/utilities/lib/amr_kdtools.pyx"],
                          libraries=["m"], depends=["yt/utilities/lib/fp_utils.pxd"])
+    config.add_extension("line_integral_convolution",
+                         ["yt/utilities/lib/line_integral_convolution.pyx"])
     config.add_subpackage("tests")
 
     if os.environ.get("GPERFTOOLS", "no").upper() != "NO":

diff -r 1e76952e015dfb14db5d97333c182edc658d4aaa -r 98b9f2c2d791d921eccc8e2edc419d863eb14650 yt/visualization/plot_modifications.py
--- a/yt/visualization/plot_modifications.py
+++ b/yt/visualization/plot_modifications.py
@@ -23,6 +23,7 @@
 
 from matplotlib.patches import Circle
 from matplotlib.colors import colorConverter
+from matplotlib import cm
 from mpl_toolkits.axes_grid1.anchored_artists import AnchoredSizeBar
 
 from yt.funcs import *
@@ -37,6 +38,8 @@
 from yt.utilities.lib.geometry_utils import triangle_plane_intersect
 from yt.analysis_modules.cosmological_observation.light_ray.light_ray \
      import periodic_ray
+from yt.utilities.lib.line_integral_convolution \
+     import line_integral_convolution_2d
 import warnings
 
 from . import _MPL
@@ -2055,3 +2058,120 @@
             lcb(plot)
 
         return plot
+
+class LineIntegralConvolutionCallback(PlotCallback):
+    """
+    annotate_line_integral_convolution(field_x, field_y, texture=None,
+                                       kernellen=50., lim=(0.5,0.6),
+                                       cmap='binary', alpha=0.8,
+                                       const_alpha=False):
+
+    Add the line integral convolution to the plot for vector fields 
+    visualization. Two component of vector fields needed to be provided
+    (i.e., velocity_x and velocity_y, magentic_field_x and magnetic_field_y).
+
+    Parameters
+    ----------
+
+    field_x, field_y : string
+        The names of two components of vector field which will be visualized
+
+    texture : 2-d array with the same shape of image, optional
+        Texture will be convolved when computing line integral convolution.
+        A white noise background will be used as default.
+
+    kernellen : float, optional
+        The lens of kernel for convolution, which is the length over which the
+        convolution will be performed. For longer kernellen, longer streamline
+        structure will appear.
+
+    lim : 2-element tuple, list, or array, optional
+        The value of line integral convolution will be clipped to the range
+        of lim, which applies upper and lower bounds to the values of line
+        integral convolution and enhance the visibility of plots. Each element
+        should be in the range of [0,1].
+
+    cmap : string, optional
+        The name of colormap for line integral convolution plot.
+
+    alpha : float, optional
+        The alpha value for line integral convolution plot.
+
+    const_alpha : boolean, optional
+        If set to False (by default), alpha will be weighted spatially by
+        the values of line integral convolution; otherwise a constant value
+        of the given alpha is used.
+
+    Example
+    -------
+
+    >>> import yt
+    >>> ds = yt.load('Enzo_64/DD0020/data0020')
+    >>> s = yt.SlicePlot(ds, 'z', 'density')
+    >>> s.annotate_line_integral_convolution('velocity_x', 'velocity_y',\
+                                             lim=(0.5,0.65))
+    """
+    _type_name = "line_integral_convolution"
+    def __init__(self, field_x, field_y, texture=None, kernellen=50.,
+                 lim=(0.5,0.6), cmap='binary', alpha=0.8, const_alpha=False):
+        PlotCallback.__init__(self)
+        self.field_x = field_x
+        self.field_y = field_y
+        self.texture = texture
+        self.kernellen = kernellen
+        self.lim = lim
+        self.cmap = cmap
+        self.alpha = alpha
+        self.const_alpha = const_alpha
+
+    def __call__(self, plot):
+        x0, x1 = plot.xlim
+        y0, y1 = plot.ylim
+        xx0, xx1 = plot._axes.get_xlim()
+        yy0, yy1 = plot._axes.get_ylim()
+        bounds = [x0,x1,y0,y1]
+        extent = [xx0,xx1,yy0,yy1]
+
+        plot._axes.hold(True)
+        nx = plot.image._A.shape[0]
+        ny = plot.image._A.shape[1]
+        pixX = plot.data.ds.coordinates.pixelize(plot.data.axis,
+                                                 plot.data,
+                                                 self.field_x,
+                                                 bounds,
+                                                 (nx,ny))
+        pixY = plot.data.ds.coordinates.pixelize(plot.data.axis,
+                                                 plot.data,
+                                                 self.field_y,
+                                                 bounds,
+                                                 (nx,ny))
+
+        vectors = np.concatenate((pixX[...,np.newaxis],
+                                  pixY[...,np.newaxis]),axis=2)
+
+        if self.texture == None:
+            self.texture = np.random.rand(nx,ny).astype(np.double)
+        elif self.texture.shape != (nx,ny):
+            raise SyntaxError("'texture' must have the same shape "
+                              "with that of output image (%d, %d)" % (nx,ny))
+
+        kernel = np.sin(np.arange(self.kernellen)*np.pi/self.kernellen)
+        kernel = kernel.astype(np.double)
+
+        lic_data = line_integral_convolution_2d(vectors,self.texture,kernel)
+        lic_data = np.flipud(lic_data / lic_data.max())
+        lic_data_clip = np.clip(lic_data,self.lim[0],self.lim[1])
+
+        if self.const_alpha:
+            plot._axes.imshow(lic_data_clip, extent=extent, cmap=self.cmap, 
+                              alpha=self.alpha)
+        else:
+            lic_data_rgba = cm.ScalarMappable(norm=None, cmap=self.cmap).\
+                            to_rgba(lic_data_clip)
+            lic_data_clip_rescale = (lic_data_clip - self.lim[0]) \
+                                    / (self.lim[1] - self.lim[0])
+            lic_data_rgba[...,3] = lic_data_clip_rescale * self.alpha
+            plot._axes.imshow(lic_data_rgba, extent=extent)
+        plot._axes.hold(False)
+
+        return plot

diff -r 1e76952e015dfb14db5d97333c182edc658d4aaa -r 98b9f2c2d791d921eccc8e2edc419d863eb14650 yt/visualization/tests/test_callbacks.py
--- a/yt/visualization/tests/test_callbacks.py
+++ b/yt/visualization/tests/test_callbacks.py
@@ -52,6 +52,7 @@
 #  X scale
 #    material_boundary
 #  X ray
+#  X line_integral_convolution
 
 @contextlib.contextmanager
 def _cleanup_fname():
@@ -326,3 +327,25 @@
             max_level=3, cmap="gist_stern")
         p.save(prefix)
 
+def test_line_integral_convolution_callback():
+    with _cleanup_fname() as prefix:
+        ds = fake_amr_ds(fields =
+            ("density", "velocity_x", "velocity_y", "velocity_z"))
+        for ax in 'xyz':
+            p = ProjectionPlot(ds, ax, "density")
+            p.annotate_line_integral_convolution("velocity_x", "velocity_y")
+            yield assert_fname, p.save(prefix)[0]
+            p = ProjectionPlot(ds, ax, "density", weight_field="density")
+            p.annotate_line_integral_convolution("velocity_x", "velocity_y")
+            yield assert_fname, p.save(prefix)[0]
+            p = SlicePlot(ds, ax, "density")
+            p.annotate_line_integral_convolution("velocity_x", "velocity_y")
+            yield assert_fname, p.save(prefix)[0]
+        # Now we'll check a few additional minor things
+        p = SlicePlot(ds, "x", "density")
+        p.annotate_line_integral_convolution("velocity_x", "velocity_y", 
+                                             kernellen=100., lim=(0.4,0.7),
+                                             cmap='algae', alpha=0.9,
+                                             const_alpha=True)
+        p.save(prefix)
+

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