Skip to Content

Wiki source code of FDSN Guide

Version 3.1 by robert on 2025/03/24 10:35
Hide last authors
robert 1.1 1 {{box cssClass="floatinginfobox" title="**Contents**"}}
2 {{toc/}}
3 {{/box}}
4
robert 1.2 5 = [[How to Install ObsPy>>url:https://github.com/obspy/obspy/wiki#installation]] =
6
7 = [[Seed-Vault>>https://github.com/AuScope/seed-vault]] =
8
robert 1.1 9 = Connecting to an FDSN Server =
10
robert 1.2 11 == How to connect to AusPass with & without authenticated access ==
12
13
14 {{code language="python"}}
15 import obspy
16 from obspy.clients.fdsn import Client
17
18 # Initialize FDSN client for AusPass
19
20 # For open access data, no username or password is required.
21 client = Client('AUSPASS')
22
23 # To access restricted data, supply your username and password
24 # Replace 'Z1' and '12345' with your actual credentials
25 client = Client('AUSPASS', user='Z1', password='12345')
26 {{/code}}
27
robert 1.1 28 = Station Metadata =
29
robert 2.2 30 Information such as site locations, sensor and data logger types, response information, etc are in the station metadata. This can be accessed directly(link) or via the obspy get_stations (link) tool.
robert 1.1 31
32
robert 2.4 33 == How to download event, station, instrument response ==
robert 1.1 34
35
robert 2.4 36 {{code language="python"}}
37 import obspy
38 from obspy.clients.fdsn import Client
39
40 # Use AusPass client for station, waveform, and earthquake information
41 client = Client("AUSPASS")
42
43
44 # Download station information for AUMTC station in S1 network at the response level
45 inv = client.get_stations(network="S1", station="AUMTC", location="*",
robert 3.1 46 channel="*", level="response")
robert 2.4 47
48 # Inventory metadata is stored in a Inventory > Network > Station > Channel hierarchy
49
50 print(inv) #inventory level
51
52 print(inv[0]) # network level (the first network in the inventory)
53
54 print(inv[0][0]) # station level (the first station of the first network in the inventory)
55
56 print(inv[0][0][0]) # channel level (the first channel of the first station of the first network in the inventoy)
57
58 # you can also select items directly
59
60 print(inv.select(station='AUMTC',channel='HHZ')[0][0][0])
61
62 # instrument response is attached to a channel object
63
64 response = inv.select(station='AUMTC',channel='HHZ')[0][0][0].response
65 {{/code}}
66
67
robert 1.1 68
69
70 = Waveform Data =
71
robert 2.2 72 Waveform data (e.g. the actual seismic data) can be accessed directly (link) or via obspy's get_waveforms (link) tool. It can also be accessed via various tools such as seed-vault, pyweed, etc (add links).
robert 1.1 73
robert 3.1 74 == Downloading and Storing data ==
robert 1.1 75
robert 3.1 76 === How to download waveform data ===
77
robert 2.3 78 {{code language="python"}}
79 from obspy import UTCDateTime
80 from obspy.clients.fdsn import Client
81
82 # Initialize the FDSN client (you can also specify other data centers)
83 client = Client("AUSPASS")
84
85 # Event information
86 network = "S1"
87 station = "AUGRF"
88 starttime = UTCDateTime("2021-09-21T23:15:53") # The time of the earthquake
89 endtime = starttime + 360 # One hour of data after the earthquake
90
91 # Download the MiniSEED data
92 st = client.get_waveforms(network=network, station=station, location="*", channel="BHZ",
93 starttime=starttime, endtime=endtime)
94 # Save the stream to a MiniSEED file
95 st.write("Woodspoint_2021.mseed", format="MSEED")
96 print("Downloaded and saved the MiniSEED file.")
97 {{/code}}
98
robert 3.1 99 === How to download a LOT of waveform data ===
robert 2.3 100
101 {{code language="python"}}
robert 3.1 102 from obspy import UTCDateTime
103 from obspy.clients.fdsn import Client
104 import datetime
105
106 # Initialize client and set parameters
107 client = Client("IRIS")
108 start_date = UTCDateTime("2023-07-01")
109 end_date = UTCDateTime("2023-07-02")
110 network = "S1"
111 station = "AUMTS"
112 channel = "BHZ"
113
114 # Loop to download data one day at a time
115 while start_date <= end_date:
116 next_date = start_date + datetime.timedelta(days=1)
117 try:
118 st = client.get_waveforms(network, station, "*", channel, start_date, next_date)
119 st.write(f"{start_date.date}.mseed", format="MSEED")
120 except:
121 print(f"Failed for {start_date} to {next_date}")
122 start_date = next_date
123 {{/code}}
124
125 === How to store and archive waveform data in SeisComP Data Structure (SDS) ===
126
127 {{code language="python"}}
128 import os
129 from obspy import UTCDateTime, read
130 from obspy.clients.fdsn import Client
131
132 # Initialize the client
133 client = Client("AUSPASS") # Replace with the correct client endpoint if different
134
135 # Define event time and time window
136 event_time = UTCDateTime("2021-09-21T23:15:53")
137 starttime = event_time - 600 # 10 minutes before the event
138 endtime = event_time + 1800 # 30 minutes after the event
139
140 # Download the waveform data
141 st = client.get_waveforms(network="S1", station="AUMTS", location="*", channel="*", starttime=starttime, endtime=endtime)
142
143 # Create SDS structure: ROOT/YEAR/NET/STA/CHAN.TYPE/NET.STA.LOC.CHAN.YEAR.DAY
144 sds_root = "." # Replace with your desired directory
145
146 for tr in st:
147 net = tr.stats.network
148 sta = tr.stats.station
149 loc = tr.stats.location
150 chan = tr.stats.channel
151 year = str(tr.stats.starttime.year)
152 jday = str(tr.stats.starttime.julday).zfill(3)
153
154 sds_path = os.path.join(sds_root, year, net, sta, f"{chan}.D", f"{net}.{sta}.{loc}.{chan}.{year}.{jday}")
155
156 # Create directories if they don't exist
157 os.makedirs(os.path.dirname(sds_path), exist_ok=True)
158
159 # Save the trace as a MiniSEED file
160 tr.write(sds_path, format="MSEED")
161 {{/code}}
162
163 ==
164 Common Data Operations ==
165
166 === How to remove instrument response ===
167
168 {{code language="python"}}
robert 2.3 169 from obspy import read
170 from obspy.core.util import AttribDict
171
172 # Load the MiniSEED file
173 st = read("Woodspoint_2021.mseed")
174
175 # Download the instrument response
176 inv = client.get_stations(network=network, station=station, location="*",
177 channel="*", starttime=starttime, endtime=endtime,
178 level="response")
179
180 # Remove the instrument response
181 output = 'VEL' # Output unit ('VEL' = velocity (default), 'DISP' = displacement, 'ACC' = acceleration)
182
183 for tr in st:
184 tr.remove_response(inventory=inv, output=output, plot=True)
185
186 # Save the corrected MiniSEED file
187 st.write("Woodspoint_2021_corrected.mseed", format="MSEED")
188 {{/code}}
189
robert 3.1 190 === How to apply a bandpass filter ===
robert 2.3 191
192 {{code language="python"}}
193 from obspy import read
194
195 # Load the MiniSEED file
196 st = read("Woodspoint_2021.mseed")
197
198 # Define the frequency band
199 freq_min = 0.1 # Minimum frequency in Hz
200 freq_max = 1.0 # Maximum frequency in Hz
201
202 # Apply the bandpass filter
203 for tr in st:
204 tr.filter(type='bandpass', freqmin=freq_min, freqmax=freq_max)
205
206 # Save the filtered MiniSEED file
207 st.write("Woodspoint_2021_filtered.mseed", format="MSEED")
208 {{/code}}
209
robert 3.1 210 === How to slice a waveform ===
robert 2.3 211
212 {{code language="python"}}
213 from obspy import read, UTCDateTime, Stream # Importing Stream here
214
215 # Load the filtered MiniSEED file
216 st = read("Woodspoint_2021_filtered.mseed")
217
218 # Define the time window for slicing
219 slice_start = UTCDateTime("2021-09-21T23:20:00")
220 slice_end = slice_start +10
221
222 # Slice the waveform for each Trace in the Stream
223 sliced_st = Stream() # Now Stream is defined
224 for tr in st:
225 sliced_tr = tr.slice(starttime=slice_start, endtime=slice_end)
226 sliced_st.append(sliced_tr)
227
228 # Save the sliced MiniSEED file
229 sliced_st.write("Woodspoint_2021_filtered_sliced.mseed", format="MSEED")
230 {{/code}}
231
robert 3.1 232 === How to save a waveform ===
robert 2.3 233
234 {{code language="python"}}
235 # Save the sliced file as MiniSEED
236 sliced_st.write("Woodspoint_2021_filtered_sliced.mseed", format="MSEED")
237
238 # Or, save the sliced SAC file
239 sliced_st.write("Woodspoint_2021_filtered_sliced.sac", format="SAC")
240 {{/code}}
241
robert 3.1 242 === How to convert miniSEED to SAC ===
robert 2.3 243
244 {{code language="python"}}
245 from obspy import read
246
247 # Read the MiniSEED file
248 st = read("Woodspoint_2021.mseed")
249
250 # Take the first Trace from the Stream
251 tr = st[0]
252
253 # Save that Trace as a SAC file
254 tr.write("Woodspoint_2021.sac", format="SAC")
255 {{/code}}
256
257
robert 1.3 258 = Earthquake Data =
robert 1.1 259
robert 3.1 260 Earthquake data can be accessed directly or via ObsPy's get_events code
robert 2.2 261
robert 1.3 262 == How to download an Earthquake Catalog ==
263
264 {{code language="python"}}
265 from obspy.clients.fdsn import Client
266 from obspy import UTCDateTime
267
268 # Initialize the AusPass FDSN client
269 client = Client("AUSPASS")
270
271 # Define the time range for the earthquake catalog
272 start_time = UTCDateTime("2021-08-01")
273 end_time = UTCDateTime("2022-01-01") # End of year
274
275 # Define the geographic region (latitude and longitude for Woodspoint, Victoria, Australia)
276 latitude = -37.47
277 longitude = 146.10
278 max_radius = 5 # in degrees
279
280 # Download the earthquake catalog
281 catalog = client.get_events(starttime=start_time, endtime=end_time,
282 minmagnitude=2, latitude=latitude, longitude=longitude,
283 maxradius=max_radius)
284
285 # Save the catalog to a file (e.g., QuakeML format)
286 catalog.write("Woodspoint_earthquakes.xml", format="QUAKEML")
287 {{/code}}
288
robert 1.5 289 == How to plot (Global) Earthquakes ==
robert 1.4 290
291 {{code language="python"}}
292 from obspy import UTCDateTime
293 from obspy.clients.fdsn import Client
294
295 # Initialize FDSN client to connect to the IRIS data center
296 client = Client("IRIS")
297
298 # Set the time range for fetching earthquake data
299 # Start time: January 1, 2023
300 # End time: Current time
301 starttime = UTCDateTime("2023-01-01")
302 endtime = UTCDateTime()
303
304 # Fetch earthquake events with a minimum magnitude of 7
305 catalog = client.get_events(starttime=starttime, endtime=endtime, minmagnitude=7)
306 #client.get_events(). This function returns a Catalog object that contains a list of Event objects.
307 #Each Event object, in turn, has an Origins attribute that contains the depth information
308
309 # Plot the fetched earthquake data using an orthographic projection
310 catalog.plot(projection="ortho", title="Global Earthquakes with Magnitude >= 7 since 2023")
311 #catalog.plot(), ObsPy automatically uses the depth information to color the events in the plot
312 {{/code}}
313
robert 2.2 314 == How to plot (Local) Earthquakes ==
315
316 {{code language="python"}}
317 from obspy import UTCDateTime
318 from obspy.clients.fdsn import Client
319
320 # Initialize FDSN client
321 client = Client("AUSPASS")
322
323 # Define time range
324 starttime = UTCDateTime("2023-01-01")
325 endtime = UTCDateTime()
326
327 # Latitude and longitude bounds for Australia
328 minlatitude = -44.0
329 maxlatitude = -10.0
330 minlongitude = 113.0
331 maxlongitude = 154.0
332
333 # Fetch event data for Australia with a minimum magnitude
334 catalog = client.get_events(starttime=starttime, endtime=endtime, minmagnitude=4,
335 minlatitude=minlatitude, maxlatitude=maxlatitude,
336 minlongitude=minlongitude, maxlongitude=maxlongitude)
337
338 # Plot the earthquakes
339 catalog.plot(projection="local", title="Australia Earthquakes", resolution="i")
340 {{/code}}
341
robert 3.1 342