Post-processing package#

seispy.core.depmodel module#

class seispy.core.depmodel.DepModel(dep_range, velmod='iasp91', elevation=0.0, layer_mod=False)[source]#

Bases: object

Class for constructing 1D velocity model and computing delay time of Ps, PpPs and PsPs Rays.

Examples

>>> model = DepModel(np.array([0, 20.1, 35.1, 100]))
>>> print(model.dz)
[ 0.  20.1 15.  64.9]
>>> print(model.vp)
[5.8        6.5        8.04001059 8.04764706]
>>> print(model.vs)
[3.36       3.75       4.47003177 4.49294118]

Methods

ccp_model([dep_range, elevation, layerd])

import ccp configure and init DepModel object for time2depth convertion if any parameters given is wrong, return a default DepModel object

plot_model([show])

plot model with matplotlib

radius_s(rayp[, phase, sphere])

calculate horizontal radius from piercing point to station postion.

raylength(rayp[, phase, sphere])

calculate ray length, P for Sp and S for Ps

read_layer_model(dep_range, h, vp, vs[, ...])

Read layer model from given parameters

save_tvel(filename)

save vel mod in tvel format for taup

tpds(rayps, raypp[, sphere])

calculate travel time of Pds

tpppds(rayps, raypp[, sphere])

calculate travel time of Ppds

tpspds(rayps[, sphere])

calculate travel time of Ppsds

classmethod ccp_model(dep_range=array([0., 20.1, 35.1, 100.]), elevation=0, layerd=False, **kwargs)[source]#

import ccp configure and init DepModel object for time2depth convertion if any parameters given is wrong, return a default DepModel object

there’s 3 types of input is allowed:

  1. mod3d, stla, stlo: for 3d model( need modification

  2. modfolder, staname: for dir

  3. mod: for single file

plot_model(show=True)[source]#

plot model with matplotlib

Parameters:

show (bool, optional) – whether to show the plot, by default True

radius_s(rayp, phase='P', sphere=True)[source]#

calculate horizontal radius from piercing point to station postion.

P for Sp and S for Ps.

Parameters:
  • rayp (float or numpy.ndarray) – ray parameter

  • phase (str, optional) – phase name, by default ‘P’

  • sphere (bool, optional) – whether to use sphere earth, by default True

Return type:

numpy.ndarray

Returns:

horizontal radius

Examples

>>> model = DepModel(np.array([0, 20.1, 35.1, 100]))
>>> model.dz
array([ 0. , 20.1, 15. , 64.9])
>>> model.R
array([6371. , 6350.9, 6335.9, 6271. ])
>>> model.radius_s(1.2,phase="S", sphere=False)*111.2
array([0.        , 0.0002478 , 0.00046823, 0.00142685])
raylength(rayp, phase='P', sphere=True)[source]#

calculate ray length, P for Sp and S for Ps

Parameters:
  • rayp (float or numpy.ndarray) – ray parameter

  • phase (str, optional) – phase name, by default ‘P’

  • sphere (bool, optional) – whether to use sphere earth, by default True

Return type:

numpy.ndarray

Returns:

ray length

classmethod read_layer_model(dep_range, h, vp, vs, rho=None, elevation=0)[source]#

Read layer model from given parameters

Parameters:
  • dep_range (numpy.ndarray) – Depth range for conversion

  • h (numpy.ndarray) – Thickness of each layer

  • vp (numpy.ndarray) – P-wave velocity of each layer

  • vs (numpy.ndarray) – S-wave velocity of each layer

  • rho (numpy.ndarray) – Density of each layer, by default None

  • elevation (float) – Elevation in km, by default 0

Return type:

DepModel

Returns:

DepModel object

Examples

>>> dep_range = np.arange(100)
>>> h = np.array([20, 15., 0])
>>> vp = np.array([5.8, 6.5, 8.04])
>>> vs = np.array([3.36, 3.75, 4.47])
>>> model = DepModel.read_layer_model(dep_range, h, vp, vs)
save_tvel(filename)[source]#

save vel mod in tvel format for taup

Parameters:

filename (str) – output file name

Examples

>>> model = DepModel(np.array([0, 20.1, 35.1, 100]))
>>> model.save_tvel("test")
    0.00     5.80     3.36     2.72
   20.10    5.800     3.36     2.72
   20.10     6.50     3.75     2.92
   35.10    6.500     3.75     2.92
   35.10     8.04     4.47     3.32
  100.00    8.040     4.47     3.32
  100.00     8.05     4.49     3.36
tpds(rayps, raypp, sphere=True)[source]#

calculate travel time of Pds

Parameters:
  • rayps (float or numpy.ndarray) – ray parameter of Ps wave

  • raypp (float or numpy.ndarray) – ray parameter of P wave

  • sphere (bool, optional) – whether to use sphere earth, by default True

Return type:

numpy.ndarray

Returns:

travel time of Pds

tpppds(rayps, raypp, sphere=True)[source]#

calculate travel time of Ppds

Parameters:
  • rayps (float or numpy.ndarray) – ray parameter of Ps wave

  • raypp (float or numpy.ndarray) – ray parameter of P wave

  • sphere (bool, optional) – whether to use sphere earth, by default True

Return type:

numpy.ndarray

Returns:

travel time of Ppds

tpspds(rayps, sphere=True)[source]#

calculate travel time of Ppsds

Parameters:
  • rayps (float or numpy.ndarray) – ray parameter of Ps wave

  • sphere (bool, optional) – whether to use sphere earth, by default True

Return type:

numpy.ndarray

Returns:

travel time of Ppsds

seispy.core.depmodel.interp_depth_model(model, lat, lon, new_dep)[source]#

Interpolate Vp and Vs from 3D velocity with a specified depth range.

Parameters:
  • mod3d (np.lib.npyio.NpzFile()) – 3D velocity loaded from a .npz file

  • lat (float) – Latitude of position in 3D velocity model

  • lon (float) – Longitude of position in 3D velocity model

  • new_dep (np.ndarray()) – 1D array of depths in km

Return type:

np.ndarray()

Returns:

Vp and Vs in new_dep

seispy.hk module#

seispy.hk.ci(allstack, h, kappa, ev_num)[source]#

Search best H and kappa from stacked matrix. Calculate error for H and kappa

Parameters:
  • allstack – stacked HK matrix

  • h – 1-D array of H

  • kappa – 1-D array of kappa

  • ev_num – event number

Returns:

seispy.hk.hk()[source]#
seispy.hk.hksta(hpara: HKPara, isplot=False, isdisplay=False)[source]#
seispy.hk.hkstack(seis, t0, dt, p, h, kappa, vp=6.3, weight=(0.7, 0.2, 0.1))[source]#
seispy.hk.plot(stack, allstack, h, kappa, besth, bestk, cvalue, cmap=<matplotlib.colors.ListedColormap object>, title=None, path=None)[source]#
seispy.hk.print_result(besth, bestk, maxhsig, maxksig, print_comment=True)[source]#
seispy.hk.time2idx(times, ti0, dt)[source]#
seispy.hk.tppps(depth, eta_p, eta_s)[source]#
seispy.hk.tps(depth, eta_p, eta_s)[source]#
seispy.hk.tpsps(depth, eta_s)[source]#
seispy.hk.transarray(array, axis=0)[source]#
seispy.hk.vslow(v, rayp)[source]#
class seispy.hkpara.HKPara[source]#

Bases: object

Attributes:
hrange
krange
property hrange#
property krange#
seispy.hkpara.hkpara(cfg_file)[source]#

seispy.slantstack module#

class seispy.slantstack.SlantStack(seis, timeaxis, dis)[source]#

Bases: object

Methods

plot([cmap, xlim, vmin, vmax, figpath, colorbar])

Imaging for slant stacking

stack([ref_dis, rayp_range, tau_range])

Slant stack for receiver function

syn_tps(phase_list[, velmodel, focal_dep])

Calculate the theoretical tau and reference rayp for the given phase list

plot(cmap='jet', xlim=None, vmin=None, vmax=None, figpath=None, colorbar=True)[source]#

Imaging for slant stacking

Parameters:
cmapstr, optional

The colormap to use, by default ‘jet’

xlim_type_, optional

Set the x limits of the current axes., by default None

vmin_type_, optional

Normalize the minium value data to the vmin, by default None

vmax_type_, optional

Normalize the maximum value data to the vmax, by default None

figpath_type_, optional

Output path to the image, by default None

colorbarbool, optional

Wether plot the colorbar, by default True

stack(ref_dis=None, rayp_range=None, tau_range=None)[source]#

Slant stack for receiver function

Parameters:
  • ref_dis (int or float, optional) – reference distance, by default None

  • rayp_range (numpy.ndarray, optional) – range of ray parameter, by default None

  • tau_range (numpy.ndarray, optional) – range of tau, by default None

syn_tps(phase_list, velmodel='iasp91', focal_dep=10)[source]#

Calculate the theoretical tau and reference rayp for the given phase list

Parameters:
  • phase_list (list) – phase list for theoretical arrival time

  • velmodel (str, optional) – velocity model, by default ‘iasp91’

  • focal_dep (int or float, optional) – focal depth, by default 10

seispy.rfani module#

class seispy.rfani.RFAni(sacdatar, tb, te, tlen=3, val=10, rayp=0.06, model='iasp91')[source]#

Bases: object

Methods

baz_stack([val])

Stack RF data with back-azimuth.

joint_ani([weight])

Joint method for crustal anisotropy estimation.

plot_correct([fvd, dt, enf, outpath])

Plot the RF data with back-azimuth and time after P corrected.

plot_polar([cmap, show, outpath])

Polar map of crustal anisotropy.

plot_stack_baz([enf, outpath])

Plot the stack of RF data with back-azimuth.

radial_energy_max()

Calculate the energy of radial component.

rotate_to_fast_slow()

Rotate RF data to fast and slow direction.

search_peak(energy[, opt])

Search the peak of energy.

search_peak_list(energy[, opt])

Search the peak of energy.

xyz2grd(energy)

Interpolate energy to grid.

baz_stack(val=10)[source]#

Stack RF data with back-azimuth.

Parameters:
valint, optional

Interval of back-azimuth, by default 10

joint_ani(weight=[0.5, 0.3, 0.2])[source]#

Joint method for crustal anisotropy estimation.

Parameters:
weightlist, optional

Weight for three energy matrix, by default [0.5, 0.3, 0.2]

plot_correct(fvd=0, dt=0.44, enf=80, outpath=None)[source]#

Plot the RF data with back-azimuth and time after P corrected.

Parameters:
fvdint, optional

Fast velocity direction, by default 0

dtfloat, optional

Time delay for correction, by default 0.44

enfint, optional

Enlarge factor for back-azimuth, by default 80

outpathstr, optional

Output path for saving the figure, by default None

plot_polar(cmap=<matplotlib.colors.ListedColormap object>, show=False, outpath='./')[source]#

Polar map of crustal anisotropy. See Liu and Niu (2012, doi: 10.1111/j.1365-246X.2011.05249.x) in detail.

Parameters:
  • cmap (optional) – Colormap of matplotlib, defaults to ‘rainbow’

  • show (bool, optional) – If show the polar map in the Matplotlib window, defaults to True

  • outpath (str, optional) – Output path to saving the figure. If show the figure in the Matplotlib window, this option will be invalid, defaults to current directory.

plot_stack_baz(enf=60, outpath='./')[source]#

Plot the stack of RF data with back-azimuth.

Parameters:
enfint, optional

Enlarge factor for back-azimuth, by default 60

outpathstr, optional

Output path for saving the figure, by default ‘./’

radial_energy_max()[source]#

Calculate the energy of radial component.

rotate_to_fast_slow()[source]#

Rotate RF data to fast and slow direction.

search_peak(energy, opt='max')[source]#

Search the peak of energy.

Parameters:
energynp.ndarray

Energy matrix

optstr, optional

Option for searching peak, by default ‘max’

search_peak_list(energy, opt='max')[source]#

Search the peak of energy.

Parameters:
energynp.ndarray

Energy matrix

optstr, optional

Option for searching peak, by default ‘max’

xyz2grd(energy)[source]#

Interpolate energy to grid.

seispy.harmonics module#

class seispy.harmonics.Harmonics(rfsta, tmin=-5, tmax=10, bin_stack=True)[source]#

Bases: object

Methods

cut_trace()

Trim RFs from tmin to tmax

harmo_trans()

Harmonic decomposition for extracting anisotropic and isotropic features from the radial and transverse RFs

plot([outpath, enf])

Plot harmonic and unmodeled components

write_constant([out_sac_path])

Write constant component to SAC file.

cut_trace()[source]#

Trim RFs from tmin to tmax

harmo_trans()[source]#

Harmonic decomposition for extracting anisotropic and isotropic features from the radial and transverse RFs

plot(outpath='./', enf=2.0)[source]#

Plot harmonic and unmodeled components

Parameters:
  • outpath (str, optional) – Output path, defaults to ‘./’

  • enf (float, optional) – Amplification factor, defaults to 2.0

write_constant(out_sac_path='./')[source]#

Write constant component to SAC file.

Parameters:

out_sac_path (str, optional) – Output path, defaults to ‘./’