1# http://pyrocko.org - GPLv3 

2# 

3# The Pyrocko Developers, 21st Century 

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

5 

6''' 

7Data model and content types handled by the Squirrel framework. 

8 

9Squirrel uses flat content types to represent waveform, station, channel, 

10response, event, and a few other objects. A common subset of the information in 

11these objects is indexed in the database, currently: kind, codes, time interval 

12and sampling rate. The :py:class:`Nut` objects encapsulate this information 

13together with information about where and how to get the associated bulk data. 

14 

15Further content types are defined here to handle waveform orders, waveform 

16promises, data coverage and sensor information. 

17''' 

18 

19from __future__ import absolute_import, print_function 

20 

21import re 

22import fnmatch 

23import hashlib 

24import numpy as num 

25from os import urandom 

26from base64 import urlsafe_b64encode 

27from collections import defaultdict, namedtuple 

28 

29from pyrocko import util 

30from pyrocko.guts import Object, SObject, String, Timestamp, Float, Int, \ 

31 Unicode, Tuple, List, StringChoice, Any 

32from pyrocko.model import squirrel_content, Location 

33from pyrocko.response import FrequencyResponse, MultiplyResponse, \ 

34 IntegrationResponse, DifferentiationResponse, simplify_responses, \ 

35 FrequencyResponseCheckpoint 

36 

37from .error import ConversionError, SquirrelError 

38 

39 

40guts_prefix = 'squirrel' 

41 

42 

43g_codes_pool = {} 

44 

45 

46class CodesError(SquirrelError): 

47 pass 

48 

49 

50class Codes(SObject): 

51 pass 

52 

53 

54def normalize_nslce(*args, **kwargs): 

55 if args and kwargs: 

56 raise ValueError('Either *args or **kwargs accepted, not both.') 

57 

58 if len(args) == 1: 

59 if isinstance(args[0], str): 

60 args = tuple(args[0].split('.')) 

61 elif isinstance(args[0], tuple): 

62 args = args[0] 

63 else: 

64 raise ValueError('Invalid argument type: %s' % type(args[0])) 

65 

66 nargs = len(args) 

67 if nargs == 5: 

68 t = args 

69 

70 elif nargs == 4: 

71 t = args + ('',) 

72 

73 elif nargs == 0: 

74 d = dict( 

75 network='', 

76 station='', 

77 location='', 

78 channel='', 

79 extra='') 

80 

81 d.update(kwargs) 

82 t = tuple(kwargs.get(k, '') for k in ( 

83 'network', 'station', 'location', 'channel', 'extra')) 

84 

85 else: 

86 raise CodesError( 

87 'Does not match NSLC or NSLCE codes pattern: %s' % '.'.join(args)) 

88 

89 if '.'.join(t).count('.') != 4: 

90 raise CodesError( 

91 'Codes may not contain a ".": "%s", "%s", "%s", "%s", "%s"' % t) 

92 

93 return t 

94 

95 

96CodesNSLCEBase = namedtuple( 

97 'CodesNSLCEBase', [ 

98 'network', 'station', 'location', 'channel', 'extra']) 

99 

100 

101class CodesNSLCE(CodesNSLCEBase, Codes): 

102 ''' 

103 Codes denominating a seismic channel (NSLC or NSLCE). 

104 

105 FDSN/SEED style NET.STA.LOC.CHA is accepted or NET.STA.LOC.CHA.EXTRA, where 

106 the EXTRA part in the latter form can be used to identify a custom 

107 processing applied to a channel. 

108 ''' 

109 

110 __slots__ = () 

111 __hash__ = CodesNSLCEBase.__hash__ 

112 

113 as_dict = CodesNSLCEBase._asdict 

114 

115 def __new__(cls, *args, safe_str=None, **kwargs): 

116 nargs = len(args) 

117 if nargs == 1 and isinstance(args[0], CodesNSLCE): 

118 return args[0] 

119 elif nargs == 1 and isinstance(args[0], CodesNSL): 

120 t = (args[0].nsl) + ('*', '*') 

121 elif nargs == 1 and isinstance(args[0], CodesX): 

122 t = ('*', '*', '*', '*', '*') 

123 elif safe_str is not None: 

124 t = safe_str.split('.') 

125 else: 

126 t = normalize_nslce(*args, **kwargs) 

127 

128 x = CodesNSLCEBase.__new__(cls, *t) 

129 return g_codes_pool.setdefault(x, x) 

130 

131 def __init__(self, *args, **kwargs): 

132 Codes.__init__(self) 

133 

134 def __str__(self): 

135 if self.extra == '': 

136 return '.'.join(self[:-1]) 

137 else: 

138 return '.'.join(self) 

139 

140 def __eq__(self, other): 

141 if not isinstance(other, CodesNSLCE): 

142 other = CodesNSLCE(other) 

143 

144 return CodesNSLCEBase.__eq__(self, other) 

145 

146 @property 

147 def safe_str(self): 

148 return '.'.join(self) 

149 

150 @property 

151 def nslce(self): 

152 return self[:4] 

153 

154 @property 

155 def nslc(self): 

156 return self[:4] 

157 

158 @property 

159 def nsl(self): 

160 return self[:3] 

161 

162 @property 

163 def ns(self): 

164 return self[:2] 

165 

166 def as_tuple(self): 

167 return tuple(self) 

168 

169 def replace(self, **kwargs): 

170 x = CodesNSLCEBase._replace(self, **kwargs) 

171 return g_codes_pool.setdefault(x, x) 

172 

173 

174def normalize_nsl(*args, **kwargs): 

175 if args and kwargs: 

176 raise ValueError('Either *args or **kwargs accepted, not both.') 

177 

178 if len(args) == 1: 

179 if isinstance(args[0], str): 

180 args = tuple(args[0].split('.')) 

181 elif isinstance(args[0], tuple): 

182 args = args[0] 

183 else: 

184 raise ValueError('Invalid argument type: %s' % type(args[0])) 

185 

186 nargs = len(args) 

187 if nargs == 3: 

188 t = args 

189 

190 elif nargs == 0: 

191 d = dict( 

192 network='', 

193 station='', 

194 location='') 

195 

196 d.update(kwargs) 

197 t = tuple(kwargs.get(k, '') for k in ( 

198 'network', 'station', 'location')) 

199 

200 else: 

201 raise CodesError( 

202 'Does not match NSL codes pattern: %s' % '.'.join(args)) 

203 

204 if '.'.join(t).count('.') != 2: 

205 raise CodesError( 

206 'Codes may not contain a ".": "%s", "%s", "%s"' % t) 

207 

208 return t 

209 

210 

211CodesNSLBase = namedtuple( 

212 'CodesNSLBase', [ 

213 'network', 'station', 'location']) 

214 

215 

216class CodesNSL(CodesNSLBase, Codes): 

217 ''' 

218 Codes denominating a seismic station (NSL). 

219 

220 NET.STA.LOC is accepted, slightly different from SEED/StationXML, where 

221 LOC is part of the channel. By setting location='*' is possible to get 

222 compatible behaviour in most cases. 

223 ''' 

224 

225 __slots__ = () 

226 __hash__ = CodesNSLBase.__hash__ 

227 

228 as_dict = CodesNSLBase._asdict 

229 

230 def __new__(cls, *args, safe_str=None, **kwargs): 

231 nargs = len(args) 

232 if nargs == 1 and isinstance(args[0], CodesNSL): 

233 return args[0] 

234 elif nargs == 1 and isinstance(args[0], CodesNSLCE): 

235 t = args[0].nsl 

236 elif nargs == 1 and isinstance(args[0], CodesX): 

237 t = ('*', '*', '*') 

238 elif safe_str is not None: 

239 t = safe_str.split('.') 

240 else: 

241 t = normalize_nsl(*args, **kwargs) 

242 

243 x = CodesNSLBase.__new__(cls, *t) 

244 return g_codes_pool.setdefault(x, x) 

245 

246 def __init__(self, *args, **kwargs): 

247 Codes.__init__(self) 

248 

249 def __str__(self): 

250 return '.'.join(self) 

251 

252 def __eq__(self, other): 

253 if not isinstance(other, CodesNSL): 

254 other = CodesNSL(other) 

255 

256 return CodesNSLBase.__eq__(self, other) 

257 

258 @property 

259 def safe_str(self): 

260 return '.'.join(self) 

261 

262 @property 

263 def ns(self): 

264 return self[:2] 

265 

266 @property 

267 def nsl(self): 

268 return self[:3] 

269 

270 def as_tuple(self): 

271 return tuple(self) 

272 

273 def replace(self, **kwargs): 

274 x = CodesNSLBase._replace(self, **kwargs) 

275 return g_codes_pool.setdefault(x, x) 

276 

277 

278CodesXBase = namedtuple( 

279 'CodesXBase', [ 

280 'name']) 

281 

282 

283class CodesX(CodesXBase, Codes): 

284 ''' 

285 General purpose codes for anything other than channels or stations. 

286 ''' 

287 

288 __slots__ = () 

289 __hash__ = CodesXBase.__hash__ 

290 __eq__ = CodesXBase.__eq__ 

291 

292 as_dict = CodesXBase._asdict 

293 

294 def __new__(cls, name='', safe_str=None): 

295 if isinstance(name, CodesX): 

296 return name 

297 elif isinstance(name, (CodesNSLCE, CodesNSL)): 

298 name = '*' 

299 elif safe_str is not None: 

300 name = safe_str 

301 else: 

302 if '.' in name: 

303 raise CodesError('Code may not contain a ".": %s' % name) 

304 

305 x = CodesXBase.__new__(cls, name) 

306 return g_codes_pool.setdefault(x, x) 

307 

308 def __init__(self, *args, **kwargs): 

309 Codes.__init__(self) 

310 

311 def __str__(self): 

312 return '.'.join(self) 

313 

314 @property 

315 def safe_str(self): 

316 return '.'.join(self) 

317 

318 @property 

319 def ns(self): 

320 return self[:2] 

321 

322 def as_tuple(self): 

323 return tuple(self) 

324 

325 def replace(self, **kwargs): 

326 x = CodesXBase._replace(self, **kwargs) 

327 return g_codes_pool.setdefault(x, x) 

328 

329 

330g_codes_patterns = {} 

331 

332 

333def match_codes(pattern, codes): 

334 spattern = pattern.safe_str 

335 scodes = codes.safe_str 

336 if spattern not in g_codes_patterns: 

337 rpattern = re.compile(fnmatch.translate(spattern), re.I) 

338 g_codes_patterns[spattern] = rpattern 

339 

340 rpattern = g_codes_patterns[spattern] 

341 return bool(rpattern.match(scodes)) 

342 

343 

344g_content_kinds = [ 

345 'undefined', 

346 'waveform', 

347 'station', 

348 'channel', 

349 'response', 

350 'event', 

351 'waveform_promise'] 

352 

353 

354g_codes_classes = [ 

355 CodesX, 

356 CodesNSLCE, 

357 CodesNSL, 

358 CodesNSLCE, 

359 CodesNSLCE, 

360 CodesX, 

361 CodesNSLCE] 

362 

363g_codes_classes_ndot = { 

364 0: CodesX, 

365 2: CodesNSL, 

366 3: CodesNSLCE, 

367 4: CodesNSLCE} 

368 

369 

370def to_codes_simple(kind_id, codes_safe_str): 

371 return g_codes_classes[kind_id](safe_str=codes_safe_str) 

372 

373 

374def to_codes(kind_id, obj): 

375 return g_codes_classes[kind_id](obj) 

376 

377 

378def to_codes_guess(s): 

379 try: 

380 return g_codes_classes_ndot[s.count('.')](s) 

381 except KeyError: 

382 raise CodesError('Cannot guess codes type: %s' % s) 

383 

384 

385g_content_kind_ids = ( 

386 UNDEFINED, WAVEFORM, STATION, CHANNEL, RESPONSE, EVENT, 

387 WAVEFORM_PROMISE) = range(len(g_content_kinds)) 

388 

389 

390g_tmin, g_tmax = util.get_working_system_time_range()[:2] 

391 

392 

393try: 

394 g_tmin_queries = max(g_tmin, util.str_to_time_fillup('1900-01-01')) 

395except Exception: 

396 g_tmin_queries = g_tmin 

397 

398 

399def to_kind(kind_id): 

400 return g_content_kinds[kind_id] 

401 

402 

403def to_kinds(kind_ids): 

404 return [g_content_kinds[kind_id] for kind_id in kind_ids] 

405 

406 

407def to_kind_id(kind): 

408 return g_content_kinds.index(kind) 

409 

410 

411def to_kind_ids(kinds): 

412 return [g_content_kinds.index(kind) for kind in kinds] 

413 

414 

415g_kind_mask_all = 0xff 

416 

417 

418def to_kind_mask(kinds): 

419 if kinds: 

420 return sum(1 << kind_id for kind_id in to_kind_ids(kinds)) 

421 else: 

422 return g_kind_mask_all 

423 

424 

425def str_or_none(x): 

426 if x is None: 

427 return None 

428 else: 

429 return str(x) 

430 

431 

432def float_or_none(x): 

433 if x is None: 

434 return None 

435 else: 

436 return float(x) 

437 

438 

439def int_or_none(x): 

440 if x is None: 

441 return None 

442 else: 

443 return int(x) 

444 

445 

446def int_or_g_tmin(x): 

447 if x is None: 

448 return g_tmin 

449 else: 

450 return int(x) 

451 

452 

453def int_or_g_tmax(x): 

454 if x is None: 

455 return g_tmax 

456 else: 

457 return int(x) 

458 

459 

460def tmin_or_none(x): 

461 if x == g_tmin: 

462 return None 

463 else: 

464 return x 

465 

466 

467def tmax_or_none(x): 

468 if x == g_tmax: 

469 return None 

470 else: 

471 return x 

472 

473 

474def time_or_none_to_str(x): 

475 if x is None: 

476 return '...'.ljust(17) 

477 else: 

478 return util.time_to_str(x) 

479 

480 

481def codes_to_str_abbreviated(codes, indent=' '): 

482 codes = [str(x) for x in codes] 

483 

484 if len(codes) > 20: 

485 scodes = '\n' + util.ewrap(codes[:10], indent=indent) \ 

486 + '\n%s[%i more]\n' % (indent, len(codes) - 20) \ 

487 + util.ewrap(codes[-10:], indent=' ') 

488 else: 

489 scodes = '\n' + util.ewrap(codes, indent=indent) \ 

490 if codes else '<none>' 

491 

492 return scodes 

493 

494 

495g_offset_time_unit_inv = 1000000000 

496g_offset_time_unit = 1.0 / g_offset_time_unit_inv 

497 

498 

499def tsplit(t): 

500 if t is None: 

501 return None, 0.0 

502 

503 t = util.to_time_float(t) 

504 if type(t) is float: 

505 t = round(t, 5) 

506 else: 

507 t = round(t, 9) 

508 

509 seconds = num.floor(t) 

510 offset = t - seconds 

511 return int(seconds), int(round(offset * g_offset_time_unit_inv)) 

512 

513 

514def tjoin(seconds, offset): 

515 if seconds is None: 

516 return None 

517 

518 return util.to_time_float(seconds) \ 

519 + util.to_time_float(offset*g_offset_time_unit) 

520 

521 

522tscale_min = 1 

523tscale_max = 365 * 24 * 3600 # last edge is one above 

524tscale_logbase = 20 

525 

526tscale_edges = [tscale_min] 

527while True: 

528 tscale_edges.append(tscale_edges[-1]*tscale_logbase) 

529 if tscale_edges[-1] >= tscale_max: 

530 break 

531 

532 

533tscale_edges = num.array(tscale_edges) 

534 

535 

536def tscale_to_kscale(tscale): 

537 

538 # 0 <= x < tscale_edges[1]: 0 

539 # tscale_edges[1] <= x < tscale_edges[2]: 1 

540 # ... 

541 # tscale_edges[len(tscale_edges)-1] <= x: len(tscale_edges) 

542 

543 return int(num.searchsorted(tscale_edges, tscale)) 

544 

545 

546@squirrel_content 

547class Station(Location): 

548 ''' 

549 A seismic station. 

550 ''' 

551 

552 codes = CodesNSL.T() 

553 

554 tmin = Timestamp.T(optional=True) 

555 tmax = Timestamp.T(optional=True) 

556 

557 description = Unicode.T(optional=True) 

558 

559 def __init__(self, **kwargs): 

560 kwargs['codes'] = CodesNSL(kwargs['codes']) 

561 Location.__init__(self, **kwargs) 

562 

563 @property 

564 def time_span(self): 

565 return (self.tmin, self.tmax) 

566 

567 def get_pyrocko_station(self): 

568 from pyrocko import model 

569 return model.Station(*self._get_pyrocko_station_args()) 

570 

571 def _get_pyrocko_station_args(self): 

572 return ( 

573 self.codes.network, 

574 self.codes.station, 

575 self.codes.location, 

576 self.lat, 

577 self.lon, 

578 self.elevation, 

579 self.depth, 

580 self.north_shift, 

581 self.east_shift) 

582 

583 

584class Sensor(Location): 

585 ''' 

586 Representation of a channel group. 

587 ''' 

588 

589 codes = CodesNSLCE.T() 

590 

591 tmin = Timestamp.T(optional=True) 

592 tmax = Timestamp.T(optional=True) 

593 

594 deltat = Float.T(optional=True) 

595 

596 @property 

597 def time_span(self): 

598 return (self.tmin, self.tmax) 

599 

600 def __init__(self, **kwargs): 

601 kwargs['codes'] = CodesNSLCE(kwargs['codes']) 

602 Location.__init__(self, **kwargs) 

603 

604 def _get_sensor_args(self): 

605 def getattr_rep(k): 

606 if k == 'codes': 

607 return self.codes.replace(self.codes[:-1]) 

608 else: 

609 return getattr(self, k) 

610 

611 return [getattr_rep(k) for k in self.T.propnames] 

612 

613 @classmethod 

614 def from_channels(cls, channels): 

615 groups = defaultdict(list) 

616 for channel in channels: 

617 groups[channel._get_sensor_args()].append(channel) 

618 

619 return [ 

620 cls(**dict((k, v) for (k, v) in zip(cls.T.propnames, args))) 

621 for args, _ in groups.items()] 

622 

623 def _get_pyrocko_station_args(self): 

624 return ( 

625 self.codes.network, 

626 self.codes.station, 

627 self.codes.location, 

628 self.lat, 

629 self.lon, 

630 self.elevation, 

631 self.depth, 

632 self.north_shift, 

633 self.east_shift) 

634 

635 

636@squirrel_content 

637class Channel(Sensor): 

638 ''' 

639 A channel of a seismic station. 

640 ''' 

641 

642 dip = Float.T(optional=True) 

643 azimuth = Float.T(optional=True) 

644 

645 @classmethod 

646 def from_channels(cls, channels): 

647 raise NotImplementedError() 

648 

649 def get_pyrocko_channel(self): 

650 from pyrocko import model 

651 return model.Channel(*self._get_pyrocko_channel_args()) 

652 

653 def _get_pyrocko_channel_args(self): 

654 return ( 

655 self.codes.channel, 

656 self.azimuth, 

657 self.dip) 

658 

659 

660observational_quantities = [ 

661 'acceleration', 'velocity', 'displacement', 'pressure', 'rotation', 

662 'temperature'] 

663 

664 

665technical_quantities = [ 

666 'voltage', 'counts'] 

667 

668 

669class QuantityType(StringChoice): 

670 ''' 

671 Choice of observational or technical quantity. 

672 

673 SI units are used for all quantities, where applicable. 

674 ''' 

675 choices = observational_quantities + technical_quantities 

676 

677 

678class ResponseStage(Object): 

679 ''' 

680 Representation of a response stage. 

681 

682 Components of a seismic recording system are represented as a sequence of 

683 response stages, e.g. sensor, pre-amplifier, digitizer, digital 

684 downsampling. 

685 ''' 

686 input_quantity = QuantityType.T(optional=True) 

687 input_sample_rate = Float.T(optional=True) 

688 output_quantity = QuantityType.T(optional=True) 

689 output_sample_rate = Float.T(optional=True) 

690 elements = List.T(FrequencyResponse.T()) 

691 log = List.T(Tuple.T(3, String.T())) 

692 

693 @property 

694 def stage_type(self): 

695 if self.input_quantity in observational_quantities \ 

696 and self.output_quantity in observational_quantities: 

697 return 'conversion' 

698 

699 if self.input_quantity in observational_quantities \ 

700 and self.output_quantity == 'voltage': 

701 return 'sensor' 

702 

703 elif self.input_quantity == 'voltage' \ 

704 and self.output_quantity == 'voltage': 

705 return 'electronics' 

706 

707 elif self.input_quantity == 'voltage' \ 

708 and self.output_quantity == 'counts': 

709 return 'digitizer' 

710 

711 elif self.input_quantity == 'counts' \ 

712 and self.output_quantity == 'counts' \ 

713 and self.input_sample_rate != self.output_sample_rate: 

714 return 'decimation' 

715 

716 elif self.input_quantity in observational_quantities \ 

717 and self.output_quantity == 'counts': 

718 return 'combination' 

719 

720 else: 

721 return 'unknown' 

722 

723 @property 

724 def summary(self): 

725 irate = self.input_sample_rate 

726 orate = self.output_sample_rate 

727 factor = None 

728 if irate and orate: 

729 factor = irate / orate 

730 return 'ResponseStage, ' + ( 

731 '%s%s => %s%s%s' % ( 

732 self.input_quantity or '?', 

733 ' @ %g Hz' % irate if irate else '', 

734 self.output_quantity or '?', 

735 ' @ %g Hz' % orate if orate else '', 

736 ' :%g' % factor if factor else '') 

737 ) 

738 

739 def get_effective(self): 

740 return MultiplyResponse(responses=list(self.elements)) 

741 

742 

743D = 'displacement' 

744V = 'velocity' 

745A = 'acceleration' 

746 

747g_converters = { 

748 (V, D): IntegrationResponse(1), 

749 (A, D): IntegrationResponse(2), 

750 (D, V): DifferentiationResponse(1), 

751 (A, V): IntegrationResponse(1), 

752 (D, A): DifferentiationResponse(2), 

753 (V, A): DifferentiationResponse(1)} 

754 

755 

756def response_converters(input_quantity, output_quantity): 

757 if input_quantity is None or input_quantity == output_quantity: 

758 return [] 

759 

760 if output_quantity is None: 

761 raise ConversionError('Unspecified target quantity.') 

762 

763 try: 

764 return [g_converters[input_quantity, output_quantity]] 

765 

766 except KeyError: 

767 raise ConversionError('No rule to convert from "%s" to "%s".' % ( 

768 input_quantity, output_quantity)) 

769 

770 

771@squirrel_content 

772class Response(Object): 

773 ''' 

774 The instrument response of a seismic station channel. 

775 ''' 

776 

777 codes = CodesNSLCE.T() 

778 tmin = Timestamp.T(optional=True) 

779 tmax = Timestamp.T(optional=True) 

780 

781 stages = List.T(ResponseStage.T()) 

782 checkpoints = List.T(FrequencyResponseCheckpoint.T()) 

783 

784 deltat = Float.T(optional=True) 

785 log = List.T(Tuple.T(3, String.T())) 

786 

787 def __init__(self, **kwargs): 

788 kwargs['codes'] = CodesNSLCE(kwargs['codes']) 

789 Object.__init__(self, **kwargs) 

790 

791 @property 

792 def time_span(self): 

793 return (self.tmin, self.tmax) 

794 

795 @property 

796 def nstages(self): 

797 return len(self.stages) 

798 

799 @property 

800 def input_quantity(self): 

801 return self.stages[0].input_quantity if self.stages else None 

802 

803 @property 

804 def output_quantity(self): 

805 return self.stages[-1].input_quantity if self.stages else None 

806 

807 @property 

808 def output_sample_rate(self): 

809 return self.stages[-1].output_sample_rate if self.stages else None 

810 

811 @property 

812 def stages_summary(self): 

813 def grouped(xs): 

814 xs = list(xs) 

815 g = [] 

816 for i in range(len(xs)): 

817 g.append(xs[i]) 

818 if i+1 < len(xs) and xs[i+1] != xs[i]: 

819 yield g 

820 g = [] 

821 

822 if g: 

823 yield g 

824 

825 return '+'.join( 

826 '%s%s' % (g[0], '(%i)' % len(g) if len(g) > 1 else '') 

827 for g in grouped(stage.stage_type for stage in self.stages)) 

828 

829 @property 

830 def summary(self): 

831 orate = self.output_sample_rate 

832 return '%s %-16s %s' % ( 

833 self.__class__.__name__, self.str_codes, self.str_time_span) \ 

834 + ', ' + ', '.join(( 

835 '%s => %s' % ( 

836 self.input_quantity or '?', self.output_quantity or '?') 

837 + (' @ %g Hz' % orate if orate else ''), 

838 self.stages_summary, 

839 )) 

840 

841 def get_effective(self, input_quantity=None): 

842 elements = response_converters(input_quantity, self.input_quantity) 

843 

844 elements.extend( 

845 stage.get_effective() for stage in self.stages) 

846 

847 return MultiplyResponse(responses=simplify_responses(elements)) 

848 

849 

850@squirrel_content 

851class Event(Object): 

852 ''' 

853 A seismic event. 

854 ''' 

855 

856 name = String.T(optional=True) 

857 time = Timestamp.T() 

858 duration = Float.T(optional=True) 

859 

860 lat = Float.T() 

861 lon = Float.T() 

862 depth = Float.T(optional=True) 

863 

864 magnitude = Float.T(optional=True) 

865 

866 def get_pyrocko_event(self): 

867 from pyrocko import model 

868 model.Event( 

869 name=self.name, 

870 time=self.time, 

871 lat=self.lat, 

872 lon=self.lon, 

873 depth=self.depth, 

874 magnitude=self.magnitude, 

875 duration=self.duration) 

876 

877 @property 

878 def time_span(self): 

879 return (self.time, self.time) 

880 

881 

882def ehash(s): 

883 return hashlib.sha1(s.encode('utf8')).hexdigest() 

884 

885 

886def random_name(n=8): 

887 return urlsafe_b64encode(urandom(n)).rstrip(b'=').decode('ascii') 

888 

889 

890@squirrel_content 

891class WaveformPromise(Object): 

892 ''' 

893 Information about a waveform potentially downloadable from a remote site. 

894 

895 In the Squirrel framework, waveform promises are used to permit download of 

896 selected waveforms from a remote site. They are typically generated by 

897 calls to 

898 :py:meth:`~pyrocko.squirrel.base.Squirrel.update_waveform_promises`. 

899 Waveform promises are inserted and indexed in the database similar to 

900 normal waveforms. When processing a waveform query, e.g. from 

901 :py:meth:`~pyrocko.squirrel.base.Squirrel.get_waveform`, and no local 

902 waveform is available for the queried time span, a matching promise can be 

903 resolved, i.e. an attempt is made to download the waveform from the remote 

904 site. The promise is removed after the download attempt (except when a 

905 network error occurs). This prevents Squirrel from making unnecessary 

906 queries for waveforms missing at the remote site. 

907 ''' 

908 

909 codes = CodesNSLCE.T() 

910 tmin = Timestamp.T() 

911 tmax = Timestamp.T() 

912 

913 deltat = Float.T(optional=True) 

914 

915 source_hash = String.T() 

916 

917 def __init__(self, **kwargs): 

918 kwargs['codes'] = CodesNSLCE(kwargs['codes']) 

919 Object.__init__(self, **kwargs) 

920 

921 @property 

922 def time_span(self): 

923 return (self.tmin, self.tmax) 

924 

925 

926class InvalidWaveform(Exception): 

927 pass 

928 

929 

930class WaveformOrder(Object): 

931 ''' 

932 Waveform request information. 

933 ''' 

934 

935 source_id = String.T() 

936 codes = CodesNSLCE.T() 

937 deltat = Float.T() 

938 tmin = Timestamp.T() 

939 tmax = Timestamp.T() 

940 gaps = List.T(Tuple.T(2, Timestamp.T())) 

941 

942 @property 

943 def client(self): 

944 return self.source_id.split(':')[1] 

945 

946 def describe(self, site='?'): 

947 return '%s:%s %s [%s - %s]' % ( 

948 self.client, site, str(self.codes), 

949 util.time_to_str(self.tmin), util.time_to_str(self.tmax)) 

950 

951 def validate(self, tr): 

952 if tr.ydata.size == 0: 

953 raise InvalidWaveform( 

954 'waveform with zero data samples') 

955 

956 if tr.deltat != self.deltat: 

957 raise InvalidWaveform( 

958 'incorrect sampling interval - waveform: %g s, ' 

959 'meta-data: %g s' % ( 

960 tr.deltat, self.deltat)) 

961 

962 if not num.all(num.isfinite(tr.ydata)): 

963 raise InvalidWaveform('waveform has NaN values') 

964 

965 

966def order_summary(orders): 

967 codes_list = sorted(set(order.codes for order in orders)) 

968 if len(codes_list) > 3: 

969 return '%i order%s: %s - %s' % ( 

970 len(orders), 

971 util.plural_s(orders), 

972 str(codes_list[0]), 

973 str(codes_list[1])) 

974 

975 else: 

976 return '%i order%s: %s' % ( 

977 len(orders), 

978 util.plural_s(orders), 

979 ', '.join(str(codes) for codes in codes_list)) 

980 

981 

982class Nut(Object): 

983 ''' 

984 Index entry referencing an elementary piece of content. 

985 

986 So-called *nuts* are used in Pyrocko's Squirrel framework to hold common 

987 meta-information about individual pieces of waveforms, stations, channels, 

988 etc. together with the information where it was found or generated. 

989 ''' 

990 

991 file_path = String.T(optional=True) 

992 file_format = String.T(optional=True) 

993 file_mtime = Timestamp.T(optional=True) 

994 file_size = Int.T(optional=True) 

995 

996 file_segment = Int.T(optional=True) 

997 file_element = Int.T(optional=True) 

998 

999 kind_id = Int.T() 

1000 codes = Codes.T() 

1001 

1002 tmin_seconds = Int.T(default=0) 

1003 tmin_offset = Int.T(default=0, optional=True) 

1004 tmax_seconds = Int.T(default=0) 

1005 tmax_offset = Int.T(default=0, optional=True) 

1006 

1007 deltat = Float.T(default=0.0) 

1008 

1009 content = Any.T(optional=True) 

1010 

1011 content_in_db = False 

1012 

1013 def __init__( 

1014 self, 

1015 file_path=None, 

1016 file_format=None, 

1017 file_mtime=None, 

1018 file_size=None, 

1019 file_segment=None, 

1020 file_element=None, 

1021 kind_id=0, 

1022 codes=CodesX(''), 

1023 tmin_seconds=None, 

1024 tmin_offset=0, 

1025 tmax_seconds=None, 

1026 tmax_offset=0, 

1027 deltat=None, 

1028 content=None, 

1029 tmin=None, 

1030 tmax=None, 

1031 values_nocheck=None): 

1032 

1033 if values_nocheck is not None: 

1034 (self.file_path, self.file_format, self.file_mtime, 

1035 self.file_size, 

1036 self.file_segment, self.file_element, 

1037 self.kind_id, codes_safe_str, 

1038 self.tmin_seconds, self.tmin_offset, 

1039 self.tmax_seconds, self.tmax_offset, 

1040 self.deltat) = values_nocheck 

1041 

1042 self.codes = to_codes_simple(self.kind_id, codes_safe_str) 

1043 self.content = None 

1044 else: 

1045 if tmin is not None: 

1046 tmin_seconds, tmin_offset = tsplit(tmin) 

1047 

1048 if tmax is not None: 

1049 tmax_seconds, tmax_offset = tsplit(tmax) 

1050 

1051 self.kind_id = int(kind_id) 

1052 self.codes = codes 

1053 self.tmin_seconds = int_or_g_tmin(tmin_seconds) 

1054 self.tmin_offset = int(tmin_offset) 

1055 self.tmax_seconds = int_or_g_tmax(tmax_seconds) 

1056 self.tmax_offset = int(tmax_offset) 

1057 self.deltat = float_or_none(deltat) 

1058 self.file_path = str_or_none(file_path) 

1059 self.file_segment = int_or_none(file_segment) 

1060 self.file_element = int_or_none(file_element) 

1061 self.file_format = str_or_none(file_format) 

1062 self.file_mtime = float_or_none(file_mtime) 

1063 self.file_size = int_or_none(file_size) 

1064 self.content = content 

1065 

1066 Object.__init__(self, init_props=False) 

1067 

1068 def __eq__(self, other): 

1069 return (isinstance(other, Nut) and 

1070 self.equality_values == other.equality_values) 

1071 

1072 def hash(self): 

1073 return ehash(','.join(str(x) for x in self.key)) 

1074 

1075 def __ne__(self, other): 

1076 return not (self == other) 

1077 

1078 def get_io_backend(self): 

1079 from . import io 

1080 return io.get_backend(self.file_format) 

1081 

1082 def file_modified(self): 

1083 return self.get_io_backend().get_stats(self.file_path) \ 

1084 != (self.file_mtime, self.file_size) 

1085 

1086 @property 

1087 def dkey(self): 

1088 return (self.kind_id, self.tmin_seconds, self.tmin_offset, self.codes) 

1089 

1090 @property 

1091 def key(self): 

1092 return ( 

1093 self.file_path, 

1094 self.file_segment, 

1095 self.file_element, 

1096 self.file_mtime) 

1097 

1098 @property 

1099 def equality_values(self): 

1100 return ( 

1101 self.file_segment, self.file_element, 

1102 self.kind_id, self.codes, 

1103 self.tmin_seconds, self.tmin_offset, 

1104 self.tmax_seconds, self.tmax_offset, self.deltat) 

1105 

1106 @property 

1107 def tmin(self): 

1108 return tjoin(self.tmin_seconds, self.tmin_offset) 

1109 

1110 @tmin.setter 

1111 def tmin(self, t): 

1112 self.tmin_seconds, self.tmin_offset = tsplit(t) 

1113 

1114 @property 

1115 def tmax(self): 

1116 return tjoin(self.tmax_seconds, self.tmax_offset) 

1117 

1118 @tmax.setter 

1119 def tmax(self, t): 

1120 self.tmax_seconds, self.tmax_offset = tsplit(t) 

1121 

1122 @property 

1123 def kscale(self): 

1124 if self.tmin_seconds is None or self.tmax_seconds is None: 

1125 return 0 

1126 return tscale_to_kscale(self.tmax_seconds - self.tmin_seconds) 

1127 

1128 @property 

1129 def waveform_kwargs(self): 

1130 network, station, location, channel, extra = self.codes 

1131 

1132 return dict( 

1133 network=network, 

1134 station=station, 

1135 location=location, 

1136 channel=channel, 

1137 extra=extra, 

1138 tmin=self.tmin, 

1139 tmax=self.tmax, 

1140 deltat=self.deltat) 

1141 

1142 @property 

1143 def waveform_promise_kwargs(self): 

1144 return dict( 

1145 codes=self.codes, 

1146 tmin=self.tmin, 

1147 tmax=self.tmax, 

1148 deltat=self.deltat) 

1149 

1150 @property 

1151 def station_kwargs(self): 

1152 network, station, location = self.codes 

1153 return dict( 

1154 codes=self.codes, 

1155 tmin=tmin_or_none(self.tmin), 

1156 tmax=tmax_or_none(self.tmax)) 

1157 

1158 @property 

1159 def channel_kwargs(self): 

1160 network, station, location, channel, extra = self.codes 

1161 return dict( 

1162 codes=self.codes, 

1163 tmin=tmin_or_none(self.tmin), 

1164 tmax=tmax_or_none(self.tmax), 

1165 deltat=self.deltat) 

1166 

1167 @property 

1168 def response_kwargs(self): 

1169 return dict( 

1170 codes=self.codes, 

1171 tmin=tmin_or_none(self.tmin), 

1172 tmax=tmax_or_none(self.tmax), 

1173 deltat=self.deltat) 

1174 

1175 @property 

1176 def event_kwargs(self): 

1177 return dict( 

1178 name=self.codes, 

1179 time=self.tmin, 

1180 duration=(self.tmax - self.tmin) or None) 

1181 

1182 @property 

1183 def trace_kwargs(self): 

1184 network, station, location, channel, extra = self.codes 

1185 

1186 return dict( 

1187 network=network, 

1188 station=station, 

1189 location=location, 

1190 channel=channel, 

1191 extra=extra, 

1192 tmin=self.tmin, 

1193 tmax=self.tmax-self.deltat, 

1194 deltat=self.deltat) 

1195 

1196 @property 

1197 def dummy_trace(self): 

1198 return DummyTrace(self) 

1199 

1200 @property 

1201 def summary(self): 

1202 if self.tmin == self.tmax: 

1203 ts = util.time_to_str(self.tmin) 

1204 else: 

1205 ts = '%s - %s' % ( 

1206 util.time_to_str(self.tmin), 

1207 util.time_to_str(self.tmax)) 

1208 

1209 return ' '.join(( 

1210 ('%s,' % to_kind(self.kind_id)).ljust(9), 

1211 ('%s,' % str(self.codes)).ljust(18), 

1212 ts)) 

1213 

1214 

1215def make_waveform_nut(**kwargs): 

1216 return Nut(kind_id=WAVEFORM, **kwargs) 

1217 

1218 

1219def make_waveform_promise_nut(**kwargs): 

1220 return Nut(kind_id=WAVEFORM_PROMISE, **kwargs) 

1221 

1222 

1223def make_station_nut(**kwargs): 

1224 return Nut(kind_id=STATION, **kwargs) 

1225 

1226 

1227def make_channel_nut(**kwargs): 

1228 return Nut(kind_id=CHANNEL, **kwargs) 

1229 

1230 

1231def make_response_nut(**kwargs): 

1232 return Nut(kind_id=RESPONSE, **kwargs) 

1233 

1234 

1235def make_event_nut(**kwargs): 

1236 return Nut(kind_id=EVENT, **kwargs) 

1237 

1238 

1239def group_channels(nuts): 

1240 by_ansl = {} 

1241 for nut in nuts: 

1242 if nut.kind_id != CHANNEL: 

1243 continue 

1244 

1245 ansl = nut.codes[:4] 

1246 

1247 if ansl not in by_ansl: 

1248 by_ansl[ansl] = {} 

1249 

1250 group = by_ansl[ansl] 

1251 

1252 k = nut.codes[4][:-1], nut.deltat, nut.tmin, nut.tmax 

1253 

1254 if k not in group: 

1255 group[k] = set() 

1256 

1257 group.add(nut.codes[4]) 

1258 

1259 return by_ansl 

1260 

1261 

1262class DummyTrace(object): 

1263 

1264 def __init__(self, nut): 

1265 self.nut = nut 

1266 self.codes = nut.codes 

1267 self.meta = {} 

1268 

1269 @property 

1270 def tmin(self): 

1271 return self.nut.tmin 

1272 

1273 @property 

1274 def tmax(self): 

1275 return self.nut.tmax 

1276 

1277 @property 

1278 def deltat(self): 

1279 return self.nut.deltat 

1280 

1281 @property 

1282 def nslc_id(self): 

1283 return self.codes.nslc 

1284 

1285 @property 

1286 def network(self): 

1287 return self.codes.network 

1288 

1289 @property 

1290 def station(self): 

1291 return self.codes.station 

1292 

1293 @property 

1294 def location(self): 

1295 return self.codes.location 

1296 

1297 @property 

1298 def channel(self): 

1299 return self.codes.channel 

1300 

1301 @property 

1302 def extra(self): 

1303 return self.codes.extra 

1304 

1305 def overlaps(self, tmin, tmax): 

1306 return not (tmax < self.nut.tmin or self.nut.tmax < tmin) 

1307 

1308 

1309def duration_to_str(t): 

1310 if t > 24*3600: 

1311 return '%gd' % (t / (24.*3600.)) 

1312 elif t > 3600: 

1313 return '%gh' % (t / 3600.) 

1314 elif t > 60: 

1315 return '%gm' % (t / 60.) 

1316 else: 

1317 return '%gs' % t 

1318 

1319 

1320class Coverage(Object): 

1321 ''' 

1322 Information about times covered by a waveform or other time series data. 

1323 ''' 

1324 kind_id = Int.T( 

1325 help='Content type.') 

1326 pattern = Codes.T( 

1327 help='The codes pattern in the request, which caused this entry to ' 

1328 'match.') 

1329 codes = Codes.T( 

1330 help='NSLCE or NSL codes identifier of the time series.') 

1331 deltat = Float.T( 

1332 help='Sampling interval [s]', 

1333 optional=True) 

1334 tmin = Timestamp.T( 

1335 help='Global start time of time series.', 

1336 optional=True) 

1337 tmax = Timestamp.T( 

1338 help='Global end time of time series.', 

1339 optional=True) 

1340 changes = List.T( 

1341 Tuple.T(2, Any.T()), 

1342 help='List of change points, with entries of the form ' 

1343 '``(time, count)``, where a ``count`` of zero indicates start of ' 

1344 'a gap, a value of 1 start of normal data coverage and a higher ' 

1345 'value duplicate or redundant data coverage.') 

1346 

1347 @classmethod 

1348 def from_values(cls, args): 

1349 pattern, codes, deltat, tmin, tmax, changes, kind_id = args 

1350 return cls( 

1351 kind_id=kind_id, 

1352 pattern=pattern, 

1353 codes=codes, 

1354 deltat=deltat, 

1355 tmin=tmin, 

1356 tmax=tmax, 

1357 changes=changes) 

1358 

1359 @property 

1360 def summary(self): 

1361 ts = '%s - %s,' % ( 

1362 util.time_to_str(self.tmin), 

1363 util.time_to_str(self.tmax)) 

1364 

1365 srate = self.sample_rate 

1366 

1367 return ' '.join(( 

1368 ('%s,' % to_kind(self.kind_id)).ljust(9), 

1369 ('%s,' % str(self.codes)).ljust(18), 

1370 ts, 

1371 '%10.3g,' % srate if srate else '', 

1372 '%4i' % len(self.changes), 

1373 '%s' % duration_to_str(self.total))) 

1374 

1375 @property 

1376 def sample_rate(self): 

1377 if self.deltat is None: 

1378 return None 

1379 elif self.deltat == 0.0: 

1380 return 0.0 

1381 else: 

1382 return 1.0 / self.deltat 

1383 

1384 @property 

1385 def labels(self): 

1386 srate = self.sample_rate 

1387 return ( 

1388 ('%s' % str(self.codes)), 

1389 '%.3g' % srate if srate else '') 

1390 

1391 @property 

1392 def total(self): 

1393 total_t = None 

1394 for tmin, tmax, _ in self.iter_spans(): 

1395 total_t = (total_t or 0.0) + (tmax - tmin) 

1396 

1397 return total_t 

1398 

1399 def iter_spans(self): 

1400 last = None 

1401 for (t, count) in self.changes: 

1402 if last is not None: 

1403 last_t, last_count = last 

1404 if last_count > 0: 

1405 yield last_t, t, last_count 

1406 

1407 last = (t, count) 

1408 

1409 

1410__all__ = [ 

1411 'to_codes', 

1412 'to_codes_guess', 

1413 'to_codes_simple', 

1414 'to_kind', 

1415 'to_kinds', 

1416 'to_kind_id', 

1417 'to_kind_ids', 

1418 'CodesError', 

1419 'Codes', 

1420 'CodesNSLCE', 

1421 'CodesNSL', 

1422 'CodesX', 

1423 'Station', 

1424 'Channel', 

1425 'Sensor', 

1426 'Response', 

1427 'Nut', 

1428 'Coverage', 

1429 'WaveformPromise', 

1430]