seisobs-webapp/scripts/mkGeolocationTable.py

#! /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 Kasper D. Fischer <kasper.fischer@rub.de>

    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/.

    $Id$
'''

def mkGeolocationTable(file=''):
    # imports
    try:
        import xml.etree.cElementTree as ET
    except ImportError:
        import xml.etree.ElementTree as ET
    from sys import stdin
    import warnings
    from time import sleep
    from geopy.geocoders import Photon
    import urllib2 as URL
    import json as JSON

    # constants
    namespaces = {'sc3': 'http://geofon.gfz-potsdam.de/ns/seiscomp3-schema/0.7',
        'qml': 'http://quakeml.org/xmlns/bed/1.2'}

    def simple_warning(message, category, filename, lineno, file=None, line=None):
        return 'Warning: %s\n' % (message)
    warnings.formatwarning = simple_warning

    # try loading the file
    geolocationTable = {}
    if file :
        try:
            jsonfile = open(file)
            jsonfileContent = jsonfile.read().split('=')[1].replace(';', '')
            geolocationTable = JSON.loads(jsonfileContent)
        except:
            geolocationTable = {}
            warnings.warn('Could not parse file %s' %file)

    # parse event.xml
    DOM = ET.parse(stdin).getroot()
    geolocator = Photon()

    # iterate over all events
    for event in DOM.iterfind('qml:eventParameters/qml:event', namespaces):
        publicID = event.attrib['publicID'].split('/')[2]
        lat = event.find('./qml:origin/qml:latitude/qml:value', namespaces).text
        lng = event.find('./qml:origin/qml:longitude/qml:value', namespaces).text
        evaluationMode = event.find('./qml:origin/qml:evaluationMode', namespaces).text

        if publicID in geolocationTable:
            warnings.warn('Skipping cached event %s' %(publicID))
        elif evaluationMode == 'automatic':
            warnings.warn('Skipping automatic event %s' %(publicID))
        else:
            try:
                location = geolocator.reverse("{lat}, {lng}".format(lat=lat, lng=lng))
            except:
                warnings.warn('Reverse Geolocation failed. Skipping event.')
                continue
            place = []
            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:
                            warnings.warn('Could not extract city for event {0} at {1} N / {2} E (URL: {3})'.format(publicID, lat, lng, url))
            geolocationTable[publicID] = place
            sleep(1.1)

            '''
            #warnings.warn('Processing event %s' %publicID)
            # send and evaluate nominatim request
            url = 'https://open.mapquestapi.com/nominatim/v1/reverse.php?key=oE512xGQxeH1n2cueacJ6xzRv7qFlCib&lat={0}&lon={1}&zoom=10&format=json'.format(lat, lng)
            response = URL.urlopen(url)
            if ( response.msg == 'OK' ):
                data = JSON.loads(response.read())
                city = []
                try:
                    try:
                        city = data['address']['city']
                    except:
                        warnings.warn('Request {3} for event {0} at {1} N / {2} E did not provide city attribute\n\t(Response: {4})'.format(publicID, lat, lng, url, data))
                        try:
                            city = data['address']['town']
                            warnings.warn('Using attribute town ({1}) for event {0}'.format(publicID, city))
                        except:
                            try:
                                city = data['address']['county']
                                warnings.warn('Using attribute county ({1}) for event {0}'.format(publicID, city))
                            except:
                                warnings.warn('Skipping event')
                                continue
                    countryCode = data['address']['country_code'].upper()
                    if ( countryCode == 'DE' ):
                        geolocationTable[publicID] = city.encode('utf-8')
                    else:
                        geolocationTable[publicID] = '{0} ({1})'.format(city.encode('utf-8'), countryCode)
                except:
                    warnings.warn('Could not extract city for event {0} at {1} N / {2} E (URL: {3})'.format(publicID, lat, lng, url))
            else:
                warnings.warn('Request {0} failed'.format(url))
            '''
    # dump json
    print('var geolocationTable = {0};'.format(JSON.dumps(geolocationTable, sort_keys=True)))

# __main__
if __name__ == "__main__":
    # parse arguments
    import argparse
    versionText = '$Revision$ ($Date$, $Author$)'.replace('$', '').replace(':','')
    parser = argparse.ArgumentParser(
        description='Reverse geocoding lookup of events in xml format (stdin).',
        epilog=versionText)
    parser.add_argument('-v', '-V', '--version', action='version',
        version=versionText)
    parser.add_argument('-f', '--file', action='store', dest='file',
        help='read in JSON file containing old output.')
    cla = parser.parse_args()

    # call mkGeolocationTable(...)
    mkGeolocationTable(file=cla.file)