diff --git a/gwinc/noise/suspensionthermal.py b/gwinc/noise/suspensionthermal.py
index 1dfdafdd64aa1b21ef60382749a54c62645b2a8e..1c16ab609a9ea3528d7e99a1225935005578988d 100644
--- a/gwinc/noise/suspensionthermal.py
+++ b/gwinc/noise/suspensionthermal.py
@@ -18,7 +18,9 @@ def suspR(f, ifo):
theta = ifo.Suspension.VHCoupling.theta
noise = zeros((1, f.size))
- if np.isscalar(Temp) or len(Temp) == 1: # if the temperature is uniform along the suspension
+
+ # if the temperature is uniform along the suspension
+ if np.isscalar(Temp) or len(Temp) == 1:
##########################################################
# Suspension TFs
##########################################################
@@ -38,7 +40,8 @@ def suspR(f, ifo):
w = 2*pi*f
noise = 4 * kB * Temp * abs(imag(dxdF)) / w
- else: # if the temperature is set for each suspension stage
+ # if the temperature is set for each suspension stage
+ else:
##########################################################
# Suspension TFs
##########################################################
@@ -61,7 +64,7 @@ def suspR(f, ifo):
# thermal noise (m^2/Hz) for one suspension
w = 2*pi*f
noise = noise + 4 * kB * Temp[ii] * abs(imag(dxdF[ii,:])) / w
-
+
if ifo.Suspension.Type == 'Quad_MB':
raise Exception('not dealing with Quad_MB suspensions currently')
#mbquadlite2lateral_20090819TM_TN;
@@ -298,7 +301,7 @@ def suspBQuad(f, ifo):
delta_h3 = delta_h3*Temp[2]/(rho_st*C_st)
# solutions to equations of motion
- B = np.array([ 0, 0, 0, 1]).T
+ B = np.array([ 0, 0, 0, 1]).T
w = 2*pi * f
# thermoelastic correction factor, silica
@@ -356,9 +359,9 @@ def suspBQuad(f, ifo):
# Equations of motion for the system
###############################################################
- m_list=np.hstack((m1, m2, m3, m4)) # array of the mass
- kh_list=np.vstack((kh1, kh2, kh3, kh4)) # array of the horiz spring constants
- kv_list=np.vstack((kv1, kv2, kv3, kv4)) # array of the vert spring constants
+ m_list = np.hstack((m1, m2, m3, m4)) # array of the mass
+ kh_list = np.vstack((kh1, kh2, kh3, kh4)) # array of the horiz spring constants
+ kv_list = np.vstack((kv1, kv2, kv3, kv4)) # array of the vert spring constants
# Calculate TFs turning on the loss of each stage one by one
hForce = mat_struct()
@@ -385,14 +388,14 @@ def suspBQuad(f, ifo):
vForce.singlylossy[ii,:], vTable = calc_transfer_functions(Av, B, k_list, f)
# horizontal
- k_list=kh_list # all of the losses are on
+ k_list = kh_list # all of the losses are on
# construct Eq of motion matrix
Ah = construct_eom_matrix(k_list, m_list, f)
# calculate TFs
hForce.fullylossy, hTable = calc_transfer_functions(Ah, B, k_list, f)
# vertical
- k_list=kv_list # all of the losses are on
+ k_list = kv_list # all of the losses are on
# construct Eq of motion matrix
Av = construct_eom_matrix(k_list, m_list, f)
# calculate TFs
@@ -474,7 +477,6 @@ def suspQuad(f, ifo):
g = scipy.constants.g
kB = scipy.constants.k
- Temp = ifo.Suspension.Temp
Temp = ifo.Suspension.Temp
if np.isscalar(Temp) or len(Temp) == 1:
Temp = [Temp, Temp, Temp, Temp]
@@ -673,9 +675,9 @@ def suspQuad(f, ifo):
# Equations of motion for the system
###############################################################
- m_list=np.hstack((m1, m2, m3, m4)) # array of the mass
- kh_list=np.vstack((kh1, kh2, kh3, kh4)) # array of the horiz spring constants
- kv_list=np.vstack((kv1, kv2, kv3, kv4)) # array of the vert spring constants
+ m_list = np.hstack((m1, m2, m3, m4)) # array of the mass
+ kh_list = np.vstack((kh1, kh2, kh3, kh4)) # array of the horiz spring constants
+ kv_list = np.vstack((kv1, kv2, kv3, kv4)) # array of the vert spring constants
# Calculate TFs turning on the loss of each stage one by one
hForce = mat_struct()
@@ -702,14 +704,14 @@ def suspQuad(f, ifo):
vForce.singlylossy[ii,:], vTable = calc_transfer_functions(Av, B, k_list, f)
# horizontal
- k_list=kh_list # all of the losses are on
+ k_list = kh_list # all of the losses are on
# construct Eq of motion matrix
Ah = construct_eom_matrix(k_list, m_list, f)
# calculate TFs
hForce.fullylossy, hTable = calc_transfer_functions(Ah, B, k_list, f)
# vertical
- k_list=kv_list # all of the losses are on
+ k_list = kv_list # all of the losses are on
# construct Eq of motion matrix
Av = construct_eom_matrix(k_list, m_list, f)
# calculate TFs
diff --git a/gwinc/util.py b/gwinc/util.py
index e6ff1b613ccbdfe5cae482f73201040ffc25190a..87bf6e178530defcd2d1425510574f832788f944 100644
--- a/gwinc/util.py
+++ b/gwinc/util.py
@@ -3,37 +3,19 @@ from numpy import pi, sqrt
from scipy.io.matlab.mio5_params import mat_struct
from scipy.io import loadmat
import scipy.special
-from noise.coatingthermal import getCoatDopt
-def SpotSizes(g1, g2, L, lambda_):
- """calculates spot sizes using FP cavity parameters
- all parameters are in SI units
-
-
- ex. [w1, w2] = SpotSizes(0.9, 0.81, 3995, 1550e-9);"""
-
- w1 = (L*lambda_/pi) * sqrt(g2/g1/(1-g1*g2))
- w1 = sqrt(w1)
-
- w2 = (L*lambda_/pi) * sqrt(g1/g2/(1-g1*g2))
- w2 = sqrt(w2)
-
- w0 = g1*g2*(1-g1*g2) / (g1 + g2 - 2*g1*g2)**2
- w0 = (L*lambda_/pi) * sqrt(w0)
- w0 = sqrt(w0)
-
- return w1, w2, w0
+from .noise.coatingthermal import getCoatDopt
-def precompIFO(ifo, PRfixed):
- """add precomputed data to the IFO model
+def precompIFO(ifo, PRfixed=0):
+ """Add precomputed data to the IFO model.
To prevent recomputation of these precomputed data, if the
ifo argument contains ifo.gwinc.PRfixed, and this matches
the argument PRfixed, no changes are made.
- (mevans June 2008)"""
-
+ """
+
# check PRfixed
#if 'gwinc' in ifo.__dict__:
# && isfield(ifo.gwinc, 'PRfixed') && ifo.gwinc.PRfixed == PRfixed
@@ -41,6 +23,11 @@ def precompIFO(ifo, PRfixed):
ifo.gwinc = mat_struct()
ifo.gwinc.PRfixed = PRfixed
+ ################################# DERIVED TEMP
+
+ if 'Temp' not in ifo.Materials.Substrate.__dict__:
+ ifo.Materials.Substrate.Temp = ifo.Constants.Temp
+
################################# DERIVED OPTICS VALES
# Calculate optics' parameters
ifo.Materials.MirrorVolume = pi*ifo.Materials.MassRadius**2 * \
@@ -86,19 +73,18 @@ def precompIFO(ifo, PRfixed):
# Seismic noise term is saved in a .mat file defined in your respective IFOModel.m
# It is loaded here and put into the ifo structure.
- darmsei = loadmat(ifo.Seismic.darmSeiSusFile)
-
- ifo.Seismic.darmseis_f = darmsei['darmseis_f'][0]
- ifo.Seismic.darmseis_x = darmsei['darmseis_x'][0]
+ if 'darmSeiSusFile' in ifo.Seismic.__dict__:
+ darmsei = loadmat(ifo.Seismic.darmSeiSusFile)
+ ifo.Seismic.darmseis_f = darmsei['darmseis_f'][0]
+ ifo.Seismic.darmseis_x = darmsei['darmseis_x'][0]
return ifo
def precompPower(ifo, PRfixed):
- """Find and return power on ifo beamsplitter limited by
- thermal lensing effects, finesse and power recycling factor of ifo"""
+ """Compute power on beamsplitter and finesse and power recycling factor.
- # Modified to use a fixed PR factor Rana, Feb 07
+ """
# constants
c = scipy.constants.c
@@ -150,10 +136,9 @@ def precompPower(ifo, PRfixed):
def precompQuantum(ifo):
+ """Compute quantum noise parameters.
- ##########################################
- # input numbers
- ##########################################
+ """
# physical constants
hbar = scipy.constants.hbar # J s
@@ -177,10 +162,6 @@ def precompQuantum(ifo):
fGammaArm = gammaArm / (2*pi)
rSR = sqrt(1 - Tsrm)
- ##########################################
- # IFO equations
- ##########################################
-
# fSQL as defined in D&D paper (eq 33 in P1400018 and/or PRD paper)
tSR = sqrt(Tsrm)
fSQL = (1/(2*pi))*(8/c)*sqrt((Parm*w0)/(m*Titm))*(tSR/(1+rSR))
@@ -189,3 +170,25 @@ def precompQuantum(ifo):
fGammaIFO = fGammaArm * ((1 + rSR) / (1 - rSR))
return fSQL, fGammaIFO, fGammaArm
+
+
+def SpotSizes(g1, g2, L, lambda_):
+ """Calculate spot sizes using FP cavity parameters.
+
+ All parameters are in SI units.
+
+ ex. w1, w2 = SpotSizes(0.9, 0.81, 3995, 1550e-9)
+
+ """
+
+ w1 = (L*lambda_/pi) * sqrt(g2/g1/(1-g1*g2))
+ w1 = sqrt(w1)
+
+ w2 = (L*lambda_/pi) * sqrt(g1/g2/(1-g1*g2))
+ w2 = sqrt(w2)
+
+ w0 = g1*g2*(1-g1*g2) / (g1 + g2 - 2*g1*g2)**2
+ w0 = (L*lambda_/pi) * sqrt(w0)
+ w0 = sqrt(w0)
+
+ return w1, w2, w0