Coverage for /usr/local/lib/python3.11/dist-packages/pyrocko/gf/fd_rotation.py: 0%

1 

2from pyrocko import gf, trace 

3 

4 

5comp_azi = { 

6 'N': 0.0, 

7 'E': 90.0, 

8 'D': 0.0} 

9 

10comp_dip = { 

11 'N': 0.0, 

12 'E': 0.0, 

13 'D': 90.0} 

14 

15 

16def fd_rot_surface_targets(delta, target): 

17 

18 assert target.interpolation == 'multilinear' 

19 assert target.quantity == 'displacement' 

20 

21 targets = [] 

22 for kn in [-1, 0, 1]: 

23 for ke in [-1, 0, 1]: 

24 for comp in 'NED': 

25 targets.append(gf.Target( 

26 codes=target.codes[:3] + (comp,), 

27 quantity='displacement', 

28 lat=target.lat, 

29 lon=target.lon, 

30 north_shift=target.north_shift + kn * delta, 

31 east_shift=target.east_shift + ke * delta, 

32 azimuth=comp_azi[comp], 

33 dip=comp_dip[comp], 

34 interpolation='multilinear', 

35 store_id=target.store_id)) 

36 

37 return targets 

38 

39 

40def fd_rot_surface_postprocess(delta, traces, rfmax=0.8): 

41 

42 fnyquist = 0.5 / traces[0].deltat 

43 fmax = rfmax * fnyquist 

44 tpad = 2.0 / fmax 

45 tmin = min(tr.tmin for tr in traces) 

46 tmax = max(tr.tmax for tr in traces) 

47 traces_prepared = [] 

48 

49 for tr in traces: 

50 tr.extend(tmin-tpad, tmax+tpad, fillmethod='repeat') 

51 tr_filt = tr.transfer( 

52 tfade=tpad, 

53 freqlimits=(-1., -1., fmax, fnyquist), 

54 demean=False) 

55 

56 traces_prepared.append(tr_filt) 

57 

58 tr_filt.set_location('filtered') 

59 

60 trace.snuffle(traces + traces_prepared) 

61 

62 d = trace.get_traces_data_as_array(traces_prepared) 

63 d = d.reshape((3, 3, 3, d.shape[1])) 

64 

65 kn, ke, kd = 0, 1, 2 

66 il, ic, ih = 0, 1, 2 

67 

68 dd_de = (d[ic, ih, kd] - d[ic, il, kd]) / (2.*delta) 

69 dd_dn = (d[ih, ic, kd] - d[il, ic, kd]) / (2.*delta) 

70 dn_de = (d[ic, ih, kn] - d[ic, il, kn]) / (2.*delta) 

71 de_dn = (d[ih, ic, ke] - d[il, ic, ke]) / (2.*delta) 

72 

73 d_rot = [ 

74 dd_de, 

75 - dd_dn, 

76 0.5 * (de_dn - dn_de)] 

77 

78 traces_rot = [] 

79 for comp, data in zip('NED', d_rot): 

80 tr = trace.Trace( 

81 network=traces[0].network, 

82 station=traces[0].station, 

83 location=traces[0].location, 

84 channel='ROT_%s' % comp, 

85 tmin=tmin, 

86 ydata=data) 

87 # tr.differentiate(1) 

88 traces_rot.append(tr) 

89 

90 traces_dis, traces_vel, traces_acc = [], [], [] 

91 for comp, data in zip('NED', d[ic, ic, :]): 

92 tr_dis = traces[0].copy() 

93 tr_dis.set_ydata(data) 

94 tr_dis.set_codes(channel='DIS_%s' % comp) 

95 tr_vel = tr_dis.differentiate(1, inplace=False) 

96 tr_vel.set_codes(channel='VEL_%s' % comp) 

97 tr_acc = tr_dis.differentiate(2, inplace=False) 

98 tr_acc.set_codes(channel='ACC_%s' % comp) 

99 traces_dis.append(tr_dis) 

100 traces_vel.append(tr_vel) 

101 traces_acc.append(tr_acc) 

102 

103 return traces_dis, traces_vel, traces_acc, traces_rot