1# http://pyrocko.org - GPLv3 

2# 

3# The Pyrocko Developers, 21st Century 

4# ---|P------/S----------~Lg---------- 

5 

6from pyrocko.gui.snuffling import Snuffling, Param, Switch, Choice 

7from pyrocko.gui.marker import PhaseMarker 

8from pyrocko import gf 

9from pyrocko import cake 

10import numpy as num 

11 

12try: 

13 newstr = unicode 

14except NameError: 

15 newstr = str 

16 

17 

18class CakePhase(Snuffling): 

19 ''' 

20 <html> 

21 <head> 

22 <style type="text/css"> 

23 body { margin-left:10px }; 

24 </style> 

25 </head> 

26 <body> 

27 <h1 align="center">Theoretical Phase Arrivals</h1> 

28 <p> 

29 This snuffling uses pyrocko's 

30 <a href="http://emolch.github.io/pyrocko/v0.3/cake_doc.html">Cake</a> 

31 module to calculate seismic rays for layered earth models. </p> 

32 <p> 

33 <b>Parameters:</b><br /> 

34 <b>&middot; Global shift</b> - Add time onset to phases. <br /> 

35 <b>&middot; Add Model</b> - Add a model to drop down menu. <br /> 

36 <b>&middot; Add Phase</b> - Add a phase definition. 

37 (GUI reset required)<br /> 

38 </p> 

39 <p> 

40 Instructions and information on Cake's syntax of seismic rays can be 

41 found in the 

42 <a href="http://emolch.github.io/pyrocko/ 

43 v0.3/cake_doc.html#cmdoption-cake--phase">Cake documentation</a>. 

44 </p> 

45 </body> 

46 </html> 

47 ''' 

48 

49 def setup(self): 

50 self.set_name('Cake Phase') 

51 

52 # self._phase_names = ('PmP ~S ~P ~P(moho)s ~P(sill-top)s' 

53 # ' ~P(sill-bottom)s Pdiff').split() 

54 self._phase_names = ('~P Pg Sg pP p P Pdiff PKP PcP PcS PKIKP pPKIKP' 

55 ' SSP PPS SPP PSP SP PS ~PS ~SP Pn s S Sn PP PPP' 

56 ' ScS Sdiff SS SSS SKS SKIKS').split() 

57 

58 for iphase, name in enumerate(self._phase_names): 

59 self.add_parameter(Switch(name, 'wantphase_%i' % iphase, 

60 iphase == 0)) 

61 

62 model_names = cake.builtin_models() 

63 model_names = [ 

64 'Cake builtin: %s' % model_name for model_name in model_names] 

65 

66 self._engine = gf.LocalEngine(use_config=True) 

67 store_ids = self._engine.get_store_ids() 

68 

69 for store_id in store_ids: 

70 model_names.append('GF Store: %s' % store_id) 

71 

72 self._models = model_names 

73 

74 self.model_choice = Choice('Model', 'chosen_model', 

75 'ak135-f-continental.m', self._models) 

76 

77 self.add_parameter(self.model_choice) 

78 

79 self.add_parameter(Param('Global shift', 'tshift', 0., -20., 20.)) 

80 self.add_parameter(Switch('Use station depth', 

81 'use_station_depth', False)) 

82 self.add_trigger('Add Phase', self.add_phase_definition) 

83 self.add_trigger('Add Model', self.add_model_to_choice) 

84 self.add_trigger('Plot Model', self.plot_model) 

85 self.add_trigger('Plot Rays', self.plot_rays) 

86 

87 self._phases = {} 

88 self._model = None 

89 

90 def panel_visibility_changed(self, bool): 

91 pass 

92 

93 def wanted_phases(self): 

94 try: 

95 wanted = [] 

96 for iphase, name in enumerate(self._phase_names): 

97 if getattr(self, 'wantphase_%i' % iphase): 

98 if name in self._phases: 

99 phases = self._phases[name] 

100 else: 

101 if name.startswith('~'): 

102 phases = [cake.PhaseDef(name[1:])] 

103 else: 

104 phases = cake.PhaseDef.classic(name) 

105 

106 self._phases[name] = phases 

107 for pha in phases: 

108 pha.name = name 

109 

110 wanted.extend(phases) 

111 except (cake.UnknownClassicPhase, cake.PhaseDefParseError) as e: 

112 self.fail(str(e)) 

113 

114 return wanted 

115 

116 def call(self, plot_rays=False): 

117 

118 self.cleanup() 

119 wanted = self.wanted_phases() 

120 

121 if not wanted: 

122 return 

123 

124 event, stations = self.get_active_event_and_stations() 

125 

126 if not stations: 

127 self.fail('No station information available.') 

128 

129 self.update_model() 

130 model = self._model[1] 

131 

132 depth = event.depth 

133 if depth is None: 

134 depth = 0.0 

135 

136 allrays = [] 

137 alldists = [] 

138 for station in stations: 

139 dist = event.distance_to(station) 

140 alldists.append(dist) 

141 

142 if self.use_station_depth: 

143 rdepth = station.depth 

144 else: 

145 rdepth = 0.0 

146 

147 multi_dists = [] 

148 nmax = 1 

149 for i in range(0, nmax): 

150 multi_dists.append(dist*cake.m2d + 360.*i) 

151 multi_dists.append((i+1)*360. - dist*cake.m2d) 

152 

153 rays = model.arrivals( 

154 phases=wanted, 

155 distances=multi_dists, 

156 zstart=depth, 

157 zstop=rdepth) 

158 

159 for ray in rays: 

160 time = ray.t 

161 name = ray.given_phase().name 

162 incidence_angle = ray.incidence_angle() 

163 takeoff_angle = ray.takeoff_angle() 

164 

165 time += event.time + self.tshift 

166 m = PhaseMarker([(station.network, station.station, '*', '*')], 

167 time, time, 2, 

168 phasename=name, 

169 event=event, 

170 incidence_angle=incidence_angle, 

171 takeoff_angle=takeoff_angle) 

172 self.add_marker(m) 

173 

174 allrays.extend(rays) 

175 

176 if plot_rays: 

177 fig = self.figure(name='Ray Paths') 

178 from pyrocko import cake_plot 

179 cake_plot.my_rays_plot(model, None, allrays, depth, 0.0, 

180 num.array(alldists)*cake.m2d, 

181 axes=fig.gca()) 

182 

183 fig.canvas.draw() 

184 

185 def update_model(self): 

186 if not self._model or self._model[0] != self.chosen_model: 

187 if self.chosen_model.startswith('Cake builtin: '): 

188 load_model = cake.load_model( 

189 self.chosen_model.split(': ', 1)[1]) 

190 

191 elif self.chosen_model.startswith('GF Store: '): 

192 store_id = self.chosen_model.split(': ', 1)[1] 

193 

194 load_model = self._engine.get_store(store_id)\ 

195 .config.earthmodel_1d 

196 else: 

197 load_model = cake.load_model(self.chosen_model) 

198 

199 self._model = (self.chosen_model, load_model) 

200 

201 def update_model_choices(self): 

202 self.set_parameter_choices('chosen_model', self._models) 

203 

204 def add_model_to_choice(self): 

205 ''' 

206 Called from trigger 'Add Model'. 

207 

208 Adds another choice to the drop down 'Model' menu. 

209 ''' 

210 

211 in_model = self.input_filename('Load Model') 

212 if in_model not in self._models: 

213 self._models.append(in_model) 

214 self.update_model_choices() 

215 

216 self.set_parameter('chosen_model', in_model) 

217 self.call() 

218 

219 def add_phase_definition(self): 

220 ''' 

221 Called from trigger 'Add Phase Definition'. 

222 

223 Adds another phase option. 

224 Requires a reset of the GUI. 

225 ''' 

226 phase_def = str(self.input_dialog('Add New Phase', 

227 'Enter Phase Definition')) 

228 self._phase_names.append(phase_def) 

229 

230 self.add_parameter( 

231 Switch(phase_def, 

232 'wantphase_%s' % str(len(self._phase_names)-1), True)) 

233 self.reset_gui(reloaded=True) 

234 self.call() 

235 

236 def plot_model(self): 

237 self.update_model() 

238 

239 from pyrocko import cake_plot 

240 

241 fig = self.figure(name='Model: %s' % self._model[0]) 

242 

243 cake_plot.my_model_plot(self._model[1], axes=fig.gca()) 

244 

245 fig.canvas.draw() 

246 

247 def plot_rays(self): 

248 self.call(plot_rays=True) 

249 

250 

251def __snufflings__(): 

252 ''' 

253 Returns a list of snufflings to be exported by this module. 

254 ''' 

255 

256 return [CakePhase()]