Skip to Content

Wiki source code of FDSN Guide

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