Post-processing package

seispy.core.depmodel module

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

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)

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)

plot model with matplotlib

Parameters:

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

radius_s(rayp, phase='P', sphere=True)

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)

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)

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)

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)

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)

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)

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)

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)

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()

seispy.hk.hksta(hpara: HKPara, isplot=False, isdisplay=False)

seispy.hk.hkstack(seis, t0, dt, p, h, kappa, vp=6.3, weight=(0.7, 0.2, 0.1))

seispy.hk.plot(stack, allstack, h, kappa, besth, bestk, cvalue, cmap=<matplotlib.colors.ListedColormap object>, title=None, path=None)

seispy.hk.print_result(besth, bestk, maxhsig, maxksig, print_comment=True)

seispy.hk.time2idx(times, ti0, dt)

seispy.hk.tppps(depth, eta_p, eta_s)

seispy.hk.tps(depth, eta_p, eta_s)

seispy.hk.tpsps(depth, eta_s)

seispy.hk.transarray(array, axis=0)

seispy.hk.vslow(v, rayp)

class seispy.hkpara.HKPara

Bases: object

Attributes:

hrange

krange

property hrange

property krange

seispy.hkpara.hkpara(cfg_file)

seispy.slantstack module

class seispy.slantstack.SlantStack(seis, timeaxis, dis)

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)

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)

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)

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')

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)

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])

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)

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='./')

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='./')

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()

Calculate the energy of radial component.

rotate_to_fast_slow()

Rotate RF data to fast and slow direction.

search_peak(energy, opt='max')

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')

Search the peak of energy.

Parameters:

energynp.ndarray

Energy matrix

optstr, optional

Option for searching peak, by default ‘max’

xyz2grd(energy)

Interpolate energy to grid.

seispy.harmonics module

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

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()

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=2.0)

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='./')

Write constant component to SAC file.

Parameters:

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