ovito.io

This module primarily provides two high-level functions for reading and writing external data files:

ovito.io.import_file(location, **params)

Imports data from an external file.

This Python function corresponds to the Load File command in OVITO’s user interface. The format of the imported file is automatically detected. However, depending on the file’s format, additional keyword parameters may need to be supplied to specify how the data should be interpreted. These keyword parameters are documented below.

Parameters:location (str) – The data file to import. This can be a local file path or a remote sftp:// URL.
Returns:The Pipeline that has been created for the imported data.

The function creates and returns a new Pipeline wired to a FileSource, which is responsible for loading the selected input data from disk. The import_file() function requests a first evaluation of the new pipeline in order to fill the pipeline’s input cache with data from the file.

Note that the newly created Pipeline is not automatically inserted into the three-dimensional scene. That means it won’t appear in rendered images or the interactive viewports of the graphical OVITO program. You need to explicitly insert the pipeline into the scene by calling add_to_scene() when desired.

Furthermore, note that the returned Pipeline may be re-used to load a different input file later on. Instead of calling import_file() again to load the next file, you can use the pipeline.source.load(...) method to replace the input file of the already existing pipeline.

File columns

When importing XYZ files or binary LAMMPS dump files, the mapping of file columns to OVITO’s particle properties must be specified using the columns keyword parameter:

pipeline = import_file("file.xyz", columns = 
  ["Particle Identifier", "Particle Type", "Position.X", "Position.Y", "Position.Z"])

The number of list entries must match the number of data columns in the input file. See the list of standard particle properties for possible target properties. If you specify a non-standard property name for a column, a user-defined ParticleProperty object is created from the values in that file column. For vector properties, the component name must be appended to the property’s base name as demonstrated for the Position property in the example above. To completely ignore a file column during import, specify None at the corresponding position in the columns list.

For text-based LAMMPS dump files, OVITO automatically determines a reasonable column-to-property mapping, but you may override it using the columns keyword. This can make sense, for example, if the file columns containing the particle coordinates do not follow the standard name scheme x, y, and z (e.g. when reading time-averaged atomic positions computed by LAMMPS).

File sequences

You can import a sequence of files by passing a filename containing a * wildcard character to import_file(). There may be only one * in the filename (and not in the directory names). The wildcard matches only to numbers in a filename.

OVITO scans the directory and imports all matching files that belong to the sequence. Note that OVITO only loads the first file into memory though.

The length of the imported time series is reported by the num_frames field of the pipeline’s FileSource and is also reflected by the global AnimationSettings object. You can step through the frames of the animation sequence as follows:

from ovito.io import import_file

# Import a sequence of files.
pipeline = import_file('simulation.*.dump')

# Loop over all frames of the sequence.
for frame in range(pipeline.source.num_frames):    

    # Evaluating the pipeline at the current animation frame loads the 
    # corresponding input file from the sequence into memory:
    data = pipeline.compute(frame)
    
    # Work with the loaded data...
    print("Particle positions at frame %i:" % frame)
    print(data.particle_properties['Position'][...])

Multi-timestep files

Some file formats store multiple frames in one file. In some cases OVITO cannot know (e.g. for XYZ and LAMMPS dump formats) that a file contains multiple frames (because, by default, reading the entire file is avoided for performance reasons). Then it is necessary to explicitly tell OVITO to scan the entire file and load a sequence of frames by supplying the multiple_frames option:

pipeline = import_file("file.dump", multiple_frames = True)

LAMMPS atom style

When loading a LAMMPS data file that is based on an atom style other than “atomic”, the style must be explicitly specified using the atom_style keyword parameter. Only the following LAMMPS atom styles are currently supported by OVITO: angle, atomic, body, bond, charge, dipole, full, molecular. LAMMPS data files may contain a hint indicating the atom style (in particular files generated by LAMMPS’s own write_data command contain that hint). Such files can be loaded without specifying the atom_style keyword.

Topology and trajectory files

Some simulation codes write a topology file and separate trajectory file. The former contains only static information like the bonding between atoms, the atom types, etc., which do not change during an MD simulation run, while the latter stores the varying data (first and foremost the atomic trajectories). To load such a topology-trajectory pair of files in OVITO, first load the topology file with the import_file() function, then insert a LoadTrajectoryModifier into the returned Pipeline to also load the trajectory data.

ovito.io.export_file(data, file, format, **params)

High-level function that exports data to a file. See the Data export section for an overview of this topic.

Parameters:
  • data – The object to be exported. See below for options.
  • file (str) – The output file path.
  • format (str) – The type of file to write. See below for options.

Data to export

Various kinds of objects are accepted as data parameter by the function:

  • Pipeline: Exports the data dynamically generated by a pipeline. Since pipelines can be evaluated at different animation times, multi-frame sequences can be produced from a single Pipeline object.
  • DataCollection: Exports the data contained in a data collection. Data objects in the collection that are not compatible with the chosen output format are ignored. Since a data collection represents a static snapshot, no frame sequences can be written.
  • DataObject: Exports the data object as if it were the only part of a DataCollection.
  • None: All pipelines that are part of the current scene (see ovito.DataSet.scene_pipelines) are exported. This option makes sense for scene descriptions such as the POV-Ray format.

Output format

The format parameter determines the type of file to write; the filename suffix is ignored. However, for filenames that end with .gz, automatic gzip compression is activated if the selected format is text-based. The following format strings are supported:

  • "txt" – Exporting global attributes to a text file (see below)
  • "lammps/dump" – LAMMPS text-based dump format
  • "lammps/data" – LAMMPS data format
  • "imd" – IMD format
  • "vasp" – POSCAR format
  • "xyz" – XYZ format
  • "fhi-aims" – FHI-aims format
  • "netcdf/amber" – Binary format for MD data following the AMBER format convention
  • "vtk/trimesh" – ParaView VTK format for exporting SurfaceMesh objects
  • "vtk/disloc" – ParaView VTK format for exporting DislocationNetwork objects
  • "ca"Text-based format for storing dislocation lines
  • "povray" – POV-Ray scene format

Depending on the selected format, additional keyword arguments must be passed to export_file():

File columns

For the output formats lammps/dump, xyz, imd and netcdf/amber, you must specify the set of particle properties to export using the columns keyword parameter:

export_file(pipeline, "output.xyz", "xyz", columns = 
  ["Particle Identifier", "Particle Type", "Position.X", "Position.Y", "Position.Z"]
)

You can export the standard particle properties and any user-defined properties present in the pipeline’s output DataCollection. For vector properties, the component name must be appended to the base name as demonstrated above for the Position property.

Exporting several simulation frames

By default, only the current animation frame (given by the current_frame global variable) is exported. To export a different frame, pass the frame keyword parameter to the export_file() function. Alternatively, you can export all frames of the current animation sequence at once by passing multiple_frames=True. Refined control of the exported frame sequence is possible using the keyword arguments start_frame, end_frame, and every_nth_frame.

The lammps/dump and xyz file formats can store multiple frames in a single output file. For other formats, or if you intentionally want to generate one file per frame, you must pass a wildcard filename to export_file(). This filename must contain exactly one * character as in the following example, which will be replaced with the animation frame number:

export_file(pipeline, "output.*.dump", "lammps/dump", multiple_frames=True)

The above call is equivalent to the following for-loop:

for i in range(pipeline.source.num_frames):
    export_file(pipeline, "output.%i.dump" % i, "lammps/dump", frame=i)

LAMMPS atom style

When writing files in the lammps/data format, the LAMMPS atom style “atomic” is used by default. If you want to create a data file that uses a different atom style, specify it with the atom_style keyword parameter:

export_file(pipeline, "output.data", "lammps/data", atom_tyle="bond")

The following LAMMPS atom styles are currently supported by OVITO: angle, atomic, bond, charge, dipole, full, molecular, sphere.

Global attributes

The txt file format allows you to export global quantities computed by the data pipeline to a text file. For example, to write out the number of FCC atoms identified by a CommonNeighborAnalysisModifier as a function of simulation time, one would do the following:

export_file(pipeline, "data.txt", "txt", 
    columns=["Timestep", "CommonNeighborAnalysis.counts.FCC"], 
    multiple_frames=True)

See the documentation of the individual analysis modifiers to find out which global quantities the compute. You can also determine which attributes are available as follows:

print(pipeline.compute().attributes)