#! /usr/bin/env python3 # -*- coding: utf-8 -*- ''' Script to lookup city names of events with Nominatim service The input should be an valid quakeML file passed to stdin. The output will will be a javascript structure to be included in the SeisObs map service. The script should be updated regularly keep the total number of all AJAX calls to the Nominatim service small, e. g. : curl -s "https://fdsnws.geophysik.ruhr-uni-bochum.de/fdsnws/event/1/query?minlat=50&maxlat=54&minlon=3&maxlon=10&minmag=1" | mkGeolocationTable.py > geolocationTable.js License Copyright 2020-2021 Kasper D. Fischer This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/. ''' def mkGeolocationTable(file=''): ## imports # XML ElementTree try: import xml.etree.cElementTree as ET except ImportError: import xml.etree.ElementTree as ET # json import json as JSON # sys stdin from sys import stdin # geopy from geopy.geocoders import Photon from geopy.extra.rate_limiter import RateLimiter from geopy.exc import GeocoderServiceError ## constants URL = 'https://photon.komoot.io/reverse?lon={lng:.3f}&lat={lat:.3f}&limit=5' NAMESPACES = {'sc3': 'http://geofon.gfz-potsdam.de/ns/seiscomp3-schema/0.7', 'qml': 'http://quakeml.org/xmlns/bed/1.2'} # try loading the file geolocationTable = {} if file : try: jsonfile = open(file) jsonfileContent = jsonfile.read().split('=')[1].replace(';', '') geolocationTable = JSON.loads(jsonfileContent) except: geolocationTable = {} logging.warning('Could not parse file %s' %file) # parse event.xml DOM = ET.parse(stdin).getroot() geolocator = Photon() reverse_geolocate = RateLimiter(geolocator.reverse, min_delay_seconds=1) # iterate over all events count = 0 for event in DOM.iterfind('qml:eventParameters/qml:event', NAMESPACES): count += 1 publicID = event.attrib['publicID'].split('/')[2] lat = float(event.find('./qml:origin/qml:latitude/qml:value', NAMESPACES).text) lng = float(event.find('./qml:origin/qml:longitude/qml:value', NAMESPACES).text) evaluationMode = event.find('./qml:origin/qml:evaluationMode', NAMESPACES).text if publicID in geolocationTable: logging.warning('Skipping cached event {id}'.format(id=publicID)) elif evaluationMode == 'automatic': logging.warning('Skipping automatic event {id}'.format(id=publicID)) else: logging.info('Processing event {id}'.format(id=publicID)) try: locations = reverse_geolocate("{lat:.3f}, {lng:.3f}".format(lat=lat, lng=lng),exactly_one=False,limit=5) except GeocoderServiceError: logging.warning('Reverse Geolocation failed. Skipping event.') continue place = [] for location in locations: try: place = location.raw['properties']['city'] except KeyError: try: place = location.raw['properties']['town'] except KeyError: try: place = location.raw['properties']['village'] except KeyError: try: place = location.raw['properties']['county'] except KeyError: logging.debug('Could not extract city for event {id} at {lat:.3f} N / {lng:.3f} E (Service: {url}), trying next returned result' .format(id=publicID, lat=lat, lng=lng, url=URL.format(lat=lat,lng=lng))) logging.debug(location.raw) if not place: logging.critical('Could not extract city for event {id} at {lat:.3f} N / {lng:.3f} E (Service: {url})' .format(id=publicID, lat=lat, lng=lng, url=URL.format(lat=lat,lng=lng))) geolocationTable[publicID] = place # dump json print('var geolocationTable = {0};'.format(JSON.dumps(geolocationTable, sort_keys=True))) logging.info("processed %d events", count) # __main__ if __name__ == "__main__": # use module logging import logging # parse arguments import argparse versionText = 'v0.1 (2021-10-16)' parser = argparse.ArgumentParser( description='Reverse geocoding lookup of events in xml format (stdin).', epilog=versionText) parser.add_argument('-V', '--version', action='version', version=versionText, help="show version") parser.add_argument('-f', '--file', action='store', dest='file', help='read in JSON file containing old output.') group = parser.add_mutually_exclusive_group() group.add_argument("-v", "--verbose", action="store_true", help="increase output verbosity") group.add_argument("-q", "--quiet", action="store_true", help="disable output") group.add_argument("-d", "--debug", action="store_true", help="show debugging output", ) cla = parser.parse_args() # set logging level if cla.quiet: logging.basicConfig(format='%(levelname)s: %(message)s', level=logging.ERROR) elif cla.verbose: logging.basicConfig(format='%(levelname)s: %(message)s', level=logging.INFO) elif cla.debug: logging.basicConfig(format='%(levelname)s: %(message)s', level=logging.DEBUG) else: logging.basicConfig(format='%(levelname)s: %(message)s', level=logging.WARN) # call mkGeolocationTable(...) mkGeolocationTable(file=cla.file)