Skip to Content

Wiki source code of FDSN Guide

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