Add vibration_analyze.breit_wigner_resonant_frequency().

Now people using the module don't have to think or read comments to
determine the resonant frequency from fitted parameters ;).

diff --git a/calibcant/vibration_analyze.py b/calibcant/vibration_analyze.py
index 1844dcf4e78be771f0e1645ed71078de68839a77..a1125b6a5c15746a5efab3771b862850a2337e49 100644 (file)
--- a/calibcant/vibration_analyze.py
+++ b/calibcant/vibration_analyze.py
@@ -454,6 +454,11 @@ def breit_wigner_area(A, B, C):
     #  <V(t)**2> = (pi*C) / (2*B*A**2)
     return (_numpy.pi * C) / (2 * B * A**2)
 
+def breit_wigner_resonant_frequency(A, B):
+    if (B**2 >= 2*A**2):
+        return 0  # over- or critically-damped
+    return _numpy.sqrt(A**2 - B**2/2)
+
 def save(filename, group='/', raw=None, config=None,
          deflection_channel_config=None, processed=None):
     specs = [
@@ -543,9 +548,8 @@ def plot(deflection=None, freq_axis=None, power=None, A=None, B=None,
             freq_axes.plot(freq_axis, fitdata, 'b-')
             noisefloor = D + 0*freq_axis;
             freq_axes.plot(freq_axis, noisefloor)
-
-            if B**2 < 2*A**2:
-                res_freq = _numpy.sqrt(A**2 - B**2/2)
+            res_freq = breit_wigner_resonant_frequency(A=A, B=B)
+            if res_freq > 0:
                 freq_axes.axvline(res_freq, color='b', zorder=-1)
 
         freq_axes.set_title('power spectral density %s' % timestamp)