Source code for openquake.commands.plot_sites

from __future__ import division
import math
import numpy
from openquake.commonlib import sap, datastore
from openquake.hazardlib.site import Tile


[docs]def make_tiles(sitecol, sites_per_tile, maximum_distance): tiles = [] hint = math.ceil(len(sitecol) / sites_per_tile) for sites in sitecol.split_in_tiles(hint): tiles.append(Tile(sites, maximum_distance)) return tiles
[docs]def plot_sites(calc_id): """ Plot the hazard sites of a calculations with one ore more a bounding boxes (a box for each tile). If point sources are present, prints them too as circles of radius `integration_distance + max rupture radius`. """ import matplotlib.pyplot as p from matplotlib.patches import Rectangle dstore = datastore.read(calc_id) sitecol = dstore['sitecol'] csm = dstore['composite_source_model'] oq = dstore['oqparam'] fig = p.figure() ax = fig.add_subplot(111) ax.grid(True) tiles = make_tiles(sitecol, oq.sites_per_tile, oq.maximum_distance) print 'There are %d tiles' % len(tiles) for tile in tiles: xs = [] ys = [] area = [] for src in csm.get_sources(): if src in tile and getattr(src, 'location', None): xs.append(src.location.x) ys.append(src.location.y) radius = src._get_max_rupture_projection_radius() r = (tile.maximum_distance[src.tectonic_region_type] + radius) / tile.KM_ONE_DEGREE a = numpy.pi * r ** 2 area.append(a) ax.add_patch(Rectangle(*tile.get_rectangle(), fill=False)) p.scatter(tile.fix_lons(xs), ys, marker='o', s=area) p.scatter(tile.fix_lons(sitecol.lons), sitecol.lats, marker='+') p.show()
parser = sap.Parser(plot_sites) parser.arg('calc_id', 'a computation id', type=int)