Source Code for Module pylal.ligolw_sicluster

  1  # Copyright (C) 2006  Duncan A. Brown 
  2  # 
  3  # This program is free software; you can redistribute it and/or modify it 
  4  # under the terms of the GNU General Public License as published by the 
  5  # Free Software Foundation; either version 2 of the License, or (at your 
  6  # option) any later version. 
  7  # 
  8  # This program is distributed in the hope that it will be useful, but 
  9  # WITHOUT ANY WARRANTY; without even the implied warranty of 
 10  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General 
 11  # Public License for more details. 
 12  # 
 13  # You should have received a copy of the GNU General Public License along 
 14  # with this program; if not, write to the Free Software Foundation, Inc., 
 15  # 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA. 
 16   
 17  # 
 18  # ============================================================================= 
 19  # 
 20  #                                   Preamble 
 21  # 
 22  # ============================================================================= 
 23  # 
 24   
 25  import sys 
 26   
 27  from lalburst import snglcluster 
 28   
 29  from glue.ligolw import table 
 30  from glue.ligolw import lsctables 
 31  from glue.ligolw.utils import process as ligolw_process 
 32  from glue.ligolw.utils import search_summary as ligolw_search_summary 
 33  from pylal import git_version 
 34  from pylal import SnglInspiralUtils 
 35   
 36  __author__ = "Duncan Brown <dbrown@ligo.caltech.edu>" 
 37   
 38   
 39  # 
 40  # ============================================================================= 
 41  # 
 42  #                                 Preparation 
 43  # 
 44  # ============================================================================= 
 45  # 
 46   
 47 -def get_tables(doc): 
 48    snglinspiraltable = lsctables.SnglInspiralTable.get_table(doc) 
 49   
 50    input_times = None 
 51    output_times = None 
 52    try: 
 53      searchsummtable = lsctables.SearchSummaryTable.get_table(doc) 
 54      input_times = searchsummtable.get_inlist().extent() 
 55      output_times = searchsummtable.get_outlist().extent() 
 56    except ValueError: 
 57      pass 
 58       
 59    return input_times, output_times, snglinspiraltable 
 60   
 61   
 62  # 
 63  # ============================================================================= 
 64  # 
 65  #                           Add Process Information 
 66  # 
 67  # ============================================================================= 
 68  # 
 69   
 70 -def append_process(doc, **kwargs): 
 71    process = ligolw_process.append_process( 
 72      doc, program = "ligolw_sicluster", version = git_version.verbose_msg, 
 73      cvs_repository = "lscsoft", cvs_entry_time = git_version.date, 
 74      comment = kwargs["comment"]) 
 75   
 76    ligolw_process.append_process_params(doc, process,  
 77      [("--cluster-window", "lstring", kwargs["cluster_window"])]) 
 78    if kwargs["snr_threshold"] > 0: 
 79      ligolw_process.append_process_params(doc, process,  
 80        [("--snr-threshold", "lstring", kwargs["snr_threshold"])]) 
 81    if kwargs["sort_descending_snr"]: 
 82      ligolw_process.append_process_params(doc, process,  
 83        [("--sort-descending-snr", "lstring", " ")]) 
 84    if kwargs["sort_ascending_snr"]: 
 85      ligolw_process.append_process_params(doc, process,  
 86        [("--sort-ascending-snr", "lstring", " ")]) 
 87   
 88    return process 
 89   
 90   
 91  # 
 92  # ============================================================================= 
 93  # 
 94  #                             Clustering Algorithm 
 95  # 
 96  # ============================================================================= 
 97  # 
 98   
 99 -def SnglInspiralCluster(a, b): 
100    """ 
101    Replace a with a cluster constructed from a and b.  
102    """ 
103    if b.snr >= a.snr: 
104      return b 
105    else: 
106      return a 
107   
108  # 
109  # ============================================================================= 
110  # 
111  #                                 Library API 
112  # 
113  # ============================================================================= 
114  # 
115   
116 -def ligolw_sicluster(doc, **kwargs): 
117    # Extract segments and tables 
118    inseg, outseg, snglinspiraltable = get_tables(doc) 
119   
120    # Add process information 
121    try: 
122      process = append_process(doc, **kwargs) 
123    except ValueError: 
124      process = None 
125   
126    # Delete all triggers below threshold 
127    if kwargs["snr_threshold"] > 0: 
128      thresh = float(kwargs["snr_threshold"]) 
129      if kwargs["verbose"]: 
130        print >>sys.stderr, "discarding triggers with snr < %f ..." % \ 
131          kwargs["snr_threshold"] 
132      for i in range(len(snglinspiraltable) - 1, -1, -1): 
133        if snglinspiraltable[i].snr <= thresh: 
134          del snglinspiraltable[i] 
135   
136    # Cluster 
137    snglcluster.cluster_events( 
138      snglinspiraltable, 
139      testfunc = lambda a, b: SnglInspiralUtils.CompareSnglInspiral(a, b, twindow = kwargs["cluster_window"]), 
140      clusterfunc = SnglInspiralCluster, 
141      sortfunc = SnglInspiralUtils.CompareSnglInspiralByEndTime, 
142      bailoutfunc = lambda a, b: SnglInspiralUtils.CompareSnglInspiral(a, b, twindow = kwargs["cluster_window"]), 
143      verbose = kwargs["verbose"] 
144    ) 
145   
146    # Sort by signal-to-noise ratio 
147    if kwargs["sort_ascending_snr"] or kwargs["sort_descending_snr"]: 
148      if kwargs["verbose"]: 
149        print >>sys.stderr, "sorting by snr ..." 
150      snglinspiraltable.sort(SnglInspiralUtils.CompareSnglInspiralBySnr) 
151      if kwargs["sort_descending_snr"]: 
152        snglinspiraltable.reverse() 
153   
154    # Add search summary information 
155    if process and inseg and outseg: 
156      ligolw_search_summary.append_search_summary(doc, process, inseg = inseg, outseg = outseg,  
157        nevents = len(snglinspiraltable)) 
158    if process: 
159      ligolw_process.set_process_end_time(process) 
160   
161    return doc 
162