Source code for openquake.engine.tools.correct_complex_sources

# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4
#
# Copyright (C) 2014-2017 GEM Foundation
#
# OpenQuake is free software: you can redistribute it and/or modify it
# under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# OpenQuake 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 Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with OpenQuake. If not, see <http://www.gnu.org/licenses/>.

"""
Correct complex fault sources to comply with Aki and Richards convention.
"""
import sys
import numpy

from openquake.hazardlib import nrml
from openquake.baselib.node import node_from_xml

from openquake.hazardlib.geo.line import Line
from openquake.hazardlib.geo.point import Point
from openquake.hazardlib.geo.surface import ComplexFaultSurface

AKI_RICH_ERR_MSG = 'Surface does not conform with Aki & Richards convention'
WRONG_ORDER_ERR_MSG = 'Edges points are not in the right order'


[docs]def make_edge(edge): ls = edge.LineString.posList.text.split() coords = numpy.array(map(float, ls)).reshape(-1, 3) return Line([Point(*coord) for coord in coords])
[docs]def reverse(edge): poslist = map(float, edge.LineString.posList.text.split()) coords = numpy.array(poslist).reshape(-1, 3)[::-1] # reversing text = '\n'.join('%s %s %s' % tuple(coord) for coord in coords) edge.LineString.posList.text = text
[docs]def fix_source_node(node): if node.tag.endswith('complexFaultSource'): geom = node.complexFaultGeometry top = geom.faultTopEdge intermediate = [edge for edge in geom.getnodes('intermediateEdge')] bottom = geom.faultBottomEdge edges = map(make_edge, [top] + intermediate + [bottom]) try: ComplexFaultSurface.from_fault_data(edges, mesh_spacing=4.) except ValueError as excp: if AKI_RICH_ERR_MSG in str(excp): print(excp) print('Reverting edges ...') reverse(geom.faultTopEdge) reverse(geom.faultBottomEdge) elif WRONG_ORDER_ERR_MSG in str(excp): print(excp) print('reverting bottom edge ...') reverse(geom.faultBottomEdge) else: raise
if __name__ == '__main__': fname = sys.argv[1] src_model = node_from_xml(fname).sourceModel for src_node in src_model: fix_source_node(src_node) with open(fname, 'wb') as f: nrml.write([src_model], f, xmlns=nrml.NAMESPACE)