From 080b9fb9d2416070a2ee7a8403efb8513022875f Mon Sep 17 00:00:00 2001
From: Jameson Graef Rollins <jrollins@finestructure.net>
Date: Thu, 19 Sep 2019 23:55:09 -0700
Subject: [PATCH] suspension thermal noise simplification

new suspension_thermal() function for calculating thermal noise of a
single suspended optic.

ifo SuspensionThermal Noise class modified to handle specific case of
Fabry-Perot Michelson.
---
 gwinc/ifo/noises.py              |  3 ++-
 gwinc/noise/suspensionthermal.py | 43 +++++++++++++++++---------------
 2 files changed, 25 insertions(+), 21 deletions(-)

diff --git a/gwinc/ifo/noises.py b/gwinc/ifo/noises.py
index f203c827..81d56cd0 100644
--- a/gwinc/ifo/noises.py
+++ b/gwinc/ifo/noises.py
@@ -103,7 +103,8 @@ class SuspensionThermal(nb.Noise):
     )
 
     def calc(self):
-        return noise.suspensionthermal.susptherm(self.freq, self.ifo)
+        n = noise.suspensionthermal.suspension_thermal(self.freq, self.ifo.Suspension)
+        return n * 4 * self.ifo.gwinc.dhdl_sqr
 
 
 class CoatingBrownian(nb.Noise):
diff --git a/gwinc/noise/suspensionthermal.py b/gwinc/noise/suspensionthermal.py
index dcc34618..c5b4910a 100644
--- a/gwinc/noise/suspensionthermal.py
+++ b/gwinc/noise/suspensionthermal.py
@@ -1,39 +1,45 @@
+'''Functions to calculate suspension thermal noise
+
+'''
 from __future__ import division
-from numpy import pi, imag, zeros
 import numpy as np
+from numpy import pi, imag, zeros
 
 from .. import const
 
 
-def susptherm(f, ifo):
-    """Suspention thermal noise.
+def suspension_thermal(f, sus):
+    """Suspention thermal noise for a single suspended test mass
 
-    'Temp' must either be set for each stage individually, or globally
-    in ifo.Suspension.Temp.  The latter will be preferred if
-    specified, so if you wish to use per-stage tempurature you must
-    remove Suspension.Temp.
+    :f: frequency array in Hz
+    :sus: gwinc suspension structure
+
+    :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:.
 
     Assumes suspension transfer functions and V-H coupling have been
-    pre-calculated and populated into the relevant `ifo` struct
-    fields.
+    pre-calculated and populated into the relevant `sus` fields.
 
     """
     # Assign Physical Constants
     kB = const.kB
 
     # and vertical to beamline coupling angle
-    theta = ifo.Suspension.VHCoupling.theta
+    theta = sus.VHCoupling.theta
 
     noise = zeros((1, f.size))
 
     # if the temperature is uniform along the suspension
-    if 'Temp' in ifo.Suspension:
+    if 'Temp' in sus:
         ##########################################################
         # Suspension TFs
         ##########################################################
 
-        hForce = ifo.Suspension.hForce
-        vForce = ifo.Suspension.vForce
+        hForce = sus.hForce
+        vForce = sus.vForce
 
         ##########################################################
         # Thermal Noise Calculation
@@ -46,7 +52,7 @@ def susptherm(f, ifo):
 
         # thermal noise (m^2/Hz) for one suspension
         w = 2*pi*f
-        noise = 4 * kB * ifo.Suspension.Temp * abs(imag(dxdF)) / w
+        noise = 4 * kB * sus.Temp * abs(imag(dxdF)) / w
 
     # if the temperature is set for each suspension stage
     else:
@@ -54,15 +60,15 @@ def susptherm(f, ifo):
         # Suspension TFs
         ##########################################################
 
-        hForce = ifo.Suspension.hForce_singlylossy[:, :]
-        vForce = ifo.Suspension.vForce_singlylossy[:, :]
+        hForce = sus.hForce_singlylossy[:, :]
+        vForce = sus.vForce_singlylossy[:, :]
 
         ##########################################################
         # Thermal Noise Calculation
         ##########################################################
 
         dxdF = zeros(hForce.shape, dtype=complex)
-        for n, stage in enumerate(reversed(ifo.Suspension.Stage)):
+        for n, stage in enumerate(reversed(sus.Stage)):
             # add up the contribution from each stage
 
             # convert to beam line motion.  theta is squared because
@@ -74,7 +80,4 @@ def susptherm(f, ifo):
             w = 2*pi*f
             noise += 4 * kB * stage.Temp * abs(imag(dxdF[n, :])) / w
 
-    # 4 masses, turn into gravitational wave strain
-    noise *= 4 * ifo.gwinc.dhdl_sqr
-
     return np.squeeze(noise)
-- 
GitLab