Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

calculates distances using MDA.distance_array #31

Open
wants to merge 4 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from 2 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 5 additions & 0 deletions CHANGELOG.rst
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,11 @@ Upon version 0.8, and before version 1, SV2 major version increments are reflect
Changelog
=========

0.8.10 (2020-xx-xx)
-------------------

* cli ``dist`` now uses ``MDA.distance_array`` function. ISSUE #30

0.8.9 (2020-03-03)
------------------

Expand Down
2 changes: 2 additions & 0 deletions docs/conf.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,8 @@
'MDAnalysis',
'MDAnalysis.analysis',
'MDAnalysis.analysis.rms',
'MDAnalysis.lib',
'MDAnalysis.lib.distances',
'mdtraj',
'simtk.openmm.app',
'simtk.openmm',
Expand Down
41 changes: 25 additions & 16 deletions src/taurenmd/cli_distances.py
Original file line number Diff line number Diff line change
@@ -1,18 +1,23 @@
"""
# Calculates distances between centers of geometry of two selections.

Distance is given in 3D XYZ coordinate space units.
Distance is given in XYZ coordinate space units.

## Algorithm

Distance between centers of geometry is calculated by::
Calculates the distance between the centers of geometry of
two different selections, and for each frame in the trajectory, using::

np.linalg.norm(np.subtract(coord1, coord2))
[MDAnalysis.lib.distances.distance_array function](https://www.mdanalysis.org/docs/documentation_pages/lib/distances.html#MDAnalysis.lib.distances.distance_array)

Where, ``coord*`` are the centers of geometry of each atom selection
``-l1`` and ``-l2``, respectively.
Read further on [np.linalg.norm](https://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.norm.html)
and [np.subtract](https://docs.scipy.org/doc/numpy/reference/generated/numpy.subtract.html?highlight=subtract#numpy-subtract).
The result is a distance value for each frame.

### Note on Periodic Boxes

The distances calculated with this client do not account for artifacts

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is the idea here that users would have to create a new trajectory of unwrapped positions before they can proceed?

Copy link
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That is how I am using it, and my original idea. I have improved the description in the CLI docs (last commit). Having these two steps separate forces the users to have a granular comprehension and understanding of their operations. So I will keep it this way in this PR but I may (almost certainly) return to this in the future. Thanks!

originated from splits across the periodic box boundaries. To correct
for this, you whould ``unwrap`` or ``image`` the molecules first. See
``trajedit`` and ``imagemol`` interfaces for this.

## Examples

Expand Down Expand Up @@ -42,10 +47,11 @@
import functools

import numpy as np
from MDAnalysis.lib.distances import distance_array

import taurenmd.core as tcore
from taurenmd import _BANNER, Path, log
from taurenmd.libs import libcli, libio, libmda, libplot # noqa: F401
from taurenmd.libs import libcli, libio, libmda, libplot
from taurenmd.logger import S, T


Expand Down Expand Up @@ -125,18 +131,21 @@ def main(
atom_sel1 = u.select_atoms(sel1)
atom_sel2 = u.select_atoms(sel2)

distances = np.ones(len(u.trajectory[frame_slice]), dtype=np.float32)
distances = np.ones((len(u.trajectory[frame_slice]), 1), dtype=np.float64)

log.info(T('Calculating distances'))
# https://www.mdanalysis.org/MDAnalysisTutorial/atomgroups.html
# https://www.mdanalysis.org/docs/documentation_pages/core/groups.html#MDAnalysis.core.groups.AtomGroup.center_of_geometry
for i, _ts in enumerate(u.trajectory[frame_slice]):
for i, ts in enumerate(u.trajectory[frame_slice]):

distances[i] = np.linalg.norm(
np.subtract(
atom_sel1.center_of_geometry(),
atom_sel2.center_of_geometry(),
)
sel1_coords = atom_sel1.center_of_geometry()
sel2_coords = atom_sel2.center_of_geometry()

distance_array(
sel1_coords,
sel2_coords,
box=ts.dimensions,
result=distances[i:i + 1, :],
)

log.info(S('calculated a total of {} distances.', len(distances)))
Expand Down Expand Up @@ -171,7 +180,7 @@ def main(

libplot.param(
list(range(len(u.trajectory))[frame_slice]),
distances,
distances[:, 0],
**plotvars,
)

Expand Down