Read waveform files into an ObsPy Stream object.
The read() function opens either one or multiple waveform files given via file name or URL using the pathname_or_url attribute.
The format of the waveform file will be automatically detected if not given. See the Supported Formats section below for available formats.
This function returns an ObsPy Stream object, a list-like object of multiple ObsPy Trace objects.
Parameters: |
|
---|---|
Returns: | An ObsPy Stream object. |
Basic Usage
In most cases a filename is specified as the only argument to read(). For a quick start you may omit all arguments and ObsPy will create and return a basic example seismogram. Further usages of the read() function can be seen in the Further Examples section underneath.
>>> from obspy import read
>>> st = read()
>>> print(st)
3 Trace(s) in Stream:
BW.RJOB..EHZ | 2009-08-24T00:20:03.000000Z - ... | 100.0 Hz, 3000 samples
BW.RJOB..EHN | 2009-08-24T00:20:03.000000Z - ... | 100.0 Hz, 3000 samples
BW.RJOB..EHE | 2009-08-24T00:20:03.000000Z - ... | 100.0 Hz, 3000 samples
Supported Formats
Additional ObsPy modules extend the functionality of the read() function. The following table summarizes all known file formats currently supported by ObsPy. The table order also reflects the order of the autodetection routine if no format option is specified.
Please refer to the Linked Function Call of each module for any extra options available at the import stage.
Format | Required Module | Linked Function Call |
---|---|---|
MSEED | obspy.mseed | obspy.mseed.core.readMSEED() |
SAC | obspy.sac | obspy.sac.core.readSAC() |
GSE2 | obspy.gse2 | obspy.gse2.core.readGSE2() |
SEISAN | obspy.seisan | obspy.seisan.core.readSEISAN() |
SACXY | obspy.sac | obspy.sac.core.readSACXY() |
GSE1 | obspy.gse2 | obspy.gse2.core.readGSE1() |
Q | obspy.sh | obspy.sh.core.readQ() |
SH_ASC | obspy.sh | obspy.sh.core.readASC() |
SLIST | obspy.core | obspy.core.ascii.readSLIST() |
TSPAIR | obspy.core | obspy.core.ascii.readTSPAIR() |
Y | obspy.y | obspy.y.core.readY() |
SEGY | obspy.segy | obspy.segy.core.readSEGY() |
SU | obspy.segy | obspy.segy.core.readSU() |
SEG2 | obspy.seg2 | obspy.seg2.seg2.readSEG2() |
WAV | obspy.wav | obspy.wav.core.readWAV() |
PICKLE | obspy.core | obspy.core.stream.readPickle() |
DATAMARK | obspy.datamark | obspy.datamark.core.readDATAMARK() |
CSS | obspy.css | obspy.css.core.readCSS() |
Next to the read() function the write() method of the returned Stream object can be used to export the data to the file system.
Further Examples
Example waveform files may be retrieved via http://examples.obspy.org.
Reading multiple local files using wildcards.
The following code uses wildcards, in this case it matches two files. Both files are then read into a single Stream object.
>>> from obspy import read
>>> st = read("/path/to/loc_R*.z")
>>> print(st)
2 Trace(s) in Stream:
.RJOB..Z | 2005-08-31T02:33:49.850000Z - ... | 200.0 Hz, 12000 samples
.RNON..Z | 2004-06-09T20:05:59.850000Z - ... | 200.0 Hz, 12000 samples
Reading a local file without format detection.
Using the format parameter disables the automatic detection and enforces reading a file in a given format.
>>> from obspy import read
>>> st = read("/path/to/loc_RJOB20050831023349.z", format="GSE2")
>>> print(st)
1 Trace(s) in Stream:
.RJOB..Z | 2005-08-31T02:33:49.850000Z - ... | 200.0 Hz, 12000 samples
Reading a remote file via HTTP protocol.
>>> from obspy import read
>>> st = read("http://examples.obspy.org/loc_RJOB20050831023349.z")
>>> print(st)
1 Trace(s) in Stream:
.RJOB..Z | 2005-08-31T02:33:49.850000Z - ... | 200.0 Hz, 12000 samples
Reading a compressed files.
>>> from obspy import read
>>> st = read("/path/to/tspair.ascii.gz")
>>> print(st)
1 Trace(s) in Stream:
XX.TEST..BHZ | 2008-01-15T00:00:00.025000Z - ... | 40.0 Hz, 635 samples
>>> st = read("http://examples.obspy.org/slist.ascii.bz2")
>>> print(st)
1 Trace(s) in Stream:
XX.TEST..BHZ | 2008-01-15T00:00:00.025000Z - ... | 40.0 Hz, 635 samples
Reading a file-like object.
>>> from StringIO import StringIO
>>> import urllib2
>>> example_url = "http://examples.obspy.org/loc_RJOB20050831023349.z"
>>> stringio_obj = StringIO(urllib2.urlopen(example_url).read())
>>> st = read(stringio_obj)
>>> print(st)
1 Trace(s) in Stream:
.RJOB..Z | 2005-08-31T02:33:49.850000Z - ... | 200.0 Hz, 12000 samples
Using ‘starttime’ and ‘endtime’ parameters
>>> from obspy import read
>>> dt = UTCDateTime("2005-08-31T02:34:00")
>>> st = read("http://examples.obspy.org/loc_RJOB20050831023349.z",
... starttime=dt, endtime=dt+10)
>>> print(st)
1 Trace(s) in Stream:
.RJOB..Z | 2005-08-31T02:34:00.000000Z - ... | 200.0 Hz, 2001 samples