Commit 8ef75e17 authored by Jameson Graef Rollins's avatar Jameson Graef Rollins

seismic: use scipy PchipInterpolator for seismic spectrum interpolation

This is the scipy equivalent of the 'pchip' interpolation algorithm that
matgwinc is using.

Also minor clean up and clarify some of the variable names
parent 57b64517
from __future__ import division
import numpy as np
from numpy import log10
from scipy.interpolate import interp1d
from scipy.interpolate import PchipInterpolator as interp1d
def seismic(f, ifo):
......@@ -16,13 +16,13 @@ def seismic(f, ifo):
theta = ifo.Suspension.VHCoupling.theta
# noise input, horizontal and vertical
nx, np_ = seisBSC(f)
nt, nr = seisBSC(f)
# horizontal noise total
nh = (abs(hTable)**2) * nx**2
nh = (abs(hTable)**2) * nt**2
# vertical noise total
nv = (abs(theta * vTable)**2) * nx**2
nv = (abs(theta * vTable)**2) * nt**2
# new total noise
n = nv + nh
......@@ -38,16 +38,17 @@ def seismic(f, ifo):
def seisBSC(f):
"""Rough ISI noise source spectra.
Returns (ISI translational DOFs, ISI rotational DOFs)
Returns ISI (translational, rotational) DOFs
"""
# translational DOFs (from Rana's bsc_seismic.m)
SEI_F = np.array([0.01, 0.03, 0.1, 0.2, 0.5, 1, 10, 30, 300])
SEI_X = np.array([3e-6, 1e-6, 2e-7, 2e-7, 8e-10, 1e-11, 3e-13, 3e-14, 3e-14])
nx = 10**(interp1d(SEI_F,log10(SEI_X), bounds_error=False, fill_value=-14)(f))
# translational DOFs
SEI_T = np.array([3e-6, 1e-6, 2e-7, 2e-7, 8e-10, 1e-11, 3e-13, 3e-14, 3e-14])
nt = 10**(interp1d(SEI_F, log10(SEI_T))(f))
# rotational DOFs
SEI_P = np.array([1e-8, 3e-8, 2e-8, 1e-8, 4e-10, 1e-11, 3e-13, 3e-14, 3e-14])
np_ = 10**(interp1d(SEI_F,log10(SEI_P), bounds_error=False, fill_value=-14)(f))
SEI_R = np.array([1e-8, 3e-8, 2e-8, 1e-8, 4e-10, 1e-11, 3e-13, 3e-14, 3e-14])
nr = 10**(interp1d(SEI_F, log10(SEI_T))(f))
return nx, np_
return nt, nr
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment