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

# http://pyrocko.org - GPLv3 

# 

# The Pyrocko Developers, 21st Century 

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

from __future__ import absolute_import, division, print_function 

 

import logging 

import os.path as op 

import numpy as num 

 

from pyrocko import gf, util 

from pyrocko.guts import Float, List 

 

from .base import TargetGenerator, NoiseGenerator 

from ..station import RandomStationGenerator, StationGenerator 

 

DEFAULT_STORE_ID = 'ak135_static' 

 

logger = logging.getLogger('pyrocko.scenario.targets.gnss_campaign') 

guts_prefix = 'pf.scenario' 

 

 

class GPSNoiseGenerator(NoiseGenerator): 

measurement_duarion_days = Float.T( 

default=2., 

help='Measurement duration in days') 

 

def add_noise(self, campaign): 

# https://www.nat-hazards-earth-syst-sci.net/15/875/2015/nhess-15-875-2015.pdf 

waterlevel = 1. - (.99 + .0015 * self.measurement_duarion_days) # noqa 

logger.warning('GPSNoiseGenerator is a work-in-progress!') 

 

for ista, sta in enumerate(campaign.stations): 

pass 

# rstate = self.get_rstate(ista) 

 

# sta.north.sigma = 8e-3 

# sta.east.sigma = 8e-3 

 

# sta.north.shift += rstate.normal(0., sta.north.sigma) 

# sta.east.shift += rstate.normal(0., sta.east.sigma) 

 

 

class GNSSCampaignGenerator(TargetGenerator): 

station_generators = List.T( 

StationGenerator.T(), 

default=[RandomStationGenerator.D( 

network_name='GN', 

channels=None)], 

help='The StationGenerator.') 

 

noise_generator = NoiseGenerator.T( 

default=GPSNoiseGenerator.D(), 

optional=True, 

help='Add Synthetic noise to the GNSS displacements.') 

 

store_id = gf.StringID.T( 

default=DEFAULT_STORE_ID, 

help='The GF store to use for forward-calculations.') 

 

def get_stations(self): 

stations = [] 

for station_generator in self.station_generators: 

stations.extend(station_generator.get_stations()) 

return stations 

 

def get_targets(self): 

stations = self.get_stations() 

lats = num.array([s.effective_lat for s in stations]) 

lons = num.array([s.effective_lon for s in stations]) 

 

target = gf.GNSSCampaignTarget( 

lats=lats, 

lons=lons, 

store_id=self.store_id) 

 

return [target] 

 

def get_gnss_campaigns(self, engine, sources, tmin=None, tmax=None): 

try: 

resp = engine.process( 

sources, 

self.get_targets(), 

nthreads=0) 

except gf.meta.OutOfBounds: 

logger.warning('Could not calculate GNSS displacements' 

' - the GF store\'s extend is too small!') 

return [] 

 

campaigns = [r.campaign for r in resp.static_results()] 

 

stacked_campaign = campaigns[0] 

stacked_campaign.name = 'Scenario Campaign' 

for camp in campaigns[1:]: 

for ista, sta in enumerate(camp.stations): 

stacked_campaign.stations[ista].north.shift += sta.north.shift 

stacked_campaign.stations[ista].east.shift += sta.east.shift 

stacked_campaign.stations[ista].up.shift += sta.up.shift 

 

for ista, sta in enumerate(stacked_campaign.stations): 

sta.code = 'SY%02d' % (ista + 1) 

 

if self.noise_generator is not None: 

self.noise_generator.add_noise(stacked_campaign) 

 

return [stacked_campaign] 

 

def ensure_data(self, engine, sources, path, tmin=None, tmax=None): 

path_gnss = op.join(path, 'gnss') 

util.ensuredir(path_gnss) 

 

networks = [] 

for sg in self.station_generators: 

try: 

networks.append(sg.network_name) 

except AttributeError: 

pass 

 

fn = op.join( 

path_gnss, 

'campaign-%s.yml' % '_'.join(networks)) 

 

if op.exists(fn): 

return 

 

campaigns = self.get_gnss_campaigns(engine, sources, tmin, tmax) 

 

with open(fn, 'w') as f: 

for camp in campaigns: 

camp.dump(stream=f) 

 

def add_map_artists(self, engine, sources, automap): 

automap.add_gnss_campaign(self.get_gnss_campaigns(engine, sources)[0])