Skip to content
GWpy 2.1.2 documentation logo
GWpy 2.1.2 5. Comparing seismic trends between LIGO sites
Type to start searching
    GWpy
    • TimeSeries examples
    GWpy
    • What is GWpy?
    • How do I install GWpy?
    • Citing GWpy
    • Data structures
    • Time Series data
    • The FrequencySeries class
    • The Spectrogram
    • State vectors
    • Data-quality segments
    • Data tables
    • Data manipulation
    • Signal processing
    • Visualising data
    • Plotting in GWpy (gwpy.plot)
    • Command line plotting with GWpy
    • Other utilities
    • The Channel class
    • Times and timestamps
    • Astrophysical modelling
    • Configuring GWpy from the environment
    • Examples
    • TimeSeries examples
    • Signal processing examples
    • FrequencySeries examples
    • Spectrogram examples
    • Segments examples
    • Tabular data examples
    • Miscellaneous examples
    • Developer notes
    • Publishing a release
    • Show Source

    5. Comparing seismic trends between LIGO sites¶

    On Feb 13 2015 there was a massive earthquake in the Atlantic Ocean, that should have had an impact on LIGO operations, I’d like to find out.

    First: we import the objects we need, one for getting the data:

    from gwpy.timeseries import TimeSeriesDict
    

    and one for plotting the data:

    from gwpy.plot import Plot
    

    Next we define the channels we want, namely the 0.03Hz-1Hz ground motion band-limited RMS channels (1-second average trends). We do this using string-replacement so we can substitute the interferometer prefix easily when we need to:

    channels = [
        '{ifo}:ISI-BS_ST1_SENSCOR_GND_STS_X_BLRMS_30M_100M.mean,s-trend',
        '{ifo}:ISI-BS_ST1_SENSCOR_GND_STS_Y_BLRMS_30M_100M.mean,s-trend',
        '{ifo}:ISI-BS_ST1_SENSCOR_GND_STS_Z_BLRMS_30M_100M.mean,s-trend',
    ]
    

    At last we can get() 12 hours of data for each interferometer:

    lho = TimeSeriesDict.get([c.format(ifo='H1') for c in channels],
                             'Feb 13 2015 16:00', 'Feb 14 2015 04:00')
    llo = TimeSeriesDict.get([c.format(ifo='L1') for c in channels],
                             'Feb 13 2015 16:00', 'Feb 14 2015 04:00')
    

    Next we can plot the data, with a separate Axes for each instrument:

    plot = Plot(lho, llo, figsize=(12, 6), sharex=True, yscale='log')
    ax1, ax2 = plot.axes
    for ifo, ax in zip(('Hanford', 'Livingston'), (ax1, ax2)):
        ax.legend(['X', 'Y', 'Z'])
        ax.text(1.01, 0.5, ifo, ha='left', va='center', transform=ax.transAxes,
                fontsize=18)
    ax1.set_ylabel(r'$1-3$\,Hz motion [nm/s]', y=-0.1)
    ax2.set_ylabel('')
    ax1.set_title('Magnitude 7.1 earthquake impact on LIGO')
    plot.show()
    

    (png)

    ../../../_images/blrms-5.png

    As we can see, the earthquake had a huge impact on the LIGO observatories, severly imparing operations for several hours.

    Previous 4. Cross-correlating two TimeSeries
    Next 6. Plotting segments for a StateVector
    © Copyright 2013, 2017-2021 Cardiff University and 2013-2017 Lousiana State University.
    Created using Sphinx 4.3.0. and Material for Sphinx