1#!/usr/bin/env python 

2# http://pyrocko.org - GPLv3 

3# 

4# The Pyrocko Developers, 21st Century 

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

6import sys 

7import re 

8import numpy as num 

9# import logger 

10from pyrocko import cake, util, orthodrome 

11from pyrocko.plot import cake_plot as plot 

12from optparse import OptionParser, OptionGroup 

13import matplotlib.pyplot as plt 

14from pyrocko.plot import mpl_init, mpl_papersize, mpl_margins 

15 

16r2d = cake.r2d 

17 

18 

19class Anon(dict): 

20 

21 def __getattr__(self, x): 

22 return self[x] 

23 

24 def getn(self, *keys): 

25 return Anon([(k, self[k]) for k in keys]) 

26 

27 def gett(self, *keys): 

28 return tuple([self[k] for k in keys]) 

29 

30 

31def process_color(s, phase_colors): 

32 m = re.match(r'^([^{]+){([^}]*)}$', s) 

33 if not m: 

34 return s 

35 

36 sphase = m.group(1) 

37 color = m.group(2) 

38 phase_colors[sphase] = plot.str_to_mpl_color(color) 

39 

40 return sphase 

41 

42 

43def optparse( 

44 required=(), 

45 optional=(), 

46 args=sys.argv, 

47 usage='%prog [options]', 

48 descr=None): 

49 

50 want = required + optional 

51 

52 parser = OptionParser( 

53 prog='cake', 

54 usage=usage, 

55 description=descr.capitalize()+'.', 

56 add_help_option=False, 

57 formatter=util.BetterHelpFormatter()) 

58 

59 parser.add_option( 

60 '-h', '--help', action='help', help='Show help message and exit.') 

61 

62 if 'phases' in want: 

63 group = OptionGroup(parser, 'Phases', ''' 

64 

65Seismic phase arrivals may be either specified as traditional phase names (e.g. 

66P, S, PP, PcP, ...) or in Cake's own syntax which is more powerful. Use the 

67--classic option, for traditional phase names. Use the --phase option if you 

68want to define phases in Cake's syntax. 

69 

70''') 

71 group.add_option( 

72 '--phase', '--phases', dest='phases', action="append", 

73 default=[], metavar='PHASE1,PHASE2,...', 

74 help='''Comma separated list of seismic phases in Cake\'s syntax. 

75 

76The definition of a seismic propagation path in Cake's phase syntax is a string 

77consisting of an alternating sequence of "legs" and "knees". 

78 

79A "leg" represents seismic wave propagation without any conversions, 

80encountering only super-critical reflections. Legs are denoted by "P", "p", 

81"S", or "s". The capital letters are used when the take-off of the "leg" is 

82in downward direction, while the lower case letters indicate a take-off in 

83upward direction. 

84 

85A "knee" is an interaction with an interface. It can be a mode conversion, a 

86reflection, or propagation as a headwave or diffracted wave. 

87 

88 * conversion is simply denoted as: "(INTERFACE)" or "DEPTH" 

89 * upperside reflection: "v(INTERFACE)" or "vDEPTH" 

90 * underside reflection: "^(INTERFACE)" or "^DEPTH" 

91 * normal kind headwave or diffracted wave: "v_(INTERFACE)" or "v_DEPTH" 

92 

93The interface may be given by name or by depth: INTERFACE is the name of an 

94interface defined in the model, DEPTH is the depth of an interface in [km] (the 

95interface closest to that depth is chosen). If two legs appear consecutively 

96without an explicit "knee", surface interaction is assumed. 

97 

98The preferred standard interface names in cake are "conrad", "moho", "cmb" 

99(core-mantle boundary), and "cb" (inner core boundary). 

100 

101The phase definition may end with a backslash "\\", to indicate that the ray 

102should arrive at the receiver from above instead of from below. It is possible 

103to restrict the maximum and minimum depth of a "leg" by appending 

104"<(INTERFACE)" or "<DEPTH" or ">(INTERFACE)" or ">DEPTH" after the leg 

105character, respectively. 

106 

107When plotting rays or travel-time curves, the color can be set by appending 

108"{COLOR}" to the phase definition, where COLOR is the name of a color or an RGB 

109or RGBA color tuple in the format "R/G/B" or "R/G/B/A", respectively. The 

110values can be normalized to the range [0, 1] or to [0, 255]. The latter is only 

111assumed when any of the values given exceeds 1.0. 

112''') 

113 

114 group.add_option( 

115 '--classic', dest='classic_phases', action='append', 

116 default=[], metavar='PHASE1,PHASE2,...', 

117 help='''Comma separated list of seismic phases in classic 

118nomenclature. Run "cake list-phase-map" for a list of available 

119phase names. When plotting, color can be specified in the same way 

120as in --phases.''') 

121 

122 parser.add_option_group(group) 

123 

124 if 'model' in want: 

125 group = OptionGroup(parser, 'Model') 

126 group.add_option( 

127 '--model', dest='model_filename', metavar='(NAME or FILENAME)', 

128 help='Use builtin model named NAME or user model from file ' 

129 'FILENAME. By default, the "ak135-f-continental.m" model is ' 

130 'used. Run "cake list-models" for a list of builtin models.') 

131 

132 group.add_option( 

133 '--format', dest='model_format', metavar='FORMAT', 

134 choices=['nd', 'hyposat'], default='nd', 

135 help='Set model file format (available: nd, hyposat; default: ' 

136 'nd).') 

137 group.add_option( 

138 '--crust2loc', dest='crust2loc', metavar='LAT,LON', 

139 help='Set model from CRUST2.0 profile at location (LAT,LON).') 

140 group.add_option( 

141 '--crust2profile', dest='crust2profile', metavar='KEY', 

142 help='Set model from CRUST2.0 profile with given KEY.') 

143 

144 parser.add_option_group(group) 

145 

146 if any(x in want for x in ( 

147 'zstart', 'zstop', 'distances', 'sloc', 'rloc')): 

148 

149 group = OptionGroup(parser, 'Source-receiver geometry') 

150 if 'zstart' in want: 

151 group.add_option( 

152 '--sdepth', dest='sdepth', type='float', default=0.0, 

153 metavar='FLOAT', 

154 help='Source depth [km] (default: 0)') 

155 if 'zstop' in want: 

156 group.add_option( 

157 '--rdepth', dest='rdepth', type='float', default=0.0, 

158 metavar='FLOAT', 

159 help='Receiver depth [km] (default: 0)') 

160 if 'distances' in want: 

161 group.add_option( 

162 '--distances', dest='sdist', metavar='DISTANCES', 

163 help='Surface distances as "start:stop:n" or ' 

164 '"dist1,dist2,..." [km]') 

165 group.add_option( 

166 '--sloc', dest='sloc', metavar='LAT,LON', 

167 help='Source location (LAT,LON).') 

168 group.add_option( 

169 '--rloc', dest='rloc', metavar='LAT,LON', 

170 help='Receiver location (LAT,LON).') 

171 parser.add_option_group(group) 

172 

173 if 'material' in want: 

174 group = OptionGroup( 

175 parser, 'Material', 

176 'An isotropic elastic material may be specified by giving ' 

177 'a combination of some of the following options. ') 

178 group.add_option( 

179 '--vp', dest='vp', default=None, type='float', metavar='FLOAT', 

180 help='P-wave velocity [km/s]') 

181 group.add_option( 

182 '--vs', dest='vs', default=None, type='float', metavar='FLOAT', 

183 help='S-wave velocity [km/s]') 

184 group.add_option( 

185 '--rho', dest='rho', default=None, type='float', metavar='FLOAT', 

186 help='density [g/cm**3]') 

187 group.add_option( 

188 '--qp', dest='qp', default=None, type='float', metavar='FLOAT', 

189 help='P-wave attenuation Qp (default: 1456)') 

190 group.add_option( 

191 '--qs', dest='qs', default=None, type='float', metavar='FLOAT', 

192 help='S-wave attenuation Qs (default: 600)') 

193 group.add_option( 

194 '--poisson', dest='poisson', default=None, type='float', 

195 metavar='FLOAT', 

196 help='Poisson ratio') 

197 group.add_option( 

198 '--lambda', dest='lame_lambda', default=None, type='float', 

199 metavar='FLOAT', 

200 help='Lame parameter lambda [GPa]') 

201 group.add_option( 

202 '--mu', dest='lame_mu', default=None, type='float', 

203 metavar='FLOAT', 

204 help='Shear modulus [GPa]') 

205 group.add_option( 

206 '--qk', dest='qk', default=None, type='float', metavar='FLOAT', 

207 help='Bulk attenuation Qk') 

208 group.add_option( 

209 '--qmu', dest='qmu', default=None, type='float', metavar='FLOAT', 

210 help='Shear attenuation Qmu') 

211 parser.add_option_group(group) 

212 

213 if any(x in want for x in ( 

214 'vred', 'as_degrees', 'accuracy', 'slowness', 'interface', 

215 'aspect', 'shade_model')): 

216 

217 group = OptionGroup(parser, 'General') 

218 if 'vred' in want: 

219 group.add_option( 

220 '--vred', dest='vred', type='float', metavar='FLOAT', 

221 help='Velocity for time reduction in plot [km/s]') 

222 

223 if 'as_degrees' in want: 

224 group.add_option( 

225 '--degrees', dest='as_degrees', action='store_true', 

226 default=False, 

227 help='Distances are in [deg] instead of [km], velocities in ' 

228 '[deg/s] instead of [km/s], slownesses in [s/deg] ' 

229 'instead of [s/km].') 

230 

231 if 'accuracy' in want: 

232 group.add_option( 

233 '--accuracy', dest='accuracy', type='float', 

234 metavar='MAXIMUM_RELATIVE_RMS', default=0.002, 

235 help='Set accuracy for model simplification.') 

236 

237 if 'slowness' in want: 

238 group.add_option( 

239 '--slowness', dest='slowness', type='float', metavar='FLOAT', 

240 default=0.0, 

241 help='Select surface slowness [s/km] (default: 0)') 

242 

243 if 'interface' in want: 

244 group.add_option( 

245 '--interface', dest='interface', metavar='(NAME or DEPTH)', 

246 help='Name or depth [km] of interface to select') 

247 

248 if 'aspect' in want: 

249 group.add_option( 

250 '--aspect', dest='aspect', type='float', metavar='FLOAT', 

251 help='Aspect ratio for plot') 

252 

253 if 'shade_model' in want: 

254 group.add_option( 

255 '--no-shade-model', dest='shade_model', action='store_false', 

256 default=True, 

257 help='Suppress shading of earth model layers') 

258 

259 parser.add_option_group(group) 

260 

261 if any(x in want for x in ('output_format', 'save', 'size', 'show')): 

262 group = OptionGroup(parser, 'Output', 'Output specifications') 

263 if 'output_format' in want: 

264 group.add_option( 

265 '--output-format', dest='output_format', metavar='FORMAT', 

266 default='textual', 

267 choices=('textual', 'nd'), 

268 help='Set model output format (available: textual, nd, ' 

269 'default: textual)') 

270 if 'save' in want: 

271 group.add_option( 

272 '-s', '--save', dest='save', metavar='PATH', 

273 help='saves plot to .png (default) or other py-supported\ 

274 endings without showing, use --show or -u for showing', 

275 default='') 

276 if 'size' in want: 

277 group.add_option( 

278 '--size', dest='size', type='string', 

279 default='a4', 

280 help='gives size of returned plot, use \'a5\' or \'a4\'') 

281 if 'show' in want: 

282 group.add_option( 

283 '-u', '--show', dest='show', action='store_true', 

284 help='shows plot when saving (-u for unhide)') 

285 

286 parser.add_option_group(group) 

287 

288 if usage == 'cake help-options': 

289 parser.print_help() 

290 

291 (options, args) = parser.parse_args(args) 

292 

293 if len(args) != 2: 

294 parser.error( 

295 'Cake arguments should look like "--option" or "--option=...".') 

296 

297 d = {} 

298 as_degrees = False 

299 if 'as_degrees' in want: 

300 as_degrees = options.as_degrees 

301 d['as_degrees'] = as_degrees 

302 

303 if 'accuracy' in want: 

304 d['accuracy'] = options.accuracy 

305 

306 if 'output_format' in want: 

307 d['output_format'] = options.output_format 

308 

309 if 'save' in want: 

310 d['save'] = options.save 

311 

312 if 'size' in want: 

313 d['size'] = options.size 

314 

315 if 'show' in want: 

316 d['show'] = options.show 

317 

318 if 'aspect' in want: 

319 d['aspect'] = options.aspect 

320 

321 if 'shade_model' in want: 

322 d['shade_model'] = options.shade_model 

323 

324 if 'phases' in want: 

325 phases = [] 

326 phase_colors = {} 

327 try: 

328 for ss in options.phases: 

329 for s in ss.split(','): 

330 s = process_color(s, phase_colors) 

331 phases.append(cake.PhaseDef(s)) 

332 

333 for pp in options.classic_phases: 

334 for p in pp.split(','): 

335 p = process_color(p, phase_colors) 

336 phases.extend(cake.PhaseDef.classic(p)) 

337 

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

339 parser.error(e) 

340 

341 if not phases and 'phases' in required: 

342 s = process_color('P', phase_colors) 

343 phases.append(cake.PhaseDef(s)) 

344 

345 if phases: 

346 d['phase_colors'] = phase_colors 

347 d['phases'] = phases 

348 

349 if 'model' in want: 

350 if options.model_filename: 

351 d['model'] = cake.load_model( 

352 options.model_filename, options.model_format) 

353 

354 if options.crust2loc or options.crust2profile: 

355 if options.crust2loc: 

356 try: 

357 args = tuple( 

358 [float(x) for x in options.crust2loc.split(',')]) 

359 except Exception: 

360 parser.error( 

361 'format for --crust2loc option is ' 

362 '"LATITUDE,LONGITUDE"') 

363 elif options.crust2profile: 

364 args = (options.crust2profile.upper(),) 

365 else: 

366 assert False 

367 

368 if 'model' in d: 

369 d['model'] = d['model'].replaced_crust(args) 

370 else: 

371 from pyrocko import crust2x2 

372 profile = crust2x2.get_profile(*args) 

373 d['model'] = cake.LayeredModel.from_scanlines( 

374 cake.from_crust2x2_profile(profile)) 

375 

376 if 'vred' in want: 

377 d['vred'] = options.vred 

378 if d['vred'] is not None: 

379 if not as_degrees: 

380 d['vred'] *= r2d * cake.km / cake.earthradius 

381 

382 if 'distances' in want: 

383 distances = None 

384 if options.sdist: 

385 if options.sdist.find(':') != -1: 

386 ssn = options.sdist.split(':') 

387 if len(ssn) != 3: 

388 parser.error( 

389 'format for distances is ' 

390 '"min_distance:max_distance:n_distances"') 

391 

392 distances = num.linspace( 

393 float(ssn[0]), float(ssn[1]), int(ssn[2])) 

394 else: 

395 distances = num.array( 

396 list(map( 

397 float, options.sdist.split(','))), dtype=float) 

398 

399 if not as_degrees: 

400 distances *= r2d * cake.km / cake.earthradius 

401 

402 if options.sloc and options.rloc: 

403 try: 

404 slat, slon = tuple([float(x) for x in options.sloc.split(',')]) 

405 rlat, rlon = tuple([float(x) for x in options.rloc.split(',')]) 

406 except Exception: 

407 parser.error( 

408 'format for --sloc and --rloc options is ' 

409 '"LATITUDE,LONGITUDE"') 

410 

411 distance_sr = orthodrome.distance_accurate50m_numpy( 

412 slat, slon, rlat, rlon) 

413 

414 distance_sr *= r2d / cake.earthradius 

415 if distances is not None: 

416 distances = num.concatenate((distances, [distance_sr])) 

417 else: 

418 distances = num.array([distance_sr], dtype=float) 

419 

420 if distances is not None: 

421 d['distances'] = distances 

422 else: 

423 if 'distances' not in required: 

424 d['distances'] = None 

425 

426 if 'slowness' in want: 

427 d['slowness'] = options.slowness/cake.d2r 

428 if not as_degrees: 

429 d['slowness'] /= cake.km*cake.m2d 

430 

431 if 'interface' in want: 

432 if options.interface: 

433 try: 

434 d['interface'] = float(options.interface)*cake.km 

435 except ValueError: 

436 d['interface'] = options.interface 

437 

438 else: 

439 d['interface'] = None 

440 

441 if 'zstart' in want: 

442 d['zstart'] = options.sdepth*cake.km 

443 

444 if 'zstop' in want: 

445 d['zstop'] = options.rdepth*cake.km 

446 

447 if 'material' in want: 

448 md = {} 

449 userfactor = dict( 

450 vp=1000., vs=1000., rho=1000., qp=1., qs=1., qmu=1., qk=1., 

451 lame_lambda=1.0e9, lame_mu=1.0e9, poisson=1.) 

452 

453 for k in userfactor.keys(): 

454 if getattr(options, k) is not None: 

455 md[k] = getattr(options, k) * userfactor[k] 

456 

457 if not (bool('lame_lambda' in md) == bool('lame_mu' in md)): 

458 parser.error('lambda and mu must be specified both.') 

459 if 'lame_lambda' in md and 'lame_mu' in md: 

460 md['lame'] = md.pop('lame_lambda'), md.pop('lame_mu') 

461 

462 if md: 

463 try: 

464 d['material'] = cake.Material(**md) 

465 except cake.InvalidArguments as e: 

466 parser.error(str(e)) 

467 

468 for k in list(d.keys()): 

469 if k not in want: 

470 del d[k] 

471 

472 for k in required: 

473 if k not in d: 

474 if k == 'model': 

475 d['model'] = cake.load_model('ak135-f-continental.m') 

476 

477 elif k == 'distances': 

478 d['distances'] = num.linspace(10*cake.km, 100*cake.km, 10) \ 

479 / cake.earthradius * r2d 

480 

481 elif k == 'phases': 

482 d['phases'] = list(map(cake.PhaseDef, 'Pp')) 

483 

484 else: 

485 parser.error('missing %s' % k) 

486 

487 return Anon(d) 

488 

489 

490def my_simplify_model(mod, accuracy): 

491 mod_simple = mod.simplify(max_rel_error=accuracy) 

492 cake.write_nd_model_fh(mod_simple, sys.stdout) 

493 

494 

495def d2u(d): 

496 return dict((k.replace('-', '_'), v) for (k, v) in d.items()) 

497 

498 

499def mini_fmt(v, d=5, try_fmts='fe'): 

500 for fmt in try_fmts: 

501 for i in range(d, -1, -1): 

502 s = ('%%%i.%i%s' % (d, i, fmt)) % v 

503 if len(s) <= d and (v == 0.0 or float(s) != 0.0): 

504 return s 

505 

506 return s 

507 

508 

509def scatter_out_fmt(d, m, v): 

510 if v < 1e-8: 

511 return ' ', ' ' 

512 else: 

513 return '%s%s' % ('\\/'[d == cake.UP], 'SP'[m == cake.P]), \ 

514 mini_fmt(v*100, 5) 

515 

516 

517def scatter_in_fmt(d, m, dwant): 

518 if d == dwant: 

519 return '%s%s' % ('\\/'[d == cake.UP], 'SP'[m == cake.P]) 

520 else: 

521 return ' ' 

522 

523 

524def print_scatter(model, p=0.0, interface=None): 

525 if interface is not None: 

526 discontinuities = [model.discontinuity(interface)] 

527 else: 

528 discontinuities = model.discontinuities() 

529 

530 for discontinuity in discontinuities: 

531 print('%s (%g km)' % (discontinuity, discontinuity.z/cake.km)) 

532 print() 

533 cols = [] 

534 for in_direction in (cake.DOWN, cake.UP): 

535 for in_mode in (cake.P, cake.S): 

536 p_critical = discontinuity.critical_ps(in_mode)[ 

537 in_direction == cake.UP] 

538 

539 if p_critical is None or p >= p_critical: 

540 continue 

541 

542 vals = [] 

543 for out_direction in (cake.UP, cake.DOWN): 

544 for out_mode in (cake.S, cake.P): 

545 vals.append( 

546 (out_direction, out_mode, discontinuity.efficiency( 

547 in_direction, out_direction, 

548 in_mode, out_mode, p))) 

549 

550 if all(x[-1] == 0.0 for x in vals): 

551 continue 

552 

553 sout = [scatter_out_fmt(d, m, v) for (d, m, v) in vals] 

554 

555 sin1 = scatter_in_fmt(in_direction, in_mode, cake.DOWN) 

556 sin2 = scatter_in_fmt(in_direction, in_mode, cake.UP) 

557 line1 = '%s %5s %5s' % ( 

558 ' '*len(sin1), sout[0][1], sout[1][1]) 

559 line2 = '%s %-5s %-5s' % ( 

560 sin1, sout[0][0], sout[1][0]) 

561 line4 = '%s %-5s %-5s' % ( 

562 sin2, sout[2][0], sout[3][0]) 

563 line5 = '%s %5s %5s' % ( 

564 ' '*len(sin2), sout[2][1], sout[3][1]) 

565 line3 = '-' * len(line1) 

566 cols.append((line1, line2, line3, line4, line5)) 

567 

568 for cols in zip(*cols): 

569 print(' ' + ' '.join(cols)) 

570 

571 print() 

572 print() 

573 

574 

575def print_arrivals( 

576 model, distances=[], phases=cake.PhaseDef('P'), zstart=0.0, zstop=0.0, 

577 as_degrees=False): 

578 

579 headers = 'slow dist time take inci effi spre phase used'.split() 

580 space = (7, 5, 6, 4, 4, 4, 4, 17, 17) 

581 

582 if as_degrees: 

583 units = 's/deg deg s deg deg % %'.split() 

584 else: 

585 units = 's/km km s deg deg % %'.split() 

586 

587 hline = ' '.join(x.ljust(s) for (x, s) in zip(headers, space)) 

588 uline = ' '.join(('%s' % x).ljust(s) for (x, s) in zip(units, space)) 

589 

590 print(hline) 

591 print(uline) 

592 print('-' * len(hline)) 

593 

594 for ray in model.arrivals( 

595 distances=distances, phases=phases, zstart=zstart, zstop=zstop): 

596 

597 if as_degrees: 

598 sd = ray.x 

599 slow = ray.p/cake.r2d 

600 else: 

601 sd = ray.x*(cake.d2r*cake.earthradius/cake.km) 

602 slow = ray.p/(r2d*cake.d2m/cake.km) 

603 

604 su = '(%s)' % ray.path.used_phase(p=ray.p, eps=1.0).used_repr() 

605 

606 print(' '.join(tuple(mini_fmt(x, s).rjust(s) for (x, s) in zip( 

607 (slow, sd, ray.t, ray.takeoff_angle(), ray.incidence_angle(), 

608 100*ray.efficiency(), 100*ray.spreading()*ray.surface_sphere()), 

609 space)) + tuple( 

610 x.ljust(17) for x in (ray.path.phase.definition(), su)))) 

611 

612 

613def plot_init(size, save, show): 

614 fontsize = 9 

615 mpl_init() 

616 fig = plt.figure(figsize=mpl_papersize(size, 'landscape')) 

617 

618 labelpos = mpl_margins(fig, w=7., h=5., units=fontsize) 

619 axes = fig.add_subplot(1, 1, 1) 

620 labelpos(axes, 2., 1.5) 

621 

622 showplt = bool(show or not save) 

623 

624 return fig, axes, showplt 

625 

626 

627class CakeError(Exception): 

628 pass 

629 

630 

631def plot_end(save, fig, show=True): 

632 if save: 

633 try: 

634 fig.savefig(save) 

635 if show: 

636 plt.show() 

637 except OSError as e: 

638 raise CakeError(str(e)) 

639 except ValueError as e: 

640 raise CakeError(str(e)) 

641 

642 

643def main(args=None): 

644 

645 if args is None: 

646 args = sys.argv[1:] 

647 

648 subcommand_descriptions = { 

649 'print': 'get information on model/phase/material properties', 

650 'arrivals': 'print list of phase arrivals', 

651 'paths': 'print ray path details', 

652 'plot-xt': 'plot traveltime vs distance curves', 

653 'plot-xp': 'plot ray parameter vs distance curves', 

654 'plot-rays': 'plot ray propagation paths', 

655 'plot': 'plot combination of ray and traveltime curves', 

656 'plot-model': 'plot velocity model', 

657 'list-models': 'list builtin velocity models', 

658 'list-phase-map': 'show translation table for classic phase names', 

659 'simplify-model': 'create a simplified version of a layered model', 

660 'scatter': 'show details about scattering at model interfaces'} 

661 

662 usage = '''cake <subcommand> [options] 

663 

664Subcommands: 

665 

666 print %(print)s 

667 arrivals %(arrivals)s 

668 paths %(paths)s 

669 plot-xt %(plot_xt)s 

670 plot-xp %(plot_xp)s 

671 plot-rays %(plot_rays)s 

672 plot %(plot)s 

673 plot-model %(plot_model)s 

674 list-models %(list_models)s 

675 list-phase-map %(list_phase_map)s 

676 simplify-model %(simplify_model)s 

677 scatter %(scatter)s 

678 

679To get further help and a list of available options for any subcommand run: 

680 

681 cake <subcommand> --help 

682 

683'''.strip() % d2u(subcommand_descriptions) 

684 

685 usage_sub = 'cake %s [options]' 

686 if len(args) < 1: 

687 sys.exit('Usage: %s' % usage) 

688 

689 command = args[0] 

690 

691 args[0:0] = ['cake'] 

692 descr = subcommand_descriptions.get(command, None) 

693 subusage = usage_sub % command 

694 

695 if command == 'print': 

696 c = optparse( 

697 (), ('model', 'phases', 'material', 'output_format'), 

698 usage=subusage, descr=descr, args=args) 

699 

700 if 'model' in c: 

701 if c.output_format == 'textual': 

702 print(c.model) 

703 print() 

704 elif c.output_format == 'nd': 

705 cake.write_nd_model_fh(c.model, sys.stdout) 

706 

707 if 'phases' in c: 

708 for phase in c.phases: 

709 print(phase) 

710 print() 

711 

712 if 'material' in c: 

713 print(c.material.describe()) 

714 print() 

715 

716 elif command == 'arrivals': 

717 c = optparse( 

718 ('model', 'phases', 'distances'), 

719 ('zstart', 'zstop', 'as_degrees'), 

720 usage=subusage, descr=descr, args=args) 

721 

722 print_arrivals( 

723 c.model, 

724 **c.getn('zstart', 'zstop', 'phases', 'distances', 'as_degrees')) 

725 

726 elif command == 'paths': 

727 c = optparse( 

728 ('model', 'phases'), 

729 ('zstart', 'zstop', 'as_degrees'), 

730 usage=subusage, descr=descr, args=args) 

731 

732 mod = c.model 

733 for path in mod.gather_paths(**c.getn('phases', 'zstart', 'zstop')): 

734 print(path.describe(path.endgaps(c.zstart, c.zstop), c.as_degrees)) 

735 

736 elif command in ('plot-xt', 'plot-xp', 'plot-rays', 'plot'): 

737 if command in ('plot-xt', 'plot'): 

738 c = optparse( 

739 ('model', 'phases',), 

740 ('zstart', 'zstop', 'distances', 'as_degrees', 'vred', 

741 'phase_colors', 'save', 'size', 'show'), 

742 usage=subusage, descr=descr, args=args) 

743 else: 

744 c = optparse( 

745 ('model', 'phases'), 

746 ('zstart', 'zstop', 'distances', 'as_degrees', 'aspect', 

747 'shade_model', 'phase_colors', 'save', 'size', 'show'), 

748 usage=subusage, descr=descr, args=args) 

749 

750 mod = c.model 

751 paths = mod.gather_paths(**c.getn('phases', 'zstart', 'zstop')) 

752 

753 if c.distances is not None: 

754 arrivals = mod.arrivals( 

755 **c.getn('phases', 'zstart', 'zstop', 'distances')) 

756 else: 

757 arrivals = None 

758 

759 fig, axes, showplt = plot_init(c.size, c.save, c.show) 

760 

761 if command == 'plot-xp': 

762 plot.my_xp_plot( 

763 paths, c.zstart, c.zstop, c.distances, 

764 c.as_degrees, show=showplt, phase_colors=c.phase_colors) 

765 

766 elif command == 'plot-xt': 

767 plot.my_xt_plot( 

768 paths, c.zstart, c.zstop, c.distances, c.as_degrees, 

769 vred=c.vred, show=showplt, 

770 phase_colors=c.phase_colors) 

771 

772 elif command == 'plot-rays': 

773 if c.as_degrees: 

774 plot.my_rays_plot_gcs( 

775 mod, paths, arrivals, c.zstart, c.zstop, c.distances, 

776 show=showplt, phase_colors=c.phase_colors) 

777 

778 else: 

779 plot.my_rays_plot( 

780 mod, paths, arrivals, c.zstart, c.zstop, c.distances, 

781 show=showplt, aspect=c.aspect, shade_model=c.shade_model, 

782 phase_colors=c.phase_colors) 

783 

784 elif command == 'plot': 

785 plot.my_combi_plot( 

786 mod, paths, arrivals, c.zstart, c.zstop, c.distances, 

787 c.as_degrees, show=showplt, vred=c.vred, 

788 phase_colors=c.phase_colors) 

789 

790 try: 

791 plot_end(save=c.save, fig=fig, show=c.show) 

792 except CakeError as e: 

793 exit('cake.py: %s' % str(e)) 

794 

795 elif command in ('plot-model',): 

796 c = optparse( 

797 ('model',), ('save', 'size', 'show'), 

798 usage=subusage, descr=descr, args=args) 

799 mod = c.model 

800 fig, axes, showplt = plot_init(c.size, c.save, c.show) 

801 plot.my_model_plot(mod, show=showplt) 

802 try: 

803 plot_end(save=c.save, fig=fig, show=c.show) 

804 except CakeError as e: 

805 exit('cake.py: %s' % str(e)) 

806 

807 elif command in ('simplify-model',): 

808 c = optparse(('model',), ('accuracy',), 

809 usage=subusage, descr=descr, args=args) 

810 my_simplify_model(c.model, c.accuracy) 

811 

812 elif command in ('list-models',): 

813 c = optparse((), (), usage=subusage, descr=descr, args=args) 

814 for x in cake.builtin_models(): 

815 print(x) 

816 

817 elif command in ('list-phase-map',): 

818 c = optparse((), (), usage=subusage, descr=descr, args=args) 

819 defs = cake.PhaseDef.classic_definitions() 

820 for k in sorted(defs.keys()): 

821 print('%-15s: %s' % (k, ', '.join(defs[k]))) 

822 

823 elif command in ('scatter',): 

824 c = optparse( 

825 ('model',), 

826 ('slowness', 'interface', 'as_degrees'), 

827 usage=subusage, descr=descr, args=args) 

828 

829 print_scatter(c.model, p=c.slowness, interface=c.interface) 

830 

831 elif command in ('help-options',): 

832 optparse( 

833 (), 

834 ('model', 'accuracy', 'slowness', 'interface', 'phases', 

835 'distances', 'zstart', 'zstop', 'distances', 'as_degrees', 

836 'material', 'vred', 'save'), 

837 usage='cake help-options', descr='list all available options', 

838 args=args) 

839 

840 elif command in ('--help', '-h', 'help'): 

841 sys.exit('Usage: %s' % usage) 

842 

843 else: 

844 sys.exit('cake: no such subcommand: %s' % command) 

845 

846 

847if __name__ == '__main__': 

848 main(sys.argv[1:])