1# http://pyrocko.org - GPLv3 

2# 

3# The Pyrocko Developers, 21st Century 

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

5from __future__ import absolute_import, division, print_function 

6 

7import logging 

8import os.path as op 

9import numpy as num 

10 

11from pyrocko import gf, util 

12from pyrocko.guts import Float, List 

13 

14from .base import TargetGenerator, NoiseGenerator 

15from ..station import RandomStationGenerator, StationGenerator 

16 

17DEFAULT_STORE_ID = 'ak135_static' 

18 

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

20guts_prefix = 'pf.scenario' 

21 

22 

23class GPSNoiseGenerator(NoiseGenerator): 

24 measurement_duarion_days = Float.T( 

25 default=2., 

26 help='Measurement duration in days') 

27 

28 def add_noise(self, campaign): 

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

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

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

32 

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

34 pass 

35 # rstate = self.get_rstate(ista) 

36 

37 # sta.north.sigma = 8e-3 

38 # sta.east.sigma = 8e-3 

39 

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

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

42 

43 

44class GNSSCampaignGenerator(TargetGenerator): 

45 station_generators = List.T( 

46 StationGenerator.T(), 

47 default=[RandomStationGenerator.D( 

48 network_name='GN', 

49 channels=None)], 

50 help='The StationGenerator.') 

51 

52 noise_generator = NoiseGenerator.T( 

53 default=GPSNoiseGenerator.D(), 

54 optional=True, 

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

56 

57 store_id = gf.StringID.T( 

58 default=DEFAULT_STORE_ID, 

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

60 

61 def get_stations(self): 

62 stations = [] 

63 for station_generator in self.station_generators: 

64 stations.extend(station_generator.get_stations()) 

65 return stations 

66 

67 def get_targets(self): 

68 stations = self.get_stations() 

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

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

71 

72 target = gf.GNSSCampaignTarget( 

73 lats=lats, 

74 lons=lons, 

75 store_id=self.store_id) 

76 

77 return [target] 

78 

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

80 try: 

81 resp = engine.process( 

82 sources, 

83 self.get_targets(), 

84 nthreads=0) 

85 except gf.meta.OutOfBounds: 

86 logger.warning('Could not calculate GNSS displacements' 

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

88 return [] 

89 

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

91 

92 stacked_campaign = campaigns[0] 

93 stacked_campaign.name = 'Scenario Campaign' 

94 for camp in campaigns[1:]: 

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

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

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

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

99 

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

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

102 

103 if self.noise_generator is not None: 

104 self.noise_generator.add_noise(stacked_campaign) 

105 

106 return [stacked_campaign] 

107 

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

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

110 util.ensuredir(path_gnss) 

111 

112 networks = [] 

113 for sg in self.station_generators: 

114 try: 

115 networks.append(sg.network_name) 

116 except AttributeError: 

117 pass 

118 

119 fn = op.join( 

120 path_gnss, 

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

122 

123 if op.exists(fn): 

124 return 

125 

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

127 

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

129 for camp in campaigns: 

130 camp.dump(stream=f) 

131 

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

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