yt.visualization.streamlines.Streamlines

class yt.visualization.streamlines.Streamlines(ds, positions, xfield='velocity_x', yfield='velocity_x', zfield='velocity_x', volume=None, dx=None, length=None, direction=1, get_magnitude=False)[source]

A collection of streamlines that flow through the volume

The Streamlines object contains a collection of streamlines defined as paths that are parallel to a specified vector field.

Parameters:

ds : ~yt.data_objects.Dataset

This is the dataset to streamline

pos : array_like

An array of initial starting positions of the streamlines.

xfield: field, optional :

The x component of the vector field to be streamlined. Default:’velocity_x’

yfield: field, optional :

The y component of the vector field to be streamlined. Default:’velocity_y’

zfield: field, optional :

The z component of the vector field to be streamlined. Default:’velocity_z’

volume : yt.extensions.volume_rendering.AMRKDTree, optional

The volume to be streamlined. Can be specified for finer-grained control, but otherwise will be automatically generated. At this point it must use the AMRKDTree. Default: None

dx : float, optional

Optionally specify the step size during the integration. Default: minimum dx

length : float, optional

Optionally specify the length of integration. Default: np.max(self.ds.domain_right_edge-self.ds.domain_left_edge)

direction : real, optional

Specifies the direction of integration. The magnitude of this value has no effect, only the sign.

get_magnitude: bool, optional :

Specifies whether the Streamlines.magnitudes array should be filled with the magnitude of the vector field at each point in the streamline. This seems to be a ~10% hit to performance. Default: False

Examples

>>> import yt
>>> import numpy as np
>>> import matplotlib.pylab as pl
>>>
>>> from yt.visualization.api import Streamlines
>>> from mpl_toolkits.mplot3d import Axes3D
>>>
>>> # Load the dataset and set some parameters
>>> ds = load('IsolatedGalaxy/galaxy0030/galaxy0030')
>>> c = np.array([0.5]*3)
>>> N = 100
>>> scale = 1.0
>>> pos_dx = np.random.random((N,3))*scale-scale/2.
>>> pos = c+pos_dx
>>>
>>> # Define and construct streamlines
>>> streamlines = Streamlines(
        ds,pos, 'velocity_x', 'velocity_y', 'velocity_z', length=1.0) 
>>> streamlines.integrate_through_volume()
>>>
>>> # Make a 3D plot of the streamlines and save it to disk
>>> fig=pl.figure() 
>>> ax = Axes3D(fig)
>>> for stream in streamlines.streamlines:
>>>     stream = stream[np.all(stream != 0.0, axis=1)]
>>>     ax.plot3D(stream[:,0], stream[:,1], stream[:,2], alpha=0.1)
>>> pl.savefig('streamlines.png')

Attributes

Methods