allow precessing BNS

bilby/gw/conversion.py

@@ -235,7 +235,7 @@ def convert_to_lal_binary_neutron_star_parameters(parameters):
 
 
     Mass: mass_1, mass_2
-    Spin: chi_1, chi_2, tilt_1, tilt_2, phi_12, phi_jl
+    Spin: a_1, a_2, tilt_1, tilt_2, phi_12, phi_jl
     Extrinsic: luminosity_distance, theta_jn, phase, ra, dec, geocent_time, psi
 
     This involves popping a lot of things from parameters.
@@ -258,19 +258,6 @@ def convert_to_lal_binary_neutron_star_parameters(parameters):
     converted_parameters, added_keys =\
         convert_to_lal_binary_black_hole_parameters(converted_parameters)
 
-    for idx in ['1', '2']:
-        mag = 'a_{}'.format(idx)
-        if mag in original_keys:
-            tilt = 'tilt_{}'.format(idx)
-            if tilt in original_keys:
-                converted_parameters['chi_{}'.format(idx)] = (
-                    converted_parameters[mag] *
-                    np.cos(converted_parameters[tilt]))
-            else:
-                converted_parameters['chi_{}'.format(idx)] = (
-                    converted_parameters[mag])
-
-    # catch if tidal parameters aren't present
     if not any([key in converted_parameters for key in
                 ['lambda_1', 'lambda_2', 'lambda_tilde', 'delta_lambda']]):
         converted_parameters['lambda_1'] = 0

bilby/gw/prior.py

@@ -336,10 +336,12 @@ class BNSPriorDict(PriorDict):
             By default this generates many additional parameters, see
             BNSPriorDict.default_conversion_function
         """
-        if not aligned_spin:
-            logger.warning('Non-aligned spins not yet supported for BNS.')
+        if aligned_spin:
+            default_file = 'binary_neutron_stars.prior'
+        else:
+            default_file = 'precessing_binary_neutron_stars.prior'
         if dictionary is None and filename is None:
-            filename = os.path.join(os.path.dirname(__file__), 'prior_files', 'binary_neutron_stars.prior')
+            filename = os.path.join(os.path.dirname(__file__), 'prior_files', default_file)
             logger.info('No prior given, using default BNS priors in {}.'.format(filename))
         elif filename is not None:
             if not os.path.isfile(filename):

bilby/gw/prior_files/precessing_binary_neutron_stars.prior

+# These are the default priors we use for BNS systems.
+# Note that you may wish to use more specific mass and distance parameters.
+# These commands are all known to bilby.gw.prior.
+# Lines beginning "#" are ignored.
+mass_1 = Uniform(name='mass_1', minimum=1, maximum=2, unit='$M_{\\odot}$', boundary=None)
+mass_2 = Uniform(name='mass_2', minimum=1, maximum=2, unit='$M_{\\odot}$', boundary=None)
+mass_ratio = Constraint(name='mass_ratio', minimum=0.125, maximum=1)
+# chirp_mass = Uniform(name='chirp_mass', minimum=0.87, maximum=1.74, unit='$M_{\\odot}$', boundary=None)
+# total_mass = Uniform(name='total_mass', minimum=2, maximum=4, unit='$M_{\\odot}$', boundary=None)
+# mass_ratio = Uniform(name='mass_ratio', minimum=0.5, maximum=1, boundary=None)
+# symmetric_mass_ratio = Uniform(name='symmetric_mass_ratio', minimum=0.22, maximum=0.25, boundary=None)
+a_1 = Uniform(name='a_1', minimum=0, maximum=0.05, boundary='reflective')
+a_2 = Uniform(name='a_2', minimum=0, maximum=0.05, boundary='reflective')
+tilt_1 = Sine(name='tilt_1', boundary='reflective')
+tilt_2 = Sine(name='tilt_2', boundary='reflective')
+# cos_tilt_1 = Uniform(name='cos_tilt_1', minimum=-1, maximum=1, boundary=None)
+# cos_tilt_2 = Uniform(name='cos_tilt_2', minimum=-1, maximum=1, boundary=None)
+phi_12 = Uniform(name='phi_12', minimum=0, maximum=2 * np.pi, boundary='periodic')
+phi_jl = Uniform(name='phi_jl', minimum=0, maximum=2 * np.pi, boundary='periodic')
+luminosity_distance = bilby.gw.prior.UniformSourceFrame(name='luminosity_distance', minimum=10, maximum=500, unit='Mpc', boundary=None)
+dec = Cosine(name='dec', boundary='reflective')
+ra = Uniform(name='ra', minimum=0, maximum=2 * np.pi, boundary='periodic')
+theta_jn = Sine(name='theta_jn', boundary='reflective')
+# cos_theta_jn = Uniform(name='cos_theta_jn', minimum=-1, maximum=1, boundary=None)
+psi = Uniform(name='psi', minimum=0, maximum=np.pi, boundary='periodic')
+phase = Uniform(name='phase', minimum=0, maximum=2 * np.pi, boundary='periodic')
+lambda_1 = Uniform(name='lambda_1', minimum=0, maximum=3000, boundary=None)
+lambda_2 = Uniform(name='lambda_2', minimum=0, maximum=3000, boundary=None)
+# lambda_tilde = Uniform(name='lambda_tilde', minimum=0, maximum=5000, boundary=None)
+# delta_lambda = Uniform(name='delta_lambda', minimum=-5000, maximum=5000, boundary=None)

bilby/gw/source.py

@@ -73,8 +73,9 @@ def lal_binary_black_hole(
 
 
 def lal_binary_neutron_star(
-        frequency_array, mass_1, mass_2, luminosity_distance, chi_1, chi_2,
-        theta_jn, phase, lambda_1, lambda_2, **kwargs):
+        frequency_array, mass_1, mass_2, luminosity_distance, a_1, tilt_1,
+        phi_12, a_2, tilt_2, phi_jl, theta_jn, phase, lambda_1, lambda_2,
+        **kwargs):
     """ A Binary Neutron Star waveform model using lalsimulation
 
     Parameters
@@ -87,27 +88,27 @@ def lal_binary_neutron_star(
         The mass of the lighter object in solar masses
     luminosity_distance: float
         The luminosity distance in megaparsec
-    chi_1: float
-        Dimensionless aligned spin
-    chi_2: float
-        Dimensionless aligned spin
+    a_1: float
+        Dimensionless primary spin magnitude
+    tilt_1: float
+        Primary tilt angle
+    phi_12: float
+        Azimuthal angle between the two component spins
+    a_2: float
+        Dimensionless secondary spin magnitude
+    tilt_2: float
+        Secondary tilt angle
+    phi_jl: float
+        Azimuthal angle between the total binary angular momentum and the
+        orbital angular momentum
     theta_jn: float
         Orbital inclination
     phase: float
         The phase at coalescence
-    ra: float
-        The right ascension of the binary
-    dec: float
-        The declination of the object
-    geocent_time: float
-        The time at coalescence
-    psi: float
-        Orbital polarisation
     lambda_1: float
         Dimensionless tidal deformability of mass_1
     lambda_2: float
         Dimensionless tidal deformability of mass_2
-
     kwargs: dict
         Optional keyword arguments
 
@@ -116,20 +117,16 @@ def lal_binary_neutron_star(
     dict: A dictionary with the plus and cross polarisation strain modes
     """
     waveform_kwargs = dict(
-        waveform_approximant='TaylorF2', reference_frequency=50.0,
+        waveform_approximant='IMRPhenomPv2_NRTidal', reference_frequency=50.0,
         minimum_frequency=20.0, maximum_frequency=frequency_array[-1],
-        pn_spin_order=-1, pn_tidal_order=-1, pn_phase_order=-1, pn_amplitude_order=0)
-
-    a_1 = abs(chi_1)
-    a_2 = abs(chi_2)
-    tilt_1 = np.arccos(np.sign(chi_1))
-    tilt_2 = np.arccos(np.sign(chi_2))
+        pn_spin_order=-1, pn_tidal_order=-1, pn_phase_order=-1,
+        pn_amplitude_order=0)
     waveform_kwargs.update(kwargs)
     return _base_lal_cbc_fd_waveform(
         frequency_array=frequency_array, mass_1=mass_1, mass_2=mass_2,
         luminosity_distance=luminosity_distance, theta_jn=theta_jn, phase=phase,
-        a_1=a_1, a_2=a_2, tilt_1=tilt_1, tilt_2=tilt_2, lambda_1=lambda_1,
-        lambda_2=lambda_2, **waveform_kwargs)
+        a_1=a_1, a_2=a_2, tilt_1=tilt_1, tilt_2=tilt_2, phi_12=phi_12,
+        phi_jl=phi_jl, lambda_1=lambda_1, lambda_2=lambda_2, **waveform_kwargs)
 
 
 def lal_eccentric_binary_black_hole_no_spins(