Core functionality for storing n-D grids — gridData.core

The core module contains classes and functions that are independent of the grid data format. In particular this module contains the Grid class that acts as a universal constructor for specific formats:

g = Grid(**kwargs)           # construct
g.export(filename, format)   # export to the desired format

Some formats can also be read:

g = Grid()                   # make an empty Grid
g.load(filename)             # populate with data from filename

Classes and functions

class gridData.core.Grid(grid=None, edges=None, origin=None, delta=None, metadata=None, interpolation_spline_order=3, file_format=None)[source]

A multidimensional grid object with origin and grid spacings.

Grid objects can be used in arithmetical calculations just like numpy arrays if they are compatible, i.e., they have the same shapes and lengths. In order to make arrays compatible, they an be resampled (resample()) on a common grid.

The attribute grid that holds the data is a standard numpy array and so the data can be directly manipulated.

Data can be read from a number of molecular volume/density formats and written out in different formats with export().

Parameters:
  • grid (numpy.ndarray or str (optional)) – Build the grid either from a histogram or density (a numpy nD array) or read data from a filename.

  • edges (list (optional)) – List of arrays, the lower and upper bin edges along the axes (same as the output by numpy.histogramdd())

  • origin (numpy.ndarray (optional)) – Cartesian coordinates of the center of grid position at index [0, 0, ..., 0].

  • delta (numpy.ndarray (optional)) – Either n x n array containing the cell lengths in each dimension, or n x 1 array for rectangular arrays.

  • metadata (dict (optional)) – A user defined dictionary of arbitrary key/value pairs associated with the density; the class does not touch metadata but stores it with save()

  • interpolation_spline_order (int (optional)) – Order of interpolation function for resampling with resample(); cubic splines = 3 and the default is 3

  • file_format (str (optional)) – Name of the file format; only necessary when grid is a filename (see load()) and autodetection of the file format fails. The default is None and normally the file format is guessed from the file extension.

Raises:
  • TypeError – If the dimensions of the various input data do not agree with each other.

  • ValueError – If some of the required data are not provided in the keyword arguments, e.g., if only the grid is supplied as an array but not the edges or only grid and one of origin and delta.

  • NotImplementedError – If triclinic (non-orthorhombic) boxes are supplied in delta .. Note:: delta can only be a 1D array of length grid.ndim

grid

This array can be any number of dimensions supported by NumPy in order to represent high-dimensional data. When used with data that represents real space densities then the axis convention in GridDataFormats is that axis 0 corresponds to the Cartesian \(x\) component, axis 1 corresponds to the \(y\) component, and axis 2 to the \(z\) component.

Type:

numpy.ndarray

delta

Length of a grid cell (spacing or voxelsize) in \(x\), \(y\), \(z\) dimensions. This is a 1D array with length Grid.grid.ndim.

Type:

numpy.ndarray

origin

Array with the Cartesian coordinates of the coordinate system origin, the center of cell Grid.grid[0, 0, .., 0].

Type:

numpy.ndarray

edges

List of arrays, one for each axis in grid. Each 1D edge array describes the edges of the grid cells along the corresponding axis. The length of an edge array for axis i is grid.shape[i] + 1 because it contains the lower boundary for the first cell, the boundaries between all grid cells, and the upper boundary for the last cell. The edges are assumed to be regular with spacing indicated in delta, namely Grid.delta[i] for axis i.

Type:

list

midpoints

List of arrays, one for each axis in grid. Each 1D midpoints array contains the midpoints of the grid cells along the corresponding axis.

Type:

list

metadata

A user-defined dictionary that can be used to annotate the data. The content is not touched by Grid. It is saved together with the other data with save().

Type:

dict

Example

Create a Grid object from data.

From numpy.histogramdd():

grid, edges = numpy.histogramdd(...)
g = Grid(grid, edges=edges)

From an arbitrary grid:

g = Grid(grid, origin=origin, delta=delta)

From a saved file:

g = Grid(filename)

or

g = Grid()
g.load(filename)

Notes

In principle, the dimension (number of axes) is arbitrary but in practice many formats only support three and almost all functionality is only tested for this special case.

The export() method with format='dx' always exports a 3D object. Other methods might work for an array of any dimension (in particular the Python pickle output).

Changed in version 0.5.0: New file_format keyword argument.

Changed in version 0.7.0: CCP4 files are now read with gridData.mrc.MRC and not anymore with the deprecated/buggy ccp4.CCP4

centers()[source]

Returns the coordinates of the centers of all grid cells as an iterator.

See also

numpy.ndindex()

check_compatible(other)[source]

Check if other can be used in an arithmetic operation.

other is compatible if

  1. other is a scalar

  2. other is a grid defined on the same edges

In order to make other compatible, resample it on the same grid as this one using resample().

Parameters:

other (Grid or float or int) – Another object to be used for standard arithmetic operations with this Grid

Raises:

TypeError – if not compatible

See also

resample()

default_format = 'DX'

Default format for exporting with export().

export(filename, file_format=None, type=None, typequote='"')[source]

export density to file using the given format.

The format can also be deduced from the suffix of the filename although the file_format keyword takes precedence.

The default format for export() is ‘dx’. Use ‘dx’ for visualization.

Implemented formats:

dx

OpenDX

pickle

pickle (use Grid.load() to restore); Grid.save() is simpler than export(format='python').

Parameters:
  • filename (str) – name of the output file

  • file_format ({'dx', 'pickle', None} (optional)) – output file format, the default is “dx”

  • type (str (optional)) –

    for DX, set the output DX array type, e.g., “double” or “float”. By default (None), the DX type is determined from the numpy dtype of the array of the grid (and this will typically result in “double”).

    New in version 0.4.0.

  • typequote (str (optional)) –

    For DX, set the character used to quote the type string; by default this is a double-quote character, ‘”’. Custom parsers like the one from NAMD-GridForces (backend for MDFF) expect no quotes, and typequote=’’ may be used to appease them.

    New in version 0.5.0.

property interpolated

B-spline function over the data grid(x,y,z).

The interpolated() function allows one to obtain data values for any values of the coordinates:

interpolated([x1,x2,...],[y1,y2,...],[z1,z2,...]) -> F[x1,y1,z1],F[x2,y2,z2],...

The interpolation order is set in Grid.interpolation_spline_order.

The interpolated function is computed once and is cached for better performance. Whenever interpolation_spline_order is modified, Grid.interpolated() is recomputed.

The value for unknown data is set in Grid.interpolation_cval (TODO: also recompute when interpolation_cval value is changed.)

Example

Example usage for resampling:

XX, YY, ZZ = numpy.mgrid[40:75:0.5, 96:150:0.5, 20:50:0.5]
FF = interpolated(XX, YY, ZZ)

Note

Values are interpolated with a spline function. It is possible that the spline will generate values that would not normally appear in the data. For example, a density is non-negative but a cubic spline interpolation can generate negative values, especially at the boundary between 0 and high values.

Internally, the function uses scipy.ndimage.map_coordinates() with mode="constant" whereby interpolated values outside the interpolated grid are determined by filling all values beyond the edge with the same constant value, defined by the interpolation_cval parameter, which when not set defaults to the minimum value in the interpolated grid.

Changed in version 0.6.0: Interpolation outside the grid is now performed with mode="constant" rather than mode="nearest", eliminating extruded volumes when interpolating beyond the grid.

property interpolation_spline_order

Order of the B-spline interpolation of the data.

3 = cubic; 4 & 5 are also supported

Only choose values that are acceptable to scipy.ndimage.spline_filter()!

See also

interpolated

load(filename, file_format=None)[source]

Load saved grid and edges from filename

The load() method calls the class’s constructor method and completely resets all values, based on the loaded data.

resample(edges)[source]

Resample data to a new grid with edges edges.

This method creates a new grid with the data from the current grid resampled to a regular grid specified by edges. The order of the interpolation is set by Grid.interpolation_spline_order: change the value before calling resample().

Parameters:

edges (tuple of arrays or Grid) – edges of the new grid or a Grid instance that provides Grid.edges

Returns:

a new Grid with the data interpolated over the new grid cells

Return type:

Grid

Examples

Providing edges (a tuple of three arrays, indicating the boundaries of each grid cell):

g = grid.resample(edges)

As a convenience, one can also supply another Grid as the argument for this method

g = grid.resample(othergrid)

and the edges are taken from Grid.edges.

resample_factor(factor)[source]

Resample to a new regular grid.

Parameters:

factor (float) – The number of grid cells are scaled with factor in each dimension, i.e., factor * N_i cells along each dimension i. Must be positive, and cannot result in fewer than 2 cells along a dimension.

Returns:

interpolated grid – The resampled data are represented on a Grid with the new grid cell sizes.

Return type:

Grid

See also

resample

Changed in version 0.6.0: Previous implementations would not alter the range of the grid edges being resampled on. As a result, values at the grid edges would creep steadily inward. The new implementation recalculates the extent of grid edges for every resampling.

save(filename)[source]

Save a grid object to filename and add “.pickle” extension.

Internally, this calls Grid.export(filename, format="python"). A grid can be regenerated from the saved data with

g = Grid(filename="grid.pickle")

Note

The pickle format depends on the Python version and therefore it is not guaranteed that a grid saved with, say, Python 2.7 can also be read with Python 3.5. The OpenDX format is a better alternative for portability.

gridData.core.ndmeshgrid(*arrs)[source]

Return a mesh grid for N dimensions.

The input are N arrays, each of which contains the values along one axis of the coordinate system. The arrays do not have to have the same number of entries. The function returns arrays that can be fed into numpy functions so that they produce values for all points spanned by the axes arrs.

Original from http://stackoverflow.com/questions/1827489/numpy-meshgrid-in-3d and fixed.