ovito.pipeline
This module contains classes that are part of OVITO’s data pipeline system.
Pipelines:
Modifier
(base class of all built-in modifiers)
ModifierInterface
(base class of user-defined modifiers)
Data sources:
- class ovito.pipeline.FileSource
This object type serves as a
Pipeline.source
and takes care of reading the input data for aPipeline
from an external file.You normally do not need to create an instance of this class yourself; the
import_file()
function does it for you and wires the fully configuredFileSource
to the newPipeline
. However, if needed, theFileSource.load()
method allows you to load a different input file later on and replace the input of the existing pipeline with a new dataset:from ovito.io import import_file # Create a new pipeline with a FileSource: pipeline = import_file('input/first_file.dump') # Get the data from the first file: data1 = pipeline.compute() # Use FileSource.load() method to replace the pipeline's input with a different file pipeline.source.load('input/second_file.dump') # Now the pipeline gets its input data from the new file: data2 = pipeline.compute()
Furthermore, you will encounter other
FileSource
objects in conjunction with certain modifiers that need secondary input data from a separate file. TheCalculateDisplacementsModifier
, for example, manages its ownFileSource
for loading reference particle positions from a separate input file. Another example is theLoadTrajectoryModifier
, which employs its own separateFileSource
instance to load the particle trajectories from disk and combine them with the topology data previously loaded by the mainFileSource
of the data pipeline.Data access
The
FileSource
class provides two ways of accessing the data that is loaded from the external input file(s). For read-only access to the data, theFileSource.compute()
method should be called. It loads the data of a specific frame from the input simulation trajectory and returns it as a newDataCollection
object:# This creates a new Pipeline with an attached FileSource. pipeline = import_file('input/simulation.dump') # Request data of animation frame 0 from the FileSource. data = pipeline.source.compute(0) print(data.particles.positions[...])
To modify or amend the
DataCollection
that was loaded by theFileSource
, you should insert a user-defined modifier function into the pipeline. A typical use case is assigning the radii and names to particle types loaded from a simulation file that doesn’t contain named atom types:pipeline = import_file('input/simulation.dump') # User-defined modifier function that assigns names and radii to numeric atom types: def setup_atom_types(frame, data): types = data.particles_.particle_types_ types.type_by_id_(1).name = "Cu" types.type_by_id_(1).radius = 1.35 types.type_by_id_(2).name = "Zr" types.type_by_id_(2).radius = 1.55 pipeline.modifiers.append(setup_atom_types)
- compute(frame=None)
Requests data from this data source. The
FileSource
will load it from the external file if needed.The optional frame parameter determines the frame to retrieve, which must be in the range 0 through (
num_frames
-1). If no frame number is specified, the current time slider position is used (will always be frame 0 for scripts not running in the context of an interactive OVITO session).The
FileSource
uses a caching mechanism to keep the data for one or more frames in memory. Thus, invokingcompute()
repeatedly to retrieve the same frame will typically be very fast.- Parameters
frame (int) – The source frame to retrieve.
- Returns
A new
DataCollection
containing the loaded data.
- property data
This field provides access to the internal
DataCollection
, i.e. the master copy of the data loaded from the input file (at frame 0).
- load(location, **params)
Sets a new input file, from which this data source will retrieve its data from.
The function accepts additional keyword arguments, which are forwarded to the format-specific file reader. For further information, please see the documentation of the
import_file()
function, which has the same function interface as this method.- Parameters
location (str) – The local file(s) or remote URL to load.
- class ovito.pipeline.Modifier
This is the base class for all modifier types in OVITO. See the
ovito.modifiers
module for a list of concrete modifier types that can be inserted into a dataPipeline
.- property enabled
Controls whether the modifier is applied to the data. Disabled modifiers are skipped during evaluation of a data pipeline.
- Default
True
- class ovito.pipeline.ModifierInterface
Base:
traits.has_traits.HasStrictTraits
Abstract base class for Python script modifiers that follow the advanced programming interface.
New in version 3.8.0.
- class InputSlot
Represents the upstream pipeline generating the input data for a custom modifier implementation.
- compute(frame)
Computes the results of the upstream pipeline connected to this input slot.
frame specifies the trajectory frame to retrieve, which must be in the range 0 to (
num_frames
-1).The slot uses a caching mechanism to keep the data for one or more frames in memory. Thus, invoking
compute()
repeatedly to retrieve the same frame will typically be very fast.- Parameters
frame (int) – The trajectory frame to retrieve from the upstream pipeline.
- Return type
- input_caching_hints(frame, *, input_slots, **kwargs)
User-defined modifiers that access multiple trajectory frames in their
modify()
method should implement this method to communicate the list of frames going to be needed. The pipeline system will keep the data of these trajectory frames in an internal cache to avoid unnecessary I/O and compute operations. See Input data caching.- Parameters
- Return type
If your modifier defines additional input slots, the function must return a dictionary that specifies for each input slot, including the standard upstream slot, which input frame(s) should be cached. For example:
extra_slot = OvitoObjectTrait(FileSource) def input_caching_hints(self, frame, **kwargs): return { 'upstream': frame, 'extra_slot': 0 }
Note
This method is supposed to be implemented as part of a user-defined modifier class but it should not be called by user code. The pipeline system will automatically invoke this method whenever necessary.
- abstract modify(data, *, frame, input_slots, data_cache, **kwargs)
The actual work function which gets called by the pipeline system to let the modifier do its thing.
- Parameters
data (DataCollection) – Data snapshot which should be modified by the modifier function in place.
frame (int) – Zero-based trajectory frame number.
input_slots (Dict[str, InputSlot]) – One or more
InputSlot
objects representing the upstream data pipeline(s) connected to this modifier.data_cache (DataCollection) – A data container (initially empty) which may be used by the modifier function to store intermediate results.
kwargs (Any) – Other arguments that may be passed by the pipeline system. This parameter should be included for forward compatibility with future versions of OVITO.
- Return type
- class ovito.pipeline.Pipeline
This class encapsulates a data pipeline, consisting of a data source and a chain of zero or more modifiers, which manipulate the data on the way through the pipeline.
Pipeline creation
Every pipeline has a data source, which loads or dynamically generates the input data entering the pipeline. This source is accessible through the
Pipeline.source
field and may be replaced with a different kind of source object if needed. For pipelines created by theimport_file()
function, the data source is automatically set to be aFileSource
object, which loads the input data from the external file and feeds it into the pipeline. Another kind of data source is theStaticSource
, which can be used if you want to programmatically specify the input data for the pipeline instead of loading it from a file.The modifiers that are part of the pipeline are accessible through the
Pipeline.modifiers
field. This list is initially empty and you can populate it with the modifier types found in theovito.modifiers
module. Note that it is possible to employ the sameModifier
instance in more than one pipeline. And it is okay to use the same data source object for several pipelines, letting them process the same input data.Pipeline evaluation
Once the pipeline is set up, its computation results can be requested by calling
compute()
, which means that the input data will be loaded/generated by thesource
and all modifiers of the pipeline are applied to the data one after the other. Thecompute()
method returns a newDataCollection
storing the data objects produced by the pipeline. Under the hood, an automatic caching system ensures that unnecessary file accesses and computations are avoided. Repeatedly callingcompute()
will not trigger a recalculation of the pipeline’s results unless you alter the pipeline’s data source, the chain of modifiers, or a modifier’s parameters.Usage example
The following code example shows how to create a new pipeline by importing an MD simulation file and inserting a
SliceModifier
to cut away some of the particles. Finally, the total number of remaining particles is printed.from ovito.io import import_file from ovito.modifiers import SliceModifier # Import a simulation file. This creates a Pipeline object. pipeline = import_file('input/simulation.dump') # Insert a modifier that operates on the data: pipeline.modifiers.append(SliceModifier(normal=(0,0,1), distance=0)) # Compute the effect of the slice modifier by evaluating the pipeline. output = pipeline.compute() print("Remaining particle count:", output.particles.count)
If you would like to access the unmodified input data of the pipeline, i.e. before it has been processed by any of the modifiers, you can call the
FileSource.compute()
method instead:# Access the pipeline's input data provided by the FileSource: input = pipeline.source.compute() print("Input particle count:", input.particles.count)
Data visualization
If you intend to produce graphical renderings of a output data produced by a pipeline, you must make the pipeline part of the current three-dimensional scene by calling the
Pipeline.add_to_scene()
method.Data export
To export the generated data of the pipeline to an output file, simply call the
ovito.io.export_file()
function with the pipeline.- add_to_scene()
Inserts the pipeline into the three-dimensional scene by appending it to the
ovito.Scene.pipelines
list. The visual representation of the pipeline’s output data will appear in rendered images and in the interactive viewports.You can remove the pipeline from the scene again using
remove_from_scene()
.
- compute(frame=None)
Computes and returns the output of this data pipeline (for one trajectory frame).
This method requests an evaluation of the pipeline and blocks until the input data has been obtained from the data
source
, e.g. a simulation file, and all modifiers have been applied to the data. If you invoke thecompute()
method repeatedly without changing the modifiers in the pipeline between calls, the pipeline may serve subsequent requests by returning cached output data.The optional frame parameter specifies the animation frame at which the pipeline should be evaluated. Frames are always consecutively numbered (0, 1, 2, …). If you don’t specify any frame number, the current time slider position is used – or frame 0 if not running in the context of an interactive OVITO Pro session.
The
compute()
method raises aRuntimeError
if the pipeline could not be successfully evaluated for some reason. This may happen due to invalid modifier settings or file I/O errors, for example.- Parameters
frame (int) – The animation frame number at which the pipeline should be evaluated.
- Returns
A
DataCollection
produced by the data pipeline.
Attention
This method returns a snapshot of the results of the current pipeline, representing an independent data copy. That means snapshot will not reflect changes you subsequently make to the pipeline or the modifiers within the pipeline. After changing the pipeline, you have to invoke
compute()
again to let the pipeline produce a new updated snapshot.Attention
The returned
DataCollection
represents a copy of the pipeline’s internal data, which means, if you subsequently make any changes to the objects in theDataCollection
, those changes will not be visible to the modifiers within the pipeline – even if you add those modifiers to the pipeline after thecompute()
call as in this example:data = pipeline.compute() data.particles_.create_property('Foo', data=...) pipeline.modifier.append(ExpressionSelectionModifier(expression='Foo > 0')) new_data = pipeline.compute() # ERROR
The second call to
compute()
will raise an error, because theExpressionSelectionModifier
references the new particle propertyFoo
, which does not exist in the original data seen by the pipeline. That’ because we’ve added the propertyFoo
only to theParticles
object that is stored in our snapshotdata
. ThisDataCollection
is independent from the transient data the pipeline operates on.To make the property
Foo
available to modifiers in the pipeline, we thus need to create the property within the pipeline. This can be accomplished by performing the modification step as a Python modifier function that is inserted into the pipeline:def add_foo(frame, data): data.particles_.create_property('Foo', data=...) pipeline.modifier.append(add_foo) pipeline.modifier.append(ExpressionSelectionModifier(expression='Foo > 0'))
Downstream modifiers now see the new particle property created by our user-defined modifier function
add_foo
, which operates on a transient data collection managed by the pipeline system.
- property modifiers
The sequence of modifiers in the pipeline.
This list contains any modifiers that are applied to the input data provided by the pipeline’s data
source
. You can add and remove modifiers as needed using standard Pythonappend()
anddel
operations. The head of the list represents the beginning of the pipeline, i.e. the first modifier receives the data from the datasource
, manipulates it and passes the results on to the second modifier in the list and so forth.Example: Adding a new modifier to the end of a data pipeline:
pipeline.modifiers.append(WrapPeriodicImagesModifier())
- remove_from_scene()
Removes the visual representation of the pipeline from the scene by deleting it from the
ovito.Scene.pipelines
list. The output data of the pipeline will disappear from viewports.
- property source
The object that provides the data entering the pipeline. This typically is a
FileSource
instance if the pipeline was created by a call toimport_file()
. You can assign a new source to the pipeline if needed. See theovito.pipeline
module for a list of available pipeline source types. Note that you can even make several pipelines share the same source object.
- class ovito.pipeline.PythonScriptSource
Example:
from ovito.pipeline import Pipeline, PythonScriptSource from ovito.io import export_file from ovito.data import Particles, DataCollection import numpy # Define a data source function, which fills an empty DataCollection with # some data objects. def create_particles(frame: int, data: DataCollection): data.particles = Particles(count = 20) coordinates = data.particles_.create_property('Position') coordinates[:,0] = numpy.linspace(0.0, 50.0, data.particles.count) coordinates[:,1] = numpy.cos(coordinates[:,0]/4.0 + frame/5.0) coordinates[:,2] = numpy.sin(coordinates[:,0]/4.0 + frame/5.0) # Create a data pipeline with a PythonScriptSource, which executes our # script function defined above. pipeline = Pipeline(source = PythonScriptSource(function = create_particles)) # Export the results of the data pipeline to an output file. # The system will invoke the Python function defined above once per animation frame. export_file(pipeline, 'output/trajectory.xyz', format='xyz', columns=['Position.X', 'Position.Y', 'Position.Z'], multiple_frames=True, start_frame=0, end_frame=10)
- property function
The Python function to be called each time the data pipeline is evaluated by the system.
The function must have a signature as shown in the example above. The frame parameter specifies the current animation frame number at which the data pipeline is being evaluated. The
DataCollection
data is initially empty should be populated with data objects by the user-defined source function.- Default
None
- class ovito.pipeline.ReferenceConfigurationModifier
Base:
ovito.pipeline.Modifier
This is the common base class of analysis modifiers that perform some kind of comparison of the current particle configuration with a reference configuration. For example, the
CalculateDisplacementsModifier
, theAtomicStrainModifier
and theWignerSeitzAnalysisModifier
are modifier types that require a reference configuration as additional input.Constant and sliding reference configurations
The
ReferenceConfigurationModifier
base class provides various fields that allow you to specify the reference particle configuration. By default, frame 0 of the currently loaded simulation sequence is used as reference. You can select any other frame with thereference_frame
field. Sometimes an incremental analysis is desired, instead of a fixed reference configuration. That means the sliding reference configuration and the current configuration are separated along the time axis by a constant period (delta t). The incremental analysis mode is activated by setting theuse_frame_offset
flag and specifying the desiredframe_offset
.External reference configuration file
By default, the reference particle positions are obtained by evaluating the same data pipeline that also provides the current particle positions, i.e. which the modifier is part of. That means any modifiers preceding this modifier in the pipeline will also act upon the reference particle configuration, but not modifiers that follow in the pipeline.
Instead of taking it from the same data pipeline, you can explicitly provide a reference configuration by loading it from a separate data file. To this end the
reference
field contains aFileSource
object and you can use itsload()
method to load the reference particle positions from a separate file.Handling of periodic boundary conditions and cell deformations
Certain analysis modifiers such as the
CalculateDisplacementsModifier
and theAtomicStrainModifier
calculate the displacements particles experienced between the reference and the current configuration. Since particle coordinates in periodic simulation cells are often stored in a wrapped form, calculating the displacement vectors is non-trivial when particles have crossed the periodic boundaries. By default, the minimum image convention is used in these cases, but you can turn if off by settingminimum_image_convention
toFalse
, for example if the input particle coordinates are given in unwrapped form.Furthermore, if the simulation cell of the reference and the current configuration are different, it makes a slight difference whether displacements are calculated in the reference or in the current frame. The
affine_mapping
property controls the type of coordinate mapping that is used.- property affine_mapping
- Selects the type of affine deformation applied to the particle coordinates of either the reference or the current configuration prior to the actual analysis computation. Must be one of the following modes:
ReferenceConfigurationModifier.AffineMapping.Off
ReferenceConfigurationModifier.AffineMapping.ToReference
ReferenceConfigurationModifier.AffineMapping.ToCurrent
When affine mapping is disabled (
AffineMapping.Off
), particle displacement vectors are simply calculated from the difference of current and reference positions, irrespective of the cell shape the reference and current configuration. Note that this can introduce a small geometric error if the shape of the periodic simulation cell changes considerably. The modeAffineMapping.ToReference
applies an affine transformation to the current configuration such that all particle positions are first mapped to the reference cell before calculating the displacement vectors. The last option,AffineMapping.ToCurrent
, does the reverse: it maps the reference particle positions to the deformed cell before calculating the displacements.- Default
ReferenceConfigurationModifier.AffineMapping.Off
- property frame_offset
The relative frame offset when using a sliding reference configuration (if
use_frame_offset
==True
). Negative frame offsets correspond to reference configurations that precede the current configuration in time.- Default
-1
- property minimum_image_convention
If
False
, then displacements are calculated from the particle coordinates in the reference and the current configuration as is. Note that in this case the calculated displacements of particles that have crossed a periodic simulation cell boundary will be wrong if their coordinates are stored in a wrapped form. IfTrue
, then the minimum image convention is applied when calculating the displacements of particles that have crossed a periodic boundary.- Default
True
- property reference
A pipeline
source
object that provides the reference particle positions. By default this field isNone
, in which case the modifier obtains the reference particle positions from data pipeline it is part of. You can explicitly assign a data source object such as aFileSource
or aStaticSource
to this field to specify an explicit reference configuration.# The modifier that requires a reference config: mod = CalculateDisplacementsModifier() # Load the reference config from a separate input file. mod.reference = FileSource() # Note: You may have to import FileSource from the ovito.pipeline module. mod.reference.load('input/simulation.0.dump')
- Default
None
- property reference_frame
The frame number to use as reference configuration. Ignored if
use_frame_offset
is set.- Default
0
- property use_frame_offset
Determines whether a sliding reference configuration is taken at a constant time offset (specified by
frame_offset
) relative to the current frame. IfFalse
, a constant reference configuration is used (set by thereference_frame
parameter) irrespective of the current frame.- Default
False
- class ovito.pipeline.StaticSource
Serves as a data
source
for aPipeline
. AStaticSource
manages aDataCollection
, which it will pass to thePipeline
as input data. One typically initializes aStaticSource
with a collection of data objects, then wiring it to aPipeline
as follows:from ovito.pipeline import StaticSource, Pipeline from ovito.data import DataCollection, SimulationCell, Particles from ovito.modifiers import CreateBondsModifier from ovito.io import export_file # Insert a new SimulationCell object into a data collection: data = DataCollection() cell = SimulationCell(pbc = (False, False, False)) cell[:,0] = (4,0,0) cell[:,1] = (0,2,0) cell[:,2] = (0,0,2) data.objects.append(cell) # Create a Particles object containing two particles: particles = Particles() particles.create_property('Position', data=[[0,0,0],[2,0,0]]) data.objects.append(particles) # Create a new Pipeline with a StaticSource as data source: pipeline = Pipeline(source = StaticSource(data = data)) # Apply a modifier: pipeline.modifiers.append(CreateBondsModifier(cutoff = 3.0)) # Write pipeline results to an output file: export_file(pipeline, 'output/structure.data', 'lammps/data', atom_style='bond')
- compute(frame=None)
Returns a copy of the
DataCollection
stored in this source’sdata
field.- Parameters
frame – This parameter is ignored, because the data of a
StaticSource
is not time-dependent.- Returns
A new
DataCollection
containing the data stored in thisStaticSource
.
- property data
The
DataCollection
managed by this object, which will be fed to the pipeline.- Default
None