1# -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; -*- 2# vim: tabstop=4 expandtab shiftwidth=4 softtabstop=4 3# 4# MDAnalysis --- https://www.mdanalysis.org 5# Copyright (c) 2006-2017 The MDAnalysis Development Team and contributors 6# (see the file AUTHORS for the full list of names) 7# 8# Released under the GNU Public Licence, v2 or any higher version 9# 10# Please cite your use of MDAnalysis in published work: 11# 12# R. J. Gowers, M. Linke, J. Barnoud, T. J. E. Reddy, M. N. Melo, S. L. Seyler, 13# D. L. Dotson, J. Domanski, S. Buchoux, I. M. Kenney, and O. Beckstein. 14# MDAnalysis: A Python package for the rapid analysis of molecular dynamics 15# simulations. In S. Benthall and S. Rostrup editors, Proceedings of the 15th 16# Python in Science Conference, pages 102-109, Austin, TX, 2016. SciPy. 17# 18# N. Michaud-Agrawal, E. J. Denning, T. B. Woolf, and O. Beckstein. 19# MDAnalysis: A Toolkit for the Analysis of Molecular Dynamics Simulations. 20# J. Comput. Chem. 32 (2011), 2319--2327, doi:10.1002/jcc.21787 21# 22# 23 24 25""" 26Distance analysis --- :mod:`MDAnalysis.analysis.distances` 27========================================================== 28 29This module provides functions to rapidly compute distances between 30atoms or groups of atoms. 31 32:func:`dist` and :func:`between` can take atom groups that do not even 33have to be from the same :class:`~MDAnalysis.core.universe.Universe`. 34 35See Also 36-------- 37:mod:`MDAnalysis.lib.distances` 38""" 39 40from __future__ import absolute_import 41 42__all__ = ['distance_array', 'self_distance_array', 43 'contact_matrix', 'dist', 'between'] 44 45import numpy as np 46import scipy.sparse 47 48from MDAnalysis.lib.distances import distance_array, self_distance_array 49from MDAnalysis.lib.c_distances import contact_matrix_no_pbc, contact_matrix_pbc 50from MDAnalysis.lib.NeighborSearch import AtomNeighborSearch 51from MDAnalysis.lib.distances import calc_bonds 52 53 54import warnings 55import logging 56logger = logging.getLogger("MDAnalysis.analysis.distances") 57 58 59def contact_matrix(coord, cutoff=15.0, returntype="numpy", box=None): 60 '''Calculates a matrix of contacts. 61 62 There is a fast, high-memory-usage version for small systems 63 (*returntype* = 'numpy'), and a slower, low-memory-usage version for 64 larger systems (*returntype* = 'sparse'). 65 66 If *box* dimensions are passed then periodic boundary conditions 67 are applied. 68 69 Parameters 70 --------- 71 coord : array 72 Array of coordinates of shape ``(N, 3)`` and dtype float32. 73 cutoff : float, optional, default 15 74 Particles within `cutoff` are considered to form a contact. 75 returntype : string, optional, default "numpy" 76 Select how the contact matrix is returned. 77 * ``"numpy"``: return as an ``(N. N)`` :class:`numpy.ndarray` 78 * ``"sparse"``: return as a :class:`scipy.sparse.lil_matrix` 79 box : array-like or ``None``, optional, default ``None`` 80 Simulation cell dimensions in the form of 81 :attr:`MDAnalysis.trajectory.base.Timestep.dimensions` when 82 periodic boundary conditions should be taken into account for 83 the calculation of contacts. 84 85 Returns 86 ------- 87 array or sparse matrix 88 The contact matrix is returned in a format determined by the `returntype` 89 keyword. 90 91 See Also 92 -------- 93 :mod:`MDAnalysis.analysis.contacts` for native contact analysis 94 95 96 .. versionchanged:: 0.11.0 97 Keyword *suppress_progmet* and *progress_meter_freq* were removed. 98 ''' 99 100 if returntype == "numpy": 101 adj = (distance_array(coord, coord, box=box) < cutoff) 102 return adj 103 elif returntype == "sparse": 104 # Initialize square List of Lists matrix of dimensions equal to number 105 # of coordinates passed 106 sparse_contacts = scipy.sparse.lil_matrix((len(coord), len(coord)), dtype='bool') 107 if box is not None: 108 # with PBC 109 contact_matrix_pbc(coord, sparse_contacts, box, cutoff) 110 else: 111 # without PBC 112 contact_matrix_no_pbc(coord, sparse_contacts, cutoff) 113 return sparse_contacts 114 115 116def dist(A, B, offset=0, box=None): 117 """Return distance between atoms in two atom groups. 118 119 The distance is calculated atom-wise. The residue ids are also 120 returned because a typical use case is to look at CA distances 121 before and after an alignment. Using the `offset` keyword one can 122 also add a constant offset to the resids which facilitates 123 comparison with PDB numbering. 124 125 Arguments 126 --------- 127 A, B : AtomGroup 128 :class:`~MDAnalysis.core.groups.AtomGroup` with the 129 same number of atoms 130 offset : integer or tuple, optional, default 0 131 An integer `offset` is added to *resids_A* and *resids_B* (see 132 below) in order to produce PDB numbers. 133 134 If `offset` is :class:`tuple` then ``offset[0]`` is added to 135 *resids_A* and ``offset[1]`` to *resids_B*. Note that one can 136 actually supply numpy arrays of the same length as the atom 137 group so that an individual offset is added to each resid. 138 139 Returns 140 ------- 141 resids_A : array 142 residue ids of the `A` group (possibly changed with `offset`) 143 resids_B : array 144 residue ids of the `B` group (possibly changed with `offset`) 145 distances : array 146 distances between the atoms 147 """ 148 149 if A.atoms.n_atoms != B.atoms.n_atoms: 150 raise ValueError("AtomGroups A and B do not have the same number of atoms") 151 try: 152 off_A, off_B = offset 153 except (TypeError, ValueError): 154 off_A = off_B = int(offset) 155 residues_A = np.array(A.resids) + off_A 156 residues_B = np.array(B.resids) + off_B 157 158 d = calc_bonds(A.positions, B.positions, box) 159 return np.array([residues_A, residues_B, d]) 160 161 162def between(group, A, B, distance): 163 """Return sub group of `group` that is within `distance` of both `A` and `B` 164 165 This function is not aware of periodic boundary conditions. 166 167 Can be used to find bridging waters or molecules in an interface. 168 169 Similar to "*group* and (AROUND *A* *distance* and AROUND *B* *distance*)". 170 171 Parameters 172 ---------- 173 group : AtomGroup 174 Find members of `group` that are between `A` and `B` 175 A : AtomGroup 176 B : AtomGroup 177 `A` and `B` are the groups of atoms between which atoms in 178 `group` are searched for. The function works is more 179 efficient if `group` is bigger than either `A` or `B`. 180 distance : float 181 maximum distance for an atom to be counted as in the vicinity of 182 `A` or `B` 183 184 Returns 185 ------- 186 AtomGroup 187 :class:`~MDAnalysis.core.groups.AtomGroup` of atoms that 188 fulfill the criterion 189 190 191 .. versionadded: 0.7.5 192 193 """ 194 ns_group = AtomNeighborSearch(group) 195 resA = set(ns_group.search(A, distance)) 196 resB = set(ns_group.search(B, distance)) 197 return sum(sorted(resB.intersection(resA))) 198