From 210e05c444373c0636456c7145494824cccc445e Mon Sep 17 00:00:00 2001 From: Ahmed Salah <72524423+AHMED-salah00@users.noreply.github.com> Date: Tue, 11 Apr 2023 16:31:51 +0200 Subject: [PATCH] Fixing some doctests. (#4118) Towards #3925 * Solved failing doctests in `core/universe.py`, `lib/util.py`, and `lib/transformations.py` --- package/MDAnalysis/core/universe.py | 30 ++++-- package/MDAnalysis/lib/transformations.py | 123 +++++++++++++++++++++- package/MDAnalysis/lib/util.py | 31 +++++- 3 files changed, 166 insertions(+), 18 deletions(-) diff --git a/package/MDAnalysis/core/universe.py b/package/MDAnalysis/core/universe.py index 572fdf628df..7de0e5aef7f 100644 --- a/package/MDAnalysis/core/universe.py +++ b/package/MDAnalysis/core/universe.py @@ -806,9 +806,12 @@ def add_TopologyAttr(self, topologyattr, values=None): ------- For example to add bfactors to a Universe: - >>> u.add_TopologyAttr('bfactors') - >>> u.atoms.bfactors - array([ 0., 0., 0., ..., 0., 0., 0.]) + >>> import MDAnalysis as mda + >>> from MDAnalysis.tests.datafiles import PSF, DCD + >>> u = mda.Universe(PSF, DCD) + >>> u.add_TopologyAttr('tempfactors') + >>> u.atoms.tempfactors + array([0., 0., 0., ..., 0., 0., 0.]) .. versionchanged:: 0.17.0 Can now also add TopologyAttrs with a string of the name of the @@ -863,8 +866,15 @@ def del_TopologyAttr(self, topologyattr): ------- For example to remove bfactors to a Universe: - >>> u.del_TopologyAttr('bfactors') - >>> hasattr(u.atoms[:3], 'bfactors') + >>> import MDAnalysis as mda + >>> from MDAnalysis.tests.datafiles import PSF, DCD + >>> u = mda.Universe(PSF, DCD) + >>> u.add_TopologyAttr('tempfactors') + >>> hasattr(u.atoms[:3], 'tempfactors') + True + >>> + >>> u.del_TopologyAttr('tempfactors') + >>> hasattr(u.atoms[:3], 'tempfactors') False @@ -987,9 +997,12 @@ def add_Residue(self, segment=None, **attrs): Adding a new GLY residue, then placing atoms within it: - >>> newres = u.add_Residue(segment=u.segments[0], resid=42, resname='GLY') + >>> import MDAnalysis as mda + >>> from MDAnalysis.tests.datafiles import PSF, DCD + >>> u = mda.Universe(PSF, DCD) + >>> newres = u.add_Residue(segment=u.segments[0], resid=42, resname='GLY', resnum=0) >>> u.atoms[[1, 2, 3]].residues = newres - >>> u.select_atoms('resname GLY and resid 42') + >>> u.select_atoms('resname GLY and resid 42 and resnum 0') """ @@ -1407,6 +1420,7 @@ def from_smiles(cls, smiles, sanitize=True, addHs=True, -------- To create a Universe with 10 conformers of ethanol: + >>> from rdkit.Chem import AllChem >>> u = mda.Universe.from_smiles('CCO', numConfs=10) >>> u @@ -1416,7 +1430,7 @@ def from_smiles(cls, smiles, sanitize=True, addHs=True, To use a different conformer generation algorithm, like ETKDGv3: >>> u = mda.Universe.from_smiles('CCO', rdkit_kwargs=dict( - params=AllChem.ETKDGv3())) + ... params=AllChem.ETKDGv3())) >>> u.trajectory diff --git a/package/MDAnalysis/lib/transformations.py b/package/MDAnalysis/lib/transformations.py index f66d1c6ea2b..826b4efcbfe 100644 --- a/package/MDAnalysis/lib/transformations.py +++ b/package/MDAnalysis/lib/transformations.py @@ -109,6 +109,8 @@ Examples -------- +>>> from MDAnalysis.lib.transformations import * +>>> import numpy as np >>> alpha, beta, gamma = 0.123, -1.234, 2.345 >>> origin, xaxis, yaxis, zaxis = (0, 0, 0), (1, 0, 0), (0, 1, 0), (0, 0, 1) >>> I = identity_matrix() @@ -117,7 +119,7 @@ >>> Rz = rotation_matrix(gamma, zaxis) >>> R = concatenate_matrices(Rx, Ry, Rz) >>> euler = euler_from_matrix(R, 'rxyz') ->>> numpy.allclose([alpha, beta, gamma], euler) +>>> np.allclose([alpha, beta, gamma], euler) True >>> Re = euler_matrix(alpha, beta, gamma, 'rxyz') >>> is_same_transform(R, Re) @@ -136,14 +138,14 @@ >>> S = scale_matrix(1.23, origin) >>> T = translation_matrix((1, 2, 3)) >>> Z = shear_matrix(beta, xaxis, origin, zaxis) ->>> R = random_rotation_matrix(numpy.random.rand(3)) +>>> R = random_rotation_matrix(np.random.rand(3)) >>> M = concatenate_matrices(T, R, Z, S) >>> scale, shear, angles, trans, persp = decompose_matrix(M) ->>> numpy.allclose(scale, 1.23) +>>> np.allclose(scale, 1.23) True ->>> numpy.allclose(trans, (1, 2, 3)) +>>> np.allclose(trans, (1, 2, 3)) True ->>> numpy.allclose(shear, (0, math.tan(beta), 0)) +>>> np.allclose(shear, (0, math.tan(beta), 0)) True >>> is_same_transform(R, euler_matrix(axes='sxyz', *angles)) True @@ -175,6 +177,8 @@ def identity_matrix(): """Return 4x4 identity/unit matrix. + >>> from MDAnalysis.lib.transformations import identity_matrix + >>> import numpy as np >>> I = identity_matrix() >>> np.allclose(I, np.dot(I, I)) True @@ -190,6 +194,8 @@ def identity_matrix(): def translation_matrix(direction): """Return matrix to translate by direction vector. + >>> from MDAnalysis.lib.transformations import translation_matrix + >>> import numpy as np >>> v = np.random.random(3) - 0.5 >>> np.allclose(v, translation_matrix(v)[:3, 3]) True @@ -203,6 +209,9 @@ def translation_matrix(direction): def translation_from_matrix(matrix): """Return translation vector from translation matrix. + >>> from MDAnalysis.lib.transformations import (translation_matrix, + ... translation_from_matrix) + >>> import numpy as np >>> v0 = np.random.random(3) - 0.5 >>> v1 = translation_from_matrix(translation_matrix(v0)) >>> np.allclose(v0, v1) @@ -215,6 +224,8 @@ def translation_from_matrix(matrix): def reflection_matrix(point, normal): """Return matrix to mirror at plane defined by point and normal vector. + >>> from MDAnalysis.lib.transformations import reflection_matrix + >>> import numpy as np >>> v0 = np.random.random(4) - 0.5 >>> v0[3] = 1.0 >>> v1 = np.random.random(3) - 0.5 @@ -241,6 +252,9 @@ def reflection_matrix(point, normal): def reflection_from_matrix(matrix): """Return mirror plane point and normal vector from reflection matrix. + >>> from MDAnalysis.lib.transformations import (reflection_matrix, + ... reflection_from_matrix, is_same_transform) + >>> import numpy as np >>> v0 = np.random.random(3) - 0.5 >>> v1 = np.random.random(3) - 0.5 >>> M0 = reflection_matrix(v0, v1) @@ -270,6 +284,10 @@ def reflection_from_matrix(matrix): def rotation_matrix(angle, direction, point=None): """Return matrix to rotate about axis defined by point and direction. + >>> from MDAnalysis.lib.transformations import (rotation_matrix, + ... is_same_transform) + >>> import random, math + >>> import numpy as np >>> R = rotation_matrix(math.pi/2.0, [0, 0, 1], [1, 0, 0]) >>> np.allclose(np.dot(R, [0, 0, 0, 1]), [ 1., -1., 0., 1.]) True @@ -319,6 +337,10 @@ def rotation_matrix(angle, direction, point=None): def rotation_from_matrix(matrix): """Return rotation angle and axis from rotation matrix. + >>> from MDAnalysis.lib.transformations import (rotation_matrix, + ... is_same_transform, rotation_from_matrix) + >>> import random, math + >>> import numpy as np >>> angle = (random.random() - 0.5) * (2*math.pi) >>> direc = np.random.random(3) - 0.5 >>> point = np.random.random(3) - 0.5 @@ -361,6 +383,9 @@ def scale_matrix(factor, origin=None, direction=None): Use factor -1 for point symmetry. + >>> from MDAnalysis.lib.transformations import scale_matrix + >>> import random + >>> import numpy as np >>> v = (np.random.rand(4, 5) - 0.5) * 20.0 >>> v[3] = 1.0 >>> S = scale_matrix(-1.234) @@ -397,6 +422,10 @@ def scale_matrix(factor, origin=None, direction=None): def scale_from_matrix(matrix): """Return scaling factor, origin and direction from scaling matrix. + >>> from MDAnalysis.lib.transformations import (scale_matrix, + ... scale_from_matrix, is_same_transform) + >>> import random + >>> import numpy as np >>> factor = random.random() * 10 - 5 >>> origin = np.random.random(3) - 0.5 >>> direct = np.random.random(3) - 0.5 @@ -444,6 +473,9 @@ def projection_matrix(point, normal, direction=None, If pseudo is True, perspective projections will preserve relative depth such that Perspective = dot(Orthogonal, PseudoPerspective). + >>> from MDAnalysis.lib.transformations import (projection_matrix, + ... is_same_transform) + >>> import numpy as np >>> P = projection_matrix((0, 0, 0), (1, 0, 0)) >>> np.allclose(P[1:, 1:], np.identity(4)[1:, 1:]) True @@ -503,6 +535,9 @@ def projection_from_matrix(matrix, pseudo=False): Return values are same as arguments for projection_matrix function: point, normal, direction, perspective, and pseudo. + >>> from MDAnalysis.lib.transformations import (projection_matrix, + ... projection_from_matrix, is_same_transform) + >>> import numpy as np >>> point = np.random.random(3) - 0.5 >>> normal = np.random.random(3) - 0.5 >>> direct = np.random.random(3) - 0.5 @@ -586,6 +621,8 @@ def clip_matrix(left, right, bottom, top, near, far, perspective=False): Homogeneous coordinates transformed by the perspective clip matrix need to be dehomogenized (devided by w coordinate). + >>> from MDAnalysis.lib.transformations import clip_matrix + >>> import numpy as np >>> frustrum = np.random.rand(6) >>> frustrum[1] += frustrum[0] >>> frustrum[3] += frustrum[2] @@ -635,6 +672,9 @@ def shear_matrix(angle, direction, point, normal): given by the angle of P-P'-P", where P' is the orthogonal projection of P onto the shear plane. + >>> from MDAnalysis.lib.transformations import shear_matrix + >>> import random, math + >>> import numpy as np >>> angle = (random.random() - 0.5) * 4*math.pi >>> direct = np.random.random(3) - 0.5 >>> point = np.random.random(3) - 0.5 @@ -658,6 +698,10 @@ def shear_matrix(angle, direction, point, normal): def shear_from_matrix(matrix): """Return shear angle, direction and plane from shear matrix. + >>> from MDAnalysis.lib.transformations import (shear_matrix, + ... shear_from_matrix, is_same_transform) + >>> import random, math + >>> import numpy as np >>> angle = (random.random() - 0.5) * 4*math.pi >>> direct = np.random.random(3) - 0.5 >>> point = np.random.random(3) - 0.5 @@ -715,6 +759,9 @@ def decompose_matrix(matrix): Raise ValueError if matrix is of wrong type or degenerative. + >>> from MDAnalysis.lib.transformations import (translation_matrix, + ... decompose_matrix, scale_matrix, euler_matrix) + >>> import numpy as np >>> T0 = translation_matrix((1, 2, 3)) >>> scale, shear, angles, trans, persp = decompose_matrix(T0) >>> T1 = translation_matrix(trans) @@ -799,6 +846,10 @@ def compose_matrix(scale=None, shear=None, angles=None, translate=None, translate : translation vector along x, y, z axes perspective : perspective partition of matrix + >>> from MDAnalysis.lib.transformations import (compose_matrix, + ... decompose_matrix, is_same_transform) + >>> import math + >>> import numpy as np >>> scale = np.random.random(3) - 0.5 >>> shear = np.random.random(3) - 0.5 >>> angles = (np.random.random(3) - 0.5) * (2*math.pi) @@ -846,6 +897,8 @@ def orthogonalization_matrix(lengths, angles): The de-orthogonalization matrix is the inverse. + >>> from MDAnalysis.lib.transformations import orthogonalization_matrix + >>> import numpy as np >>> O = orthogonalization_matrix((10., 10., 10.), (90., 90., 90.)) >>> np.allclose(O[:3, :3], np.identity(3, float) * 10) True @@ -881,6 +934,11 @@ def superimposition_matrix(v0, v1, scaling=False, usesvd=True): The returned matrix performs rotation, translation and uniform scaling (if specified). + >>> from MDAnalysis.lib.transformations import (superimposition_matrix, + ... random_rotation_matrix, scale_matrix, translation_matrix, + ... concatenate_matrices) + >>> import random + >>> import numpy as np >>> v0 = np.random.rand(3, 10) >>> M = superimposition_matrix(v0, v0) >>> np.allclose(M, np.identity(4)) @@ -976,6 +1034,10 @@ def euler_matrix(ai, aj, ak, axes='sxyz'): ai, aj, ak : Euler's roll, pitch and yaw angles axes : One of 24 axis sequences as string or encoded tuple + >>> from MDAnalysis.lib.transformations import (euler_matrix + ... _AXES2TUPLE, _TUPLE2AXES) + >>> import math + >>> import numpy as np >>> R = euler_matrix(1, 2, 3, 'syxz') >>> np.allclose(np.sum(R[0]), -1.34786452) True @@ -1040,6 +1102,10 @@ def euler_from_matrix(matrix, axes='sxyz'): Note that many Euler angle triplets can describe one matrix. + >>> from MDAnalysis.lib.transformations import (euler_matrix, + ... euler_from_matrix, _AXES2TUPLE) + >>> import math + >>> import numpy as np >>> R0 = euler_matrix(1, 2, 3, 'syxz') >>> al, be, ga = euler_from_matrix(R0, 'syxz') >>> R1 = euler_matrix(al, be, ga, 'syxz') @@ -1094,6 +1160,8 @@ def euler_from_matrix(matrix, axes='sxyz'): def euler_from_quaternion(quaternion, axes='sxyz'): """Return Euler angles from quaternion for specified axis sequence. + >>> from MDAnalysis.lib.transformations import euler_from_quaternion + >>> import numpy as np >>> angles = euler_from_quaternion([0.99810947, 0.06146124, 0, 0]) >>> np.allclose(angles, [0.123, 0, 0]) True @@ -1108,6 +1176,7 @@ def quaternion_from_euler(ai, aj, ak, axes='sxyz'): ai, aj, ak : Euler's roll, pitch and yaw angles axes : One of 24 axis sequences as string or encoded tuple + >>> from MDAnalysis.lib.transformations import quaternion_from_euler >>> q = quaternion_from_euler(1, 2, 3, 'ryxz') >>> np.allclose(q, [0.435953, 0.310622, -0.718287, 0.444435]) True @@ -1162,6 +1231,8 @@ def quaternion_from_euler(ai, aj, ak, axes='sxyz'): def quaternion_about_axis(angle, axis): """Return quaternion for rotation about axis. + >>> from MDAnalysis.lib.transformations import quaternion_about_axis + >>> import numpy as np >>> q = quaternion_about_axis(0.123, (1, 0, 0)) >>> np.allclose(q, [0.99810947, 0.06146124, 0, 0]) True @@ -1181,6 +1252,9 @@ def quaternion_about_axis(angle, axis): def quaternion_matrix(quaternion): """Return homogeneous rotation matrix from quaternion. + >>> from MDAnalysis.lib.transformations import (identity_matrix, + ... quaternion_matrix, rotation_matrix) + >>> import numpy as np >>> M = quaternion_matrix([0.99810947, 0.06146124, 0, 0]) >>> np.allclose(M, rotation_matrix(0.123, (1, 0, 0))) True @@ -1213,6 +1287,10 @@ def quaternion_from_matrix(matrix, isprecise=False): If isprecise=True, the input matrix is assumed to be a precise rotation matrix and a faster algorithm is used. + >>> from MDAnalysis.lib.transformations import (identity_matrix, + ... quaternion_from_matrix, rotation_matrix, random_rotation_matrix, + ... is_same_transform, quaternion_matrix) + >>> import numpy as np >>> q = quaternion_from_matrix(identity_matrix(), True) >>> np.allclose(q, [1., 0., 0., 0.]) True @@ -1289,6 +1367,8 @@ def quaternion_from_matrix(matrix, isprecise=False): def quaternion_multiply(quaternion1, quaternion0): """Return multiplication of two quaternions. + >>> from MDAnalysis.lib.transformations import quaternion_multiply + >>> import numpy as np >>> q = quaternion_multiply([4, 1, -2, 3], [8, -5, 6, 7]) >>> np.allclose(q, [28, -44, -14, 48]) True @@ -1307,6 +1387,9 @@ def quaternion_multiply(quaternion1, quaternion0): def quaternion_conjugate(quaternion): """Return conjugate of quaternion. + >>> from MDAnalysis.lib.transformations import (random_quaternion, + ... quaternion_conjugate) + >>> import numpy as np >>> q0 = random_quaternion() >>> q1 = quaternion_conjugate(q0) >>> q1[0] == q0[0] and all(q1[1:] == -q0[1:]) @@ -1322,6 +1405,9 @@ def quaternion_conjugate(quaternion): def quaternion_inverse(quaternion): """Return inverse of quaternion. + >>> from MDAnalysis.lib.transformations import (random_quaternion, + ... quaternion_inverse) + >>> import numpy as np >>> q0 = random_quaternion() >>> q1 = quaternion_inverse(q0) >>> np.allclose(quaternion_multiply(q0, q1), [1, 0, 0, 0]) @@ -1334,6 +1420,7 @@ def quaternion_inverse(quaternion): def quaternion_real(quaternion): """Return real part of quaternion. + >>> from MDAnalysis.lib.transformations import quaternion_real >>> quaternion_real([3.0, 0.0, 1.0, 2.0]) 3.0 @@ -1344,6 +1431,7 @@ def quaternion_real(quaternion): def quaternion_imag(quaternion): """Return imaginary part of quaternion. + >>> from MDAnalysis.lib.transformations import quaternion_imag >>> quaternion_imag([3.0, 0.0, 1.0, 2.0]) [0.0, 1.0, 2.0] @@ -1354,6 +1442,10 @@ def quaternion_imag(quaternion): def quaternion_slerp(quat0, quat1, fraction, spin=0, shortestpath=True): r"""Return spherical linear interpolation between two quaternions. + >>> from MDAnalysis.lib.transformations import (random_quaternion, + ... quaternion_slerp) + >>> import math + >>> import numpy as np >>> q0 = random_quaternion() >>> q1 = random_quaternion() >>> q = quaternion_slerp(q0, q1, 0.0) @@ -1399,6 +1491,9 @@ def random_quaternion(rand=None): Three independent random variables that are uniformly distributed between 0 and 1. + >>> from MDAnalysis.lib.transformations import (random_quaternion, + ... vector_norm) + >>> import numpy as np >>> q = random_quaternion() >>> np.allclose(1.0, vector_norm(q)) True @@ -1431,6 +1526,8 @@ def random_rotation_matrix(rand=None): Three independent random variables that are uniformly distributed between 0 and 1 for each returned quaternion. + >>> from MDAnalysis.lib.transformations import random_rotation_matrix + >>> import numpy as np >>> R = random_rotation_matrix() >>> np.allclose(np.dot(R.T, R), np.identity(4)) True @@ -1442,6 +1539,8 @@ def random_rotation_matrix(rand=None): class Arcball(object): """Virtual Trackball Control. + >>> from MDAnalysis.lib.transformations import Arcball + >>> import numpy as np >>> ball = Arcball() >>> ball = Arcball(initial=np.identity(4)) >>> ball.place([320, 320], 320) @@ -1624,6 +1723,8 @@ def arcball_nearest_axis(point, axes): def vector_norm(data, axis=None, out=None): """Return length, i.e. eucledian norm, of ndarray along axis. + >>> from MDAnalysis.lib.transformations import vector_norm + >>> import numpy as np >>> v = np.random.random(3) >>> n = vector_norm(v) >>> np.allclose(n, np.linalg.norm(v)) @@ -1663,6 +1764,8 @@ def vector_norm(data, axis=None, out=None): def unit_vector(data, axis=None, out=None): """Return ndarray normalized by length, i.e. eucledian norm, along axis. + >>> from MDAnalysis.lib.transformations import unit_vector + >>> import numpy as np >>> v0 = np.random.random(3) >>> v1 = unit_vector(v0) >>> np.allclose(v1, v0 / np.linalg.norm(v0)) @@ -1707,6 +1810,8 @@ def unit_vector(data, axis=None, out=None): def random_vector(size): """Return array of random doubles in the half-open interval [0.0, 1.0). + >>> from MDAnalysis.lib.transformations import random_vector + >>> import numpy as np >>> v = random_vector(10000) >>> np.all(v >= 0.0) and np.all(v < 1.0) True @@ -1722,6 +1827,9 @@ def random_vector(size): def inverse_matrix(matrix): """Return inverse of square transformation matrix. + >>> from MDAnalysis.lib.transformations import (random_rotation_matrix, + ... inverse_matrix) + >>> import numpy as np >>> M0 = random_rotation_matrix() >>> M1 = inverse_matrix(M0.T) >>> np.allclose(M1, np.linalg.inv(M0.T)) @@ -1738,6 +1846,8 @@ def inverse_matrix(matrix): def concatenate_matrices(*matrices): """Return concatenation of series of transformation matrices. + >>> from MDAnalysis.lib.transformations import concatenate_matrices + >>> import numpy as np >>> M = np.random.rand(16).reshape((4, 4)) - 0.5 >>> np.allclose(M, concatenate_matrices(M)) True @@ -1754,6 +1864,9 @@ def concatenate_matrices(*matrices): def is_same_transform(matrix0, matrix1): """Return True if two matrices perform same transformation. + >>> from MDAnalysis.lib.transformations import (is_same_transform, + ... random_rotation_matrix) + >>> import numpy as np >>> is_same_transform(np.identity(4), np.identity(4)) True >>> is_same_transform(np.identity(4), random_rotation_matrix()) diff --git a/package/MDAnalysis/lib/util.py b/package/MDAnalysis/lib/util.py index ab36932dbf1..4bf21361a45 100644 --- a/package/MDAnalysis/lib/util.py +++ b/package/MDAnalysis/lib/util.py @@ -474,6 +474,8 @@ def greedy_splitext(p): Example ------- + + >>> from MDAnalysis.lib.util import greedy_splitext >>> greedy_splitext("/home/joe/protein.pdb.bz2") ('/home/joe/protein', '.pdb.bz2') @@ -1620,10 +1622,14 @@ def unique_rows(arr, return_index=False): -------- Remove dupicate rows from an array: + >>> import numpy as np + >>> from MDAnalysis.lib.util import unique_rows >>> a = np.array([[0, 1], [1, 2], [1, 2], [0, 1], [2, 3]]) >>> b = unique_rows(a) >>> b - array([[0, 1], [1, 2], [2, 3]]) + array([[0, 1], + [1, 2], + [2, 3]]) See Also -------- @@ -1680,6 +1686,8 @@ def blocks_of(a, n, m): Examples -------- + >>> import numpy as np + >>> from MDAnalysis.lib.util import blocks_of >>> arr = np.arange(16).reshape(4, 4) >>> view = blocks_of(arr, 2, 2) >>> view[:] = 100 @@ -1737,6 +1745,7 @@ def group_same_or_consecutive_integers(arr): list of :class:`numpy.ndarray` Examples + >>> import numpy as np >>> arr = np.array([ 2, 3, 4, 7, 8, 9, 10, 11, 15, 16]) >>> group_same_or_consecutive_integers(arr) [array([2, 3, 4]), array([ 7, 8, 9, 10, 11]), array([15, 16])] @@ -1791,6 +1800,7 @@ def ltruncate_int(value, ndigits): Examples -------- + >>> from MDAnalysis.lib.util import ltruncate_int >>> ltruncate_int(123, 2) 23 >>> ltruncate_int(1234, 5) @@ -1840,6 +1850,7 @@ def static_variables(**kwargs): Example ------- + >>> from MDAnalysis.lib.util import static_variables >>> @static_variables(msg='foo calls', calls=0) ... def foo(): ... foo.calls += 1 @@ -2011,6 +2022,10 @@ def check_coords(*coord_names, **options): Example ------- + >>> import numpy as np + >>> import MDAnalysis as mda + >>> from MDAnalysis.tests.datafiles import PSF, DCD + >>> from MDAnalysis.lib.util import check_coords >>> @check_coords('coords1', 'coords2', allow_atomgroup=True) ... def coordsum(coords1, coords2): ... assert coords1.dtype == np.float32 @@ -2027,12 +2042,18 @@ def check_coords(*coord_names, **options): >>> >>> # automatic handling of AtomGroups >>> u = mda.Universe(PSF, DCD) - >>> coordsum(u.atoms, u.select_atoms("index 1 to 10")) - ... + >>> try: + ... coordsum(u.atoms, u.select_atoms("index 1 to 10")) + ... except ValueError as err: + ... err + ValueError('coordsum(): coords1, coords2 must contain the same number of coordinates, got [3341, 10].') >>> >>> # automatic shape checking: - >>> coordsum(np.zeros(3), np.ones(6)) - ValueError: coordsum(): coords2.shape must be (3,) or (n, 3), got (6,). + >>> try: + ... coordsum(np.zeros(3), np.ones(6)) + ... except ValueError as err: + ... err + ValueError('coordsum(): coords2.shape must be (3,) or (n, 3), got (6,)') .. versionadded:: 0.19.0