Fixed merge conflicts and updated CHANGELOG

package/CHANGELOG

@@ -28,7 +28,7 @@ Metadata
 API Changes
 
   * rmsd doesn't superimpose by default anymore. The superposition
-    is controlled by the 'superposition' keyword now. (see issue #562)
+    is controlled by the 'superposition' keyword now. (see issue #562, #822)
 
 Enhancements
 
@@ -37,11 +37,13 @@ Enhancements
 
 Fixes
 
+  * rmsd now returns proper value when given array of weights (Issue #814)
   * change_release now finds number and dev (Issue #776)
   * test_shear_from_matrix doesn't fail for MKL builds anymore (Issue #757)
   * HEADER and TITLE now appear just once in the PDB. (Issue #741) (PR #761)
   * MOL2 files without substructure section can now be read (Issue #816)
   * MOL2 files can be written without substructure section (Issue #816)
+  * GRO files with an incomplete set of velocities can now be read (Issue #820)
 
 Changes
 

package/MDAnalysis/__init__.py

@@ -84,14 +84,14 @@
 
 Calculate the CA end-to-end distance (in angstroem)::
   >>> import numpy as np
-  >>> coord = ca.coordinates()
+  >>> coord = ca.positions
   >>> v = coord[-1] - coord[0]   # last Ca minus first one
   >>> np.sqrt(np.dot(v, v,))
   10.938133
 
 Define a function eedist():
   >>> def eedist(atoms):
-  ...     coord = atoms.coordinates()
+  ...     coord = atoms.positions
   ...     v = coord[-1] - coord[0]
   ...     return sqrt(dot(v, v,))
   ...
@@ -162,7 +162,7 @@
 del logging
 
 # DeprecationWarnings are loud by default
-warnings.simplefilter('always', DeprecationWarning)
+warnings.simplefilter('once', DeprecationWarning)
 
 
 from . import units

package/MDAnalysis/analysis/align.py

@@ -217,8 +217,8 @@ def rotation_matrix(a, b, weights=None):
     :meth:`MDAnalysis.core.AtomGroup.AtomGroup.rotate` to generate a rotated
     selection, e.g. ::
 
-    >>> R = rotation_matrix(A.select_atoms('backbone').coordinates(),
-    >>>                     B.select_atoms('backbone').coordinates())[0]
+    >>> R = rotation_matrix(A.select_atoms('backbone').positions,
+    >>>                     B.select_atoms('backbone').positions)[0]
     >>> A.atoms.rotate(R)
     >>> A.atoms.write("rotated.pdb")
 
@@ -356,8 +356,8 @@ def alignto(mobile, reference, select="all", mass_weighted=False,
         ref_com = ref_atoms.center_of_geometry()
         mobile_com = mobile_atoms.center_of_geometry()
 
-    ref_coordinates = ref_atoms.coordinates() - ref_com
-    mobile_coordinates = mobile_atoms.coordinates() - mobile_com
+    ref_coordinates = ref_atoms.positions - ref_com
+    mobile_coordinates = mobile_atoms.positions - mobile_com
 
     old_rmsd = rms.rmsd(mobile_coordinates, ref_coordinates)
 
@@ -494,10 +494,10 @@ def rms_fit_trj(traj, reference, select='all', filename=None, rmsdfile=None, pre
 
     # reference centre of mass system
     ref_com = ref_atoms.center_of_mass()
-    ref_coordinates = ref_atoms.coordinates() - ref_com
+    ref_coordinates = ref_atoms.positions - ref_com
 
     # allocate the array for selection atom coords
-    traj_coordinates = traj_atoms.coordinates().copy()
+    traj_coordinates = traj_atoms.positions.copy()
 
     # RMSD timeseries
     nframes = len(frames)
@@ -516,7 +516,7 @@ def rms_fit_trj(traj, reference, select='all', filename=None, rmsdfile=None, pre
         # shift coordinates for rotation fitting
         # selection is updated with the time frame
         x_com = traj_atoms.center_of_mass().astype(np.float32)
-        traj_coordinates[:] = traj_atoms.coordinates() - x_com
+        traj_coordinates[:] = traj_atoms.positions - x_com
 
         # Need to transpose coordinates such that the coordinate array is
         # 3xN instead of Nx3. Also qcp requires that the dtype be float64

package/MDAnalysis/analysis/contacts.py

@@ -606,7 +606,7 @@ def __init__(self, *args, **kwargs):
 
         # compute reference contacts
         dref = MDAnalysis.lib.distances.distance_array(
-            self.references[0].coordinates(), self.references[1].coordinates())
+            self.references[0].positions, self.references[1].positions)
         self.qref = self.qarray(dref)
         self.nref = self.qref.sum()
 
@@ -682,7 +682,7 @@ def run(self, store=True, force=False, start=0, stop=None, step=1, **kwargs):
             for ts in self.universe.trajectory[start:stop:step]:
                 frame = ts.frame
                 # use pre-allocated distance array to save a little bit of time
-                MDAnalysis.lib.distances.distance_array(A.coordinates(), B.coordinates(), result=self.d)
+                MDAnalysis.lib.distances.distance_array(A.positions, B.positions, result=self.d)
                 self.qarray(self.d, out=self.q)
                 n1, q1 = self.qN(self.q, out=self._qtmp)
                 self.qavg += self.q

package/MDAnalysis/analysis/density.py

@@ -470,7 +470,7 @@ def current_coordinates():
         group = u.select_atoms(atomselection)
 
         def current_coordinates():
-            return group.coordinates()
+            return group.positions
 
     coord = current_coordinates()
     logger.info("Selected {0:d} atoms out of {1:d} atoms ({2!s}) from {3:d} total.".format(coord.shape[0], len(u.select_atoms(atomselection)), atomselection, len(u.atoms)))
@@ -600,14 +600,14 @@ def notwithin_coordinates(cutoff=cutoff):
                 ns_w = NS.AtomNeighborSearch(solvent)  # build kd-tree on solvent (N_w > N_protein)
                 solvation_shell = ns_w.search_list(protein, cutoff)  # solvent within CUTOFF of protein
                 group = MDAnalysis.core.AtomGroup.AtomGroup(set_solvent - set(solvation_shell))  # bulk
-                return group.coordinates()
+                return group.positions
         else:
             def notwithin_coordinates(cutoff=cutoff):
                 # acts as '<solvent> WITHIN <cutoff> OF <protein>'
                 # must update every time step
                 ns_w = NS.AtomNeighborSearch(solvent)  # build kd-tree on solvent (N_w > N_protein)
                 group = ns_w.search_list(protein, cutoff)  # solvent within CUTOFF of protein
-                return group.coordinates()
+                return group.positions
     else:
         # slower distance matrix based (calculate all with all distances first)
         dist = np.zeros((len(solvent), len(protein)), dtype=np.float64)
@@ -620,8 +620,8 @@ def notwithin_coordinates(cutoff=cutoff):
             aggregatefunc = np.any
 
         def notwithin_coordinates(cutoff=cutoff):
-            s_coor = solvent.coordinates()
-            p_coor = protein.coordinates()
+            s_coor = solvent.positions
+            p_coor = protein.positions
             # Does water i satisfy d[i,j] > r for ALL j?
             d = MDAnalysis.analysis.distances.distance_array(s_coor, p_coor, box=box, result=dist)
             return s_coor[aggregatefunc(compare(d, cutoff), axis=1)]
@@ -751,7 +751,7 @@ def __init__(self, pdb, delta=1.0, atomselection='resname HOH and name O',
         """
         u = MDAnalysis.as_Universe(pdb)
         group = u.select_atoms(atomselection)
-        coord = group.coordinates()
+        coord = group.positions
         logger.info("Selected {0:d} atoms ({1!s}) out of {2:d} total.".format(coord.shape[0], atomselection, len(u.atoms)))
         smin = np.min(coord, axis=0) - padding
         smax = np.max(coord, axis=0) + padding

package/MDAnalysis/analysis/gnm.py

@@ -210,7 +210,7 @@ def generate_kirchoff(self):
         the cutoff. Returns the resulting matrix
         '''
         #ca = self.u.select_atoms(self.selection)
-        positions = self.ca.coordinates()
+        positions = self.ca.positions
 
         natoms = len(positions)
 
@@ -320,7 +320,7 @@ def __init__(self, universe, selection='protein', cutoff=4.5, ReportVector=None,
     def generate_kirchoff(self):
         natoms = len(self.ca.atoms)
         nresidues = len(self.ca.residues)
-        positions = self.ca.coordinates()
+        positions = self.ca.positions
         [res_positions, grid, low_x, low_y, low_z] = generate_grid(positions, self.cutoff)
         residue_index_map = [resnum for [resnum, residue] in enumerate(self.ca.residues) for atom in residue]
         matrix = np.zeros((nresidues, nresidues), "float")

package/MDAnalysis/analysis/hbonds/hbond_analysis.py

@@ -916,15 +916,15 @@ def _get_timestep():
 
     def calc_angle(self, d, h, a):
         """Calculate the angle (in degrees) between two atoms with H at apex."""
-        v1 = h.pos - d.pos
-        v2 = h.pos - a.pos
+        v1 = h.position - d.position
+        v2 = h.position - a.position
         if np.all(v1 == v2):
             return 0.0
         return np.rad2deg(angle(v1, v2))
 
     def calc_eucl_distance(self, a1, a2):
         """Calculate the Euclidean distance between two atoms. """
-        return norm(a2.pos - a1.pos)
+        return norm(a2.position - a1.position)
 
     def generate_table(self):
         """Generate a normalised table of the results.

package/MDAnalysis/analysis/helanal.py

@@ -290,7 +290,7 @@ def helanal_trajectory(universe, selection="name CA", start=None, end=None, begi
             if trajectory.time > finish:
                 break
 
-        ca_positions = ca.coordinates()
+        ca_positions = ca.positions
         twist, bending_angles, height, rnou, origins, local_helix_axes, local_screw_angles = \
             main_loop(ca_positions, ref_axis=ref_axis)
 
@@ -515,7 +515,7 @@ def helanal_main(pdbfile, selection="name CA", start=None, end=None, ref_axis=No
     logger.info("Analysing %d/%d residues", ca.n_atoms, universe.atoms.n_residues)
 
     twist, bending_angles, height, rnou, origins, local_helix_axes, local_screw_angles = \
-        main_loop(ca.coordinates(), ref_axis=ref_axis)
+        main_loop(ca.positions, ref_axis=ref_axis)
 
     #TESTED- origins are correct
     #print current_origin

package/MDAnalysis/analysis/nuclinfo.py

@@ -135,7 +135,7 @@ def wc_pair(universe, i, bp, seg1="SYSTEM", seg2="SYSTEM"):
     wc_dist = universe.select_atoms("(segid {0!s} and resid {1!s} and name {2!s}) "
                                     "or (segid {3!s} and resid {4!s} and name {5!s}) "
                                     .format(seg1, i, a1, seg2, bp, a2))
-    wc = mdamath.norm(wc_dist[0].pos - wc_dist[1].pos)
+    wc = mdamath.norm(wc_dist[0].position - wc_dist[1].position)
     return wc
 
 
@@ -168,7 +168,7 @@ def minor_pair(universe, i, bp, seg1="SYSTEM", seg2="SYSTEM"):
     c2o2_dist = universe.select_atoms("(segid {0!s} and resid {1!s} and name {2!s}) "
                                       "or (segid {3!s} and resid {4!s} and name {5!s})"
                                       .format(seg1, i, a1, seg2, bp, a2))
-    c2o2 = mdamath.norm(c2o2_dist[0].pos - c2o2_dist[1].pos)
+    c2o2 = mdamath.norm(c2o2_dist[0].position - c2o2_dist[1].position)
     return c2o2
 
 
@@ -208,7 +208,7 @@ def major_pair(universe, i, bp, seg1="SYSTEM", seg2="SYSTEM"):
     no_dist = universe.select_atoms("(segid {0!s} and resid {1!s} and name {2!s}) "
                                     "or (segid {3!s} and resid {4!s} and name {5!s}) "
                                     .format(seg1, i, a1, seg2, bp, a2))
-    major = mdamath.norm(no_dist[0].pos - no_dist[1].pos)
+    major = mdamath.norm(no_dist[0].position - no_dist[1].position)
     return major
 
 
@@ -234,11 +234,11 @@ def phase_cp(universe, seg, i):
     atom4 = universe.select_atoms(" atom {0!s} {1!s} C3\' ".format(seg, i))
     atom5 = universe.select_atoms(" atom {0!s} {1!s} C4\' ".format(seg, i))
 
-    data1 = atom1.coordinates()
-    data2 = atom2.coordinates()
-    data3 = atom3.coordinates()
-    data4 = atom4.coordinates()
-    data5 = atom5.coordinates()
+    data1 = atom1.positions
+    data2 = atom2.positions
+    data3 = atom3.positions
+    data4 = atom4.positions
+    data5 = atom5.positions
 
     r0 = (data1 + data2 + data3 + data4 + data5) * (1.0 / 5.0)
     r1 = data1 - r0

package/MDAnalysis/analysis/rms.py

@@ -204,7 +204,11 @@ def rmsd(a, b, weights=None, center=False, superposition=False):
         return qcp.CalcRMSDRotationalMatrix(a.T, b.T, N, None,
                                             relative_weights)
     else:
-        return np.sqrt(np.sum((a - b) ** 2) / a.size)
+        if weights is not None:
+            return np.sqrt(np.sum(relative_weights[:, np.newaxis]
+                * (( a - b ) ** 2)) / N)
+        else:
+            return np.sqrt(np.sum((a - b) ** 2) / N)
 
 
 def _process_selection(select):
@@ -443,7 +447,7 @@ def run(self, start=None, stop=None, step=None,
         ref_coordinates_T_64 = ref_coordinates.T.astype(np.float64)
 
         # allocate the array for selection atom coords
-        traj_coordinates = traj_atoms.coordinates().copy()
+        traj_coordinates = traj_atoms.positions.copy()
 
         if self.groupselections_atoms:
             # Only carry out a rotation if we want to calculate secondary
@@ -468,7 +472,7 @@ def run(self, start=None, stop=None, step=None,
             # shift coordinates for rotation fitting
             # selection is updated with the time frame
             x_com = traj_atoms.center_of_mass().astype(np.float32)
-            traj_coordinates[:] = traj_atoms.coordinates() - x_com
+            traj_coordinates[:] = traj_atoms.positions - x_com
 
             rmsd[k, :2] = ts.frame, trajectory.time
 

package/MDAnalysis/coordinates/CRD.py

@@ -145,7 +145,7 @@ def write(self, selection, frame=None):
                 frame = 0  # should catch cases when we are analyzing a single PDB (?)
 
         atoms = selection.atoms  # make sure to use atoms (Issue 46)
-        coor = atoms.coordinates()  # can write from selection == Universe (Issue 49)
+        coor = atoms.positions  # can write from selection == Universe (Issue 49)
         with util.openany(self.filename, 'w') as self.crd:
             self._TITLE("FRAME " + str(frame) + " FROM " + str(u.trajectory.filename))
             self._TITLE("")

package/MDAnalysis/coordinates/GRO.py

@@ -130,12 +130,15 @@ def _read_first_frame(self):
             # (dependent upon the GRO file precision)
             first_atomline = grofile.readline()
             cs = first_atomline[25:].find('.') + 1
-            has_velocities = first_atomline[20:].count('.') > 3
+
+            # Always try, and maybe add them later
+            velocities = np.zeros((n_atoms, 3), dtype=np.float32)
 
             self.ts = ts = self._Timestep(n_atoms,
-                                          velocities=has_velocities,
                                           **self._ts_kwargs)
 
+            missed_vel = False
+
             grofile.seek(0)
             for pos, line in enumerate(grofile, start=-2):
                 # 2 header lines, 1 box line at end
@@ -144,13 +147,19 @@ def _read_first_frame(self):
                     continue
                 if pos < 0:
                     continue
-                for i in range(3):
-                    ts._pos[pos, i] = float(line[20 + cs*i: 20 + cs*(i+1)])
 
-                if not has_velocities:
-                    continue
-                for i, j in enumerate(range(3, 6)):
-                    ts._velocities[pos, i] = float(line[20+cs*j:20+cs*(j+1)])
+                ts._pos[pos] = [line[20 + cs*i:20 + cs*(i+1)] for i in range(3)]
+                try:
+                    velocities[pos] = [line[20 + cs*i:20 + cs*(i+1)] for i in range(3, 6)]
+                except ValueError:
+                    # Remember that we got this error
+                    missed_vel = True
+
+        if np.any(velocities):
+            ts.velocities = velocities
+            if missed_vel:
+                warnings.warn("Not all velocities were present.  "
+                              "Unset velocities set to zero.")
 
         self.ts.frame = 0  # 0-based frame number
 

package/MDAnalysis/coordinates/PDBQT.py

@@ -293,7 +293,7 @@ def write(self, selection, frame=None):
         self.TITLE("FRAME " + str(frame) + " FROM " + str(u.trajectory.filename))
         self.CRYST1(self.convert_dimensions_to_unitcell(u.trajectory.ts))
         atoms = selection.atoms  # make sure to use atoms (Issue 46)
-        coor = atoms.coordinates()  # can write from selection == Universe (Issue 49)
+        coor = atoms.positions  # can write from selection == Universe (Issue 49)
 
         # check if any coordinates are illegal (coordinates are already in Angstroem per package default)
         if not self.has_valid_coordinates(self.pdb_coor_limits, coor):

package/MDAnalysis/coordinates/PQR.py

@@ -231,7 +231,7 @@ def write(self, selection, frame=None):
                 frame = 0  # should catch cases when we are analyzing a single frame(?)
 
         atoms = selection.atoms  # make sure to use atoms (Issue 46)
-        coordinates = atoms.coordinates()  # can write from selection == Universe (Issue 49)
+        coordinates = atoms.positions  # can write from selection == Universe (Issue 49)
         if self.convert_units:
             self.convert_pos_to_native(coordinates)  # inplace because coordinates is already a copy