Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
Parameters
----------
xyz : np.ndarray, shape=(n_frames, n_atoms, 3)
The cartesian coordinates of the atoms to write. By convention for
this trajectory format, the lengths should be in units of angstroms.
types : np.ndarray, shape(3, )
The type of each particle.
"""
if not self._mode == 'w':
raise ValueError('write() is only available when file is opened '
'in mode="w"')
if not types:
# Make all particles the same type.
types = ['X' for _ in xrange(xyz.shape[1])]
xyz = ensure_type(xyz, np.float32, 3, 'xyz', can_be_none=False,
shape=(None, None, 3), warn_on_cast=False,
add_newaxis_on_deficient_ndim=True)
for i in range(xyz.shape[0]):
self._fh.write('{0}\n'.format(xyz.shape[1]))
self._fh.write("Created with MDTraj {0}, {1}\n".format(version, str(date.today())))
for j, coord in enumerate(xyz[i]):
self._fh.write('{0} {1:8.3f} {2:8.3f} {3:8.3f}\n'.format(
types[j], coord[0], coord[1], coord[2]))
ibonds : np.ndarray, shape=[n_bonds, 2], dtype=int
n_bonds x 2 array of indices, where each row is the index of two
atom who participate in a bond.
Returns
-------
iangles : np.ndarray, shape[n_angles, 3], dtype=int
n_angles x 3 array of indices, where each row is the index of three
atoms m,n,o such that n is bonded to both m and o.
"""
nx = import_('networkx')
graph = nx.from_edgelist(ibonds)
n_atoms = graph.number_of_nodes()
iangles = []
for i in xrange(n_atoms):
for (m, n) in combinations(graph.neighbors(i), 2):
# so now the there is a bond angle m-i-n
iangles.append((m, i, n))
return np.array(iangles)
elements = np.zeros(n_atoms, dtype='S1')
atom_names = [a.name for a in conf.top.atoms()]
for i in xrange(n_atoms):
# name of the element that is atom[i]
# take the first character of the AtomNames string,
# after stripping off any digits
elements[i] = atom_names[i].strip('123456789 ')[0]
if not elements[i] in COVALENT_RADII.keys():
raise ValueError("I don't know about this AtomName: {}".format(
atom_names[i]))
distance_mtx = squareform(pdist(xyz))
connectivity = []
for i in xrange(n_atoms):
for j in xrange(i + 1, n_atoms):
# Regular bonds are assigned to all pairs of atoms where
# the interatomic distance is less than or equal to 1.3 times the
# sum of their respective covalent radii.
d = distance_mtx[i, j]
if d < 1.3 * (COVALENT_RADII[elements[i]] + COVALENT_RADII[elements[j]]):
connectivity.append((i, j))
return np.array(connectivity)
the interatomic distance is less than or equal to 1.3 times the
sum of their respective covalent radii.
References
----------
Bakken and Helgaker, JCP Vol. 117, Num. 20 22 Nov. 2002
http://folk.uio.no/helgaker/reprints/2002/JCP117b_GeoOpt.pdf
"""
from scipy.spatial.distance import squareform, pdist
xyz = conf.xyz[0, :, :]
n_atoms = xyz.shape[0]
elements = np.zeros(n_atoms, dtype='S1')
atom_names = [a.name for a in conf.top.atoms()]
for i in xrange(n_atoms):
# name of the element that is atom[i]
# take the first character of the AtomNames string,
# after stripping off any digits
elements[i] = atom_names[i].strip('123456789 ')[0]
if not elements[i] in COVALENT_RADII.keys():
raise ValueError("I don't know about this AtomName: {}".format(
atom_names[i]))
distance_mtx = squareform(pdist(xyz))
connectivity = []
for i in xrange(n_atoms):
for j in xrange(i + 1, n_atoms):
# Regular bonds are assigned to all pairs of atoms where
# the interatomic distance is less than or equal to 1.3 times the
-------
idihedrals : np.ndarray, shape[n_dihedrals, 4], dtype=int
All sets of 4 atoms A,B,C,D such that A is bonded to B, B is bonded
to C, and C is bonded to D
"""
nx = import_('networkx')
graph = nx.from_edgelist(ibonds)
n_atoms = graph.number_of_nodes()
idihedrals = []
# TODO: CHECK FOR DIHEDRAL ANGLES THAT ARE 180 and recover
# conf : msmbuilder.Trajectory
# An msmbuilder trajectory, only the first frame will be used. This
# is used purely to make the check for angle(ABC) != 180.
for a in xrange(n_atoms):
for b in graph.neighbors(a):
for c in filter(lambda c: c not in [a, b], graph.neighbors(b)):
for d in filter(lambda d: d not in [a, b, c], graph.neighbors(c)):
idihedrals.append((a, b, c, d))
return np.array(idihedrals)
def _compute_traces(self):
if not self._centered:
self._center()
self._g = np.zeros((self.n_frames,), dtype=np.float32)
for i in xrange(self.n_frames):
for j in xrange(self.n_dims):
if self.axis_major:
self._g[i] += np.dot(self.cords[i, j, :],
self.cords[i, j, :])
else:
self._g[i] += np.dot(self.cords[i, :, j],
self.cords[i, :, j])
return
def _compute_traces(self):
if not self._centered:
self._center()
self._g = np.zeros((self.n_frames,), dtype=np.float32)
for i in xrange(self.n_frames):
for j in xrange(self.n_dims):
if self.axis_major:
self._g[i] += np.dot(self.cords[i, j, :],
self.cords[i, j, :])
else:
self._g[i] += np.dot(self.cords[i, :, j],
self.cords[i, :, j])
return
def save_pdb(self, filename, force_overwrite=True):
"""Save trajectory to RCSB PDB format
Parameters
----------
filename : str
filesystem path in which to save the trajectory
force_overwrite : bool, default=True
Overwrite anything that exists at filename, if its already there
"""
self._check_valid_unitcell()
with PDBTrajectoryFile(filename, 'w', force_overwrite=force_overwrite) as f:
for i in xrange(self.n_frames):
if self._have_unitcell:
f.write(convert(self._xyz[i], Trajectory._distance_unit, f.distance_unit),
self.topology,
modelIndex=i,
unitcell_lengths=convert(self.unitcell_lengths[i], Trajectory._distance_unit, f.distance_unit),
unitcell_angles=self.unitcell_angles[i])
else:
f.write(convert(self._xyz[i], Trajectory._distance_unit, f.distance_unit),
self.topology,
modelIndex=i)
xyz : np.ndarray, shape=(n_frames, n_atoms, 3), dtype=np.float32
cell_lengths : np.ndarray, None
The lengths (a,b,c) of the unit cell for each frame, or None if
the information is not present in the file.
cell_angles : np.ndarray, None
The angles (\alpha, \beta, \gamma) defining the unit cell for
each frame, or None if the information is not present in the file.
"""
if not self._mode == 'r':
raise ValueError('read() is only available when file is opened '
'in mode="r"')
if n_frames is None:
frame_counter = itertools.count()
else:
frame_counter = xrange(n_frames)
if stride is None:
stride = 1
all_coords, all_lengths, all_angles = [], [], []
for _ in frame_counter:
try:
frame_coords, frame_lengths, frame_angles = self._read()
if atom_indices is not None:
frame_coords = frame_coords[atom_indices, :]
except _EOF:
break
all_coords.append(frame_coords)
all_lengths.append(frame_lengths)
all_angles.append(frame_angles)