Figure 4: pandas / tabular data - scatter plot - seismicity around Indonesia

Fig4_PyGMT_pandas.ipynb

@@ -0,0 +1,35 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "id": "533bc728-db2c-4f3a-ab3f-1f4f2c9ecd69",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Copy finale version of script from normal Python file"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "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.12.3"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}

Fig4_PyGMT_pandas.py

@@ -0,0 +1,63 @@
+import io
+
+import pandas as pd
+import pygmt
+import requests
+from pygmt.params import Box
+
+params = {
+    "format": "csv",
+    "starttime": "2000-01-01",
+    "endtime": "2025-10-30",
+    "mindepth": 70,
+    "minmagnitude": 5,
+    "orderby": "magnitude",
+}
+r = requests.get("https://earthquake.usgs.gov/fdsnws/event/1/query", params=params)
+df_eqs = pd.read_csv(io.StringIO(r.text))
+df_eqs = df_eqs[
+    (df_eqs["longitude"] >= 91)
+    & (df_eqs["longitude"] <= 134)
+    & (df_eqs["latitude"] >= -22)
+    & (df_eqs["latitude"] <= 19)
+    & (df_eqs["magType"] != "mb")
+]
+
+fig = pygmt.Figure()
+fig.basemap(region=[91, 134, -22, 19], projection="M15c", frame=True)
+fig.coast(land="gray95", shorelines="gray50")
+
+# Plot epicenters with color (hypocentral depth) or size (moment magnitude)
+pygmt.makecpt(cmap="SCM/navia", series=[0, 700], reverse=True, transparency=30)
+fig.colorbar(frame=["xaf+lHypocentral depth", "y+lkm"])
+fig.plot(
+    x=df_eqs.longitude,
+    y=df_eqs.latitude,
+    size=0.005 * 2**df_eqs.mag,
+    fill=df_eqs.depth,
+    cmap=True,
+    style="c",
+    pen="gray10",
+)
+# Add legend for size-coding
+legend = io.StringIO(
+    "\n".join(f"S 0.4 c {0.005 * 2**mag:.2f} - 1p 1 Mw {mag}" for mag in [5, 6, 7])
+)
+fig.legend(spec=legend, position="jBR+o0.2c+l2", box=Box(fill="white", pen="black"))
+
+# Add histogram for moment magnitude distribution
+with fig.inset(position="jBL+w7c/4c+o0.2c", margin=(1.2, 0.2, 0.9, 0.2), box=True):
+    with pygmt.config(FONT="8p"):
+        fig.histogram(
+            region=[4.8, 10.2, 0, 0],
+            projection="X?/?",
+            frame=["WSrt", "xa1f0.2+lMw", "yaf+lCounts"],
+            data=df_eqs.mag,
+            series=0.2,
+            fill="darkgray",
+            pen="lightgray",
+            histtype=0,
+        )
+
+fig.show()
+fig.savefig("Fig4_PyGMT_pandas.png")