diff --git a/gwinc/ifo/noises.py b/gwinc/ifo/noises.py index 16958042c267542ad965b3570f934cf046b041e9..e1768bf2b5c6844dee42b7c729d028663dbe122b 100644 --- a/gwinc/ifo/noises.py +++ b/gwinc/ifo/noises.py @@ -1,3 +1,4 @@ +import copy import numpy as np from numpy import pi, sin, exp, sqrt @@ -10,14 +11,33 @@ from .. import suspension ################################################## -def coating_thickness(ifo, optic): - optic = ifo.Optics.get(optic) + +def coating_thickness(materials, optic): if 'CoatLayerOpticalThickness' in optic: return optic['CoatLayerOpticalThickness'] T = optic.Transmittance dL = optic.CoatingThicknessLown dCap = optic.CoatingThicknessCap - return noise.coatingthermal.getCoatDopt(ifo.Materials, T, dL, dCap=dCap) + return noise.coatingthermal.getCoatDopt(materials, T, dL, dCap=dCap) + + +def mirror_struct(ifo, tm): + """Create a "mirror" Struct for a LIGO core optic + + This is a copy of the ifo.Materials Struct, containing Substrate + and Coating sub-Structs, as well as some basic geometrical + properties of the optic. + + """ + # NOTE: we deepcopy this struct since we'll be modifying it (and + # it would otherwise be stored by reference) + mirror = copy.deepcopy(ifo.Materials) + optic = ifo.Optics.get(tm) + mirror.Coating.dOpt = coating_thickness(ifo.Materials, optic) + mirror.update(optic) + mirror.MassVolume = pi * mirror.MassRadius**2 * mirror.MassThickness + mirror.MirrorMass = mirror.MassVolume * mirror.Substrate.MassDensity + return mirror def arm_cavity(ifo): @@ -118,15 +138,6 @@ def dhdl(f, armlen): ################################################## -def precomp_mirror(f, ifo): - ifo.Materials.MirrorVolume = \ - pi * ifo.Materials.MassRadius**2 \ - * ifo.Materials.MassThickness - ifo.Materials.MirrorMass = \ - ifo.Materials.MirrorVolume \ - * ifo.Materials.Substrate.MassDensity - - def precomp_suspension(f, ifo): pc = Struct() pc.VHCoupling = Struct() @@ -142,10 +153,10 @@ def precomp_suspension(f, ifo): return pc -def precomp_coating(f, ifo): +def precomp_mirror(f, ifo): pc = Struct() - pc.dOpt_ITM = coating_thickness(ifo, 'ITM') - pc.dOpt_ETM = coating_thickness(ifo, 'ETM') + pc.ITM = mirror_struct(ifo, 'ITM') + pc.ETM = mirror_struct(ifo, 'ETM') return pc @@ -186,10 +197,10 @@ class QuantumVacuum(nb.Noise): color='#ad03de', ) - @nb.precomp(precomp_mirror) + @nb.precomp(mirror=precomp_mirror) @nb.precomp(power=precomp_power) - def calc(self, power): - return noise.quantum.shotrad(self.freq, self.ifo, power=power) + def calc(self, mirror, power): + return noise.quantum.shotrad(self.freq, self.ifo, mirror.ETM.MirrorMass, power) class StandardQuantumLimit(nb.Noise): @@ -202,9 +213,9 @@ class StandardQuantumLimit(nb.Noise): linestyle=":", ) - @nb.precomp(precomp_mirror) - def calc(self): - return 8 * const.hbar / (self.ifo.Materials.MirrorMass * (2 * np.pi * self.freq) ** 2) + @nb.precomp(mirror=precomp_mirror) + def calc(self, mirror): + return 8 * const.hbar / (mirror.ETM.MirrorMass * (2 * np.pi * self.freq) ** 2) class Seismic(nb.Noise): @@ -314,20 +325,16 @@ class CoatingBrownian(nb.Noise): color='#fe0002', ) - @nb.precomp(precomp_mirror) - @nb.precomp(power=precomp_power) + @nb.precomp(mirror=precomp_mirror) @nb.precomp(cavity=precomp_cavity) - @nb.precomp(coat=precomp_coating) - def calc(self, power, cavity, coat): + @nb.precomp(power=precomp_power) + def calc(self, mirror, cavity, power): wavelength = self.ifo.Laser.Wavelength - materials = self.ifo.Materials nITM = noise.coatingthermal.coating_brownian( - self.freq, materials, wavelength, - cavity.wBeam_ITM, coat.dOpt_ITM, power.parm, + self.freq, mirror.ITM, wavelength, cavity.wBeam_ITM, power.parm, ) nETM = noise.coatingthermal.coating_brownian( - self.freq, materials, wavelength, - cavity.wBeam_ETM, coat.dOpt_ETM, power.parm, + self.freq, mirror.ETM, wavelength, cavity.wBeam_ETM, power.parm, ) return (nITM + nETM) * 2 @@ -342,18 +349,16 @@ class CoatingThermoOptic(nb.Noise): linestyle='--', ) + @nb.precomp(mirror=precomp_mirror) @nb.precomp(cavity=precomp_cavity) - @nb.precomp(coat=precomp_coating) - def calc(self, cavity, coat): + def calc(self, mirror, cavity): wavelength = self.ifo.Laser.Wavelength materials = self.ifo.Materials nITM, junk1, junk2, junk3 = noise.coatingthermal.coating_thermooptic( - self.freq, materials, wavelength, - cavity.wBeam_ITM, coat.dOpt_ITM, + self.freq, mirror.ITM, wavelength, cavity.wBeam_ITM, ) nETM, junk1, junk2, junk3 = noise.coatingthermal.coating_thermooptic( - self.freq, materials, wavelength, - cavity.wBeam_ETM, coat.dOpt_ETM, + self.freq, mirror.ETM, wavelength, cavity.wBeam_ETM, ) return (nITM + nETM) * 2 diff --git a/gwinc/noise/coatingthermal.py b/gwinc/noise/coatingthermal.py index c22951b2ed64519938d4634dfa78896026dce395..4ecd346f7b9db7981754ab97c3984ba5319510a4 100644 --- a/gwinc/noise/coatingthermal.py +++ b/gwinc/noise/coatingthermal.py @@ -10,7 +10,7 @@ from ..const import BESSEL_ZEROS as zeta from ..const import J0M as j0m -def coating_brownian(f, materials, wavelength, wBeam, dOpt, power): +def coating_brownian(f, mirror, wavelength, wBeam, power): """Coating brownian noise for a given collection of coating layers This function calculates Coating Brownian noise using @@ -26,11 +26,9 @@ def coating_brownian(f, materials, wavelength, wBeam, dOpt, power): low-n first. Inputs: f = frequency vector in Hz - materials = gwinc optic materials structure + mirror = mirror properties Struct wavelength = laser wavelength wBeam = beam radius (at 1 / e**2 power) - dOpt = the optical thickness, normalized by lambda, of each - coating layer power = Circulating Laser Power falling on the Mirror (W). If the material.Coating.IncCoatBrAmpNoise parameter is present and evaluates to True the amplitude noise due to coating brownian @@ -64,8 +62,8 @@ def coating_brownian(f, materials, wavelength, wBeam, dOpt, power): """ # extract substructures - sub = materials.Substrate - coat = materials.Coating + sub = mirror.Substrate + coat = mirror.Coating # Constants kBT = const.kB * sub.Temp @@ -77,6 +75,7 @@ def coating_brownian(f, materials, wavelength, wBeam, dOpt, power): nsub = sub.RefractiveIndex # Refractive Index # coating properties + dOpt = mirror.Coating.dOpt Yhighn = coat.Yhighn sigmahighn = coat.Sigmahighn nH = coat.Indexhighn @@ -218,7 +217,7 @@ def coating_brownian(f, materials, wavelength, wBeam, dOpt, power): if coat.get('IncCoatBrAmpNoise'): # get/calculate optic transmittance - mTi = optic.get('Transmittance', 1-np.abs(rho)**2) + mTi = mirror.get('Transmittance', 1-np.abs(rho)**2) # Define Re(epsilon)/2 Rp1 = np.real(Ep1) / 2 # First part of Re(epsilon)/2 @@ -242,7 +241,7 @@ def coating_brownian(f, materials, wavelength, wBeam, dOpt, power): + np.tensordot(p_Sk, S_Sk, axes=1)) AmpToDispConvFac = ((32 * power) - / (materials.MirrorMass * wavelength + / (mirror.MirrorMass * wavelength * f**2 * c * 2 * pi * sqrt(mTi))) # Adding the two pathways of noise contribution as correlated noise SbrZ = (sqrt(S_Xi) + AmpToDispConvFac * sqrt(S_Zeta))**2 @@ -252,14 +251,13 @@ def coating_brownian(f, materials, wavelength, wBeam, dOpt, power): return SbrZ -def coating_thermooptic(f, materials, wavelength, wBeam, dOpt): +def coating_thermooptic(f, mirror, wavelength, wBeam): """Optical coating thermo-optic displacement noise spectrum :f: frequency array in Hz - :materials: gwinc optic materials structure + :mirror: mirror parameter Struct :wavelength: laser wavelength :wBeam: beam radius (at 1 / e**2 power) - :dOpt: coating layer thickness array (Nlayer x 1) :returns: tuple of: StoZ = displacement noise power spectrum at :f: @@ -269,15 +267,15 @@ def coating_thermooptic(f, materials, wavelength, wBeam, dOpt): """ # compute coefficients - dTO, dTR, dTE, T, junk = getCoatTOPos(materials, wavelength, wBeam, dOpt) + dTO, dTR, dTE, T, junk = getCoatTOPos(mirror, wavelength, wBeam) # compute correction factors - gTO = getCoatThickCorr(f, materials, wavelength, dOpt, dTE, dTR) - gTE = getCoatThickCorr(f, materials, wavelength, dOpt, dTE, 0) - gTR = getCoatThickCorr(f, materials, wavelength, dOpt, 0, dTR) + gTO = getCoatThickCorr(f, mirror, wavelength, dTE, dTR) + gTE = getCoatThickCorr(f, mirror, wavelength, dTE, 0) + gTR = getCoatThickCorr(f, mirror, wavelength, 0, dTR) # compute thermal source spectrum - SsurfT, junk = getCoatThermal(f, materials, wBeam) + SsurfT, junk = getCoatThermal(f, mirror, wBeam) StoZ = SsurfT * gTO * dTO**2 SteZ = SsurfT * gTE * dTE**2 @@ -286,13 +284,12 @@ def coating_thermooptic(f, materials, wavelength, wBeam, dOpt): return (StoZ, SteZ, StrZ, T) -def getCoatTOPos(materials, wavelength, wBeam, dOpt): +def getCoatTOPos(mirror, wavelength, wBeam): """Mirror position derivative wrt thermal fluctuations - :materials: gwinc optic materials structure + :mirror: mirror parameter Struct :wavelength: laser wavelength :wBeam: beam radius (at 1 / e**2 power) - :dOpt: coating layer thickness array (Nlayer x 1) :returns: tuple of: dTO = total thermo-optic dz/dT @@ -307,13 +304,14 @@ def getCoatTOPos(materials, wavelength, wBeam, dOpt): """ # parameters - nS = materials.Substrate.RefractiveIndex + nS = mirror.Substrate.RefractiveIndex + dOpt = mirror.Coating.dOpt # compute refractive index, effective alpha and beta - nLayer, aLayer, bLayer, dLayer, sLayer = getCoatLayers(materials, wavelength, dOpt) + nLayer, aLayer, bLayer, dLayer, sLayer = getCoatLayers(mirror, wavelength) # compute coating average parameters - dc, Cc, Kc, aSub = getCoatAvg(materials, wavelength, dOpt) + dc, Cc, Kc, aSub = getCoatAvg(mirror, wavelength) # compute reflectivity and parameters dphi_dT, dphi_TE, dphi_TR, rCoat = getCoatTOPhase(1, nS, nLayer, dOpt, aLayer, bLayer, sLayer) @@ -328,7 +326,7 @@ def getCoatTOPos(materials, wavelength, wBeam, dOpt): dTE = dphi_TE * wavelength / (4 * pi) - aSub * dc # mirror finite size correction - Cfsm = getCoatFiniteCorr(materials, wavelength, wBeam, dOpt) + Cfsm = getCoatFiniteCorr(mirror, wavelength, wBeam) dTE = dTE * Cfsm # add TE and TR effects (sign is already included) @@ -337,14 +335,13 @@ def getCoatTOPos(materials, wavelength, wBeam, dOpt): return dTO, dTR, dTE, T, R -def getCoatThickCorr(f, materials, wavelength, dOpt, dTE, dTR): +def getCoatThickCorr(f, mirror, wavelength, dTE, dTR): """Finite coating thickness correction :f: frequency array in Hz - :materials: gwinc optic materials structure + :mirror: gwinc optic mirror structure :wavelength: laser wavelength :wBeam: beam radius (at 1 / e**2 power) - :dOpt: coating layer thickness array (Nlayer x 1) Uses correction factor from LIGO-T080101, "Thick Coating Correction" (Evans). @@ -362,12 +359,12 @@ def getCoatThickCorr(f, materials, wavelength, dOpt, dTE, dTR): # gTC = 1 - xi * (R - 1 / (3 * R)); # parameter extraction - pS = materials.Substrate + pS = mirror.Substrate Cs = pS.MassCM * pS.MassDensity Ks = pS.MassKappa # compute coating average parameters - dc, Cc, Kc, junk = getCoatAvg(materials, wavelength, dOpt) + dc, Cc, Kc, junk = getCoatAvg(mirror, wavelength) # R and xi (from T080101, Thick Coating Correction) w = 2 * pi * f @@ -396,11 +393,11 @@ def getCoatThickCorr(f, materials, wavelength, dOpt, dTE, dTR): return gTC -def getCoatThermal(f, materials, wBeam): +def getCoatThermal(f, mirror, wBeam): """Thermal noise spectra for a surface layer :f: frequency array in Hz - :materials: gwinc optic material structure + :mirror: mirror parameter Struct :wBeam: beam radius (at 1 / e**2 power) :returns: tuple of: @@ -408,7 +405,7 @@ def getCoatThermal(f, materials, wBeam): rdel = thermal diffusion length at each frequency in m """ - pS = materials.Substrate + pS = mirror.Substrate C_S = pS.MassCM * pS.MassDensity K_S = pS.MassKappa kBT2 = const.kB * pS.Temp**2 @@ -425,12 +422,11 @@ def getCoatThermal(f, materials, wBeam): return SsurfT, rdel -def getCoatLayers(materials, wavelength, dOpt): +def getCoatLayers(mirror, wavelength): """Layer vectors for refractive index, effective alpha and beta and geometrical thickness - :materials: gwinc optic materials structure + :mirror: mirror parameter Struct :wavelength: laser wavelength - :dOpt: coating layer thickness array (Nlayer x 1) :returns: tuple of: nLayer = refractive index of each layer, ordered input to output (N x 1) @@ -444,8 +440,9 @@ def getCoatLayers(materials, wavelength, dOpt): """ # coating parameters - pS = materials.Substrate - pC = materials.Coating + pS = mirror.Substrate + pC = mirror.Coating + dOpt = mirror.Coating.dOpt Y_S = pS.MirrorY sigS = pS.MirrorSigma @@ -499,12 +496,11 @@ def getCoatLayers(materials, wavelength, dOpt): return nLayer, aLayer, bLayer, dLayer, sLayer -def getCoatAvg(materials, wavelength, dOpt): +def getCoatAvg(mirror, wavelength): """Coating average properties - :materials: gwinc optic materials structure + :mirror: gwinc optic mirror structure :wavelength: laser wavelength - :dOpt: coating layer thickness array (Nlayer x 1) :returns: tuple of: dc = total thickness (meters) @@ -514,8 +510,9 @@ def getCoatAvg(materials, wavelength, dOpt): """ # coating parameters - pS = materials.Substrate - pC = materials.Coating + pS = mirror.Substrate + pC = mirror.Coating + dOpt = mirror.Coating.dOpt alphaS = pS.MassAlpha C_S = pS.MassCM * pS.MassDensity @@ -528,7 +525,7 @@ def getCoatAvg(materials, wavelength, dOpt): K_H = pC.ThermalDiffusivityhighn # compute refractive index, effective alpha and beta - junk1, junk2, junk3, dLayer, junk4 = getCoatLayers(materials, wavelength, dOpt) + junk1, junk2, junk3, dLayer, junk4 = getCoatLayers(mirror, wavelength) # heat capacity dc = np.sum(dLayer) @@ -598,13 +595,12 @@ def getCoatTOPhase(nIn, nOut, nLayer, dOpt, aLayer, bLayer, sLayer): return dphi_dT, dphi_TE, dphi_TR, rCoat -def getCoatFiniteCorr(materials, wavelength, wBeam, dOpt): +def getCoatFiniteCorr(mirror, wavelength, wBeam): """Finite mirror size correction - :materials: gwinc optic materials structure + :mirror: mirror parameter Struct :wavelength: laser wavelength :wBeam: beam radius (at 1 / e**2 power) - :dOpt: coating layer thickness array (Nlayer x 1) Uses correction factor from PLA 2003 vol 312 pg 244-255 "Thermodynamical fluctuations in optical mirror coatings" @@ -617,25 +613,26 @@ def getCoatFiniteCorr(materials, wavelength, wBeam, dOpt): """ # parameter extraction - R = materials.MassRadius - H = materials.MassThickness - - alphaS = materials.Substrate.MassAlpha - C_S = materials.Substrate.MassCM * materials.Substrate.MassDensity - Y_S = materials.Substrate.MirrorY - sigS = materials.Substrate.MirrorSigma - - alphaL = materials.Coating.Alphalown - C_L = materials.Coating.CVlown - Y_L = materials.Coating.Ylown - sigL = materials.Coating.Sigmalown - nL = materials.Coating.Indexlown - - alphaH = materials.Coating.Alphahighn - C_H = materials.Coating.CVhighn - Y_H = materials.Coating.Yhighn - sigH = materials.Coating.Sigmahighn - nH = materials.Coating.Indexhighn + R = mirror.MassRadius + H = mirror.MassThickness + dOpt = mirror.Coating.dOpt + + alphaS = mirror.Substrate.MassAlpha + C_S = mirror.Substrate.MassCM * mirror.Substrate.MassDensity + Y_S = mirror.Substrate.MirrorY + sigS = mirror.Substrate.MirrorSigma + + alphaL = mirror.Coating.Alphalown + C_L = mirror.Coating.CVlown + Y_L = mirror.Coating.Ylown + sigL = mirror.Coating.Sigmalown + nL = mirror.Coating.Indexlown + + alphaH = mirror.Coating.Alphahighn + C_H = mirror.Coating.CVhighn + Y_H = mirror.Coating.Yhighn + sigH = mirror.Coating.Sigmahighn + nH = mirror.Coating.Indexhighn # coating sums dL = wavelength * np.sum(dOpt[::2]) / nL diff --git a/gwinc/noise/quantum.py b/gwinc/noise/quantum.py index 01702040a0548fb9fd48a2071127f5240b95cd1c..5adf510ac04ee27ed46700841987c2782ec5b6aa 100644 --- a/gwinc/noise/quantum.py +++ b/gwinc/noise/quantum.py @@ -17,7 +17,7 @@ def sqzOptimalSqueezeAngle(Mifo, eta): return alpha -def shotrad(f, ifo, power): +def shotrad(f, ifo, mirror_mass, power): """Quantum noise noise strain spectrum :f: frequency array in Hz @@ -48,7 +48,7 @@ def shotrad(f, ifo, power): else: namespace = globals() fname = namespace['shotrad' + ifo.Optics.Type] - coeff, Mifo, Msig, Mn = fname(f, ifo, power=power) + coeff, Mifo, Msig, Mn = fname(f, ifo, mirror_mass, power) # check for consistent dimensions Nfield = Msig.shape[0] @@ -291,7 +291,7 @@ def compile_SEC_RES_TF(): print('RES', '=', str(SEC_RES_expr[0]).replace('Matrix', 'np.array')) -def shotradSignalRecycled(f, ifo, power): +def shotradSignalRecycled(f, ifo, mirror_mass, power): """Quantum noise model for signal recycled IFO (see shotrad for more info) New version July 2016 by JH based on transfer function formalism @@ -312,7 +312,7 @@ def shotradSignalRecycled(f, ifo, power): L = ifo.Infrastructure.Length # Length of arm cavities [m] l = ifo.Optics.SRM.CavityLength # SRC Length [m] T = ifo.Optics.ITM.Transmittance # ITM Transmittance [Power] - m = ifo.Materials.MirrorMass # Mirror mass [kg] + m = mirror_mass # Mirror mass [kg] # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - bsloss = ifo.Optics.BSLoss # BS Loss [Power] @@ -461,7 +461,7 @@ def shotradSignalRecycled(f, ifo, power): return coeff, Mifo, Msig, Mnoise -def shotradSignalRecycledBnC(f, ifo, power): +def shotradSignalRecycledBnC(f, ifo, mirror_mass, power): """Quantum noise model for signal recycled IFO See shotrad for more info. @@ -493,7 +493,7 @@ def shotradSignalRecycledBnC(f, ifo, power): L = ifo.Infrastructure.Length # Length of arm cavities [m] l = ifo.Optics.SRM.CavityLength # SRC Length [m] T = ifo.Optics.ITM.Transmittance # ITM Transmittance [Power] - m = ifo.Materials.MirrorMass # Mirror mass [kg] + m = mirror_mass # Mirror mass [kg] # - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - bsloss = ifo.Optics.BSLoss # BS Loss [Power]