Verified Commit ba0d18b6 authored by Anonymous's avatar Anonymous Committed by Daniel Williams
Browse files

Added more work towards ringdowns, and moving xml reading to glue rather than...

Added more work towards ringdowns, and moving xml reading to glue rather than weird bits of lalsuite.
parent b506a0e6
Loading
Loading
Loading
Loading
+42 −22
Original line number Diff line number Diff line
@@ -31,6 +31,9 @@ from pylal.xlal.datatypes.ligotimegps import LIGOTimeGPS
from glue.ligolw.utils import process
import glue

import glue.ligolw
import gzip 

import lal, lalframe
from pylal import Fr

@@ -42,6 +45,7 @@ import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import re
import random


table_types = {
@@ -57,7 +61,8 @@ table_types = {
    # Long Duration
    "adi" : lsctables.SimBurstTable,
    # Ringdown
    "rng" : lsctables.SimRingdownTable
    "rng" : lsctables.SimRingdownTable,
    "gng" : lsctables.SimRingdownTable,
    }

tables = {
@@ -99,7 +104,8 @@ class MDCSet():
                          # Long-duration
                          'adi' : 'ADI',
                          # Ringdown
                          'rng' : "Ringdown",
                          'rng' : "BBHRingdown",
                          'gng' : "Ringdown",
                          }

    inj_families_abb = dict((v,k) for k,v in inj_families_names.iteritems())
@@ -182,13 +188,15 @@ class MDCSet():
                waveform_row = waveform._row(sim)
                waveform_row.process_id = procrow.process_id
            except:

                row = sim.RowType()

                for a in lsctables.SimBurstTable.validcolumns.keys():
                for a in self.table_type.validcolumns.keys():
                    if not hasattr(waveform, a):
                        setattr(row, a, 0)
                    else:
                        setattr(row, a, getattr(waveform, a))

                row.waveform = waveform.waveform
                if self.table_type == lsctables.SimBurstTable:
                    # Fill in the time
                    row.set_time_geocent(GPS(float(waveform.time)))
                    # Get the sky locations
@@ -235,24 +243,30 @@ class MDCSet():
        file. This should be fixed so that the object works symmetrically.
        """
        i = 0
        sim_burst_table = lalburst.SimBurstTableFromLIGOLw(filename, start, stop)
        while True:
        #sim_burst_table = lalburst.SimBurstTableFromLIGOLw(filename, start, stop)

        xml = glue.ligolw.utils.load_filename("test_simbursttable.xml.gz",
                                              contenthandler = glue.ligolw.ligolw.LIGOLWContentHandler,
                                              verbose = True)
        sim_burst_table = glue.ligolw.table.get_table(xml, self.table_type.tableName)
        
        for i,simrow in enumerate(sim_burst_table):
            # This is an ugly kludge to get around the poor choice of wavform name in the xmls, and
            if sim_burst_table.waveform[:3]=="s15": 
            if simrow.waveform[:3]=="s15": 
                self.numrel_file = str(sim_burst_table.waveform)
                sim_burst_table.waveform = "Dimmelmeier+08"

            self.waveforms.append(sim_burst_table)
            self.waveforms.append(simrow)#_burst_table)
            if full:
                self._measure_hrss(i)
                self._measure_egw_rsq(i)
            
            self.times = np.append(self.times, float(sim_burst_table.time_geocent_gps))
            self.times = np.append(self.times, float(simrow.time_geocent_gps))
            #np.insert(self.times, len(self.times), sim_burst_table.time_geocent_gps)
            if sim_burst_table.next is None: break
            sim_burst_table = sim_burst_table.next
            i += 1
        del(sim_burst_table)
            #if sim_burst_table.next is None: break
            #sim_burst_table = sim_burst_table.next
            #i += 1
        #del(sim_burst_table)
            
    def _generate_burst(self,row,rate=16384.0):
        """
@@ -280,7 +294,7 @@ class MDCSet():
        # This is a temporary kludge to allow LALSimulation to
        # be bypassed for pre-calculated waveforms.
        # A more robust solution should be considered.
        exceptions = ["Ott+13", "Mueller+12", "Scheidegger+10"]
        exceptions = ["Ott+13", "Mueller+12", "Scheidegger+10", "ADI"]
        row = self.waveforms[row]
        swig_row = lalburst.CreateSimBurst()
        for a in lsctables.SimBurstTable.validcolumns.keys():
@@ -291,7 +305,13 @@ class MDCSet():
                print a, getattr(row,a)
                continue # the structure is different than the TableRow
        theta, phi = np.cos(swig_row.incl), swig_row.phi
        try:
            swig_row.numrel_data = row.numrel_data
        except AttributeError:
            # Don't fail badly if there's a problem with numerical relativity data unless this is a supernova injection
            swig_row.numrel_data = ""
            if row.waveform in exceptions:
                raise AttributeError
        hp, hx = lalburst.GenerateSimBurst(swig_row, 1.0/rate)
        hp0, hx0 = lalburst.GenerateSimBurst(swig_row, 1.0/rate)
        return hp, hx, hp0, hx0
+27 −2
Original line number Diff line number Diff line
@@ -44,6 +44,7 @@ class Waveform(object):
    """

    sim = lsctables.New(lsctables.SimBurstTable)
    table_type = lsctables.SimBurstTable
    
    numrel_data = []
    waveform = "Generic"
@@ -205,7 +206,7 @@ class Waveform(object):
        if not sim: sim = self.sim
        row = sim.RowType()

        for a in lsctables.SimBurstTable.validcolumns.keys():
        for a in self.table_type.validcolumns.keys():
            setattr(row, a, self.params[a])

        row.waveform = self.waveform
@@ -276,6 +277,7 @@ class SineGaussian(Waveform):
        self.params['pol_ellipse_e'], self.params['pol_ellipse_angle'] = self.parse_polarisation(self.polarisation)    



class Gaussian(Waveform):
    """
    A class to represent a Gaussian injection.
@@ -960,6 +962,7 @@ class Ringdown(Waveform):
    """
    A class to handle Rindown waveforms.
    """
    table_type = lsctables.SimRingdownTable
    waveform = "Ringdown"
    pass

@@ -1091,3 +1094,25 @@ class ADI(LongDuration):
        output[:,2] = strainc_new

        return output

class BBHRingdown(Ringdown):
    """
    A class to represent BBH ringdowns.
    """
    #lalsimfunction = SimBlackHoleRingdown

    def __init__(self, time, hrss,  phi0, deltaT, mass, spin, massloss, distance, inclination, l, m, sky_dist=uniform_sky):
        self._clear_params()
        self.time = self.v_start_time_ns = time
        self.sky_dist = sky_dist
        self.params['simulation_id'] = self.simulation_id =  self.sim.get_next_id()
        self.params['phi0'] = phi0
        self.params['deltaT'] = deltaT
        self.params['mass'] = mass
        self.params['spin'] = spin
        self.params['massloss'] = massloss
        self.params['eff_dist_l'] = self.eff_dist_l = distance
        self.params['hrss'] = self.hrss = hrss
        self.params['inclination'] = inclination
        self.params['l'] = l
        self.params['m'] = m
+68 −1
Original line number Diff line number Diff line
@@ -11,7 +11,7 @@ Tests for `minke` module.
import unittest

import minke

from minke import mdctools, sources

class TestMinke(unittest.TestCase):

@@ -21,6 +21,73 @@ class TestMinke(unittest.TestCase):
    def tearDown(self):
        pass

    # Check that failure occurs if the wrong waveforms are inserted into the wrong tables

    def test_insert_wrong_waveform_to_table(self):
        """Test whether adding a ringdown waveform to a SimBurstTable throws
        an error.
        """
        mdcset = mdctools.MDCSet(["L1"], table_type = "burst")
        ring = sources.Ringdown()

        with self.assertRaises(mdctools.TableTypeError):
            mdcset + ring

    def test_insert_burst_waveform_to_burst_table(self):
        """
        Test whether inserting the correct type of waveform into a table works
        """
        mdcset = mdctools.MDCSet(["L1"], table_type = "burst")
        waveform = sources.Gaussian(0,0,0)

        mdcset + waveform

        self.assertEqual(len(mdcset.waveforms), 1)

    # Check that the correct XML files are produced for different table types

    def test_write_simbursttable(self):
        """
        Write out a simburst table xml file
        """
        mdcset = mdctools.MDCSet(["L1"], table_type = "burst")
        waveform = sources.Gaussian(0.1,1e-23,1000)

        mdcset + waveform

        mdcset.save_xml("test_simbursttable.xml.gz")

    def test_verify_simbursttable(self):
        """
        Read-in the xml simburst table.
        """
        mdcset = mdctools.MDCSet(["L1"])
        mdcset.load_xml("test_simbursttable.xml.gz", full = False)

        self.assertEqual(len(mdcset.waveforms), 1)

    def test_write_simringdowntable(self):
        """
        Write out a simburst table xml file
        """
        mdcset = mdctools.MDCSet(["L1"], table_type = "ringdown")
        waveform = sources.BBHRingdown(1000, 1e-22, 0 ,0.1, 10, 0.1, 0.01, 10, 0, 2, 2,)

        mdcset + waveform

        mdcset.save_xml("test_simringdowntable.xml.gz")

    def test_verify_simringdowntable(self):
        """
        Read-in the xml simburst table.
        """
        mdcset = mdctools.MDCSet(["L1"])
        mdcset.load_xml("test_simringdowntable.xml.gz", full = False)

        self.assertEqual(len(mdcset.waveforms), 1)

        
        
if __name__ == '__main__':
    import sys
    sys.exit(unittest.main())