Adds unwrap keyword to Asphericity method

package/CHANGELOG

@@ -60,6 +60,7 @@ Enhancements
   * added unwrap keyword to center_of_geometry (PR #2279)
   * added unwrap keyword to center_of_mass (PR #2286)
   * added unwrap keyword to moment_of_inertia (PR #2287)
+  * added unwrap keyword to asphericity (PR #2290)
 
 Changes
   * added official support for Python 3.7 (PR #1963)

package/MDAnalysis/core/topologyattrs.py

@@ -1028,7 +1028,8 @@ def shape_parameter(group, **kwargs):
         ('shape_parameter', shape_parameter))
 
     @warn_if_not_unique
-    def asphericity(group, pbc=None):
+    @check_pbc_and_unwrap
+    def asphericity(group, pbc=None, unwrap=None, compound='group'):
         """Asphericity.
 
         See [Dima2004b]_ for background information.
@@ -1039,6 +1040,10 @@ def asphericity(group, pbc=None):
             If ``True``, move all atoms within the primary unit cell before
             calculation. If ``None`` use value defined in
             MDAnalysis.core.flags['use_pbc']
+        unwrap : bool, optional
+            If ``True``, compounds will be unwrapped before computing their centers.
+        compound : {'group', 'segments', 'residues', 'molecules', 'fragments'}, optional
+            Which type of component to keep together during unwrapping.
 
         Note
         ----
@@ -1059,19 +1064,24 @@ def asphericity(group, pbc=None):
 
         .. versionadded:: 0.7.7
         .. versionchanged:: 0.8 Added *pbc* keyword
+        .. versionchanged:: 0.20.0 Added *unwrap* and *compound* parameter
 
         """
         atomgroup = group.atoms
         if pbc is None:
             pbc = flags['use_pbc']
         masses = atomgroup.masses
 
+        com = atomgroup.center_of_mass(pbc=pbc, unwrap=unwrap, compound=compound)
+        if compound is not 'group':
+            com = (com * group.masses[:, None]).sum(axis=0) / group.masses.sum()
+
         if pbc:
-            recenteredpos = (atomgroup.pack_into_box(inplace=False) -
-                             atomgroup.center_of_mass(pbc=True))
+            recenteredpos = (atomgroup.pack_into_box(inplace=False) - com)
+        elif unwrap:
+            recenteredpos = (atomgroup.unwrap(inplace=False) - com)
         else:
-            recenteredpos = (atomgroup.positions -
-                             atomgroup.center_of_mass(pbc=False))
+            recenteredpos = (atomgroup.positions - com)
 
         tensor = np.zeros((3, 3))
         for x in range(recenteredpos.shape[0]):

testsuite/MDAnalysisTests/core/test_atomgroup.py

@@ -848,17 +848,17 @@ def ag(self):
     @pytest.fixture()
     def ref_noUnwrap_residues(self):
         return {
-            'COG': np.array([[21.356, 28.52 , 36.762],
-                             [32.062, 36.16 , 27.679],
+            'COG': np.array([[21.356, 28.52, 36.762],
+                             [32.062, 36.16, 27.679],
                              [27.071, 29.997, 28.506]], dtype=np.float32),
             'COM': np.array([[21.286, 28.407, 36.629],
                              [31.931, 35.814, 27.916],
-                             [26.817, 29.41 , 29.05 ]]),
+                             [26.817, 29.41, 29.05]]),
             'MOI': np.array([
                 [7333.79167791, -211.8997285, -721.50785456],
                 [-211.8997285, 7059.07470427, -91.32156884],
                 [-721.50785456, -91.32156884, 6509.31735029]]),
-
+            'Asph': 0.02060121,
         }
 
     @pytest.fixture()
@@ -872,7 +872,8 @@ def ref_Unwrap_residues(self):
                              [2.045, 28.243, 54.127]], dtype=np.float32),
             'MOI': np.array([[16747.486, -1330.489,  2938.243],
                              [-1330.489, 19315.253,  3306.212],
-                             [ 2938.243,  3306.212,  8990.481]])
+                             [ 2938.243,  3306.212,  8990.481]]),
+            'Asph': 0.2969491080,
         }
 
     @pytest.fixture()
@@ -884,6 +885,7 @@ def ref_noUnwrap(self):
                 [0.0, 0.0, 0.0],
                 [0.0, 98.6542, 0.0],
                 [0.0, 0.0, 98.65421327]]),
+            'Asph': 1.0,
         }
 
     @pytest.fixture()
@@ -895,17 +897,20 @@ def ref_Unwrap(self):
                 [0.0, 0.0, 0.0],
                 [0.0, 132.673, 0.0],
                 [0.0, 0.0, 132.673]]),
+            'Asph': 1.0,
         }
 
     def test_default_residues(self, ag, ref_noUnwrap_residues):
         assert_almost_equal(ag.center_of_geometry(compound='residues'), ref_noUnwrap_residues['COG'], self.prec)
         assert_almost_equal(ag.center_of_mass(compound='residues'), ref_noUnwrap_residues['COM'], self.prec)
         assert_almost_equal(ag.moment_of_inertia(compound='residues'), ref_noUnwrap_residues['MOI'], self.prec)
+        assert_almost_equal(ag.asphericity(compound='residues'), ref_noUnwrap_residues['Asph'], self.prec)
 
     def test_UnWrapFlag_residues(self, ag, ref_Unwrap_residues):
         assert_almost_equal(ag.center_of_geometry(unwrap=True, compound='residues'), ref_Unwrap_residues['COG'], self.prec)
         assert_almost_equal(ag.center_of_mass(unwrap=True, compound='residues'), ref_Unwrap_residues['COM'], self.prec)
         assert_almost_equal(ag.moment_of_inertia(unwrap=True, compound='residues'), ref_Unwrap_residues['MOI'], self.prec)
+        assert_almost_equal(ag.asphericity(unwrap=True, compound='residues'), ref_Unwrap_residues['Asph'], self.prec)
 
     def test_default(self, ref_noUnwrap):
         u = UnWrapUniverse(is_triclinic=False)
@@ -916,6 +921,7 @@ def test_default(self, ref_noUnwrap):
         assert_almost_equal(group.center_of_geometry(), ref_noUnwrap['COG'], self.prec)
         assert_almost_equal(group.center_of_mass(), ref_noUnwrap['COM'], self.prec)
         assert_almost_equal(group.moment_of_inertia(), ref_noUnwrap['MOI'], self.prec)
+        assert_almost_equal(group.asphericity(), ref_noUnwrap['Asph'], self.prec)
 
     def test_UnWrapFlag(self, ref_Unwrap):
         u = UnWrapUniverse(is_triclinic=False)
@@ -925,6 +931,8 @@ def test_UnWrapFlag(self, ref_Unwrap):
         assert_almost_equal(group.center_of_geometry(unwrap=True), ref_Unwrap['COG'], self.prec)
         assert_almost_equal(group.center_of_mass(unwrap=True), ref_Unwrap['COM'], self.prec)
         assert_almost_equal(group.moment_of_inertia(unwrap=True), ref_Unwrap['MOI'], self.prec)
+        assert_almost_equal(group.asphericity(unwrap=True), ref_Unwrap['Asph'], self.prec)
+
 
 class TestPBCFlag(object):