Use interpolant caching for performance improvements

In the ROQLikelihood and calibration, interp1d instances are repeatedly
instantiated. Minor improvements are made by made by caching the
interpolant and only updating what is required.

bilby/gw/detector/calibration.py

@@ -94,6 +94,34 @@ class CubicSpline(Recalibrate):
         return self.__class__.__name__ + '(prefix=\'{}\', minimum_frequency={}, maximum_frequency={}, n_points={})'\
             .format(self.prefix, self.minimum_frequency, self.maximum_frequency, self.n_points)
 
+    @property
+    def delta_amplitude_interp(self):
+        """ Return the updated interpolant or generate a new one if required """
+        amplitude_parameters = [self.params['amplitude_{}'.format(ii)]
+                                for ii in range(self.n_points)]
+
+        if hasattr(self, '_delta_amplitude_interp'):
+            self._delta_amplitude_interp.y = amplitude_parameters
+        else:
+            self._delta_amplitude_interp = interp1d(
+                self.log_spline_points, amplitude_parameters, kind='cubic',
+                bounds_error=False, fill_value=0)
+        return self._delta_amplitude_interp
+
+    @property
+    def delta_phase_interp(self):
+        """ Return the updated interpolant or generate a new one if required """
+        phase_parameters = [
+            self.params['phase_{}'.format(ii)] for ii in range(self.n_points)]
+
+        if hasattr(self, '_delta_phase_interp'):
+            self._delta_phase_interp.y = phase_parameters
+        else:
+            self._delta_phase_interp = interp1d(
+                self.log_spline_points, phase_parameters, kind='cubic',
+                bounds_error=False, fill_value=0)
+        return self._delta_phase_interp
+
     def get_calibration_factor(self, frequency_array, **params):
         """Apply calibration model
 
@@ -113,17 +141,10 @@ class CubicSpline(Recalibrate):
             The factor to multiply the strain by.
         """
         self.set_calibration_parameters(**params)
-        amplitude_parameters = [self.params['amplitude_{}'.format(ii)]
-                                for ii in range(self.n_points)]
-        delta_amplitude = interp1d(
-            self.log_spline_points, amplitude_parameters, kind='cubic',
-            bounds_error=False, fill_value=0)(np.log10(frequency_array))
+        log10_frequency_array = np.log10(frequency_array)
 
-        phase_parameters = [
-            self.params['phase_{}'.format(ii)] for ii in range(self.n_points)]
-        delta_phase = interp1d(
-            self.log_spline_points, phase_parameters, kind='cubic',
-            bounds_error=False, fill_value=0)(np.log10(frequency_array))
+        delta_amplitude = self.delta_amplitude_interp(log10_frequency_array)
+        delta_phase = self.delta_phase_interp(log10_frequency_array)
 
         calibration_factor = (1 + delta_amplitude) * (2 + 1j * delta_phase) / (2 - 1j * delta_phase)
 

bilby/gw/likelihood.py

@@ -865,6 +865,18 @@ class ROQGravitationalWaveTransient(GravitationalWaveTransient):
         self.frequency_nodes_quadratic = \
             waveform_generator.waveform_arguments['frequency_nodes_quadratic']
 
+    def get_d_inner_h_interp(self, d_inner_h_tc_array, indices):
+        """ Return the updated interpolant or generate a new one if required """
+        if hasattr(self, '_d_inner_h_interp'):
+            self._d_inner_h_interp.x = self.weights['time_samples'][indices]
+            self._d_inner_h_interp.y = d_inner_h_tc_array
+        else:
+            self._d_inner_h_interp = interp1d(
+                self.weights['time_samples'][indices],
+                d_inner_h_tc_array, kind='cubic', assume_sorted=True)
+
+        return self._d_inner_h_interp
+
     def calculate_snrs(self, waveform_polarizations, interferometer):
         """
         Compute the snrs for ROQ
@@ -916,9 +928,7 @@ class ROQGravitationalWaveTransient(GravitationalWaveTransient):
             'i,ji->j', np.conjugate(h_plus_linear + h_cross_linear),
             self.weights[interferometer.name + '_linear'][indices])
 
-        d_inner_h = interp1d(
-            self.weights['time_samples'][indices],
-            d_inner_h_tc_array, kind='cubic', assume_sorted=True)(ifo_time)
+        d_inner_h = self.get_d_inner_h_interp(d_inner_h_tc_array, indices)(ifo_time)
 
         optimal_snr_squared = \
             np.vdot(np.abs(h_plus_quadratic + h_cross_quadratic)**2,