Updated the PP catalogue to allow building GPR models from it.

parent 74c7f59f
Pipeline #55380 failed with stage
in 19 seconds
......@@ -37,7 +37,7 @@ class PPCatalogue(Catalogue):
solMass = 1.988e30 # kg
parsec = 3.0857e16 # m
def __init__(self, approximant, total_mass=20, fmin=30.):
def __init__(self, approximant, total_mass=20, fmin=30., waveforms=[]):
"""
Assemble the catalogue.
......@@ -52,61 +52,25 @@ class PPCatalogue(Catalogue):
self.total_mass = total_mass
self.fmin = fmin
def waveform(self, p, time_range, distance=1.0, coa_phase=0, t0=0, f_ref=100):
"""
Generate a single waveform from the catalogue.
"""
delta_t = (time_range[1] - time_range[0]) / time_range[2]
mass1, mass2 = self.components_from_total(self.total_mass,
p['mass ratio'])
self.waveforms = waveforms
hp, hx = get_td_waveform(approximant=self.approximant,
mass1=mass1,
mass2=mass2,
spin1x=p['spin 1x'],
spin1y=p['spin 1y'],
spin1z=p['spin 1z'],
spin2x=p['spin 2x'],
spin2y=p['spin 2y'],
spin2z=p['spin 2z'],
distance=distance,
delta_t=delta_t/1e4,
coa_phase=coa_phase,
f_ref=f_ref,
f_lower=self.fmin)
self.parameters = ["mass_ratio", "spin_1x", "spin_1y", "spin_1z",
"spin_2x", "spin_2y", "spin_2z"]
hp = Timeseries(hp)
hx = Timeseries(hx)
self.data_parameters = {0: "time",
1: "mass ratio",
2: "spin 1x",
3: "spin 1y",
4: "spin 1z",
5: "spin 2x",
6: "spin 2y",
7: "spin 2z",
8: "h+",
9: "hx"}
self.c_ind = {j: i for i, j in self.data_parameters.items()}
# Recenter the waveforms on the maximum strain
hp.times -= hp.times[np.argmax(np.abs(hp.data - 1j * hx.data))]
hx.times -= hx.times[np.argmax(np.abs(hp.data - 1j * hx.data))]
# Recenter the waveforms now to some arbitrary time
hp.times -= t0
hx.times -= t0
tix = (time_range[0] < hp.times*1e4) & (hp.times*1e4 < time_range[1])
hp.times = hp.times[tix]
hx.times = hx.times[tix]
hp.data = hp.data[tix]
hx.data = hx.data[tix]
return hp, hx
def create_training_data(self,
total_mass,
sampling="factorial",
sample_space = {"mass ratio": [0, 4, 0.5]},
f_min=None,
ma=None,
sample_rate=4096,
distance=1,
tmax=0.005, tmin=-0.010):
def create_training_data(self, total_mass, f_min=None, ma=None,
sample_rate=4096., distance=1, tmax=0.005, tmin=-0.010):
"""
Produce an array of data suitable for use as training data
using the waveforms in this catalogue.
......@@ -137,14 +101,12 @@ class PPCatalogue(Catalogue):
gen_sample = 4096
skip = int(gen_sample/sample_rate) # this is nasty
# Unlike in NR catalogues where there are a defined number of waveforms available,
# we need to construct a set of locations
for waveform in self.waveforms:
time_range = [1e4*tmin, 1e4*tmax, (tmax-tmin)*sample_rate]
for waveform_p in self.waveforms:
try:
hp, hx = waveform.timeseries(total_mass, gen_sample,
f_min, distance, ma=ma)
hp, hx = self.waveform(waveform_p, time_range)
#hp, hx = waveform.timeseries(total_mass, gen_sample,
# f_min, distance, ma=ma)
except LalsuiteError:
print(
"There was an error producing a waveform for {}"
......@@ -159,14 +121,59 @@ class PPCatalogue(Catalogue):
export[:, 0] = hp.times[ixs][::skip]
export[:, 1] = waveform.mass_ratio
export[:, [2, 3, 4, 5, 6, 7]] *= waveform.spins
export[:, 1] = waveform_p['mass ratio']
# TODO Fix this
#export[:, [2, 3, 4, 5, 6, 7]] *= waveform_p['spin 1x']
export[:, 8] = hp.data[ixs][::skip]
export[:, 9] = hx.data[ixs][::skip]
big_export = np.vstack([big_export, export])
return big_export[1:, :]
def waveform(self, p, time_range, distance=1.0, coa_phase=0, t0=0, f_ref=100):
"""
Generate a single waveform from the catalogue.
"""
delta_t = (time_range[1] - time_range[0]) / time_range[2]
mass1, mass2 = self.components_from_total(self.total_mass,
p['mass ratio'])
hp, hx = get_td_waveform(approximant=self.approximant,
mass1=mass1,
mass2=mass2,
spin1x=p['spin 1x'],
spin1y=p['spin 1y'],
spin1z=p['spin 1z'],
spin2x=p['spin 2x'],
spin2y=p['spin 2y'],
spin2z=p['spin 2z'],
distance=distance,
delta_t=delta_t/1e4,
coa_phase=coa_phase,
f_ref=f_ref,
f_lower=self.fmin)
hp = Timeseries(hp)
hx = Timeseries(hx)
# Recenter the waveforms on the maximum strain
hp.times -= hp.times[np.argmax(np.abs(hp.data - 1j * hx.data))]
hx.times -= hx.times[np.argmax(np.abs(hp.data - 1j * hx.data))]
# Recenter the waveforms now to some arbitrary time
hp.times -= t0
hx.times -= t0
tix = (time_range[0] < hp.times*1e4) & (hp.times*1e4 < time_range[1])
hp.times = hp.times[tix]
hx.times = hx.times[tix]
hp.data = hp.data[tix]
hx.data = hx.data[tix]
return hp, hx
......@@ -191,6 +198,7 @@ class NRCatalogue(Catalogue):
The file type containing the waveforms.
"""
self.origin = origin
if ftype == "hdf5":
self.file_suffix = "h5"
else:
......@@ -288,14 +296,16 @@ class NRCatalogue(Catalogue):
A query expression, for example 's1x == 0' to find all
the waveforms where the s1x component is zero.
"""
query_table = self.table.query(expression, inplace=False)
query_waveforms = [waveform for waveform
query_table = self.table.query(expression)
query_waveforms = np.array([waveform for waveform
in self.waveforms if waveform.tag
in list(query_table['tag'])]
in query_table.tag.values])
new_cat = NRCatalogue(origin="GeorgiaTech", ftype="hdf5",
new_cat = NRCatalogue(origin=self.origin, ftype="hdf5",
table=query_table,
waveforms=query_waveforms)
return new_cat
def distances(self, point, metric=np.ones(7)):
......
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