Skip to content
Snippets Groups Projects
Commit 426864c1 authored by Colm Talbot's avatar Colm Talbot
Browse files

make noise_realisation_tests.py run

parent 22874e1b
No related branches found
No related tags found
No related merge requests found
Pipeline #
......@@ -8,40 +8,40 @@ class TestNoiseRealisation(unittest.TestCase):
time_duration = 1.
sampling_frequency = 4096.
factor = np.sqrt(2./time_duration)
navg = 1000
n_avg = 1000
psd_avg = 0
H1 = peyote.detector.H1
for x in range(0, navg):
H1_hf_noise, frequencies = H1.power_spectral_density.get_noise_realisation(sampling_frequency,
time_duration)
H1.set_data(sampling_frequency, time_duration, frequency_domain_strain=H1_hf_noise)
hf_tmp = H1.data
psd_avg += abs(hf_tmp)**2
psd_avg = psd_avg/navg
interferometer = peyote.detector.H1
for x in range(0, n_avg):
interferometer.set_data(sampling_frequency, time_duration, from_power_spectral_density=True)
psd_avg += abs(interferometer.data)**2
psd_avg = psd_avg/n_avg
asd_avg = np.sqrt(abs(psd_avg))
a = H1.amplitude_spectral_density_array/factor
a = interferometer.amplitude_spectral_density_array/factor
b = asd_avg
self.assertTrue(np.isclose(a[2]/b[2], 1.00, atol=1e-1))
def test_noise_normalisation(self):
time_duration = 1.
sampling_frequency = 4096.
number_of_samples = int(np.round(time_duration*sampling_frequency))
time_array = (1./sampling_frequency) * np.linspace(0, number_of_samples, number_of_samples)
time_array = peyote.utils.create_time_series(sampling_frequency=sampling_frequency, duration=time_duration)
# generate some toy-model signal for matched filtering SNR testing
navg = 10000
snr = np.zeros(navg)
n_avg = 10000
snr = np.zeros(n_avg)
mu = np.exp(-(time_array-time_duration/2.)**2 / (2.*0.1**2)) * np.sin(2*np.pi*100*time_array)
muf, frequency_array = peyote.utils.nfft(mu, sampling_frequency)
for x in range(0, navg):
H1 = peyote.detector.H1
H1_hf_noise, frequencies = H1.power_spectral_density.get_noise_realisation(sampling_frequency, time_duration)
H1.set_data(sampling_frequency, time_duration,frequency_domain_strain=H1_hf_noise)
hf_tmp = H1.data
Sh = H1.power_spectral_density
snr[x] = peyote.utils.inner_product(hf_tmp, muf, frequency_array, Sh) / np.sqrt(peyote.utils.inner_product(muf, muf, frequency_array, Sh))
for x in range(0, n_avg):
interferometer = peyote.detector.H1
interferometer.set_data(sampling_frequency, time_duration, from_power_spectral_density=True)
hf_tmp = interferometer.data
psd = interferometer.power_spectral_density
snr[x] = peyote.utils.inner_product(hf_tmp, muf, frequency_array, psd) \
/ np.sqrt(peyote.utils.inner_product(muf, muf, frequency_array, psd))
self.assertTrue(np.isclose(np.std(snr), 1.00, atol=1e-2))
if __name__ == '__main__':
unittest.main()
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment