CCP package

seispy.rfcorrect module

class seispy.rfcorrect.RFStation(data_path, only_r=False, prime_comp='R')

Bases: object

Attributes

stel

Methods

harmonic([tb, te])	Harmonic decomposition for extracting anisotropic and isotropic features from the radial and transverse RFs
jointani(tb, te[, tlen, stack_baz_val, ...])	Eastimate crustal anisotropy with a joint method.
moveoutcorrect([ref_rayp, dep_range, ...])	Moveout correction with specified reference ray-parameter and depth
normalize([method])	Normalize amplitude of each RFs. :param method: Method of normalization with single and average avaliable. - single for normalization with max amplitude of current RF. - average for normalization with average amplitude of current station. :type method: str, optional.
psrf2depth([dep_range])	Time-to-depth conversion with specified depth series.
psrf_1D_raytracing([dep_range])	1D back ray tracing to obtained Ps conversion points at discret depths
resample(dt)	Resample RFs with specified dt
sort([key])	Sort RFs by keys in given event, evla, evlo, evdp, dis, bazi, rayp, mag, f0

psrf_3D_moveoutcorrect
psrf_3D_raytracing
psrf_3D_timecorrect
read_sample
slantstack

harmonic(tb=- 5, te=10)

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

Parameters

• tb (float, optional) – Start time relative to P, defaults to -5

• te (float, optional) – End time relative to P, defaults to 10

Returns

harmonic_trans: numpy.ndarray, float

Harmonic components with shape of (5, nsamp), nsamp = (te-tb)/RFStation.sampling

unmodel_trans: numpy.ndarray, float

Unmodel components with shape same as harmonic_trans.

jointani(tb, te, tlen=3.0, stack_baz_val=10, rayp=0.06, velmodel='iasp91', weight=[0.4, 0.4, 0.2])

Eastimate crustal anisotropy with a joint method. See Liu and Niu (2012, doi: 10.1111/j.1365-246X.2011.05249.x) in detail.

Parameters

• tb (float) – Time before Pms for search Ps peak

• te (float) – Time after Pms for search Ps peak

• tlen (float, optional) – Half time length for cut out Ps phase, defaults to 3.0

• stack_baz_val (float, optional) – The interval for stacking binned by back-azimuth, defaults to 10

• rayp (float, optional) – Reference ray-parameter for moveout correction, defaults to 0.06

• velmodel (str, optional) – velocity model for moveout correction. ‘iasp91’, ‘prem’ and ‘ak135’ is valid for internal model. Specify path to velocity model for the customized model. The format is the same as in Taup, but the depth should be monotonically increasing, defaults to ‘iasp91’

• weight (list, optional) – Weight for three different method, defaults to [0.4, 0.4, 0.2]

Returns

Dominant fast velocity direction and time delay

Return type

list, list

moveoutcorrect(ref_rayp=0.06, dep_range=array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149]), velmod='iasp91', replace=False)

Moveout correction with specified reference ray-parameter and depth

Parameters

• ref_rayp (float, optional) – reference ray-parameter in s/km, defaults to 0.06

• dep_range (np.ndarray(), optional) – Depth range used for extracting velocity in velocity model, defaults to np.arange(0, 150)

• velmod (str, optional) – Velocity model for moveout correction. ‘iasp91’, ‘prem’ and ‘ak135’ is valid for internal model. Specify path to velocity model for the customized model. The format is the same as in Taup, but the depth should be monotonically increasing, defaults to ‘iasp91’

• replace (bool, optional) – whether replace original data, False to return new array, defaults to False

Returns

rf_corr:meth:np.ndarray

Corrected RFs with component of RFStation.comp()

t_corr:meth:np.ndarray or None

Corrected RFs in transverse component. If only_r is True, this variable is None

normalize(method='single')

Normalize amplitude of each RFs. :param method: Method of normalization with single and average avaliable.

• single for normalization with max amplitude of current RF.

• average for normalization with average amplitude of current station.

psrf2depth(dep_range=array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149]), **kwargs)

Time-to-depth conversion with specified depth series.

Parameters

• dep_range (np.ndarray() , optional) – Discret conversion depth, defaults to np.arange(0, 150)

• velmod (str, optional) – Velocity model for time-to-depth conversion. ‘iasp91’, ‘prem’ and ‘ak135’ is valid for internal model. Specify path to velocity model for the customized model. The format is the same as in Taup, but the depth should be monotonically increasing, defaults to ‘iasp91’

• srayp (numpy.lib.npyio.NpzFile, optional) – Ray-parameter lib for Ps phases, If set up to None the rayp of direct is used, defaults to None

Returns

2D array of RFs in depth

Return type

np.ndarray()

psrf_1D_raytracing(dep_range=array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149]), **kwargs)

1D back ray tracing to obtained Ps conversion points at discret depths

Parameters

• dep_range (numpy.ndarray, optional) – Discret conversion depth, defaults to np.arange(0, 150)

• velmod (str, optional) – Velocity model for time-to-depth conversion. ‘iasp91’, ‘prem’ and ‘ak135’ is valid for internal model. Specify path to velocity model for the customized model. The format is the same as in Taup, but the depth should be monotonically increasing, defaults to ‘iasp91’

• srayp (numpy.lib.npyio.NpzFile, optional) – Ray-parameter lib for Ps phases, If set up to None the rayp of direct is used, defaults to None

Return pplat_s

Latitude of conversion points

Return pplon_s

Longitude of conversion points

Return tps

Time difference of Ps at each depth

Return type

list

psrf_3D_moveoutcorrect(mod3dpath, **kwargs)

psrf_3D_raytracing(mod3dpath, dep_range=array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149]), srayp=None)

psrf_3D_timecorrect(mod3dpath, dep_range=array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149]), normalize='single', **kwargs)

read_sample(data_path)

resample(dt)

Resample RFs with specified dt

Parameters

dtfloat

Target sampling interval in sec

slantstack(ref_dis=None, rayp_range=None, tau_range=None)

sort(key='bazi')

Sort RFs by keys in given event, evla, evlo, evdp, dis, bazi, rayp, mag, f0

Parameters

key (str, optional) – key to sort, defaults to bazi

property stel

class seispy.rfcorrect.SACStation(data_path, only_r=False)

Bases: seispy.rfcorrect.RFStation

Attributes

stel

Methods

harmonic([tb, te])	Harmonic decomposition for extracting anisotropic and isotropic features from the radial and transverse RFs
jointani(tb, te[, tlen, stack_baz_val, ...])	Eastimate crustal anisotropy with a joint method.
moveoutcorrect([ref_rayp, dep_range, ...])	Moveout correction with specified reference ray-parameter and depth
normalize([method])	Normalize amplitude of each RFs. :param method: Method of normalization with single and average avaliable. - single for normalization with max amplitude of current RF. - average for normalization with average amplitude of current station. :type method: str, optional.
psrf2depth([dep_range])	Time-to-depth conversion with specified depth series.
psrf_1D_raytracing([dep_range])	1D back ray tracing to obtained Ps conversion points at discret depths
resample(dt)	Resample RFs with specified dt
sort([key])	Sort RFs by keys in given event, evla, evlo, evdp, dis, bazi, rayp, mag, f0

psrf_3D_moveoutcorrect
psrf_3D_raytracing
psrf_3D_timecorrect
read_sample
slantstack

seispy.rfcorrect.interp_depth_model(model, lat, lon, new_dep)

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

Parameters

mod3dnp.lib.npyio.NpzFile()

3D velocity loaded from a .npz file

latfloat

Latitude of position in 3D velocity model

lonfloat

Longitude of position in 3D velocity model

new_depnp.ndarray()

1D array of depths in km

Returns

Vpnp.ndarray()

Vp in new_dep

Vsnp.ndarray()

Vs in new_dep

seispy.rfcorrect.moveoutcorrect_ref(stadatar, raypref, YAxisRange, chan='r', velmod='iasp91', sphere=True, phase=1)

Moveout correction refer to a specified ray-parameter

Parameters

• stadatar – data class of RFStation

• raypref – referred ray parameter in rad

• YAxisRange – Depth range in nd.array type

• velmod – Path to velocity model

• chan – channel name for correction, ‘r’, ‘t’…

Returns

Newdatar, EndIndex

seispy.rfcorrect.psrf2depth(stadatar, YAxisRange, velmod='iasp91', srayp=None, normalize='single', sphere=True, phase=1)

Time-to-depth conversion with S-wave backprojection.

Parameters

• stadatar (RFStation()) – Data class of RFStation

• YAxisRange (numpy.ndarray) – Depth range for converison

• velmod (str, optional) – Velocity for conversion, whcih can be a path to velocity file, defaults to ‘iasp91’

• srayp (str or seispy.psrayp.PsRayp(), optional) – ray-parameter library of conversion phases. See seispy.psrayp() in detail, defaults to None

• normalize (str, optional) – method of normalization, defaults to ‘single’. Please refer to RFStation.normalize()

• sphere (bool, optional) – Wether do earth-flattening transformation, defaults to True

Returns

ps_rfdepth: 2-D numpy.ndarray, float

RFs in depth with shape of (stadatar.ev_num, YAxisRange.size), stadatar.ev_num is the number of RFs in current station. YAxisRange.size is the size of depth axis.

endindex: numpy.ndarray, int

End index of each RF in depth

x_s: 2-D numpy.ndarray, float

Horizontal distance between station and S-wave conversion points with shape of (stadatar.ev_num, YAxisRange.size)

x_p: 2-D numpy.ndarray, float

Horizontal distance between station and P-wave conversion points with shape of (stadatar.ev_num, YAxisRange.size)

seispy.rfcorrect.psrf_1D_raytracing(stadatar, YAxisRange, velmod='iasp91', srayp=None, sphere=True, phase=1)

seispy.rfcorrect.psrf_3D_migration(pplat_s, pplon_s, pplat_p, pplon_p, raylength_s, raylength_p, Tpds, dep_range, mod3d)

3D time difference correction with specified ray path and 3D velocity model.

The input parameters can be generated with psrf_1D_raytracing().

Parameters

pplat_snp.ndarray()

2D array of latitude of S-wave in dep_range, (RFStation.ev_num(), dep_range.size)

pplon_snp.ndarray()

2D array of longitude of S-wave in dep_range, (RFStation.ev_num(), dep_range.size)

pplat_pnp.ndarray()

2D array of latitude of P-wave in dep_range, (RFStation.ev_num(), dep_range.size)

pplon_pnp.ndarray()

2D array of longitude of P-wave in dep_range, (RFStation.ev_num(), dep_range.size)

raylength_snp.ndarray()

2D array of ray path length of S-wave in dep_range, (RFStation.ev_num(), dep_range.size)

raylength_pnp.ndarray()

2D array of ray path length of P-wave in dep_range, (RFStation.ev_num(), dep_range.size)

Tpdsnp.ndarray()

1D array of time difference in dep_range (dep_range.size)

dep_rangenp.ndarray()

1D array of depths in km, (dep_range.size)

mod3dnp.lib.npyio.NpzFile()

3D velocity loaded from a .npz file

Returns

np.ndarray()

Corrected time difference in dep_range

seispy.rfcorrect.psrf_3D_raytracing(stadatar, YAxisRange, mod3d, srayp=None, elevation=0, sphere=True)

Back ray trace the S wavs with a assumed ray parameter of P.

Parameters

• stadatar (object RFStation) – The data class including PRFs and more parameters

• YAxisRange (numpy.ndarray) – The depth array with the same intervals

• mod3d ('Mod3DPerturbation' object) – The 3D velocity model with fields of dep, lat, lon, vp and vs.

• elevation (float) – Elevation of this station relative to sea level

Returns

pplat_s, pplon_s, pplat_p, pplon_p, tps

Type

numpy.ndarray * 5

seispy.rfcorrect.time2depth(stadatar, dep_range, Tpds, normalize='single')

Interpolate RF amplitude with specified time difference and depth range

Parameters

stadatarRFStation()

Data class of RFStation()

dep_rangenp.ndarray()

1D array of depths in km, (dep_range.size)

Tpdsnp.ndarray()

1D array of time difference in dep_range (dep_range.size)

normalizestr, optional

Normlization option, 'sinlge' and 'average' are available , by default ‘single’See RFStation.normalize() in detail.

Returns

PS_RFdepth:meth:np.ndarray

2D array of RFs in depth, (RFStation.ev_num(), dep_range.size)

end_index:meth:np.ndarray

1D array of the last digit of effective depth (RFStation.ev_num())

seispy.rfcorrect.xps_tps_map(dep_mod, srayp, prayp, is_raylen=False, sphere=True, phase=1)

Calculate horizontal distance and time difference at depths

Parameters

• dep_mod (seispy.util.DepModel()) – 1D velocity model class

• srayp (float) – conversion phase ray-parameters

• prayp (float) – S-wave ray-parameters

• is_raylen (bool, optional) – Wether calculate ray length at depths, defaults to False

• sphere (bool, optional) – Wether do earth-flattening transformation, defaults to True, defaults to True

• phase (int, optional) – Phases to calculate 1 for Ps, 2 for PpPs, 3 for PsPs+PpSs, defaults to 1

Returns

If is_raylen = False

tps: 2-D numpy.ndarray, float

RFs in depth with shape of (stadatar.ev_num, YAxisRange.size), stadatar.ev_num is the number of RFs in current station. YAxisRange.size is the size of depth axis.

x_s: 2-D numpy.ndarray, float

Horizontal distance between station and S-wave conversion points with shape of (stadatar.ev_num, YAxisRange.size)

x_p: 2-D numpy.ndarray, float

Horizontal distance between station and P-wave conversion points with shape of (stadatar.ev_num, YAxisRange.size)

otherwise, two more variables will be returned

raylength_s: 2-D numpy.ndarray, float

raylength_p: 2-D numpy.ndarray, float

seispy.rf2depth_makedata module

class seispy.rf2depth_makedata.Station(sta_lst)

Bases: object

seispy.rf2depth_makedata.init_mat(sta_num)

seispy.rf2depth_makedata.makedata(cpara, velmod3d=None, modfolder1d=None, log=<seispy.setuplog.setuplog object>)

seispy.rf2depth_makedata.makedata3d(cpara, velmod3d, log=<seispy.setuplog.setuplog object>, raytracing3d=True)

seispy.rf2depth_makedata.rf2depth()

seispy.ccpprofile module

class seispy.ccpprofile.CCPProfile(cfg_file=None, log=None)

Bases: object

Methods

initial_profile()	Initialize bins of profile
read_rfdep()	Read RFdepth file
save_stack_data([format])	If format is 'npz', saving stacked data and parameters to local as a npz file.
stack()	Stack RFs in bins

load_para

initial_profile()

Initialize bins of profile

load_para(cfg_file)

read_rfdep()

Read RFdepth file

Raises

FileNotFoundError – Not Found RFdepth file

save_stack_data(format='npz')

If format is ‘npz’, saving stacked data and parameters to local as a npz file. To load the file, please use data = np.load(fname, allow_pickle=True). data[‘cpara’] is the parameters when CCP stacking. data[‘stack_data’] is the result of stacked data.

If format is ‘dat’ the stacked data will be save into a txt file with 8 columns, including bin_lat, bin_lon, profile_dis, depth, amp, ci_low, ci_high and count. where bin_lat and bin_lon represent the position of each bin; profile_dis represents the distance in km between each bin and the start point of the profile; depth represents depth of each bin; amp means the stacked amplitude; ci_low and ci_high mean confidence interval with bootstrap method; count represents stacking number of each bin.

Parameters

format (str) – Format for stacked data

stack()

Stack RFs in bins

seispy.ccpprofile.bin_shape(cpara)

seispy.ccpprofile.fix_filename(filename, typ='dat')

seispy.ccpprofile.init_profile(lat1, lon1, lat2, lon2, val)

Initial bins along a profile with given position of two points.

Parameters

• lat1 (float) – The latitude of the start point

• lon1 (float) – The lontitude of the start point

• lat2 (float) – The latitude of the end point

• lon2 (float) – The lontitude of the end point

• val (float) – The interval between two points in km

Returns

The location of bins (bin_loca), and length between each bin and the start point (profile_range)

The bin_loca is positions of bins with a numpy.array with two column. The profile_range is distance between bin center and the start point with an 1D numpy.array. :rtype: (numpy.array, numpy.array)

seispy.ccpprofile.line_proj(lat1, lon1, lat2, lon2)

seispy.cc3d module

class seispy.ccp3d.CCP3D(cfg_file=None, log=None)

Bases: object

Methods

save_stack_data(fname)	Save stacked data and parameters to local as a npz file.
stack()	Search conversion points falling within a bin and stack them with bootstrap method.

get_depth_err
initial_grid
load_para
read_rfdep
read_stack_data
save_good_410_660
search_good_410_660

get_depth_err(type='std')

initial_grid()

load_para(cfg_file)

read_rfdep()

classmethod read_stack_data(stack_data_path, cfg_file=None, good_depth_path=None, ismtz=False)

save_good_410_660(fname)

save_stack_data(fname)

Save stacked data and parameters to local as a npz file. To load the file, please use data = np.load(fname, allow_pickle=True). data[‘cpara’] is the parameters when CCP stacking. data[‘stack_data’] is the result of stacked data.

Parameters

fname (str) – file name of stacked data

search_good_410_660(peak_410_min=380, peak_410_max=440, peak_660_min=630, peak_660_max=690)

stack()

Search conversion points falling within a bin and stack them with bootstrap method.

seispy.ccp3d.bin_shape(cpara)

seispy.ccp3d.boot_bin_stack(data_bin, n_samples=3000)

seispy.ccp3d.gen_center_bin(center_lat, center_lon, len_lat, len_lon, val)

Create spaced grid point with coordinates of the center point in the area in spherical coordinates.

Parameters

• center_lat (float) – Latitude of the center point.

• center_lon (float) – Longitude of the center point.

• len_lat (float) – Half length in degree along latitude axis.

• len_lon (float) – Half length in degree along longitude axis.

• val (float) – Interval in degree between adjacent grid point.

Returns

Coordinates of Grid points.

Return type

2-D ndarray of floats with shape (n, 2), where n is the number of grid points.