{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Miscellaneous Utilities\n", "The follow demonstrates some miscellaneous utilities included in ObsPlus.\n", "\n", "## Geodetics\n", "It is often necessary to calculate geometric parameters (distance, azimuth, etc.) for pairs of entities in two different groups. For example, distance from each event in a catalog to each receiver in an inventory. ObsPlus provides a simple class for completing this task." ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:57.826116Z", "iopub.status.busy": "2026-03-18T12:02:57.825952Z", "iopub.status.idle": "2026-03-18T12:02:59.502447Z", "shell.execute_reply": "2026-03-18T12:02:59.501635Z" } }, "outputs": [], "source": [ "# Load the catalog and inventory from the crandall dataset\n", "import obsplus\n", "from obsplus.utils.geodetics import SpatialCalculator\n", "\n", "crandall = obsplus.load_dataset(\"crandall_test\")\n", "cat = crandall.event_client.get_events()\n", "inv = crandall.station_client.get_stations()" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.504183Z", "iopub.status.busy": "2026-03-18T12:02:59.503997Z", "iopub.status.idle": "2026-03-18T12:02:59.506713Z", "shell.execute_reply": "2026-03-18T12:02:59.506020Z" } }, "outputs": [], "source": [ "# init a SpatialCalculator instance (defaults to Earth's params)\n", "spatial_calc = SpatialCalculator()" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.508299Z", "iopub.status.busy": "2026-03-18T12:02:59.508127Z", "iopub.status.idle": "2026-03-18T12:02:59.632348Z", "shell.execute_reply": "2026-03-18T12:02:59.631538Z" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/home/runner/micromamba/envs/test/lib/python3.13/site-packages/obspy/taup/taup_geo.py:105: UserWarning: Assuming spherical planet when calculating epicentral distance. Install the Python module 'geographiclib' to solve this.\n", " warnings.warn(msg)\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
distance_mazimuthback_azimuthdistance_degreesvertical_distance_m
id1id2
smi:local/248925TA.P17A..BHE41321.405319268.56672988.2616170.3711966307.0
smi:local/248883TA.P17A..BHE41577.826856268.65399188.3469710.3735005867.0
smi:local/248882TA.P17A..BHE41788.158988268.90596788.5973520.3753893477.0
smi:local/248891TA.P17A..BHE41690.761121268.24628387.9385080.3745144927.0
smi:local/248887TA.P17A..BHE40909.197252268.51810988.2160480.3674945847.0
\n", "
" ], "text/plain": [ " distance_m azimuth back_azimuth \\\n", "id1 id2 \n", "smi:local/248925 TA.P17A..BHE 41321.405319 268.566729 88.261617 \n", "smi:local/248883 TA.P17A..BHE 41577.826856 268.653991 88.346971 \n", "smi:local/248882 TA.P17A..BHE 41788.158988 268.905967 88.597352 \n", "smi:local/248891 TA.P17A..BHE 41690.761121 268.246283 87.938508 \n", "smi:local/248887 TA.P17A..BHE 40909.197252 268.518109 88.216048 \n", "\n", " distance_degrees vertical_distance_m \n", "id1 id2 \n", "smi:local/248925 TA.P17A..BHE 0.371196 6307.0 \n", "smi:local/248883 TA.P17A..BHE 0.373500 5867.0 \n", "smi:local/248882 TA.P17A..BHE 0.375389 3477.0 \n", "smi:local/248891 TA.P17A..BHE 0.374514 4927.0 \n", "smi:local/248887 TA.P17A..BHE 0.367494 5847.0 " ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# create distance dataframe\n", "df = spatial_calc(entity_1=cat, entity_2=inv)\n", "df.head()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Since a multi-index is used it provides a fairly intuitive way to look up particular event-channel pairs using a tuple of (event_id, seed_id) in conjunction with the `.loc` DataFrame attribute like so:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.661055Z", "iopub.status.busy": "2026-03-18T12:02:59.660840Z", "iopub.status.idle": "2026-03-18T12:02:59.663692Z", "shell.execute_reply": "2026-03-18T12:02:59.662946Z" } }, "outputs": [], "source": [ "event_id = str(cat[0].resource_id)\n", "seed_id = \"UU.MPU..HHZ\"\n" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.665363Z", "iopub.status.busy": "2026-03-18T12:02:59.665181Z", "iopub.status.idle": "2026-03-18T12:02:59.667751Z", "shell.execute_reply": "2026-03-18T12:02:59.666977Z" } }, "outputs": [], "source": [ "# or just get a particular parameter" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The distances can be converted to km and the distrtibutions can be described." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.669308Z", "iopub.status.busy": "2026-03-18T12:02:59.669126Z", "iopub.status.idle": "2026-03-18T12:02:59.684696Z", "shell.execute_reply": "2026-03-18T12:02:59.684068Z" } }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
azimuthback_azimuthdistance_degreesdistance_kmvertical_distance_km
count408.00408.00408.00408.00408.00
mean201.06211.690.8291.245.22
std104.3798.060.3336.551.85
min9.1648.810.1618.261.70
25%126.47122.450.6370.463.69
50%163.87189.790.8290.855.50
75%293.44307.271.12124.366.25
max358.12344.771.29144.009.31
\n", "
" ], "text/plain": [ " azimuth back_azimuth distance_degrees distance_km \\\n", "count 408.00 408.00 408.00 408.00 \n", "mean 201.06 211.69 0.82 91.24 \n", "std 104.37 98.06 0.33 36.55 \n", "min 9.16 48.81 0.16 18.26 \n", "25% 126.47 122.45 0.63 70.46 \n", "50% 163.87 189.79 0.82 90.85 \n", "75% 293.44 307.27 1.12 124.36 \n", "max 358.12 344.77 1.29 144.00 \n", "\n", " vertical_distance_km \n", "count 408.00 \n", "mean 5.22 \n", "std 1.85 \n", "min 1.70 \n", "25% 3.69 \n", "50% 5.50 \n", "75% 6.25 \n", "max 9.31 " ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Convert add km columns and delete m columns\n", "m_columns = [x for x in df.columns if x.endswith(\"_m\")]\n", "km_columns = [x.replace(\"_m\", \"_km\") for x in m_columns]\n", "\n", "df_km = df.assign(**{x: df[y] / 1000.0 for x, y in zip(km_columns, m_columns)}).drop(\n", " columns=m_columns\n", ")\n", "\n", "# Calculate stats for source reseiver distances\n", "df_km.describe().round(decimals=2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Time\n", "\n", "Working with dates/times can be frustrating, especially since (as of 2020) ObsPy, numpy, and pandas all use slightly different methods for working with time. ObsPlus provides some utilities to make things a little easier." ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.686455Z", "iopub.status.busy": "2026-03-18T12:02:59.686275Z", "iopub.status.idle": "2026-03-18T12:02:59.688962Z", "shell.execute_reply": "2026-03-18T12:02:59.688322Z" } }, "outputs": [], "source": [ "\n", "import numpy as np\n", "import obspy\n", "\n", "import obsplus\n", "from obsplus.utils import to_datetime64, to_timedelta64, to_utc" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Note:** All ObsPlus datafames use numpy/pandas datatypes. " ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.690726Z", "iopub.status.busy": "2026-03-18T12:02:59.690558Z", "iopub.status.idle": "2026-03-18T12:02:59.735160Z", "shell.execute_reply": "2026-03-18T12:02:59.734315Z" } }, "outputs": [], "source": [ "df = obsplus.events_to_df(obspy.read_events())" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.736777Z", "iopub.status.busy": "2026-03-18T12:02:59.736620Z", "iopub.status.idle": "2026-03-18T12:02:59.740451Z", "shell.execute_reply": "2026-03-18T12:02:59.739796Z" } }, "outputs": [ { "data": { "text/plain": [ "0 2012-04-04 14:21:42.300\n", "1 2012-04-04 14:18:37.000\n", "2 2012-04-04 14:08:46.000\n", "Name: time, dtype: datetime64[ns]" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "df[\"time\"]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Numpy and ObsPy time differences" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "One difference between numpy's `datatime64` and ObsPy's `UTCDateTime` is how offsets are applied. For ObsPy, numbers are simply taken as seconds, but numpy requires explicitly using `timedeta64` instances." ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.742140Z", "iopub.status.busy": "2026-03-18T12:02:59.741987Z", "iopub.status.idle": "2026-03-18T12:02:59.744871Z", "shell.execute_reply": "2026-03-18T12:02:59.744086Z" } }, "outputs": [], "source": [ "time_str = \"2020-01-03T11:00:00\"\n", "utc = obspy.UTCDateTime(time_str)\n", "dt64 = np.datetime64(time_str, \"ns\")" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.746417Z", "iopub.status.busy": "2026-03-18T12:02:59.746225Z", "iopub.status.idle": "2026-03-18T12:02:59.748849Z", "shell.execute_reply": "2026-03-18T12:02:59.748133Z" } }, "outputs": [], "source": [ "# add 1/2 second to UTCDateTime\n", "utc2 = utc + 0.5" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.750546Z", "iopub.status.busy": "2026-03-18T12:02:59.750395Z", "iopub.status.idle": "2026-03-18T12:02:59.753000Z", "shell.execute_reply": "2026-03-18T12:02:59.752350Z" } }, "outputs": [], "source": [ "# however doing the same thing with datetime64 raises a TypeError\n", "try:\n", " dt64 + 0.5\n", "except TypeError:\n", " pass" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.754713Z", "iopub.status.busy": "2026-03-18T12:02:59.754562Z", "iopub.status.idle": "2026-03-18T12:02:59.757069Z", "shell.execute_reply": "2026-03-18T12:02:59.756457Z" } }, "outputs": [], "source": [ "# so you need to use a timedelta64\n", "dt64_2 = dt64 + np.timedelta64(500, \"ms\")" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.758721Z", "iopub.status.busy": "2026-03-18T12:02:59.758572Z", "iopub.status.idle": "2026-03-18T12:02:59.762603Z", "shell.execute_reply": "2026-03-18T12:02:59.761849Z" } }, "outputs": [ { "data": { "text/plain": [ "0 2012-04-04 14:21:42.800\n", "1 2012-04-04 14:18:37.500\n", "2 2012-04-04 14:08:46.500\n", "Name: time, dtype: datetime64[ns]" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# This is, of course, also the case with datetime64 columns\n", "df[\"time\"] + np.timedelta64(500, \"ms\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "ObsPlus provides some common converters between ObsPy and Pandas/Numpy time objects." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Converting between ObsPy and Numpy time datatypes\n", "Conversion between ObsPy `UTCDateTime` and numpy `datetime64` objects when needed." ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.764163Z", "iopub.status.busy": "2026-03-18T12:02:59.763999Z", "iopub.status.idle": "2026-03-18T12:02:59.766959Z", "shell.execute_reply": "2026-03-18T12:02:59.766256Z" } }, "outputs": [], "source": [ "# Convert a time column to an array of ObsPy UTCDateTime objects\n", "utc_array = to_utc(df[\"time\"])" ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.768544Z", "iopub.status.busy": "2026-03-18T12:02:59.768381Z", "iopub.status.idle": "2026-03-18T12:02:59.771256Z", "shell.execute_reply": "2026-03-18T12:02:59.770612Z" } }, "outputs": [], "source": [ "# Convert back to datetime64\n", "dt64_array = to_datetime64(utc_array)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "
\n", "\n", "**Note**: `datetime64` arrays are much more efficient in terms of memory usage and computational efficiency than arrays of `UTCDateTime` *objects*.\n", "\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If you prefer not to manually define `timedelta64` to perform offsets, `to_timedelta64` simply converts a real number to an offset in seconds." ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "execution": { "iopub.execute_input": "2026-03-18T12:02:59.773099Z", "iopub.status.busy": "2026-03-18T12:02:59.772813Z", "iopub.status.idle": "2026-03-18T12:02:59.776248Z", "shell.execute_reply": "2026-03-18T12:02:59.775647Z" } }, "outputs": [ { "data": { "text/plain": [ "np.timedelta64(3255000000,'ns')" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "to_timedelta64(3.255)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.12" } }, "nbformat": 4, "nbformat_minor": 4 }