source: Tools/ScXML2evt/scxml2evt.py @ 14

Revision 14, 6.5 KB checked in by marcus, 14 years ago (diff)
  • improved: error output redirected to stderr
  • bugfix: clean exit if no stations left for processing after filtering
  • Property svn:executable set to *
Line 
1#! /usr/bin/env python
2# -*- coding: UTF8 -*-
3
4"""
5use lxml parsing to reformat seiscomp3's xml dump into evt format for seismic handler
6maybe there the possibility to use xslt
7"""
8
9__licence__ = "GPLv3"
10__author__ = "Marcus Walther, walther@szgrf.bgr.de"
11__url__ = "http://www.seismic-handler.org/portal/wiki/ShSoftScpToEvt"
12__version__ = 0.1
13__revision__  = __id__ = "$Id$"
14
15import sys
16from time import strftime, strptime
17
18try:
19    from lxml import etree
20except ImportError:
21    print >> sys.stderr, "fatal error: lxml module not found!"
22    sys.exit(1)
23
24# list of stations to be put in evt file
25# note: stations = [] will allow all
26stations = [
27    "WET", "TNS",
28]
29def stationFilter(data):
30    if not len(stations):
31        return True
32
33    if data["station"] in stations:
34        return True
35    else:
36        return False
37
38class sc2evt(object):
39    """convert seiscomp3 xmldump to evt file"""
40
41    def __init__(self, src, trgt):
42        self.target = trgt
43
44        # older version don't support XMLParser
45        try:
46            parser = etree.XMLParser(remove_blank_text=True)
47        except AttributeError:
48            parser = None
49
50        try:
51            x = etree.parse(src, parser)
52        except Exception, e:
53            print >> sys.stderr, e
54            sys.exit(3)
55
56        self.xml = x
57        self._output()
58
59    def _output(self):
60        if "arrivals" not in self.__dict__:
61            self._read()
62
63        # width of first column
64        WIDTH = 23
65        FINAL = "--- End of Phase ---\n\n"
66
67        fields_std = {
68            "station": "Station code",
69            "time": "Onset time",
70            "phase": "Phase name",
71            "channel": "Component",
72            "distance": "Distance (deg)",
73            "azimuth": "Epi-Azimuth (deg)",
74        }
75
76        fields_event = {
77            "latitude": "Latitude",
78            "longitude": "Longitude",
79            "depth": "Depth (km)",
80            "time": "Origin time"
81        }
82
83        fields_mag = {
84            "Mw(mB)": "Broadband Magnitude",
85        }
86
87        fields_info = {
88            "Pick Type": "automatic",
89        }
90
91        formatted = {}
92        for a in filter(stationFilter, self.arrivals):
93            for field in a:
94                try:
95                    formfield = fields_std[field]
96                except KeyError:
97                    continue
98
99                if field == "channel":
100                    value = a[field][-1]
101                else:
102                    value = a[field]
103
104                try:
105                    formatted[a["station"]][formfield] = value
106                except KeyError:
107                    formatted[a["station"]] = {formfield: value}
108
109            # add info field(s)
110            formatted[a["station"]].update(fields_info)
111
112        # add event info to last pick
113        try:
114            last = formatted[a["station"]]
115        except NameError:
116            # only if no stations left after filtering
117            print >> sys.stderr, "no picks on filtered station list"
118            sys.exit(4)
119
120        for e in self.event:
121            try:
122                formfield = fields_event[e]
123            except KeyError:
124                continue
125
126            last[formfield] = self.event[e]
127
128        # add magnitude information
129        try:
130            last[fields_mag[self.event["magnitude_type"]]] = self.event["magnitude"]
131        except:
132            pass
133
134        for entry in formatted.values():
135            for line in entry:
136                print >> self.target, "%s%s: %s" % (line, " "*(WIDTH-len(line)), entry[line])
137            print >> self.target, FINAL
138
139    def _read(self):
140        x = self.xml
141
142        # get information by xpath usage
143        picks_raw = x.xpath("/seiscomp/EventParameters/pick")
144        arrivals_raw = x.xpath("/seiscomp/EventParameters/origin/arrival")
145
146        # no picks found
147        if not len(picks_raw) or not len(arrivals_raw):
148            print >> sys.stderr, "error: no pick or arrival information found"
149            print >> sys.stderr, "error: maybe you have a very old version of libxml/libxslt"
150            sys.exit(4)
151
152        # event information
153        self.event = {
154            "time": self.sctime2sh(x.xpath("/seiscomp/EventParameters/origin/time/value")[0].text),
155            "latitude": x.xpath("/seiscomp/EventParameters/origin/latitude/value")[0].text,
156            "longitude": x.xpath("/seiscomp/EventParameters/origin/longitude/value")[0].text,
157            "depth": x.xpath("/seiscomp/EventParameters/origin/depth/value")[0].text,
158            "magnitude": x.xpath("/seiscomp/EventParameters/origin/networkMagnitude/magnitude/value")[0].text,
159            "magnitude_type": x.xpath("/seiscomp/EventParameters/origin/networkMagnitude/type")[0].text
160        }
161
162        # build hash of picks
163        picks = {}
164        for pick in picks_raw:
165            wf = pick.xpath("./waveformID")[0]
166            picks[pick.get("publicID")] = {
167                "network": wf.get("networkCode"),
168                "station": wf.get("stationCode"),
169                "channel": wf.get("channelCode"),
170                "time": self.sctime2sh(pick.xpath("./time/value")[0].text)
171            }
172
173        arrivals = []
174        # connect arrivals and picks
175        for arrival in arrivals_raw:
176            # skip bad arrivals
177            if arrival.get("weight") == "0":
178                continue
179
180            try:
181                pinfo = picks[arrival.xpath("./pickID")[0].text]
182            except KeyError:
183                continue
184
185            pinfo["phase"] = arrival.xpath("./phase")[0].text
186            pinfo["distance"] = arrival.xpath("./distance")[0].text
187            pinfo["azimuth"] = arrival.xpath("./azimuth")[0].text
188
189            arrivals.append(pinfo)
190
191        self.arrivals = arrivals
192
193    def sctime2sh(self, dtstr):
194        """
195        reformat seiscomp time string into seismic handler format
196
197        input format: 2008-06-14T00:03:49.76134Z
198        output format: 14-JUN-2008_00:03:49.761
199        """
200
201        # we need to reformat date part only
202        dt = strftime("%d-%b-%Y", strptime(dtstr[:10], "%Y-%m-%d")).upper()
203
204        # shorten milliseconds part if necessary
205        tt = dtstr[11:-1]
206        if tt.index(".") < len(tt)-4:
207            tt = tt[:tt.index(".")+4]
208
209        return "_".join((dt, tt))
210
211if __name__ == "__main__":
212    try:
213        i = sys.argv[1]
214    except IndexError:
215        print >> sys.stderr, "error: please give source of seiscomp3 xmldump as command argument"
216        sys.exit(2)
217
218    try:
219        o = open(sys.argv[2], "w")
220    except IndexError:
221        o = sys.stdout
222
223    a = sc2evt(i, o)
Note: See TracBrowser for help on using the repository browser.