Give all measured values to FunctionalForm

- no longer limited to using mean, var

pycgtool/bondset.py

@@ -63,22 +63,21 @@ class Bond:
         if not self.values:
             raise ValueError("No bonds were measured between atoms {0}".format(self.atoms))
 
-        mean, var = stat_moments(self.values)
-
+        values = np.array(self.values)
         if len(self.atoms) > 2:
-            mean_tmp = math.radians(mean)
-            var *= (math.pi * math.pi) / (180. * 180.)
-        else:
-            mean_tmp = mean
+            values = np.radians(values)
 
-        self.eqm = mean
         with np.errstate(divide="raise"):
+            self.eqm = self._func_form.eqm(values, temp)
             try:
-                self.fconst = self._func_form(mean_tmp, var, temp)
+                self.fconst = self._func_form.fconst(values, temp)
             except FloatingPointError:
                 # Happens when variance is 0, i.e. we only have one value
                 self.fconst = float("inf")
 
+        if len(self.atoms) > 2:
+            self.eqm = math.degrees(self.eqm)
+
     def __repr__(self):
         try:
             return "<Bond containing atoms {0} with r_0 {1:.3f} and force constant {2:.3e}>".format(

pycgtool/functionalforms.py

 import abc
 import math
 
+import numpy as np
+
 from pycgtool.util import SimpleEnum
 
 
@@ -53,15 +55,26 @@ class FunctionalForm(object, metaclass=abc.ABCMeta):
 
     New functional forms must define a static __call__ method.
     """
+    @staticmethod
+    def eqm(values, temp):
+        """
+        Calculate equilibrium value.
+        May be overridden by functional forms.
+
+        :param values: Measured internal coordinate values from which to calculate equilibrium value
+        :param temp: Temperature of simulation
+        :return: Calculated equilibrium value
+        """
+        return np.nanmean(values)
+
     @staticmethod
     @abc.abstractstaticmethod
-    def __call__(mean, var, temp):
+    def fconst(values, temp):
         """
         Calculate force constant.
         Abstract static method to be defined by all functional forms.
 
-        :param mean: Mean of internal coordinate distribution
-        :param var: Variance of internal coordinate distribution
+        :param values: Measured internal coordinate values from which to calculate force constant
         :param temp: Temperature of simulation
         :return: Calculated force constant
         """
@@ -70,31 +83,34 @@ class FunctionalForm(object, metaclass=abc.ABCMeta):
 
 class Harmonic(FunctionalForm):
     @staticmethod
-    def __call__(mean, var, temp):
+    def fconst(values, temp):
         rt = 8.314 * temp / 1000.
+        var = np.nanvar(values)
         return rt / var
 
 
 class CosHarmonic(FunctionalForm):
     @staticmethod
-    def __call__(mean, var, temp):
+    def fconst(values, temp):
         rt = 8.314 * temp / 1000.
+        mean = CosHarmonic.eqm(values, temp)
+        var = np.nanvar(values)
         return rt / (math.sin(mean)**2 * var)
 
 
 class MartiniDefaultLength(FunctionalForm):
     @staticmethod
-    def __call__(mean, var, temp):
+    def fconst(values, temp):
         return 1250.
 
 
 class MartiniDefaultAngle(FunctionalForm):
     @staticmethod
-    def __call__(mean, var, temp):
+    def fconst(values, temp):
         return 25.
 
 
 class MartiniDefaultDihedral(FunctionalForm):
     @staticmethod
-    def __call__(mean, var, temp):
+    def fconst(values, temp):
         return 50.

test/test_functionalforms.py

@@ -14,12 +14,17 @@ class FunctionalFormTest(unittest.TestCase):
     def test_functional_form_new(self):
         class TestFunc(FunctionalForm):
             @staticmethod
-            def __call__(mean, var, temp):
-                return "TestResult"
+            def eqm(values, temp):
+                return "TestResultEqm"
+
+            @staticmethod
+            def fconst(values, temp):
+                return "TestResultFconst"
 
         funcs = FunctionalForms()
         self.assertIn("TestFunc", funcs)
-        self.assertEqual("TestResult", funcs.TestFunc(None, None, None))
+        self.assertEqual("TestResultEqm", funcs.TestFunc.eqm(None, None))
+        self.assertEqual("TestResultFconst", funcs.TestFunc.fconst(None, None))
 
 if __name__ == '__main__':
     unittest.main()