diff --git a/gracedb/serialize/__init__.py b/gracedb/serialize/__init__.py new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/gracedb/serialize/utils.py b/gracedb/serialize/utils.py new file mode 100644 index 0000000000000000000000000000000000000000..95e14ffd710b116ecace8e055cee9d058c35523c --- /dev/null +++ b/gracedb/serialize/utils.py @@ -0,0 +1,437 @@ +#!/usr/bin/python + +from math import log +from time import gmtime, strftime + +from pylal import Fr +from glue.lal import LIGOTimeGPS +from glue.ligolw import ligolw +from glue.ligolw import table +from glue.ligolw import lsctables + +############################################################################## +# +# useful variables +# +############################################################################## + +#Need these for each search: inspiral, burst, etc. +InspiralCoincDef = lsctables.CoincDef(search = u"inspiral", \ + search_coinc_type = 0, \ + description = \ + u"sngl_inspiral<-->sngl_inspiral coincidences") +#these should work for both Omega and CWB +BurstCoincDef = lsctables.CoincDef(search = u"burst", \ + search_coinc_type = 0, \ + description = \ + u"coherent burst coincidences") + +#list of detectors participating in the coinc +#MBTA only sends triples to gracedb at the time being so this list is +#simply for convenience. burst and future inspiral searches should +#construct this list on the fly +H1L1V1_detlist = ['H1', 'L1', 'V1'] +H1L1_detlist = ['H1', 'L1'] +H1V1_detlist = ['H1', 'V1'] +L1V1_detlist = ['L1', 'V1'] +H1_detlist = ['H1'] +L1_detlist = ['L1'] +V1_detlist = ['V1'] + +#this is the subset of SnglInspiralTable.validcolumn.keys() that +#are assigned from MBTA coinc triggers +MBTA_set_keys = ['ifo', 'search', 'end_time', 'end_time_ns', 'mass1', 'mass2',\ + 'mchirp', 'mtotal', 'eta', 'snr', 'eff_distance', 'event_id',\ + 'process_id', 'channel'] +#Omega +Omega_set_keys = ['process_id', 'ifos', 'start_time', 'start_time_ns',\ + 'duration', 'confidence', 'coinc_event_id'] +#CWB +CWB_set_keys = ['process_id', 'ifos', 'start_time', 'start_time_ns',\ + 'coinc_event_id'] + +#this dictionary is the simplest way to assign event_id's +#collisions are are taken care of in the process of conversion to sqlite +insp_event_id_dict = {'H1': 'sngl_inspiral:event_id:0',\ + 'L1': 'sngl_inspiral:event_id:1',\ + 'V1': 'sngl_inspiral:event_id:2'} +#this one is designed for coherent searches, which don't have event_ids +coherent_event_id_dict = None + +#the names of the variables we're going to get from omega +omega_vars = ['time', 'frequency', 'duration', 'bandwidth', 'modeTheta',\ + 'modePhi', 'probSignal', 'probGlitch', 'logSignal','logGlitch',\ + 'network', 'URL_web', 'URL_file'] + +############################################################################## +# +# convenience functions +# +############################################################################## + +def compute_mchirp_eta(m1,m2): + """ + compute and return mchirp and eta for a given pair of masses + """ + + mtot = m1 + m2 + mu = m1*m2/mtot + eta = mu/mtot + mchirp = pow(eta,3.0/5.0)*mtot + + return float(mchirp), float(eta) + +def write_output_files(root_dir, xmldoc, log_content, \ + xml_fname = 'coinc.xml', log_fname = 'event.log'): + """ + write the xml-format coinc tables and log file + """ + + f = open(root_dir+'/'+xml_fname,'w') + xmldoc.write(f) + f.close() + + f = open(root_dir+'/'+log_fname,'w') + f.write(log_content) + f.close() + +def get_ifos_for_cwb(cwb_ifos): + """ + get human-readable things from CWB detector labels + """ + + ifos = [] + for i in cwb_ifos: + if i == '1': ifos.append('L1') + if i == '2' : ifos.append('H1') + if i == '3' : ifos.append('H2') + if i == '4' : ifos.append('G1') + if i == '5' : ifos.append('T1') + if i == '6' : ifos.append('V1') + if i == '7' : ifos.append('A1') + + return ifos + +############################################################################## +# +# table populators +# +############################################################################## + +def populate_inspiral_tables(MBTA_frame, set_keys = MBTA_set_keys, \ + event_id_dict = insp_event_id_dict): + """ + create xml file and populate the SnglInspiral and CoincInspiral tables from a + coinc .gwf file from MBTA + xmldoc: xml file to append the tables to + MBTA_frame: frame file to get info about triggers from + set_keys: columns in the SnglInspiral Table to set + process_id: process_id + event_id_dict: {ifo:event_id} dictionary to assign event_id's + coinc_event_id: coinc_event_id + detectors: detectors participating in the coinc + + returns xmldoc and contents of the comment field + """ + #initialize xml document + xmldoc = ligolw.Document() + xmldoc.appendChild(ligolw.LIGO_LW()) + #dictionaries to store about individual triggers + end_time = {} + snr = {} + mass1 = {} + mass2 = {} + Deff = {} + mchirp = {} + eta = {} + + #extract the information from the frame file + events = Fr.frgetevent(MBTA_frame) + #get the ifos from the event name + for event in events: + if 'MbtaHLV' in event['name']: + detectors = H1L1V1_detlist + elif 'MbtaHL' in event['name']: + detectors = H1L1_detlist + elif 'MbtaHV' in event['name']: + detectors = H1V1_detlist + elif 'MbtaH' in event['name']: + detectors = H1_detlist + elif 'MbtaLV' in event['name']: + detectors = L1V1_detlist + elif 'MbtaL' in event['name']: + detectors = L1_detlist + elif 'MbtaV' in event['name']: + detectors = V1_detlist + else: + raise ValueError, "Invalid FrEvent name" + + log_data = event['comment'] + '\n' + try: + far = 1/(float(event['IFAR_year'])*365.0*24*60*60) + except KeyError: + far = None + for ifo in detectors: + end_time[ifo] = LIGOTimeGPS(event[ifo+':end_time']) + snr[ifo] = float(event[ifo+':SNR']) + mass1[ifo] = float(event[ifo+':mass1']) + mass2[ifo] = float(event[ifo+':mass2']) + mchirp[ifo], eta[ifo] = compute_mchirp_eta(mass1[ifo],mass2[ifo]) + Deff[ifo] = float(event[ifo+':eff_distance']) + + #fill the SnglInspiralTable + sin_table = lsctables.New(lsctables.SnglInspiralTable) + xmldoc.childNodes[0].appendChild(sin_table) + process_id = lsctables.ProcessTable.get_next_id() + for ifo in detectors: + row = sin_table.RowType() + row.ifo = ifo + row.search = 'MBTA' + row.end_time = end_time[ifo].seconds + row.end_time_ns = end_time[ifo].nanoseconds + row.mass1 = mass1[ifo] + row.mass2 = mass2[ifo] + row.mchirp = mchirp[ifo] + row.mtotal = mass1[ifo] + mass2[ifo] + row.eta = eta[ifo] + row.snr = snr[ifo] + row.eff_distance = Deff[ifo] + row.event_id = event_id_dict[ifo] + row.process_id = process_id + row.channel = '' + #zero out the rest of the columns + #should work in chi2 and chi2cut + for key in sin_table.validcolumns.keys(): + if key not in set_keys: + setattr(row,key,None) + sin_table.append(row) + + #CoincInspiralTable + #using the conventions found in: + #https://www.lsc-group.phys.uwm.edu/ligovirgo/cbcnote/S6Plan/ + #090505160219S6PlanningNotebookCoinc_and_Experiment_Tables_ihope_implementation? + #highlight=%28coinc%29|%28table%29 + + temp_data_loc = None + + if len(detectors) < 2: + return xmldoc, log_data, temp_data_loc + + #coinc_event_id = coinc_event_id_base + str(UID) + + #this next line is to guard against future violence to the table definitions + cintcols = ['coinc_event_id','ifos','end_time','end_time_ns','mass','mchirp','snr','false_alarm_rate','combined_far'] + cin_table = lsctables.New(lsctables.CoincInspiralTable,columns=cintcols) + xmldoc.childNodes[0].appendChild(cin_table) + row = cin_table.RowType() + row.set_ifos(detectors) + cid = lsctables.CoincTable.get_next_id() + row.coinc_event_id = cid + row.end_time = end_time['H1'].seconds + row.end_time_ns = end_time['H1'].nanoseconds + row.mass = (sum(mass1.values()) + sum(mass2.values()))/3 + row.mchirp = sum(mchirp.values())/3 + #the snr here is really the snr NOT effective snr + row.snr = pow(sum([x*x for x in snr.values()]),0.5) + if far is not None: + #far is triggers in Hz + row.combined_far = float(far) + else: + row.combined_far = None + row.false_alarm_rate = None + + cin_table.append(row) + + + xmldoc = populate_coinc_tables(xmldoc,cid,insp_event_id_dict,\ + InspiralCoincDef,detectors) + + return xmldoc, log_data, temp_data_loc + +def populate_omega_tables(datafile, set_keys = Omega_set_keys): + """ + """ + #initialize xml document + xmldoc = ligolw.Document() + xmldoc.appendChild(ligolw.LIGO_LW()) + + #extract the data from the intial Omega file + f = open(datafile, 'r') + omega_list = [] + for line in f.readlines(): + if not line.strip(): continue # ignore blank lines + elif '#' in line.strip()[0]: continue # ignore comments + elif '=' not in line: raise ValueError, "Improperly formatted line" + else: + omega_list.extend([dat.strip() for dat in line.split('=',1)]) + f.close() + omega_data = dict(zip(omega_list[::2],omega_list[1::2])) + # basic error checking +# for key in omega_data: +# if not (key in omega_vars): +# raise ValueError, "Unknown variable" + + #create the content for the event.log file + log_data = '\nLog File created '\ + +strftime("%a, %d %b %Y %H:%M:%S", gmtime())\ + +'\n' + + for var in omega_vars: + log_data += var + ': ' + omega_data[var] + '\n' + + #pull out the ifos + detectors = [ifo for ifo in omega_data['network'].split(',')] + + #fill the MutliBurstTable + mb_table = lsctables.New(lsctables.MultiBurstTable) + xmldoc.childNodes[0].appendChild(mb_table) + row = mb_table.RowType() + row.process_id = lsctables.ProcessTable.get_next_id() + row.set_ifos(detectors) + st = LIGOTimeGPS(omega_data['time']) + row.start_time = st.seconds + row.start_time_ns = st.nanoseconds + row.duration = None + row.confidence = -log(float(omega_data['probGlitch'])) + cid = lsctables.CoincTable.get_next_id() + row.coinc_event_id = cid + for key in mb_table.validcolumns.keys(): + if key not in set_keys: + setattr(row,key,None) + mb_table.append(row) + + xmldoc = populate_coinc_tables(xmldoc,cid, coherent_event_id_dict,\ + BurstCoincDef, detectors) + + return xmldoc, log_data, omega_data['URL_file'] + +def populate_cwb_tables(datafile, set_keys=CWB_set_keys): + """ + """ + #initialize xml document + xmldoc = ligolw.Document() + xmldoc.appendChild(ligolw.LIGO_LW()) + + #extract the data from the file + f = open(datafile,'r') + cwb_list = [] + for line in f.readlines(): + if not line.strip(): continue #ignore blanks + elif '#' in line.strip()[0]: continue #skip comments + elif 'H1:' in line.strip() or 'L1:' in line.strip() or 'V1:' in line.strip(): continue #skip DQ stuff + elif ':' in line.strip(): cwb_list.extend([dat.strip() for dat in line.split(':',1)]) + + f.close() + cwb_data = dict(zip(cwb_list[::2],cwb_list[1::2])) + + #create the content for the event.log file + log_data = '\nLog File created '\ + +strftime("%a, %d %b %Y %H:%M:%S", gmtime())\ + +'\n' + detectors = get_ifos_for_cwb(cwb_data['ifo'].split()) + + for var in cwb_data: + log_data += var + ': ' + cwb_data[var] + '\n' + + #fill the MutliBurstTable + mb_table = lsctables.New(lsctables.MultiBurstTable) + xmldoc.childNodes[0].appendChild(mb_table) + row = mb_table.RowType() + row.process_id = lsctables.ProcessTable.get_next_id() + row.set_ifos(detectors) + st = LIGOTimeGPS(cwb_data['start'][0]) + row.start_time = st.seconds + row.start_time_ns = st.nanoseconds + cid = lsctables.CoincTable.get_next_id() + row.coinc_event_id = cid + for key in mb_table.validcolumns.keys(): + if key not in set_keys: + setattr(row,key,None) + mb_table.append(row) + + xmldoc = populate_coinc_tables(xmldoc,cid, coherent_event_id_dict,\ + BurstCoincDef, detectors) + + return xmldoc, log_data, None + + + +def populate_coinc_tables(xmldoc, coinc_event_id, event_id_dict,\ + CoincDef, detectors, \ + time_slide_id = None, likelihood = None): + """ + populate a set of coinc tables + xmldoc: xml file to append the tables to + CoincDef: pre-initialized CoincDef table row + detectors: detectors participating in the coinc + """ + #make sure there's actually a coinc there to write + if len(detectors) < 2: + return xmldoc + else: + #CoincTable + coinc_table = lsctables.New(lsctables.CoincTable) + xmldoc.childNodes[0].appendChild(coinc_table) + row = coinc_table.RowType() + row.process_id = lsctables.ProcessTable.get_next_id() + row.coinc_event_id = coinc_event_id + coinc_def_id = lsctables.CoincDefTable.get_next_id() + row.coinc_def_id = coinc_def_id + row.time_slide_id = time_slide_id + row.set_instruments(detectors) + if 'inspiral' in CoincDef.search: + row.nevents = len(detectors) + elif 'burst' in CoincDef.search: + row.nevents = 1 + else: + raise ValueError, "Unrecognize CoincDef.search" + row.likelihood = likelihood + coinc_table.append(row) + + #CoincMapTable + coinc_map_table = lsctables.New(lsctables.CoincMapTable) + xmldoc.childNodes[0].appendChild(coinc_map_table) + for ifo in detectors: + row = coinc_map_table.RowType() + row.coinc_event_id = coinc_event_id + if 'inspiral' in CoincDef.search: + row.table_name = lsctables.SnglInspiralTable.tableName.split(':')[0] + elif 'burst' in CoincDef.search: + row.table_name = lsctables.MultiBurstTable.tableName.split(':')[0] + else: + raise ValueError, "Unrecognize CoincDef.search" + if event_id_dict: + row.event_id = event_id_dict[ifo] + coinc_map_table.append(row) + if not event_id_dict: + row.event_id = coinc_event_id + coinc_map_table.append(row) + + #CoincDefTable + coinc_def_table = lsctables.New(lsctables.CoincDefTable) + xmldoc.childNodes[0].appendChild(coinc_def_table) + row = coinc_def_table.RowType() + row.coinc_def_id = coinc_def_id + row.search = CoincDef.search + row.search_coinc_type = CoincDef.search_coinc_type + row.description = CoincDef.description + coinc_def_table.append(row) + + return xmldoc + + +############################################################################## +# +# usage example +# +############################################################################## + +#here's how it works for inspirals +#populate the tables +#xmldoc, log_data, temp_data_loc = populate_inspiral_tables("MbtaFake-930909680-16.gwf") +#write the output +#write_output_files('.', xmldoc, log_data) + +#here's how it works for bursts +#xmldoc, log_data, temp_data_loc = populate_burst_tables("initial.data") +#write_output_files('.', final_xmldoc, log_data) diff --git a/gracedb/translator.py b/gracedb/translator.py index 9aad66cfada04abee2c86855d39c40c7c3e9c8ec..8f998ac05d4b7a97a7bec0df7c9c0893da49d4d4 100644 --- a/gracedb/translator.py +++ b/gracedb/translator.py @@ -6,7 +6,7 @@ from models import EventLog import glue import glue.ligolw.utils -from ligo.gracedb.utils import populate_inspiral_tables, \ +from gracedb.serialize.utils import populate_inspiral_tables, \ populate_omega_tables, \ write_output_files