Coverage for C: \ Users \ peaco \ OneDrive \ Documents \ GitHub \ mt_metadata \ mt_metadata \ timeseries \ filters \ obspy_stages.py: 72%

1""" 

2Idea here is to add logic to interrogate stage filters received from StationXML 

3 

4""" 

5 

6# ============================================================================= 

7# Imports 

8# ============================================================================= 

9import numpy as np 

10 

11 

12try: 

13 import obspy 

14except ImportError: 

15 raise ImportError("obspy_stages requires obspy to be installed.") 

16from loguru import logger 

17 

18from mt_metadata.timeseries.filters import ( 

19 CoefficientFilter, 

20 FIRFilter, 

21 FrequencyResponseTableFilter, 

22 PoleZeroFilter, 

23 TimeDelayFilter, 

24) 

25 

26 

27# ============================================================================= 

28 

29 

30def create_time_delay_filter_from_stage(stage): 

31 time_delay_filter = TimeDelayFilter() 

32 time_delay_filter = time_delay_filter.from_obspy_stage(stage) 

33 return time_delay_filter 

34 

35 

36def create_coefficent_filter_from_stage(stage): 

37 coeff_filter = CoefficientFilter() 

38 coeff_filter = coeff_filter.from_obspy_stage(stage) 

39 return coeff_filter 

40 

41 

42def create_pole_zero_filter_from_stage(stage): 

43 pz_filter = PoleZeroFilter() 

44 pz_filter = pz_filter.from_obspy_stage(stage) 

45 return pz_filter 

46 

47 

48def create_fir_filter_from_stage(stage): 

49 fir_filter = FIRFilter() 

50 fir_filter = fir_filter.from_obspy_stage(stage) 

51 return fir_filter 

52 

53 

54def create_frequency_response_table_filter_from_stage(stage): 

55 """ 

56 Notes: Issue here is that the stage has a list of response_list_elements, it does not actually 

57 have attrs frequencies, amplitudes, phases... will try assigning them on the fly here 

58 

59 Parameters 

60 ---------- 

61 stage 

62 

63 Returns 

64 ------- 

65 

66 """ 

67 n_freq = len(stage.response_list_elements) 

68 frequencies = np.full(n_freq, np.nan) 

69 amplitudes = np.full(n_freq, np.nan) 

70 phases = np.full(n_freq, np.nan) 

71 for i, response_list_element in enumerate(stage.response_list_elements): 

72 frequencies[i] = response_list_element.frequency 

73 amplitudes[i] = response_list_element.amplitude 

74 phases[i] = response_list_element.phase 

75 stage.frequencies = frequencies 

76 stage.amplitudes = amplitudes 

77 stage.phases = phases 

78 

79 fap_filter = FrequencyResponseTableFilter() 

80 fap_filter = fap_filter.from_obspy_stage(stage) 

81 return fap_filter 

82 

83 

84def check_if_coefficient_filter_is_delay_only(stage): 

85 """ 

86 stage: obspy_stage type in obspy.core.inventory.response 

87 This function may wind up being a method of the CoefficientFilter class, but leaving it 

88 separate for now. 

89 Conditions to check: 

90 1. Gain is unity 

91 2. Decimation Factor is 1 

92 3. delay is non-zero 

93 4. delay has not been corrected on the data-center side. 

94 If all four of these are true we assume this is a pathological CoefficienFilter 

95 and would actually be a TimeDelay() filter instance if one existed in StationXML. 

96 Returns 

97 ------- 

98 

99 """ 

100 cond1 = stage.stage_gain == 1.0 

101 cond2 = stage.decimation_factor == 1 

102 cond3 = stage.decimation_delay != 0.0 

103 cond4 = stage.decimation_correction == 0.0 

104 

105 if cond1 & cond2 & cond3 & cond4: 

106 return True 

107 else: 

108 return False 

109 

110 

111def create_filter_from_stage(stage): 

112 """ 

113 This works on a single stage, we need a catalog of stage classes 

114 obspy.core.inventory.response.CoefficientsTypeResponseStage 

115 obspy.core.inventory.response.FIRResponseStage 

116 obspy.core.inventory.response.PolesZerosResponseStage 

117 

118 CoefficientsTypeResponseStage: cases include 

119 -numerator=[],denominator=[] 

120 

121 Sometimes filter stages in obspy are used to kluge-represent filters of other types 

122 # Encountered Cases This Far: 

123 # CoefficientTypeResponseStage Used to package a time-delay filter 

124 

125 Parameters 

126 ---------- 

127 stage 

128 

129 Returns 

130 ------- 

131 

132 """ 

133 

134 if isinstance(stage, obspy.core.inventory.response.PolesZerosResponseStage): 

135 if stage.poles == [] and stage.zeros == []: 

136 if ( 

137 "counts" not in stage.input_units.lower() 

138 and "counts" not in stage.output_units.lower() 

139 ): 

140 logger.info( 

141 f"Converting PoleZerosResponseStage {stage.name} to a " 

142 "CoefficientFilter." 

143 ) 

144 return create_coefficent_filter_from_stage(stage) 

145 

146 return create_pole_zero_filter_from_stage(stage) 

147 

148 elif isinstance(stage, obspy.core.inventory.response.CoefficientsTypeResponseStage): 

149 is_a_delay_filter = check_if_coefficient_filter_is_delay_only(stage) 

150 if is_a_delay_filter: 

151 obspy_filter = create_time_delay_filter_from_stage(stage) 

152 else: 

153 obspy_filter = create_coefficent_filter_from_stage(stage) 

154 return obspy_filter 

155 

156 elif isinstance(stage, obspy.core.inventory.response.FIRResponseStage): 

157 try: 

158 if isinstance(stage.coefficients, list): 

159 pass 

160 else: 

161 msg = f"Expected list of coefficients, got {type(stage.coefficients)}" 

162 logger.error(msg) 

163 raise TypeError(msg) 

164 except TypeError: 

165 msg = "Something seems off with this FIR" 

166 logger.info(msg) 

167 raise ValueError(msg) 

168 obspy_filter = create_fir_filter_from_stage(stage) 

169 return obspy_filter 

170 

171 elif isinstance(stage, obspy.core.inventory.response.ResponseListResponseStage): 

172 obspy_filter = create_frequency_response_table_filter_from_stage(stage) 

173 return obspy_filter 

174 

175 elif isinstance(stage, obspy.core.inventory.response.ResponseStage): 

176 obspy_filter = create_coefficent_filter_from_stage(stage) 

177 return obspy_filter 

178 

179 else: 

180 msg = "Filter Stage of type {type(stage)} not known, or supported" 

181 logger.info(msg) 

182 raise TypeError(msg)