gwinc/ifo/noises.py

@@ -141,21 +141,9 @@ def precomp_suspension(f, ifo):
         ifo.Suspension.VHCoupling = Struct()
         ifo.Suspension.VHCoupling.theta = ifo.Infrastructure.Length / const.R_earth
 
-    if ifo.Suspension.Type == 'BQuad':
-        material = 'Silicon'
-    else:
-        material = 'Silica'
-    hForce, vForce, hTable, vTable = suspension.suspQuad(f, ifo.Suspension, material)
-    try:
-        # full TF (conventional)
-        ifo.Suspension.hForce = hForce.fullylossy
-        ifo.Suspension.vForce = vForce.fullylossy
-        # TFs with each stage lossy
-        ifo.Suspension.hForce_singlylossy = hForce.singlylossy
-        ifo.Suspension.vForce_singlylossy = vForce.singlylossy
-    except:
-        ifo.Suspension.hForce = hForce
-        ifo.Suspension.vForce = vForce
+    hForce, vForce, hTable, vTable = suspension.suspQuad(f, ifo.Suspension)
+    ifo.Suspension.hForce = hForce
+    ifo.Suspension.vForce = vForce
     ifo.Suspension.hTable = hTable
     ifo.Suspension.vTable = vTable
 

gwinc/noise/suspensionthermal.py

@@ -5,11 +5,12 @@ from __future__ import division
 import numpy as np
 from numpy import pi, imag
 
-from .. import const
+from ..const import kB
+from ..suspension import getJointParams
 
 
 def suspension_thermal(f, sus):
-    """Suspention thermal noise for a single suspended test mass
+    """Suspension thermal noise.
 
     :f: frequency array in Hz
     :sus: gwinc suspension structure
@@ -17,67 +18,45 @@ def suspension_thermal(f, sus):
     :returns: displacement noise power spectrum at :f:
 
     :Temp: must either be set for each stage individually, or globally
-    in :sus.Temp:.  The latter will be preferred if specified, so if
-    you wish to use per-stage tempurature you must remove :sus.Temp:.
+    in :sus.Temp:.  If both are specified, :sus.Temp: is interpreted as
+    the temperature of the upper joint of the top stage.
 
     Assumes suspension transfer functions and V-H coupling have been
     pre-calculated and populated into the relevant `sus` fields.
 
     """
-    # Assign Physical Constants
-    kB = const.kB
 
     # and vertical to beamline coupling angle
     theta = sus.VHCoupling.theta
 
-    noise = np.zeros((1, f.size))
+    noise = np.zeros_like(f)
 
-    # if the temperature is uniform along the suspension
-    if 'Temp' in sus:
-        ##########################################################
-        # Suspension TFs
-        ##########################################################
+    ##########################################################
+    # Suspension TFs
+    ##########################################################
 
-        hForce = sus.hForce
-        vForce = sus.vForce
+    hForce = sus.hForce
+    vForce = sus.vForce
 
-        ##########################################################
-        # Thermal Noise Calculation
-        ##########################################################
-
-        # convert to beam line motion
-        #  theta is squared because we rotate by theta into the suspension
-        #  basis, and by theta to rotate back to the beam line basis
-        dxdF = hForce + theta**2 * vForce
-
-        # thermal noise (m^2/Hz) for one suspension
-        w = 2*pi*f
-        noise = 4 * kB * sus.Temp * abs(imag(dxdF)) / w
-
-    # if the temperature is set for each suspension stage
-    else:
-        ##########################################################
-        # Suspension TFs
-        ##########################################################
-
-        hForce = sus.hForce_singlylossy[:, :]
-        vForce = sus.vForce_singlylossy[:, :]
-
-        ##########################################################
-        # Thermal Noise Calculation
-        ##########################################################
-
-        dxdF = np.zeros(hForce.shape, dtype=complex)
-        for n, stage in enumerate(reversed(sus.Stage)):
-            # add up the contribution from each stage
+    ##########################################################
+    # Thermal Noise Calculation
+    ##########################################################
 
+    # Here the suspension stage list is reversed so that the last stage
+    # in the list is the test mass
+    stages = sus.Stage[::-1]
+    w = 2*pi*f
+    dxdF = np.zeros_like(hForce)  # complex valued
+    for n, stage in enumerate(stages):
+        # add up the contribution from each joint
+        jointParams = getJointParams(sus, n)
+        for (isLower, Temp, wireMat, bladeMat) in jointParams:
             # convert to beam line motion.  theta is squared because
             # we rotate by theta into the suspension basis, and by
             # theta to rotate back to the beam line basis
-            dxdF[n, :] = hForce[n, :] + theta**2 * vForce[n, :]
-
-            # thermal noise (m^2/Hz) for one suspension
-            w = 2*pi*f
-            noise += 4 * kB * stage.Temp * abs(imag(dxdF[n, :])) / w
+            m = 2*n + isLower
+            dxdF[m, :] = hForce[m, :] + theta**2 * vForce[m, :]
+            noise += 4 * kB * Temp * abs(imag(dxdF[m, :])) / w
 
-    return np.squeeze(noise)
+    # thermal noise (m^2/Hz) for one suspension
+    return noise