CCP package¶
seispy.rfcorrect module¶
- class seispy.rfcorrect.RFStation(data_path, only_r=False, prime_comp='R')[source]¶
Bases:
object
- Attributes:
- stel
Methods
bin_stack
([key, lim, val])Stack RFs by bins of
key
with interval ofval
harmonic
([tb, te, is_stack])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.
plotr
([outformat])Plot radial RFs
plotrt
([outformat])Plot radial and transverse RFs
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
psrf_3D_moveoutcorrect
(mod3dpath, **kwargs)3D moveout correction with 3D velocity model
psrf_3D_raytracing
(mod3dpath[, dep_range, srayp])3D back ray tracing to obtained Ps conversion points at discret depths
psrf_3D_timecorrect
(mod3dpath[, dep_range, ...])3D time-to-depth conversion with 3D velocity model
read_stream
(stream, rayp, baz[, prime_comp, ...])Create RFStation instance from
obspy.Stream
resample
(dt)Resample RFs with specified dt
slantstack
([ref_dis, rayp_range, tau_range])Slant stack for receiver function
sort
([key])Sort RFs by keys in given
event
,evla
,evlo
,evdp
,dis
,bazi
,rayp
,mag
,f0
- bin_stack(key='bazi', lim=[0, 360], val=10)[source]¶
Stack RFs by bins of
key
with interval ofval
- Parameters:
key (str, optional) – Key to stack, valid in [‘bazi’, ‘rayp’ (in s/km)], defaults to ‘bazi’
val (int, optional) – Interval of bins, defaults to 10 degree for bazi
- Returns:
Stacked RFs
- Return type:
dict
with keys ofdata_prime
,datat
andcount
- harmonic(tb=-5, te=10, is_stack=True)[source]¶
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 -5te (
float
, optional) – End time relative to P, defaults to 10is_stack (bool, optional) – Wether stack the result, defaults to True
- Returns:
Harmonic components and unmodel components
- Return type:
numpy.ndarray
,numpy.ndarray
- jointani(tb, te, tlen=3.0, stack_baz_val=10, rayp=0.06, velmodel='iasp91', weight=[0.4, 0.4, 0.2])[source]¶
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, **kwargs)[source]¶
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:
np.ndarray
Corrected RFs with component of
RFStation.comp
- t_corr:
np.ndarray
orNone
Corrected RFs in transverse component. If
only_r
isTrue
, this variable isNone
- rf_corr:
- normalize(method='single')[source]¶
Normalize amplitude of each RFs.
- Parameters:
method (str, optional) – Method of normalization with
single
andaverage
avaliable. -single
for normalization with max amplitude of current RF. -average
for normalization with average amplitude of current station.
- plotr(outformat=None, **kwargs)[source]¶
Plot radial RFs
- Parameters:
outformat (str, optional) – Output format for plot, defaults to None
- plotrt(outformat=None, **kwargs)[source]¶
Plot radial and transverse RFs
- Parameters:
outformat (str, optional) – Output format for plot, defaults to None
- 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)[source]¶
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)[source]¶
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)[source]¶
3D moveout correction with 3D velocity model
- Parameters:
mod3dpath (str) – Path to 3D velocity model
dep_range (numpy.ndarray, optional) – Discret conversion depth, defaults to np.arange(0, 150)
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:
numpy.ndarray
- 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)[source]¶
3D back ray tracing to obtained Ps conversion points at discret depths
- Parameters:
mod3dpath (str) – Path to 3D velocity model
dep_range (numpy.ndarray, optional) – Discret conversion depth, defaults to np.arange(0, 150)
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_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)[source]¶
3D time-to-depth conversion with 3D velocity model
- Parameters:
mod3dpath (str) – Path to 3D velocity model
dep_range (numpy.ndarray, optional) – Discret conversion depth, defaults to np.arange(0, 150)
normalize (str, optional) – Normalization method, defaults to ‘single’, see RFStation.normalize for detail
- Returns:
2D array of RFs in depth
- Return type:
numpy.ndarray
- classmethod read_stream(stream, rayp, baz, prime_comp='R', stream_t=None)[source]¶
Create RFStation instance from
obspy.Stream
- Parameters:
stream (
obspy.Stream()
) – Stream of RFsrayp (float or
np.ndarray()
) – Ray-parameter of RFsbaz (float or
np.ndarray()
) – Back-azimuth of RFsprime_comp (str, optional) – Prime component in RF filename.
R
orQ
for PRF andL
orZ
for SRF, defaults to ‘R’stream_t (
obspy.Stream()
, optional) – Stream of transverse RFs, defaults to None
- Returns:
RFStation instance
- Return type:
- slantstack(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
- sort(key='bazi')[source]¶
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)[source]¶
Bases:
RFStation
- Attributes:
- stel
Methods
bin_stack
([key, lim, val])Stack RFs by bins of
key
with interval ofval
harmonic
([tb, te, is_stack])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.
plotr
([outformat])Plot radial RFs
plotrt
([outformat])Plot radial and transverse RFs
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
psrf_3D_moveoutcorrect
(mod3dpath, **kwargs)3D moveout correction with 3D velocity model
psrf_3D_raytracing
(mod3dpath[, dep_range, srayp])3D back ray tracing to obtained Ps conversion points at discret depths
psrf_3D_timecorrect
(mod3dpath[, dep_range, ...])3D time-to-depth conversion with 3D velocity model
read_stream
(stream, rayp, baz[, prime_comp, ...])Create RFStation instance from
obspy.Stream
resample
(dt)Resample RFs with specified dt
slantstack
([ref_dis, rayp_range, tau_range])Slant stack for receiver function
sort
([key])Sort RFs by keys in given
event
,evla
,evlo
,evdp
,dis
,bazi
,rayp
,mag
,f0
- seispy.rfcorrect.moveoutcorrect_ref(stadatar, raypref, YAxisRange, chan='', velmod='iasp91', sphere=True, phase=1)[source]¶
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)[source]¶
Time-to-depth conversion with S-wave backprojection.
- Parameters:
stadatar (
RFStation()
) – Data class of RFStationYAxisRange (
numpy.ndarray
) – Depth range for conversionvelmod (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. Seeseispy.psrayp()
in detail, defaults to Nonenormalize (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)[source]¶
- seispy.rfcorrect.psrf_3D_migration(pplat_s, pplon_s, pplat_p, pplon_p, raylength_s, raylength_p, Tpds, dep_range, mod3d)[source]¶
- 3D time difference correction with specified ray path and 3D velocity model.
The input parameters can be generated with
psrf_1D_raytracing()
.
- Parameters:
- pplat_s
np.ndarray()
2D array of latitude of S-wave in
dep_range
, (RFStation.ev_num()
,dep_range.size
)- pplon_s
np.ndarray()
2D array of longitude of S-wave in
dep_range
, (RFStation.ev_num()
,dep_range.size
)- pplat_p
np.ndarray()
2D array of latitude of P-wave in
dep_range
, (RFStation.ev_num()
,dep_range.size
)- pplon_p
np.ndarray()
2D array of longitude of P-wave in
dep_range
, (RFStation.ev_num()
,dep_range.size
)- raylength_s
np.ndarray()
2D array of ray path length of S-wave in
dep_range
, (RFStation.ev_num()
,dep_range.size
)- raylength_p
np.ndarray()
2D array of ray path length of P-wave in
dep_range
, (RFStation.ev_num()
,dep_range.size
)- Tpds
np.ndarray()
1D array of time difference in
dep_range
(dep_range.size
)- dep_range
np.ndarray()
1D array of depths in km, (
dep_range.size
)- mod3d
np.lib.npyio.NpzFile()
3D velocity loaded from a
.npz
file
- pplat_s
- Returns:
np.ndarray()
Corrected time difference in dep_range
- seispy.rfcorrect.psrf_3D_raytracing(stadatar, YAxisRange, mod3d, srayp=None, elevation=0, sphere=True)[source]¶
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
andvs
.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')[source]¶
Interpolate RF amplitude with specified time difference and depth range
- Parameters:
- stadatar
RFStation()
Data class of
RFStation()
- dep_range
np.ndarray()
1D array of depths in km, (
dep_range.size
)- Tpds
np.ndarray()
1D array of time difference in
dep_range
(dep_range.size
)- normalizestr, optional
Normlization option,
'sinlge'
and'average'
are available , by default ‘single’SeeRFStation.normalize()
in detail.
- stadatar
- Returns:
- PS_RFdepth:
np.ndarray()
2D array of RFs in depth, (
RFStation.ev_num()
,dep_range.size
)- end_index:
np.ndarray()
1D array of the last digit of effective depth (
RFStation.ev_num()
)
- PS_RFdepth:
- seispy.rfcorrect.xps_tps_map(dep_mod: DepModel, srayp, prayp, is_raylen=False, sphere=True, phase=1)[source]¶
Calculate horizontal distance and time difference at depths
- Parameters:
dep_mod (
seispy.util.DepModel()
) – 1D velocity model classsrayp (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 forPpPs
, 3 forPsPs+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
- If
seispy.rf2depth_makedata module¶
class RF2depth : process cal RF2depth class sta_part, sta_full, _RFInd
- class seispy.rf2depth_makedata.RFDepth(cpara: ~seispy.ccppara.CCPPara, log: ~logging.Logger = <Logger RF2depth (INFO)>, raytracing3d=False, velmod3d=None, modfolder1d=None)[source]¶
Bases:
object
Convert receiver function to depth axis
Methods
makedata
([psphase])Convert receiver function to depth axis
- seispy.rf2depth_makedata.rf2depth()[source]¶
CLI for Convert receiver function to depth axis There’s 4 branch provided to do RF 2 depth conversion
only -d :do moveout correction
only -r : do raytracing but no moveout correction
-d and -r : do moveout correction and raytracing
-m : use {staname}.vel file for RF2depth conversion
- class seispy.rf2depth_makedata.sta_full(station, stla, stlo, stel)¶
Bases:
tuple
Methods
count
(value, /)Return number of occurrences of value.
index
(value[, start, stop])Return first index of value.
- station¶
Alias for field number 0
- stel¶
Alias for field number 3
- stla¶
Alias for field number 1
- stlo¶
Alias for field number 2
seispy.ccpprofile module¶
- class seispy.ccpprofile.CCPProfile(cfg_file=None, log: ~seispy.setuplog.SetupLog = <seispy.setuplog.SetupLog object>)[source]¶
Bases:
object
Methods
Initialize bins of profile
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
- save_stack_data(format='npz')[source]¶
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, includingbin_lat
,bin_lon
,profile_dis
,depth
,amp
,ci_low
,ci_high
andcount
.bin_lat
andbin_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
andci_high
mean confidence interval with bootstrap method;count
represents stacking number of each bin.
- Parameters:
format (str) – Format for stacked data
- seispy.ccpprofile.create_center_bin_profile(stations, val=5, method='linear')[source]¶
Create bins along a profile with given stations
- Parameters:
stations (
seispy.rf2depth_makedata.Station
) – Stations along a profileval (float) – The interval between two points in km
method (str) – Method for interpolation
- Returns:
The location of bins (bin_loca), and length between each bin and the start point (profile_range)
- Return type:
(numpy.array, numpy.array)
- seispy.ccpprofile.init_profile(lat1, lon1, lat2, lon2, val)[source]¶
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.cc3d module¶
- class seispy.ccp3d.CCP3D(cfg_file=None, log=None)[source]¶
Bases:
object
Class for 3-D CCP stacking, Usually used to study mantle transition zone structure.
- Parameters:
cfg_file (str, optional) – Path to configure file. If not defined a instance of CCP3D.cpara will be initialed, defaults to None
log (seispy.sutuplog.logger , optional) – A logger instance. If not defined, seispy.sutuplog.logger will be initialed, defaults to None
Methods
Initial grid points and search stations within a distance.
load_para
(cfg_file)Load parameters from configure file.
Read RFdepth data from npz file.
read_stack_data
(stack_data_path[, cfg_file, ...])Read stacked data from local.
save_good_410_660
(fname)Save good 410 and 660 km discontinuities to local.
save_stack_data
(fname)Save stacked data and parameters to local as a npz file.
search_good_410_660
([peak_410_min, ...])Search good 410 and 660 km discontinuities from stacked data.
stack
()Search conversion points falling within a bin and stack them with bootstrap method.
get_depth_err
- load_para(cfg_file)[source]¶
Load parameters from configure file.
- Parameters:
cfg_file (str) – Path to configure file.
- classmethod read_stack_data(stack_data_path, cfg_file=None, good_depth_path=None, ismtz=False)[source]¶
Read stacked data from local.
- Parameters:
stack_data_path (str) – Path to stacked data.
cfg_file (str, optional) – Path to configure file, defaults to None
good_depth_path (str, optional) – Path to good depth file, defaults to None
ismtz (bool, optional) – Whether the good depth file is in mtz format, defaults to False
- Returns:
A instance of CCP3D.
- Return type:
- save_good_410_660(fname)[source]¶
Save good 410 and 660 km discontinuities to local.
- Parameters:
fname (str) – file name of good 410 and 660 km discontinuities
- save_stack_data(fname)[source]¶
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)[source]¶
Search good 410 and 660 km discontinuities from stacked data.
- Parameters:
peak_410_min (float, optional) – Minimum depth of 410 km discontinuity, defaults to 380
peak_410_max (float, optional) – Maximum depth of 410 km discontinuity, defaults to 440
peak_660_min (float, optional) – Minimum depth of 660 km discontinuity, defaults to 630
peak_660_max (float, optional) – Maximum depth of 660 km discontinuity, defaults to 690
- seispy.ccp3d.bin_shape(cpara)[source]¶
Compute the radius of bins in degree.
- Parameters:
cpara (CCPPara) – Parameters of CCP stacking.
- Returns:
Radius of bins in degree.
- Return type:
1-D ndarray of floats with shape (n, ), where n is the number of bins.
- seispy.ccp3d.boot_bin_stack(data_bin, n_samples=3000)[source]¶
Stack data with bootstrap method.
- Parameters:
data_bin (1-D ndarray of floats) – Data falling within a bin.
n_samples (int, optional) – Number of bootstrap samples, defaults to 3000
- Returns:
Mean, confidence interval and number of data falling within a bin.
- Return type:
tuple of floats and int
- seispy.ccp3d.gen_center_bin(center_lat, center_lon, len_lat, len_lon, val)[source]¶
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.