Compute the H/V spectral ratio from the ratios of PSDs

import os

if "SPHINX_DOC_BUILD" in os.environ:
    if "MSNOISE_DOC" in os.environ:
        os.chdir(os.environ["MSNOISE_DOC"])

import matplotlib
# matplotlib.use("agg")

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from pandas.plotting import register_matplotlib_converters
import datetime

register_matplotlib_converters()

plt.style.use("ggplot")

from msnoise.api import *


# connect to the database
db = connect()

# Obtain a list of dates between ``start_date`` and ``enddate``
start, end, datelist = build_movstack_datelist(db)

# Get the list of parameters from the DB:
params = get_params(db)

# Get the time axis for plotting the CCF:
taxis = get_t_axis(db)

# Get the first mov_stack configured:
mov_stack = params.mov_stack[0]



def convert_to_velocity(df):
    df = df.resample("30Min").mean()
    df.columns = 1. / df.columns
    df = df.sort_index(axis=1)
    df = df.dropna(axis=1, how='all')

    w2f = (2.0 * np.pi * df.columns)
    # The acceleration amplitude spectrum and velocity spectral amplitude (not power)
    vamp = np.sqrt(10.0 ** (df / 10.) / w2f ** 2)
    return vamp



Z = hdf_open_store("PF.FJS.00.HHZ", mode="r").PSD
Z = convert_to_velocity(Z)
E = hdf_open_store("PF.FJS.00.HHE", mode="r").PSD
E = convert_to_velocity(E)
N = hdf_open_store("PF.FJS.00.HHN", mode="r").PSD
N = convert_to_velocity(N)



r = "6h"
rZ = Z.resample(r).mean()
rE = E.resample(r).mean()
rN = N.resample(r).mean()

merged = pd.concat({'Z': rZ, 'E': rE, 'N': rN}, axis=0)
merged.index.names = ['channel', 'date']

# Swap the levels of the MultiIndex
result = merged.swaplevel('channel', 'date').sort_index(level='date')  # .iloc[:500]
result



HVSRs = {}
for date, group in result.groupby(level='date'):
    print(f"Date: {date}")
    group = group.droplevel(0)
    try:
        iZ = group.loc["Z"].values
        iE = group.loc["E"].values
        iN = group.loc["N"].values
    except:
        continue
    hvsr = np.sqrt((iE + iN) / iZ)
    HVSRs[date] = hvsr



HVSR = pd.DataFrame(HVSRs, index=result.columns)



X = HVSR.copy().fillna(0.0)
X = X.T.resample(r).mean().T
Date: 2014-03-20 00:00:00
Date: 2014-03-20 06:00:00
Date: 2014-03-20 12:00:00
Date: 2014-03-20 18:00:00
Date: 2014-03-21 00:00:00
Date: 2014-03-21 06:00:00
Date: 2014-03-21 12:00:00
Date: 2014-03-21 18:00:00
Date: 2014-03-22 00:00:00
Date: 2014-03-22 06:00:00
Date: 2014-03-22 12:00:00
Date: 2014-03-22 18:00:00
Date: 2014-03-23 00:00:00
Date: 2014-03-23 06:00:00
Date: 2014-03-23 12:00:00
Date: 2014-03-23 18:00:00
Date: 2014-03-24 00:00:00
Date: 2014-03-24 06:00:00
Date: 2014-03-24 12:00:00
Date: 2014-03-24 18:00:00
Date: 2014-03-25 00:00:00
Date: 2014-03-25 06:00:00
Date: 2014-03-25 12:00:00
Date: 2014-03-25 18:00:00
Date: 2014-03-26 00:00:00
Date: 2014-03-26 06:00:00
Date: 2014-03-26 12:00:00
Date: 2014-03-26 18:00:00
Date: 2014-03-27 00:00:00
Date: 2014-03-27 06:00:00
Date: 2014-03-27 12:00:00
Date: 2014-03-27 18:00:00
Date: 2014-03-28 00:00:00
Date: 2014-03-28 06:00:00
Date: 2014-03-28 12:00:00
Date: 2014-03-28 18:00:00
Date: 2014-03-29 00:00:00
Date: 2014-03-29 06:00:00
Date: 2014-03-29 12:00:00
Date: 2014-03-29 18:00:00
Date: 2014-03-30 00:00:00
Date: 2014-03-30 06:00:00
Date: 2014-03-30 12:00:00
Date: 2014-03-30 18:00:00
Date: 2014-03-31 00:00:00
Date: 2014-03-31 06:00:00
Date: 2014-03-31 12:00:00
Date: 2014-03-31 18:00:00
Date: 2014-04-01 00:00:00
Date: 2014-04-01 06:00:00
Date: 2014-04-01 12:00:00
Date: 2014-04-01 18:00:00
Date: 2014-04-02 00:00:00
Date: 2014-04-02 06:00:00
Date: 2014-04-02 12:00:00
Date: 2014-04-02 18:00:00
Date: 2014-04-03 00:00:00
Date: 2014-04-03 06:00:00
Date: 2014-04-03 12:00:00
Date: 2014-04-03 18:00:00
Date: 2014-04-04 00:00:00
Date: 2014-04-04 06:00:00
Date: 2014-04-04 12:00:00
Date: 2014-04-04 18:00:00
Date: 2014-04-05 00:00:00
Date: 2014-04-05 06:00:00
Date: 2014-04-05 12:00:00
Date: 2014-04-05 18:00:00
Date: 2014-04-06 00:00:00
Date: 2014-04-06 06:00:00
Date: 2014-04-06 12:00:00
Date: 2014-04-06 18:00:00
Date: 2014-04-07 00:00:00
Date: 2014-04-07 06:00:00
Date: 2014-04-07 12:00:00
Date: 2014-04-07 18:00:00
Date: 2014-04-08 00:00:00
Date: 2014-04-08 06:00:00
Date: 2014-04-08 12:00:00
Date: 2014-04-08 18:00:00
Date: 2014-04-09 00:00:00
Date: 2014-04-09 06:00:00
Date: 2014-04-09 12:00:00
Date: 2014-04-09 18:00:00
Date: 2014-04-10 00:00:00
Date: 2014-04-10 06:00:00
Date: 2014-04-10 12:00:00
Date: 2014-04-10 18:00:00
Date: 2014-04-11 00:00:00
Date: 2014-04-11 06:00:00
Date: 2014-04-11 12:00:00
Date: 2014-04-11 18:00:00
Date: 2014-04-12 00:00:00
Date: 2014-04-12 06:00:00
Date: 2014-04-12 12:00:00
Date: 2014-04-12 18:00:00
Date: 2014-04-13 00:00:00
Date: 2014-04-13 06:00:00
Date: 2014-04-13 12:00:00
Date: 2014-04-13 18:00:00
Date: 2014-04-14 00:00:00
Date: 2014-04-14 06:00:00
Date: 2014-04-14 12:00:00
Date: 2014-04-14 18:00:00
Date: 2014-04-15 00:00:00
Date: 2014-04-15 06:00:00
Date: 2014-04-15 12:00:00
Date: 2014-04-15 18:00:00
Date: 2014-04-16 00:00:00
Date: 2014-04-16 06:00:00
Date: 2014-04-16 12:00:00
Date: 2014-04-16 18:00:00
Date: 2014-04-17 00:00:00
Date: 2014-04-17 06:00:00
Date: 2014-04-17 12:00:00
Date: 2014-04-17 18:00:00
Date: 2014-04-18 00:00:00
Date: 2014-04-18 06:00:00
Date: 2014-04-18 12:00:00
Date: 2014-04-18 18:00:00
Date: 2014-04-19 00:00:00
Date: 2014-04-19 06:00:00
Date: 2014-04-19 12:00:00
Date: 2014-04-19 18:00:00
Date: 2014-04-20 00:00:00
Date: 2014-04-20 06:00:00
Date: 2014-04-20 12:00:00
Date: 2014-04-20 18:00:00
Date: 2014-04-21 00:00:00
Date: 2014-04-21 06:00:00
Date: 2014-04-21 12:00:00
Date: 2014-04-21 18:00:00
Date: 2014-04-22 00:00:00
Date: 2014-04-22 06:00:00
Date: 2014-04-22 12:00:00
Date: 2014-04-22 18:00:00
Date: 2014-04-23 00:00:00
Date: 2014-04-23 06:00:00
Date: 2014-04-23 12:00:00
Date: 2014-04-23 18:00:00
Date: 2014-04-24 00:00:00
Date: 2014-04-24 06:00:00
Date: 2014-04-24 12:00:00
Date: 2014-04-24 18:00:00
Date: 2014-04-25 00:00:00
Date: 2014-04-25 06:00:00
Date: 2014-04-25 12:00:00
Date: 2014-04-25 18:00:00
Date: 2014-04-26 00:00:00
Date: 2014-04-26 06:00:00
Date: 2014-04-26 12:00:00
Date: 2014-04-26 18:00:00
Date: 2014-04-27 00:00:00
Date: 2014-04-27 06:00:00
Date: 2014-04-27 12:00:00
Date: 2014-04-27 18:00:00
Date: 2014-04-28 00:00:00
Date: 2014-04-28 06:00:00
Date: 2014-04-28 12:00:00
Date: 2014-04-28 18:00:00
Date: 2014-04-29 00:00:00
Date: 2014-04-29 06:00:00
Date: 2014-04-29 12:00:00
Date: 2014-04-29 18:00:00
Date: 2014-04-30 00:00:00
Date: 2014-04-30 06:00:00
Date: 2014-04-30 12:00:00
Date: 2014-04-30 18:00:00
Date: 2014-05-01 00:00:00
Date: 2014-05-01 06:00:00
Date: 2014-05-01 12:00:00
Date: 2014-05-01 18:00:00
Date: 2014-05-02 00:00:00
Date: 2014-05-02 06:00:00
Date: 2014-05-02 12:00:00
Date: 2014-05-02 18:00:00
Date: 2014-05-03 00:00:00
Date: 2014-05-03 06:00:00
Date: 2014-05-03 12:00:00
Date: 2014-05-03 18:00:00
Date: 2014-05-04 00:00:00
Date: 2014-05-04 06:00:00
Date: 2014-05-04 12:00:00
Date: 2014-05-04 18:00:00
Date: 2014-05-05 00:00:00
Date: 2014-05-05 06:00:00
Date: 2014-05-05 12:00:00
Date: 2014-05-05 18:00:00
Date: 2014-05-06 00:00:00
Date: 2014-05-06 06:00:00
Date: 2014-05-06 12:00:00
Date: 2014-05-06 18:00:00
Date: 2014-05-07 00:00:00
Date: 2014-05-07 06:00:00
Date: 2014-05-07 12:00:00
Date: 2014-05-07 18:00:00
Date: 2014-05-08 00:00:00
Date: 2014-05-08 06:00:00
Date: 2014-05-08 12:00:00
Date: 2014-05-08 18:00:00
Date: 2014-05-09 00:00:00
Date: 2014-05-09 06:00:00
Date: 2014-05-09 12:00:00
Date: 2014-05-09 18:00:00
Date: 2014-05-10 00:00:00
Date: 2014-05-10 06:00:00
Date: 2014-05-10 12:00:00
Date: 2014-05-10 18:00:00
Date: 2014-05-11 00:00:00
Date: 2014-05-11 06:00:00
Date: 2014-05-11 12:00:00
Date: 2014-05-11 18:00:00
Date: 2014-05-12 00:00:00
Date: 2014-05-12 06:00:00
Date: 2014-05-12 12:00:00
Date: 2014-05-12 18:00:00
Date: 2014-05-13 00:00:00
Date: 2014-05-13 06:00:00
Date: 2014-05-13 12:00:00
Date: 2014-05-13 18:00:00
Date: 2014-05-14 00:00:00
Date: 2014-05-14 06:00:00
Date: 2014-05-14 12:00:00
Date: 2014-05-14 18:00:00
Date: 2014-05-15 00:00:00
Date: 2014-05-15 06:00:00
Date: 2014-05-15 12:00:00
Date: 2014-05-15 18:00:00
Date: 2014-05-16 00:00:00
Date: 2014-05-16 06:00:00
Date: 2014-05-16 12:00:00
Date: 2014-05-16 18:00:00
Date: 2014-05-17 00:00:00
Date: 2014-05-17 06:00:00
Date: 2014-05-17 12:00:00
Date: 2014-05-17 18:00:00
Date: 2014-05-18 00:00:00
Date: 2014-05-18 06:00:00
Date: 2014-05-18 12:00:00
Date: 2014-05-18 18:00:00
Date: 2014-05-19 00:00:00
Date: 2014-05-19 06:00:00
Date: 2014-05-19 12:00:00
Date: 2014-05-19 18:00:00
Date: 2014-05-20 00:00:00
Date: 2014-05-20 06:00:00
Date: 2014-05-20 12:00:00
Date: 2014-05-20 18:00:00
Date: 2014-05-21 00:00:00
Date: 2014-05-21 06:00:00
Date: 2014-05-21 12:00:00
Date: 2014-05-21 18:00:00
Date: 2014-05-22 00:00:00
Date: 2014-05-22 06:00:00
Date: 2014-05-22 12:00:00
Date: 2014-05-22 18:00:00
Date: 2014-05-23 00:00:00
Date: 2014-05-23 06:00:00
Date: 2014-05-23 12:00:00
Date: 2014-05-23 18:00:00
Date: 2014-05-24 00:00:00
Date: 2014-05-24 06:00:00
Date: 2014-05-24 12:00:00
Date: 2014-05-24 18:00:00
Date: 2014-05-25 00:00:00
Date: 2014-05-25 06:00:00
Date: 2014-05-25 12:00:00
Date: 2014-05-25 18:00:00
Date: 2014-05-26 00:00:00
Date: 2014-05-26 06:00:00
Date: 2014-05-26 12:00:00
Date: 2014-05-26 18:00:00
Date: 2014-05-27 00:00:00
Date: 2014-05-27 06:00:00
Date: 2014-05-27 12:00:00
Date: 2014-05-27 18:00:00
Date: 2014-05-28 00:00:00
Date: 2014-05-28 06:00:00
Date: 2014-05-28 12:00:00
Date: 2014-05-28 18:00:00
Date: 2014-05-29 00:00:00
Date: 2014-05-29 06:00:00
Date: 2014-05-29 12:00:00
Date: 2014-05-29 18:00:00
Date: 2014-05-30 00:00:00
Date: 2014-05-30 06:00:00
Date: 2014-05-30 12:00:00
Date: 2014-05-30 18:00:00
Date: 2014-05-31 00:00:00
Date: 2014-05-31 06:00:00
Date: 2014-05-31 12:00:00
Date: 2014-05-31 18:00:00
Date: 2014-06-01 00:00:00
Date: 2014-06-01 06:00:00
Date: 2014-06-01 12:00:00
Date: 2014-06-01 18:00:00
Date: 2014-06-02 00:00:00
Date: 2014-06-02 06:00:00
Date: 2014-06-02 12:00:00
Date: 2014-06-02 18:00:00
Date: 2014-06-03 00:00:00
Date: 2014-06-03 06:00:00
Date: 2014-06-03 12:00:00
Date: 2014-06-03 18:00:00
Date: 2014-06-04 00:00:00
Date: 2014-06-04 06:00:00
Date: 2014-06-04 12:00:00
Date: 2014-06-04 18:00:00
Date: 2014-06-05 00:00:00
Date: 2014-06-05 06:00:00
Date: 2014-06-05 12:00:00
Date: 2014-06-05 18:00:00
Date: 2014-06-06 00:00:00
Date: 2014-06-06 06:00:00
Date: 2014-06-06 12:00:00
Date: 2014-06-06 18:00:00
Date: 2014-06-07 00:00:00
Date: 2014-06-07 06:00:00
Date: 2014-06-07 12:00:00
Date: 2014-06-07 18:00:00
Date: 2014-06-08 00:00:00
Date: 2014-06-08 06:00:00
Date: 2014-06-08 12:00:00
Date: 2014-06-08 18:00:00
Date: 2014-06-09 00:00:00
Date: 2014-06-09 06:00:00
Date: 2014-06-09 12:00:00
Date: 2014-06-09 18:00:00
Date: 2014-06-10 00:00:00
Date: 2014-06-10 06:00:00
Date: 2014-06-10 12:00:00
Date: 2014-06-10 18:00:00
Date: 2014-06-11 00:00:00
Date: 2014-06-11 06:00:00
Date: 2014-06-11 12:00:00
Date: 2014-06-11 18:00:00
Date: 2014-06-12 00:00:00
Date: 2014-06-12 06:00:00
Date: 2014-06-12 12:00:00
Date: 2014-06-12 18:00:00
Date: 2014-06-13 00:00:00
Date: 2014-06-13 06:00:00
Date: 2014-06-13 12:00:00
Date: 2014-06-13 18:00:00
Date: 2014-06-14 00:00:00
Date: 2014-06-14 06:00:00
Date: 2014-06-14 12:00:00
Date: 2014-06-14 18:00:00
Date: 2014-06-15 00:00:00
Date: 2014-06-15 06:00:00
Date: 2014-06-15 12:00:00
Date: 2014-06-15 18:00:00
Date: 2014-06-16 00:00:00
Date: 2014-06-16 06:00:00
Date: 2014-06-16 12:00:00
Date: 2014-06-16 18:00:00
Date: 2014-06-17 00:00:00
Date: 2014-06-17 06:00:00
Date: 2014-06-17 12:00:00
Date: 2014-06-17 18:00:00
Date: 2014-06-18 00:00:00
Date: 2014-06-18 06:00:00
Date: 2014-06-18 12:00:00
Date: 2014-06-18 18:00:00
Date: 2014-06-19 00:00:00
Date: 2014-06-19 06:00:00
Date: 2014-06-19 12:00:00
Date: 2014-06-19 18:00:00
Date: 2014-06-20 00:00:00
Date: 2014-06-20 06:00:00
Date: 2014-06-20 12:00:00
Date: 2014-06-20 18:00:00
Date: 2014-06-21 00:00:00
Date: 2014-06-21 06:00:00
Date: 2014-06-21 12:00:00
Date: 2014-06-21 18:00:00
Date: 2014-06-22 00:00:00
Date: 2014-06-22 06:00:00
Date: 2014-06-22 12:00:00
Date: 2014-06-22 18:00:00
Date: 2014-06-23 00:00:00
Date: 2014-06-23 06:00:00
Date: 2014-06-23 12:00:00
Date: 2014-06-23 18:00:00
Date: 2014-06-24 00:00:00
Date: 2014-06-24 06:00:00
Date: 2014-06-24 12:00:00
Date: 2014-06-24 18:00:00
Date: 2014-06-25 00:00:00
Date: 2014-06-25 06:00:00
Date: 2014-06-25 12:00:00
Date: 2014-06-25 18:00:00
Date: 2014-06-26 00:00:00
Date: 2014-06-26 06:00:00
Date: 2014-06-26 12:00:00
Date: 2014-06-26 18:00:00
Date: 2014-06-27 00:00:00
Date: 2014-06-27 06:00:00
Date: 2014-06-27 12:00:00
Date: 2014-06-27 18:00:00
Date: 2014-06-28 00:00:00
Date: 2014-06-28 06:00:00
Date: 2014-06-28 12:00:00
Date: 2014-06-28 18:00:00
Date: 2014-06-29 00:00:00
Date: 2014-06-29 06:00:00
Date: 2014-06-29 12:00:00
Date: 2014-06-29 18:00:00
Date: 2014-06-30 00:00:00
Date: 2014-06-30 06:00:00
Date: 2014-06-30 12:00:00
Date: 2014-06-30 18:00:00
Date: 2014-07-01 00:00:00
Date: 2014-07-01 06:00:00
Date: 2014-07-01 12:00:00
Date: 2014-07-01 18:00:00
Date: 2014-07-02 00:00:00
Date: 2014-07-02 06:00:00
Date: 2014-07-02 12:00:00
Date: 2014-07-02 18:00:00
Date: 2014-07-03 00:00:00
Date: 2014-07-03 06:00:00
Date: 2014-07-03 12:00:00
Date: 2014-07-03 18:00:00
Date: 2014-07-04 00:00:00
Date: 2014-07-04 06:00:00
Date: 2014-07-04 12:00:00
Date: 2014-07-04 18:00:00
Date: 2014-07-05 00:00:00
Date: 2014-07-05 06:00:00
Date: 2014-07-05 12:00:00
Date: 2014-07-05 18:00:00
Date: 2014-07-06 00:00:00
Date: 2014-07-06 06:00:00
Date: 2014-07-06 12:00:00
Date: 2014-07-06 18:00:00
Date: 2014-07-07 00:00:00
Date: 2014-07-07 06:00:00
Date: 2014-07-07 12:00:00
Date: 2014-07-07 18:00:00
Date: 2014-07-08 00:00:00
Date: 2014-07-08 06:00:00
Date: 2014-07-08 12:00:00
Date: 2014-07-08 18:00:00
Date: 2014-07-09 00:00:00
Date: 2014-07-09 06:00:00
Date: 2014-07-09 12:00:00
Date: 2014-07-09 18:00:00
Date: 2014-07-10 00:00:00
Date: 2014-07-10 06:00:00
Date: 2014-07-10 12:00:00
Date: 2014-07-10 18:00:00
Date: 2014-07-11 00:00:00
Date: 2014-07-11 06:00:00
Date: 2014-07-11 12:00:00
Date: 2014-07-11 18:00:00
Date: 2014-07-12 00:00:00
Date: 2014-07-12 06:00:00
Date: 2014-07-12 12:00:00
Date: 2014-07-12 18:00:00
Date: 2014-07-13 00:00:00
Date: 2014-07-13 06:00:00
Date: 2014-07-13 12:00:00
Date: 2014-07-13 18:00:00
Date: 2014-07-14 00:00:00
Date: 2014-07-14 06:00:00
Date: 2014-07-14 12:00:00
Date: 2014-07-14 18:00:00
Date: 2014-07-15 00:00:00
Date: 2014-07-15 06:00:00
Date: 2014-07-15 12:00:00
Date: 2014-07-15 18:00:00
Date: 2014-07-16 00:00:00
Date: 2014-07-16 06:00:00
Date: 2014-07-16 12:00:00
Date: 2014-07-16 18:00:00
Date: 2014-07-17 00:00:00
Date: 2014-07-17 06:00:00
Date: 2014-07-17 12:00:00
Date: 2014-07-17 18:00:00
Date: 2014-07-18 00:00:00
Date: 2014-07-18 06:00:00
Date: 2014-07-18 12:00:00
Date: 2014-07-18 18:00:00
Date: 2014-07-19 00:00:00
Date: 2014-07-19 06:00:00
Date: 2014-07-19 12:00:00
Date: 2014-07-19 18:00:00

Plotting

plt.subplots(figsize=(18, 18))
plt.pcolormesh(X.index, X.columns, X.T, rasterized=True, vmax=4)
plt.colorbar()
plt.grid()
plt.xlabel("Frequency (Hz)")
plt.show()
plot compute hvsr psd

Plotting & zooming around 0.5-4 Hz

plt.subplots(figsize=(18, 18))
plt.pcolormesh(X.index, X.columns, X.T, rasterized=True, vmax=4)
plt.colorbar()
plt.xlim(0.5, 4)
plt.grid()
plt.xlabel("Frequency (Hz)")
plt.show()
plot compute hvsr psd

Plotting the HVSR curve (truncating the low frequency here):

X.loc[0.2:20].mean(axis=1).plot()
plt.xlabel("Frequency (Hz)")
plt.ylabel("Amplitude")
plt.show()
plot compute hvsr psd

Total running time of the script: (0 minutes 2.856 seconds)

Gallery generated by Sphinx-Gallery