-
Notifications
You must be signed in to change notification settings - Fork 9
/
user_functions.f90
140 lines (118 loc) · 4.02 KB
/
user_functions.f90
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
! -------------------------------------------------------------
! edit here to change T0 and T1 on some condition
! Note, this function is called AFTER the seismogram has been
! read but before it is filtered.
! -------------------------------------------------------------
subroutine modify_T0_T1_on_condition
use seismo_variables
implicit none
! do nothing
! adjust fstart and fend accordingly
FSTART=1./WIN_MAX_PERIOD
FEND=1./WIN_MIN_PERIOD
end subroutine
! -------------------------------------------------------------
! edit here to change the time dependent properties of the
! selection criteria
! Note, this function is called AFTER the seismogram has been
! read and filtered.
! -------------------------------------------------------------
subroutine set_up_criteria_arrays
use seismo_variables
implicit none
integer :: i
double precision :: time
double precision :: R_vel, R_time
double precision :: Q_vel, Q_time
! -----------------------------------------------------------------
! This is the basic version of the subroutine - no variation with time
! -----------------------------------------------------------------
do i = 1, npts
time=b+(i-1)*dt
DLNA_LIMIT(i)=DLNA_BASE
CC_LIMIT(i)=CC_BASE
TSHIFT_LIMIT(i)=TSHIFT_BASE
STALTA_W_LEVEL(i)=STALTA_BASE
S2N_LIMIT(i)=WINDOW_S2N_BASE
enddo
! these values will be used for signal2noise calculations
! if DATA_QUALITY=.true.
if (DATA_QUALITY) then
noise_start=b
noise_end=max(ph_times(1)-WIN_MIN_PERIOD*1.5,b+dt)
signal_start=noise_end
signal_end=b+(npts-1)*dt
endif
! -----------------------------------------------------------------
! Start of user-dependent portion
! This is where you reset the signal_end and noise_end values
! if you want to modify them from the defaults.
! This is also where you modulate the time dependence of the selection
! criteria. You have access to the following parameters from the
! seismogram itself:
!
! dt, b, kstnm, knetwk, kcmpnm
! evla, evlo, stla, stlo, evdp, azimuth, backazimuth, dist_deg, dist_km
! num_phases, ph_names, ph_times
!
! Example of modulation:
!-----------------------
! To increase s2n limit after arrival of R1 try
!
! R_vel=3.2
! R_time=dist_km/R_vel
! do i = 1, npts
! time=b+(i-1)*dt
! if (time.gt.R_time) then
! S2N_LIMIT(i)=2*WINDOW_S2N_BASE
! endif
! enddo
!
! --------------------------------
! Set approximate end of rayleigh wave arrival
R_vel=3.2
R_time=dist_km/R_vel
! --------------------------------
! Set approximate start of love wave arrival
Q_vel=4.2
Q_time=dist_km/Q_vel
! reset the signal_end time to be the end of the Rayleigh waves
if (DATA_QUALITY) then
signal_end=R_time
endif
! --------------------------------
! modulate criteria in time
do i = 1, npts
time=b+(i-1)*dt
! --------------------------------
! if we are beyond the Rayleigh wave, then make the all criteria stronger
! ratio criterion stronger
if (time.gt.R_time) then
S2N_LIMIT(i)=10*WINDOW_S2N_BASE ! only pick big signals
CC_LIMIT(i)= 0.95 ! only pick very similar signals
TSHIFT_LIMIT(i)=TSHIFT_BASE/3.0 ! only pick small timeshifts
DLNA_LIMIT(i)=DLNA_BASE/3.0 ! only pick small amplitude anomalies
STALTA_W_LEVEL(i)=STALTA_BASE*2.0 ! pick only distinctive arrivals
endif
! --------------------------------
! if we are in the surface wave times, then make the cross-correlation
! criterion less severe
if (time.gt.Q_time .and. time.lt.R_time) then
CC_LIMIT(i)=0.9*CC_LIMIT(i)
endif
! --------------------------------
! modulate criteria according to event depth
!
! if an intermediate depth event
if (evdp.ge.70 .and. evdp.lt.300) then
TSHIFT_LIMIT(i)=TSHIFT_BASE*1.4
! if a deep event
elseif (evdp.ge.300) then
TSHIFT_LIMIT(i)=TSHIFT_BASE*1.7
endif
enddo
!
! End of user-dependent portion
! -----------------------------------------------------------------
end subroutine
! -------------------------------------------------------------