simplified bns example and deleted not used spins

examples/injection_examples/binary_neutron_star_example.py

@@ -22,11 +22,9 @@ np.random.seed(88170235)
 
 # We are going to inject a binary neutron star waveform.  We first establish a dictionary of parameters that
 # includes all of the different waveform parameters, including masses of the two black holes (mass_1, mass_2),
-# spins of both black holes (a_1,a_2) , etc. Take into account the the waveform approximants TaylorF2 and 
-# IMRPHenomD_NRTidal can only handle aligned spins, so the parameters tilt_1, tilt_2, phi_12, phi_jl must be
-# set to 0.
-injection_parameters = dict(mass_1=1.5, mass_2=1.3, a_1=0.0, a_2=0.0, tilt_1=0.0, tilt_2=0.0, phi_12=0.0, phi_jl=0,
-                            luminosity_distance=50., iota=0.4, psi=2.659, phase=1.3, geocent_time=1126259642.413,
+# spins of both black holes (a_1,a_2) , etc. 
+injection_parameters = dict(mass_1=1.5, mass_2=1.3, a_1=0.0, a_2=0.0, luminosity_distance=50., 
+                            iota=0.4, psi=2.659, phase=1.3, geocent_time=1126259642.413,
                             ra=1.375, dec=-1.2108, lambda1=400, lambda2=450)
 
 # Set the duration and sampling frequency of the data segment that we're going to inject the signal into. For the 
@@ -49,60 +47,36 @@ hf_signal = waveform_generator.frequency_domain_strain()
 
 # Set up interferometers.  In this case we'll use three interferometers (LIGO-Hanford (H1), LIGO-Livingston (L1),
 # and Virgo (V1)).  These default to their design sensitivity and start at 40 Hz.
-H1 = tupak.gw.detector.get_empty_interferometer('H1')
-H1.minimum_frequency = 40
-H1.set_strain_data_from_power_spectral_density(sampling_frequency=sampling_frequency, duration=duration,
-                                               start_time = start_time)
-H1.inject_signal(parameters=injection_parameters,
-        injection_polarizations=hf_signal,
-        waveform_generator=waveform_generator)
-
-H1.save_data(outdir, label=label)
-H1.plot_data(signal=H1.get_detector_response(hf_signal,injection_parameters), outdir=outdir, label=label)
-
-
-#second interferometer
-L1 = tupak.gw.detector.get_empty_interferometer('L1')
-L1.minimum_frequency = 40
-L1.set_strain_data_from_power_spectral_density(sampling_frequency=sampling_frequency, duration=duration,
-                                               start_time = start_time)
-L1.inject_signal(parameters=injection_parameters,
-        injection_polarizations=hf_signal,
-        waveform_generator=waveform_generator)
-
-L1.save_data(outdir, label=label)
-L1.plot_data(signal=L1.get_detector_response(hf_signal,injection_parameters), outdir=outdir, label=label)
-
-#third interferometer
-V1 = tupak.gw.detector.get_empty_interferometer('V1')
-V1.minimum_frequency = 40
-V1.set_strain_data_from_power_spectral_density(sampling_frequency=sampling_frequency, duration=duration,
-                                               start_time = start_time)
-V1.inject_signal(parameters=injection_parameters,
-        injection_polarizations=hf_signal,
-        waveform_generator=waveform_generator)
-
-V1.save_data(outdir, label=label)
-V1.plot_data(signal=V1.get_detector_response(hf_signal,injection_parameters), outdir=outdir, label=label)
-IFOs = np.array([H1,L1,V1])
+interferometers = tupak.gw.detector.InterferometerList(['H1', 'L1', 'V1'])
+for interferometer in interferometers:
+    interferometer.minimum_frequency = 40
+interferometers.set_strain_data_from_power_spectral_densities(
+    sampling_frequency=sampling_frequency, duration=duration,
+    start_time = start_time)
+interferometers.inject_signal(parameters=injection_parameters, waveform_generator=waveform_generator)
+
 
 #priors
 priors = tupak.gw.prior.BBHPriorSet()
+priors.pop('tilt_1')
+priors.pop('tilt_2')
+priors.pop('phi_12')
+priors.pop('phi_jl')
 priors['lambda1'] = tupak.prior.Uniform(0, 3000, '$\\Lambda_1$')
 priors['lambda2'] = tupak.prior.Uniform(0, 3000, '$\\Lambda_2$')
 priors['mass_1'] = tupak.prior.Uniform(1, 2, '$m_1$')
 priors['mass_2'] = tupak.prior.Uniform(1, 2, '$m_2$')
-for key in ['a_1', 'a_2', 'tilt_1', 'tilt_2', 'phi_12', 'phi_jl', 'psi',
-           'geocent_time','ra','dec','iota','luminosity_distance','phase']:
+for key in ['a_1', 'a_2', 'psi','geocent_time','ra','dec',
+            'iota','luminosity_distance','phase']:
     priors[key] = injection_parameters[key]
     
 # Initialise the likelihood by passing in the interferometer data (IFOs) and the waveoform generator
-likelihood = tupak.gw.GravitationalWaveTransient(interferometers=IFOs, waveform_generator=waveform_generator,
+likelihood = tupak.gw.GravitationalWaveTransient(interferometers=interferometers, waveform_generator=waveform_generator,
                                               time_marginalization=False, phase_marginalization=False,
                                               distance_marginalization=False, prior=priors)
 
 # Run sampler.  In this case we're going to use the `nestle` sampler
-result = tupak.run_sampler(likelihood=likelihood, priors=priors, sampler='nestle', npoints=500,
+result = tupak.run_sampler(likelihood=likelihood, priors=priors, sampler='nestle', npoints=1000,
                            injection_parameters=injection_parameters, outdir=outdir, label=label)
 
 result.plot_corner()

tupak/gw/source.py

@@ -254,7 +254,7 @@ def supernova_pca_model(
     return {'plus': h_plus, 'cross': h_cross}
 
 def lal_binary_neutron_star(
-        frequency_array, mass_1, mass_2, luminosity_distance, a_1, tilt_1, phi_12, a_2, tilt_2, phi_jl,
+        frequency_array, mass_1, mass_2, luminosity_distance, a_1, a_2,
         iota, phase, ra, dec, geocent_time, psi, lambda1, lambda2, **kwargs):
     """ A Binary Black Hole waveform model using lalsimulation
 
@@ -269,17 +269,9 @@ def lal_binary_neutron_star(
     luminosity_distance: float
         The luminosity distance in megaparsec
     a_1: float
-        Dimensionless primary spin magnitude
-    tilt_1: float
-        Primary tilt angle. TaylorF2 and IMRPhenomD_NRTidal only handle aligned spin, set this value to 0
-    phi_12: float
-        TaylorF2 and IMRPhenomD_NRTidal only handle aligned spin, set this value to 0
+        Dimensionless spin magnitude
     a_2: float
         Dimensionless secondary spin magnitude. 
-    tilt_2: float
-        Secondary tilt angle. TaylorF2 and IMRPhenomD_NRTidal only handle aligned spin, set this value to 0
-    phi_jl: float
-        TaylorF2 and IMRPhenomD_NRTidal only handle aligned spin, set this value to 0
     iota: float
         Orbital inclination
     phase: float
@@ -319,15 +311,12 @@ def lal_binary_neutron_star(
     mass_1 = mass_1 * utils.solar_mass
     mass_2 = mass_2 * utils.solar_mass
 
-    if tilt_1 == 0 and tilt_2 == 0 and phi_12 == 0 and phi_jl==0:
-        spin_1x = 0
-        spin_1y = 0
-        spin_1z = a_1
-        spin_2x = 0
-        spin_2y = 0
-        spin_2z = a_2
-    else:
-        raise ValueError('The waveform approximants TaylorF2 and IMRPhenomD_NRTidal only support aligned spins')
+    spin_1x = 0
+    spin_1y = 0
+    spin_1z = a_1
+    spin_2x = 0
+    spin_2y = 0
+    spin_2z = a_2
 
     longitude_ascending_nodes = 0.0
     eccentricity = 0.0