Skip to content
Snippets Groups Projects
Commit e00a1702 authored by frcojimenez's avatar frcojimenez
Browse files

rd snr

parent 95dfc9c2
Branches
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
code_new/ringdown_overtones_example.ipynb 0 → 100644
+ 407
− 0
View file @ e00a1702
%% Cell type:code id: tags:
``` python
import numpy
from pycbc.waveform import generator
from matplotlib import pyplot as plt
from pycbc import psd, filter
import numpy as np
from random import seed
from random import gauss
import pycbc
import pylab
```
%% Cell type:code id: tags:
``` python
genclass = generator.TDomainFreqTauRingdownGenerator
detectors = ['H1', 'L1']
delta_f = 1./8
delta_t = 1./2048
f_lower = 20.
f_final = 1024.
params = {}
params['amp220'] = 5.e-21
params['amp221'] = 2
params['f_220'] = 256.607043264397
params['f_221'] = 250.9439959543337
params['tau_220'] = 0.004047225085452862
params['tau_221'] = 0.001338559280517331
params['dec'] = -1.25
params['distance'] = 410.0
#params['f_lower'] = 18.0
#params['f_ref'] = 20.0
#params['final_mass'] = 63.0
#params['final_spin'] = 0.64
params['inclination'] = 2.5
params['injtype'] = 'ringdown'
params['lmns'] = '222'
params['phi220'] = 2.0
params['phi221'] = 1.1
params['polarization'] = 1.75
params['ra'] = 2.2
params['t_final'] = 1.
params['tc'] = 0.
params['approximant'] = 'TdQNMfromFreqTau'
tc = params['tc']
gen = generator.FDomainDetFrameGenerator(genclass,
epoch=tc,
detectors=detectors,
variable_args=params.keys(),
f_lower = f_lower,
f_final = f_final,
delta_f=delta_f,
delta_t=delta_t)
signal = gen.generate(**params)
```
%% Cell type:code id: tags:
``` python
fig = plt.figure(0)
plt.plot(signal['H1'].to_timeseries().sample_times, signal['H1'].to_timeseries())
plt.xlim(-0.01,0.05)
fig.set_dpi(150)
plt.show()
```
%% Output
%% Cell type:code id: tags:
``` python
adv_psd = {}
for det in detectors:
adv_psd[det] = psd.analytical.aLIGOZeroDetHighPower(length=1024*8+1, delta_f=1./8, low_freq_cutoff=20)
```
%% Cell type:code id: tags:
``` python
help(psd.analytical.aLIGOZeroDetHighPower)
```
%% Output
Help on function aLIGOZeroDetHighPower in module pycbc.psd.analytical:
aLIGOZeroDetHighPower(length, delta_f, low_freq_cutoff)
Return a FrequencySeries containing the aLIGOZeroDetHighPower PSD from LALSimulation.
%% Cell type:code id: tags:
``` python
delta_t = 1.0 / 4096
tsamples = int(32 / delta_t)
ts = pycbc.noise.noise_from_psd(tsamples, delta_t, psd, seed=127)
```
%% Output
---------------------------------------------------------------------------
AttributeError Traceback (most recent call last)
<ipython-input-58-e64cb88395a2> in <module>
1 delta_t = 1.0 / 4096
2 tsamples = int(32 / delta_t)
----> 3 ts = pycbc.noise.noise_from_psd(tsamples, delta_t, psd, seed=127)
/work/francisco.jimenez/sio/pyenv/versions/3.6.6/lib/python3.6/site-packages/pycbc/noise/gaussian.py in noise_from_psd(length, delta_t, psd, seed)
102 randomness = lal.gsl_rng("ranlux", seed)
103
--> 104 N = int (1.0 / delta_t / psd.delta_f)
105 n = N//2+1
106 stride = N//2
AttributeError: module 'pycbc.psd' has no attribute 'delta_f'
%% Cell type:code id: tags:
``` python
optsnr = {}
snrsq = {}
addsnr = 0
for det in signal.keys():
snrsq[det] = filter.sigmasq(signal[det], psd=adv_psd[det],
low_frequency_cutoff=f_lower,
high_frequency_cutoff=f_final)
addsnr += snrsq[det]
optsnr = numpy.sqrt(addsnr)
print(optsnr)
```
%% Output
72.31060179760523
%% Cell type:code id: tags:
``` python
for _ in range(10):
value = gauss(0, np.sqrt(adv_psd[det][i]))
print(value)
```
%% Output
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
<ipython-input-56-8af46ede6ee5> in <module>
1 for _ in range(10):
----> 2 value = gauss(0, np.sqrt(adv_psd[det][i]))
3 print(value)
NameError: name 'i' is not defined
%% Cell type:code id: tags:
``` python
# The color of the noise matches a PSD which you provide
flen = int(2048 / delta_f) + 1
psd = pycbc.psd.aLIGOZeroDetHighPower(len(signal['H1'].to_timeseries()), delta_f, f_lower)
# Generate 4 seconds of noise at 4096 Hz
tsamples = int(4 / delta_t)
ts = pycbc.noise.noise_from_psd(tsamples, delta_t, psd, seed=127)
```
%% Cell type:code id: tags:
``` python
pylab.plot(ts.sample_times, ts+signal['H1'].to_timeseries())
pylab.ylabel('Strain')
pylab.xlabel('Time (s)')
pylab.show()
```
%% Output
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-78-d70d9b4b69c6> in <module>
----> 1 pylab.plot(ts.sample_times, ts+signal['H1'].to_timeseries())
2 pylab.ylabel('Strain')
3 pylab.xlabel('Time (s)')
4 pylab.show()
<decorator-gen-154> in __add__(self, other)
/work/francisco.jimenez/sio/pyenv/versions/3.6.6/lib/python3.6/site-packages/pycbc/types/array.py in _returntype(fn, self, *args)
230 @decorator
231 def _returntype(fn, self, *args):
--> 232 ary = fn(self,*args) # pylint:disable=not-callable
233 if ary is NotImplemented:
234 return NotImplemented
<decorator-gen-153> in __add__(self, other)
/work/francisco.jimenez/sio/pyenv/versions/3.6.6/lib/python3.6/site-packages/pycbc/types/array.py in _convert(fn, self, *args)
64 # Convert this array to the current processing scheme
65 _convert_to_scheme(self)
---> 66 return fn(self, *args)
67
68 @decorator
<decorator-gen-152> in __add__(self, other)
/work/francisco.jimenez/sio/pyenv/versions/3.6.6/lib/python3.6/site-packages/pycbc/types/array.py in _checkother(fn, self, *args)
251 elif isinstance(other, type(self)) or type(other) is Array:
252 if len(other) != len(self):
--> 253 raise ValueError('lengths do not match')
254 if other.precision == self.precision:
255 _convert_to_scheme(other)
ValueError: lengths do not match
%% Cell type:code id: tags:
``` python
```
%% Output
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
<ipython-input-6-32b2ef174fcb> in <module>
----> 1 get_lm_f0tau(mass, spin, l, m, nmodes)
NameError: name 'get_lm_f0tau' is not defined
%% Cell type:code id: tags:
``` python
from pycbc.conversions import get_lm_f0tau
import qnm
```
%% Cell type:code id: tags:
``` python
def wRD_to_f_Phys(f,M):
c=2.99792458*10**8;G=6.67259*10**(-11);MS=1.9885*10**30;
return (c**3/(M*MS*G*2*np.pi))*f
def tauRD_to_t_Phys(tau,M):
c=2.99792458*10**8;G=6.67259*10**(-11);MS=1.9885*10**30;
return ((M*MS*G)/c**3)*tau
```
%% Cell type:code id: tags:
``` python
omegas = [qnm.modes_cache(s=-2,l=2,m=2,n=i)(a=0.7)[0] for i in range (0,2)]
w = (np.real(omegas))/1
tau=-1/(np.imag(omegas))*1
```
%% Cell type:code id: tags:
``` python
tauRD_to_t_Phys(tau[0],69)
```
%% Output
0.004205658110224193
%% Cell type:code id: tags:
``` python
frompycbc=get_lm_f0tau(69, 0.7, 2, 2, 2)
```
%% Cell type:code id: tags:
``` python
frompycbc[1]
```
%% Output
array([0.00420654, 0.00139151])
%% Cell type:code id: tags:
``` python
0.00139151*2.5
```
%% Output
0.0034787750000000004
%% Cell type:code id: tags:
``` python
```
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment

AEI Privacy Policy Imprint