diff --git a/minke/mdctools.py b/minke/mdctools.py index 79f64aa054baf0062f9f7120a684b99a52cf75c5..1b1996fde3aea066f1f8a53abde74ba3fa714daa 100644 --- a/minke/mdctools.py +++ b/minke/mdctools.py @@ -415,6 +415,7 @@ class MDCSet(): The name of the injection in the cWB format """ row = self.waveforms[row] + name = '' numberspart = '' if row.waveform in ("Dimmelmeier+08", "Scheidegger+10", "Mueller+12", "Ott+13", "Yakunin+10"): @@ -618,11 +619,19 @@ class Frame(): """ Represents a frame, in order to prepare the injection frames """ - def __init__(self, start, duration, ifo): + def __init__(self, start, duration, ifo, number = -1): + """ + + Parameters + ---------- + number : int + The frame's number within the project. Defaults to -1. + """ self.start = start self.duration = duration self.end = self.start + duration self.ifos = ifo + self.number = -1 def __repr__(self): out = '' @@ -653,7 +662,7 @@ class Frame(): log += "\n" return log - def generate_gwf(self, mdc, directory, channel="SCIENCE", force=False, rate=16384.0): + def generate_gwf(self, mdc, directory, project = "Minke", channel="SCIENCE", force=False, rate=16384.0): """ Produce the gwf file which corresponds to the MDC set over the period of this frame. @@ -666,6 +675,8 @@ class Frame(): "/home/albert.einstein/data/mdc/frames/" would cause the SineGaussian injections to be made in the directories under "/home/albert.einstein/data/mdc/frames/sg" + project : str + The name of the project which this frame is a part of. Defaults to 'Minke'. channel : str The name of the channel which the injections should be made into. This is prepended by the initials for each interferometer, so there will be a channel for each interferometer in the gwf. @@ -681,29 +692,37 @@ class Frame(): family = mdc.waveforms[0].waveform epoch = lal.LIGOTimeGPS(self.start) filename = "{}-{}-{}-{}.gwf".format(ifosstr, family, self.start, self.duration) + self.frame = lalframe.FrameNew(epoch = epoch, duration = self.duration, project='', run=1, frnum=1, detectorFlags=lal.LALDETECTORTYPE_ABSENT) + + + ifobits = np.array([getattr(lal,"{}_DETECTOR_BIT".format(lal.cached_detector_by_prefix[ifo].frDetector.name.upper())) + for ifo in self.ifos]) + ifoflag = numpy.bitwise_or.reduce(ifobits) + RUN_NUM = -1 # Simulated data should have a negative run number + head_date = str(self.start)[:5] frameloc = directory+"/"+mdc.directory_path()+"/"+head_date+"/" mkdir(frameloc) if not os.path.isfile(frameloc + filename) or force: + epoch = lal.LIGOTimeGPS(self.start) + frame = lalframe.FrameNew(epoch, self.duration, project, RUN_NUM, self.number, ifoflag) data = [] - # Define the start point of the time series top be generated for the injection - # Loop through each interferometer print self.ifos for ifo in self.ifos: # Calculate the number of samples in the timeseries nsamp = int((self.end-self.start)*rate) # Make the timeseries - h_resp = lal.CreateREAL8TimeSeries("inj time series", epoch, 0, 1.0/rate, lal.StrainUnit, nsamp) + h_resp = lal.CreateREAL8TimeSeries("{}:{}".format(ifo, channel), epoch, 0, 1.0/rate, lal.StrainUnit, nsamp) # Loop over all of the injections corresponding to this frame rowlist = self.get_rowlist(mdc) print rowlist if len(rowlist)==0: return for row in rowlist: - sim_burst = mdc.waveforms[row] + sim_burst = mdc.waveforms[row]._row() if sim_burst.hrss > 1: distance = sim_burst.amplitude @@ -715,6 +734,7 @@ class Frame(): # Apply detector response det = lalsimulation.DetectorPrefixToLALDetector(ifo) # Produce the total strains + print hp, hx, sim_burst.ra h_tot = lalsimulation.SimDetectorStrainREAL8TimeSeries(hp, hx, sim_burst.ra, sim_burst.dec, sim_burst.psi, det) # Inject the waveform into the overall timeseries @@ -729,19 +749,14 @@ class Frame(): len(h_resp.data.data)) data.data.data = h_resp.data.data lalframe.FrameAddREAL8TimeSeriesSimData(self.frame, data) - # data.append({"name": "%s:%s" % (ifo, channel), - # "data": h_resp.data.data, - # "start": float(epoch), - # "dx": h_resp.deltaT, - # "kind": "SIM"}) + + lalframe.FrameAddREAL8TimeSeriesSimData(frame, h_resp) # Make the directory in which to store the files # if it doesn't exist already # Write out the frame file - print "Frame Writing" - lalframe.FrameWrite(self.frame, frameloc+filename) - #Fr.frputvect(frameloc+filename, data) + lalframe.FrameWrite(frame, frameloc+filename) class HWInj(Frame): diff --git a/minke/sources.py b/minke/sources.py index bfb44a9849cb4d1dc343b9154b6e3fb1514d1b9a..625f71d0a1800ac1be043b3b7a5e535e614a4166 100644 --- a/minke/sources.py +++ b/minke/sources.py @@ -37,6 +37,10 @@ except ImportError: import matplotlib.pyplot as plt +if "numrel_data" in lsctables.SimBurstTable.validcolumns.keys(): + NROK = True +else: + NROK = False class Waveform(object): @@ -58,7 +62,9 @@ class Waveform(object): for a in self.table_type.validcolumns.keys(): self.params[a] = None - + def __getattr__(self, name): + return self.params[name] + def parse_polarisation(self, polarisation): """ Convert a string description of a polarisation to an ellipse eccentricity and an ellipse angle. @@ -156,11 +162,12 @@ class Waveform(object): continue burstobj.waveform = str(self.waveform) - - if swig_row.numrel_data: - burstobj.numrel_data = str(swig_row.numrel_data) - else: - burstobj.numrel_data = str("") + + if NROK: + if swig_row.numrel_data: + burstobj.numrel_data = str(swig_row.numrel_data) + else: + burstobj.numrel_data = str("") return burstobj @@ -200,10 +207,11 @@ class Waveform(object): for a in self.table_type.validcolumns.keys(): setattr(row, a, self.params[a]) - if self.numrel_data: - row.numrel_data = str(self.numrel_data) - else: - row.numrel_data = self.params['numrel_data'] + if NROK: + if self.numrel_data: + row.numrel_data = str(self.numrel_data) + else: + row.numrel_data = self.params['numrel_data'] row.waveform = self.waveform # Fill in the time @@ -211,6 +219,9 @@ class Waveform(object): # Get the sky locations if not row.ra: row.ra, row.dec, row.psi = self.sky_dist() + row.ra = row.ra[0] + row.dec = row.dec[0] + row.psi = row.psi[0] row.simulation_id = sim.get_next_id() row.waveform_number = random.randint(0,int(2**32)-1) ### !! This needs to be updated. @@ -294,13 +305,13 @@ class SineGaussian(Waveform): """ self._clear_params() self.sky_dist = sky_dist - self.params['hrss'] = hrss - self.params['seed'] = seed - self.params['frequency'] = frequency - self.params['q'] = q + self.hrss = self.params['hrss'] = hrss + self.seed = self.params['seed'] = seed + self.frequency = self.params['frequency'] = frequency + self.q = self.params['q'] = q self.time = time self.polarisation = polarisation - self.params['pol_ellipse_e'], self.params['pol_ellipse_angle'] = self.parse_polarisation(self.polarisation) + self.pol_ellipse_e, self.ellipse_angle = self.params['pol_ellipse_e'], self.params['pol_ellipse_angle'] = self.parse_polarisation(self.polarisation)