Cpuopt

bin/pycbc_inspiral

@@ -124,7 +124,9 @@ parser.add_argument("--user-tag", type=str, metavar="TAG", help="""
                     This is used to identify FULL_DATA jobs for 
                     compatibility with pipedown post-processing. 
                     Option will be removed when no longer needed.""")
-parser.add_argument("--psd-recalculate-segments", type=int, help="Number of segments to use before recalculating the PSD", default=0)
+parser.add_argument("--psd-recalculate-segments", type=int, 
+                    help="Number of segments to use before recalculating the PSD", default=0)
+parser.add_argument("--gpu-callback-method", default='none')
 
 # Add options groups
 psd.insert_psd_option_group(parser)
@@ -185,11 +187,16 @@ with ctx:
     associate_psd(strain_segments, gwstrain, segments, 
                   opt.psd_recalculate_segments, flen, delta_f, flow)
 
+    template_mem = zeros(tlen, dtype = complex64)
+    cluster_window = int(opt.cluster_window * gwstrain.sample_rate)
+
     matched_filter = MatchedFilterControl(opt.low_frequency_cutoff, None,
                                    opt.snr_threshold, tlen, delta_f, complex64,
+                                   segments, template_mem, cluster_window,
                                    downsample_factor=opt.downsample_factor,
                                    upsample_threshold=opt.upsample_threshold,
-                                   upsample_method=opt.upsample_method)
+                                   upsample_method=opt.upsample_method,
+                                   gpu_callback_method=opt.gpu_callback_method)
 
     bank_chisq = vetoes.SingleDetBankVeto(opt.bank_veto_bank_file,
                                           flen, delta_f, flow, complex64, 
@@ -233,24 +240,25 @@ with ctx:
 
     logging.info("Read in template bank")
     bank = waveform.FilterBank(opt.bank_file, flen, delta_f,
-                    flow, dtype=complex64, phase_order=opt.order,
-                    taper=opt.taper_template, approximant=opt.approximant,
-                    out=zeros(tlen, dtype=complex64))
+                    flow, dtype = complex64, phase_order = opt.order,
+                    taper = opt.taper_template, approximant = opt.approximant,
+                    out = template_mem)
 
     for t_num, template in enumerate(bank):
         event_mgr.new_template(tmplt=template.params, sigmasq=template.sigmasq(segments[0].psd))
-        if opt.cluster_method == "window":
-            cluster_window = int(opt.cluster_window * gwstrain.sample_rate)
-        elif opt.cluster_method == "template":
-            cluster_window = int(template.chirp_length * gwstrain.sample_rate)
+        #if opt.cluster_method == "window":
+            # At present, moved above with class-based approach
+            #cluster_window = int(opt.cluster_window * gwstrain.sample_rate)
+        if opt.cluster_method == "template":
+            #cluster_window = int(template.chirp_length * gwstrain.sample_rate)
+            raise ValueError("FIXME: template clustering not supported yet with optimizations!")
 
         for s_num, stilde in enumerate(segments):
             logging.info("Filtering template %d/%d segment %d/%d" % 
                          (t_num + 1, len(bank), s_num + 1, len(segments)))
 
             snr, norm, corr, idx, snrv = \
-               matched_filter.matched_filter_and_cluster(template, 
-                                  template.sigmasq(stilde.psd), stilde, cluster_window)
+               matched_filter.matched_filter_and_cluster(s_num, template.sigmasq(stilde.psd))
 
             if not len(idx):
                 continue
@@ -274,6 +282,7 @@ with ctx:
         event_mgr.cluster_template_events("time_index", "snr", cluster_window)
         event_mgr.finalize_template_events()
 
+event_mgr.finalize_template_events()
 logging.info("Found %s triggers" % str(len(event_mgr.events)))
 
 if opt.chisq_threshold and opt.chisq_bins:

pycbc/events/events.py

@@ -44,6 +44,59 @@ def threshold_and_cluster(series, threshold, window):
     """Return list of values and indices values over threshold in series.
     """
 
+@schemed("pycbc.events.threshold_")
+def batched_threshold_and_cluster(series, threshold, window):
+    """Return list of values and indices values over threshold in series.
+    """
+@schemed("pycbc.events.threshold_")
+def _threshold_cluster_factory(series, window):
+    pass
+
+class ThresholdCluster(object):
+    """Create a threshold and cluster engine
+
+    Parameters
+    ---------
+    series : complex64
+      Input pycbc.types.Array (or subclass); it will be searched for
+      points above threshold that are then clustered
+    window : int
+      Size in samples of the window over which to cluster
+    """
+    def __new__(cls, *args, **kwargs):
+        real_cls = _threshold_cluster_factory(*args, **kwargs)
+        return real_cls(*args, **kwargs)
+
+# The class below should serve as the parent for all schemed classes.
+# The intention is that this class serves simply as the location for
+# all documentation of the class and its methods, though that is not
+# yet implemented.  Perhaps something along the lines of:
+#
+#    http://stackoverflow.com/questions/2025562/inherit-docstrings-in-python-class-inheritance
+#
+# will work? Is there a better way?
+class _BaseThresholdCluster(object):
+    def threshold_and_cluster(self, threshold):
+        """
+        Threshold and cluster the memory specified at instantiation with the
+        threshold specified at creation and the window size specified at creation.
+
+        Parameters:
+        -----------
+        threshold : float32
+          The minimum absolute value of the series given at object initialization
+          to return when thresholding and clustering.
+
+        Returns:
+        --------
+        event_vals : complex64
+          Numpy array, complex values of the clustered events
+        event_locs : uint32
+          Numpy array, indices into series of location of events          
+        """
+        pass
+
+
 def fc_cluster_over_window_fast(times, values, window_length):
     """ Reduce the events by clustering over a window using
     the FindChirp clustering algorithm
@@ -933,7 +986,7 @@ def get_weighted_snr(self, fac):
         coinc_def_row.search_coinc_type = 0
         coinc_def_table.append(coinc_def_row)
 
-__all__ = ['threshold_and_cluster', 'newsnr',
+__all__ = ['threshold_and_cluster', 'newsnr', 'batched_threshold_and_cluster',
            'findchirp_cluster_over_window', 'fc_cluster_over_window_fast',
-           'threshold', 'cluster_reduce',
+           'threshold', 'cluster_reduce', 'ThresholdCluster',
            'EventManager', 'EventManagerMultiDet']