Source code for seismic.xcorqc.analytic_plot_utils
#!/usr/bin/env python
"""
Utility functions supporting plotting for cross-correlation visualizations.
"""
import math
import datetime
import numpy as np
from shapely.geometry import Polygon
from descartes import PolygonPatch
[docs]def distance(origin, destination):
"""
Compute the distance in km between origin coordinates and destination coordinates.
The coordinates are (latitude, longitude) couplets in units of degrees.
:param origin: Coordinates of origin point
:type origin: tuple(float, float)
:param destination: Coordinates of destination point
:type destination: tuple(float, float)
:return: Epicentral distance between origin and destination in kilometres
:rtype: float
"""
lat1, lon1 = origin
lat2, lon2 = destination
radius = 6371 # km
dlat = math.radians(lat2 - lat1)
dlon = math.radians(lon2 - lon1)
a = math.sin(dlat / 2) * math.sin(dlat / 2) + math.cos(math.radians(lat1)) \
* math.cos(math.radians(lat2)) * math.sin(dlon / 2) * math.sin(dlon / 2)
c = 2 * math.atan2(math.sqrt(a), math.sqrt(1 - a))
d = radius * c
return d
[docs]def drawBBox(min_lon, min_lat, max_lon, max_lat, base_map, **kwargs):
"""
Draw bounding box on a basemap
:param min_lon: Minimum longitude
:type min_lon: float
:param min_lat: Minimum latitude
:type min_lat: float
:param max_lon: Maximum longitude
:type max_lon: float
:param max_lat: Maximum latitude
:type max_lat: float
:param base_map: Basemap on which to draw the bounding box
:type base_map: mpl_toolkits.basemap.Basemap
"""
bblons = np.array([min_lon, max_lon, max_lon, min_lon, min_lon])
bblats = np.array([min_lat, min_lat, max_lat, max_lat, min_lat])
x, y = base_map(bblons, bblats)
xy = zip(x, y)
poly = Polygon(xy)
base_map.ax.add_patch(PolygonPatch(poly, **kwargs))
[docs]def timestamps_to_plottable_datetimes(time_series):
"""
Convert a series of float (or equivalent) timestamp values to matplotlib plottable datetimes.
:param time_series: Series of timestamps
:type time_series: iterable container
:return: Equivalent series of plottable timestamps
:rtype: numpy.array('datetime64[ms]') with millisecond resolution
"""
plt_times = np.array([datetime.datetime.utcfromtimestamp(v)
for v in time_series]).astype('datetime64[ms]')
return plt_times