source: Tools/ScXML2evt/scxml2evt.py @ 6

#! /usr/bin/env python
# -*- coding: UTF8 -*-

"""
use lxml parsing to reformat seiscomp3's xml dump into evt format for seismic handler
maybe there the possibility to use xslt
"""

__licence__ = "GPLv3"
__author__ = "Marcus Walther, walther@szgrf.bgr.de"
__url__ = "http://www.seismic-handler.org/portal/wiki/ShSoftScpToEvt"
__version__ = 0.1
__revision__  = __id__ = "$Id$"

import sys
from time import strftime, strptime

try:
    from lxml import etree
except ImportError:
    print "fatal error: lxml module not found!"
    sys.exit(1)

# list of stations to be put in evt file
# note: stations = [] will allow all
stations = [
    "WET", "TNS",
]
def stationFilter(data):
    if not len(stations):
        return True

    if data["station"] in stations:
        return True
    else:
        return False

class sc2evt(object):
    """convert seiscomp3 xmldump to evt file"""

    def __init__(self, src, trgt):
        self.target = trgt

        # older version don't support XMLParser
        try:
            parser = etree.XMLParser(remove_blank_text=True)
        except AttributeError:
            parser = None

        try:
            x = etree.parse(src, parser)
        except Exception, e:
            print e
            sys.exit(3)

        self.xml = x
        self._output()

    def _output(self):
        if "arrivals" not in self.__dict__:
            self._read()

        # width of first column
        WIDTH = 23
        FINAL = "--- End of Phase ---\n\n"

        fields_std = {
            "station": "Station code",
            "time": "Onset time",
            "phase": "Phase name",
            "channel": "Component",
            "distance": "Distance (deg)",
            "azimuth": "Epi-Azimuth (deg)",
        }

        fields_event = {
            "latitude": "Latitude",
            "longitude": "Longitude",
            "depth": "Depth (km)",
            "time": "Origin time"
        }

        fields_mag = {
            "Mw(mB)": "Broadband Magnitude",
        }

        fields_info = {
            "Pick Type": "automatic",
        }

        formatted = {}
        for a in filter(stationFilter, self.arrivals):
            for field in a:
                try:
                    formfield = fields_std[field]
                except KeyError:
                    continue

                if field == "channel":
                    value = a[field][-1]
                else:
                    value = a[field]

                try:
                    formatted[a["station"]][formfield] = value
                except KeyError:
                    formatted[a["station"]] = {formfield: value}

            # add info field(s)
            formatted[a["station"]].update(fields_info)

        # add event info to last pick
        last = formatted[a["station"]]
        for e in self.event:
            try:
                formfield = fields_event[e]
            except KeyError:
                continue

            last[formfield] = self.event[e]

        # add magnitude information
        try:
            last[fields_mag[self.event["magnitude_type"]]] = self.event["magnitude"]
        except:
            pass

        for entry in formatted.values():
            for line in entry:
                print >> self.target, "%s%s: %s" % (line, " "*(WIDTH-len(line)), entry[line])
            print >> self.target, FINAL

    def _read(self):
        x = self.xml

        # get information by xpath usage
        picks_raw = x.xpath("/seiscomp/EventParameters/pick")
        arrivals_raw = x.xpath("/seiscomp/EventParameters/origin/arrival")

        # no picks found
        if not len(picks_raw) or not len(arrivals_raw):
            print "error: no pick or arrival information found"
            print "error: maybe you have a very old version of libxml/libxslt"
            sys.exit(4)

        # event information
        self.event = {
            "time": self.sctime2sh(x.xpath("/seiscomp/EventParameters/origin/time/value")[0].text),
            "latitude": x.xpath("/seiscomp/EventParameters/origin/latitude/value")[0].text,
            "longitude": x.xpath("/seiscomp/EventParameters/origin/longitude/value")[0].text,
            "depth": x.xpath("/seiscomp/EventParameters/origin/depth/value")[0].text,
            "magnitude": x.xpath("/seiscomp/EventParameters/origin/networkMagnitude/magnitude/value")[0].text,
            "magnitude_type": x.xpath("/seiscomp/EventParameters/origin/networkMagnitude/type")[0].text
        }

        # build hash of picks
        picks = {}
        for pick in picks_raw:
            wf = pick.xpath("./waveformID")[0]
            picks[pick.get("publicID")] = {
                "network": wf.get("networkCode"),
                "station": wf.get("stationCode"),
                "channel": wf.get("channelCode"),
                "time": self.sctime2sh(pick.xpath("./time/value")[0].text)
            }

        arrivals = []
        # connect arrivals and picks
        for arrival in arrivals_raw:
            # skip bad arrivals
            if arrival.get("weight") == "0":
                continue

            try:
                pinfo = picks[arrival.xpath("./pickID")[0].text]
            except KeyError:
                continue

            pinfo["phase"] = arrival.xpath("./phase")[0].text
            pinfo["distance"] = arrival.xpath("./distance")[0].text
            pinfo["azimuth"] = arrival.xpath("./azimuth")[0].text

            arrivals.append(pinfo)

        self.arrivals = arrivals

    def sctime2sh(self, dtstr):
        """
        reformat seiscomp time string into seismic handler format

        input format: 2008-06-14T00:03:49.76134Z
        output format: 14-JUN-2008_00:03:49.761
        """

        # we need to reformat date part only
        dt = strftime("%d-%b-%Y", strptime(dtstr[:10], "%Y-%m-%d")).upper()

        # shorten milliseconds part if necessary
        tt = dtstr[11:-1]
        if tt.index(".") < len(tt)-4:
            tt = tt[:tt.index(".")+4]

        return "_".join((dt, tt))

if __name__ == "__main__":
    try:
        i = sys.argv[1]
    except IndexError:
        print "error: please give source of seiscomp3 xmldump as command argument"
        sys.exit(2)

    try:
        o = open(sys.argv[2], "w")
    except IndexError:
        o = sys.stdout

    a = sc2evt(i, o)