diff --git a/README.md b/README.md index 454ac0c..430b599 100644 --- a/README.md +++ b/README.md @@ -1,267 +1,187 @@
-
+
- 
+ English · 中文简介 · 简体中文完整文档 · 日本語
+
+ Gallery · Quick start · Journal toolkit · API · Research-figure guide · Cite
+
+ 
+ 
+
python examples/journal_figures.py and python examples/journal_gallery.py.
-```python
-from plot_function import plot_map
+## Features
-result = plot_map(
- "data/ERA5temp_1978_monthly.nc",
- variable="t2m",
- isel={"time": 0},
- offset=-273.15, units="°C", # This file stores temperature in Kelvin.
- cmap="RdYlBu_r",
- title="January 1978 · 2 m air temperature",
- output="january.png",
-)
-```
+- **One-call maps from NetCDF** — `plot_map("file.nc", "t2m", isel={"time": 0})` handles variable/dimension selection, coordinate normalisation, projection, colourbar and export.
+- **Journal typography** — `journal_style()` applies Arial/Helvetica (with Liberation/Nimbus fallbacks), 7–9 pt text, 0.5 pt lines, editable text in PDF/SVG; `figsize("single" | "double")` gives 89/183 mm column widths.
+- **Discrete colour scales** — `discrete_cmap()` returns a `(cmap, norm)` pair whose extend triangles take the darkest colours, exactly like printed atlases.
+- **Significance & agreement** — `add_stippling()` (staggered dot lattice) and `add_hatching()` from a mask, or `Significance(p_values, correction="fdr_bh")` in `plot_map`.
+- **Lower-left insets** — `add_inset_bars()` with hatched low-agreement fractions; `add_inset_histogram()` with stacked groups, a total outline, log bins and a cumulative curve on an offset twin axis. Labels get a white halo so they stay legible over the map.
+- **Aligned marginal profiles** — `add_latitude_profile()` / `add_longitude_profile()` push latitudes through the map projection, so 40°N on the profile is level with 40°N on Robinson or Equal Earth maps; IQR/ensemble bands, min–max envelopes or individual members.
+- **Cartographic polish** — `geo_ticks()` (degree labels, dashed light graticule), `add_land()` (cream land, thin grey borders), `add_size_legend()`, `add_panel_label()`, `save_figure(dpi=600)`.
+- **Batteries included** — CLI (`plot-function plot …`), offline tests, CI, and the original notebook helpers (`from utils import plot`) kept for backwards compatibility.
-First use of coastlines may download Natural Earth data. Add `coastlines=False` to avoid that download. In a notebook, display `result.figure`; in a Python GUI session, call `matplotlib.pyplot.show()`.
-
-Want a tiny example without external downloads or the bundled climate dataset?
+## Quick start
```bash
-python examples/quickstart.py
+git clone https://github.com/GISWLH/plot-function.git
+cd plot-function
+python -m pip install -e . # Python 3.10+, not yet on PyPI
+python examples/journal_quickstart.py
```
-This creates a clearly labeled **synthetic** NetCDF field and a PNG in `examples/output/`, using only the core dependencies.
-
-## Research figures: statistics, context, and significance
-
-Build a complete composition with reusable option objects. The new layers work with the same NetCDF inputs and editable Matplotlib results.
+**A complete map from the bundled NetCDF file:**
```python
-from plot_function import plot_map, Profile, Distribution, MapFeatures
+import numpy as np
+from plot_function import plot_map, Profile, Distribution
result = plot_map(
"data/ERA5temp_1978_monthly.nc", "t2m",
- isel={"time": 6}, offset=-273.15, units="°C",
- projection="platecarree", extent=[92, 142, 8, 53], cmap="RdYlBu_r",
- profiles=[Profile(position="right", style="band"),
- Profile(position="top", style="line", weights="coslat")],
+ reduce="time", offset=-273.15, units="°C",
+ cmap="RdYlBu_r", levels=np.arange(-40, 41, 5),
+ title="Annual-mean 2 m air temperature",
+ profiles=Profile(style="band", reference=0, label="Zonal mean ± 1 s.d."),
distribution=Distribution(style="bars", weights="coslat", color="map"),
- features=MapFeatures(rivers=True, lakes=True, borders=True, resolution="50m"),
- title="East Asia · July 1978", figsize=(12, 9),
- output="research-figure.png",
+ panel_label="a",
+ output="annual_mean.png", dpi=600,
)
```
-| Layer | Styles and controls |
-| :--- | :--- |
-| Marginal profiles | Right/top; `line`, `band`, `bars`; mean or median; spatial SD or IQR; reference lines |
-| Distribution inset | Bottom-left by default; `bars`, `step`, `line`, `ecdf`; custom bins and placement; mean marker |
-| Weighting | Equal cells, cosine latitude, or explicitly aligned cell weights |
-| Significance | `stipple`, `hatch`, `contour`; supplied p-values or boolean masks; optional BH FDR correction |
-| Geographic detail | Rivers, lakes, borders, land/ocean backgrounds; 110m, 50m, or 10m Natural Earth layers |
-| Composition | Panel labels, shared color scales, editable axes, PNG/PDF/SVG export |
-
-
-
-
@@ -13,7 +11,18 @@
新工作流无需准备 Shapefile 或 GeoTIFF，也不依赖 Salem。可选海岸线由 Cartopy 的 Natural Earth 缓存提供；设置 `coastlines=False` 即可完全离线绘图。
-[快速开始](#快速开始) · [示例画廊](#示例画廊) · [使用自己的数据](#使用自己的数据) · [命令行](#命令行) · [API 参考](docs/API.md)
+[顶刊图库](#顶刊图库04-新增) · [快速开始](#快速开始) · [示例画廊](#示例画廊) · [使用自己的数据](#使用自己的数据) · [命令行](#命令行) · [API 参考](docs/API.md)
+
+## 顶刊图库（0.4 新增）
+
+新增 `plot_function.journal` 模块，参照 Nature / Science 版式：Arial 7–9 pt 字体、带三角端点的离散色标、低一致性打点与纹理、左下角统计小图（分类柱状图带低一致性纹理 / 对数直方图 + 累积曲线双轴）、右侧**与地图纬度严格对齐**的纬向剖面（Robinson 等任意投影均对齐，含 IQR / 集合范围阴影）、`40°N` 式刻度与浅色虚线经纬网、600 dpi 导出。
+
+| | |
+| :---: | :---: |
+| | |
+| | |
+
+图 1、2、4 使用**合成示例数据**（图中已注明），图 3 使用仓库自带的 ERA5 1978 年数据。运行 `python examples/journal_figures.py` 可全部复现；最小示例见 `python examples/journal_quickstart.py`，用法见英文 [README](README.md#journal-toolkit) 与 [API 参考](docs/API.md#plot_functionjournal)。
## 快速开始
diff --git a/data/countries.shp.DESKTOP-1H76CSN.10272.3740.sr.lock b/data/countries.shp.DESKTOP-1H76CSN.10272.3740.sr.lock
deleted file mode 100644
index e69de29..0000000
diff --git a/docs/API.md b/docs/API.md
index acc93a2..5374629 100644
--- a/docs/API.md
+++ b/docs/API.md
@@ -73,9 +73,9 @@ Accepts all `open_field` options, plus:
| `colorbar` | `True` | Add a horizontal colorbar. Disable for shared custom colorbars. |
| `plotfunc` | `"pcolormesh"` | `pcolormesh` or `contourf`. |
| `ax` | `None` | An existing Cartopy GeoAxes; overrides `projection`. Its figure background remains caller-controlled. |
-| `figsize` | `(10, 5.8)` | Inches, used only when creating a new figure. |
+| `figsize` | `(7.2, 4.2)` | Inches (≈183 mm double column), used only when creating a new figure. |
| `output` | `None` | Optional path to save immediately. Parent directories are created. |
-| `dpi` | `180` | Resolution for immediate export. |
+| `dpi` | `300` | Resolution for immediate export; use 600 for print. |
Regional extents control the viewport; they do not crop, mask, or resample the
returned data. For very large arrays, select/crop the data before plotting.
@@ -89,7 +89,7 @@ returned data. For very large arrays, select/crop the data before plotting.
| `.artist` | The plotted mappable, usable for a shared colorbar. |
| `.colorbar` | Colorbar or `None`. |
| `.data` | The normalized 2D field, including explicitly converted values. |
-| `.save(path, dpi=180, **kwargs)` | Save with a tight bounding box and the figure's background. Additional options go to `Figure.savefig`. Returns a `Path`. |
+| `.save(path, dpi=300, **kwargs)` | Save with a tight bounding box and the figure's background. Additional options go to `Figure.savefig`. Returns a `Path`. |
```python
result.axes.text(0.02, 0.02, "ERA5 · 1978", transform=result.axes.transAxes)
@@ -137,3 +137,45 @@ preserving the original API's behavior.
`add_profile`, `add_distribution`, `add_significance`, and `add_features` accept
an option object or its keyword arguments. See [Research figures](RESEARCH_FIGURES.md)
for styles, weights, scope, algorithms, defaults, and examples.
+
+## Data contract
+
+- **Rectilinear geographic grids:** latitude and longitude must each be one-dimensional, finite, and contain at least two distinct values. CF `standard_name` / geographic `units` and common `lat`/`lon` or `latitude`/`longitude` names are recognized.
+- **Explicit scientific choices:** multiple fields require `variable=`; nonspatial dimensions with more than one value require selection or reduction. Reductions are unweighted and skip missing values. Supported statistics: `mean`, `median`, `min`, `max`, `sum`, `std`.
+- **Coordinate normalization:** latitude is sorted, longitude is wrapped to `[-180, 180)` and sorted; a redundant global 0°/360° endpoint is removed. Other duplicates are rejected.
+- **Missing values:** xarray decodes NetCDF fill values; nonfinite values are masked. An entirely missing field produces an actionable error.
+- **Bounded scope:** curvilinear grids, unstructured meshes, projected x/y grids in metres, and antimeridian-crossing regional grids/extents need preprocessing. The API does not reproject source grids or infer units, vertical levels, time weights, or area weights.
+- **Memory:** the selected/reduced 2D field is loaded before the file is closed. For large workflows, preprocess with xarray and pass a prepared `DataArray`.
+
+## `plot_function.journal`
+
+Publication helpers that decorate *any* Matplotlib/Cartopy axes. They never
+modify global `rcParams`; wrap figure code in `journal_style()`. Gallery:
+[`examples/journal_figures.py`](../examples/journal_figures.py); minimal example:
+[`examples/journal_quickstart.py`](../examples/journal_quickstart.py).
+
+| Function | Key arguments | Returns |
+| :--- | :--- | :--- |
+| `journal_style(**overrides)` | any rcParams | context manager (Arial → Helvetica → Liberation Sans → Nimbus Sans → DejaVu Sans, 7 pt base, 0.5 pt lines, `pdf.fonttype=42`) |
+| `figsize(width="double", ratio=0.55)` | `"single"` (89 mm), `"onehalf"` (120 mm), `"double"` (183 mm) or mm | `(w, h)` inches |
+| `discrete_cmap(cmap, levels, extend="both")` | colormap name/object or explicit colour list (`len(levels)-1` + number of extend ends) | `(ListedColormap, BoundaryNorm)` |
+| `add_colorbar(mappable, ax, label=, title=, orientation=, bounds=, tick_every=)` | `bounds=[x, y, w, h]` in figure fraction, or placement relative to `ax` | `Colorbar` |
+| `geo_ticks(ax, xticks=, yticks=, lat_labels=, lon_labels=, gridlines=)` | degree tick positions | gridliner list |
+| `add_land(ax, color=LAND, borders=True, resolution="110m")` | Natural Earth (downloads once) | feature artists |
+| `add_stippling(ax, lon, lat, mask, stride=2, size=1.2, offset=True)` | 1-D or 2-D lon/lat, boolean mask | `PathCollection` or `None` |
+| `add_hatching(ax, lon, lat, mask, hatch="....")` | | contour set or `None` |
+| `hatch_patch(hatch)` | | legend handle |
+| `add_inset_bars(ax, labels, values, hatched=, colors=, bounds=, frame="box")` | `hatched` = uncertain part of each bar | inset axes |
+| `add_inset_histogram(ax, groups, bins=, log=, cumulative=, bounds=, frame="open")` | array or `{label: values}`; `cumulative=True` or weights | `(hist_axes, twin_axes or None)` |
+| `add_latitude_profile(ax, lat, center, lower=, upper=, outer=, members=, width=, pad=)` | values per latitude | profile axes, y axis locked to the map |
+| `add_longitude_profile(ax, lon, center, lower=, upper=)` | values per longitude | profile axes above the map |
+| `add_panel_label(ax, "a", style="({})")` | | `Text` |
+| `add_size_legend(ax, values, labels, scale=)` | `scale(value) -> marker area (pt²)` | `Legend` |
+| `save_figure(fig, path, dpi=600, formats=None)` | e.g. `formats=("png", "pdf")` | path(s) |
+
+**Profile alignment.** Profile coordinates are transformed with
+`ax.projection.transform_points` at the central longitude (or the centre of
+the extent), and the profile's limits are synchronised with the map's at draw
+time. On pseudo-cylindrical projections (Robinson, Equal Earth, Mollweide)
+latitude maps to a single y value, so alignment is exact; on conic projections
+it is exact only along the reference meridian.
diff --git a/docs/assets/README.md b/docs/assets/README.md
index 2d271ca..8de7de2 100644
--- a/docs/assets/README.md
+++ b/docs/assets/README.md
@@ -28,12 +28,25 @@ plots, not image-generated depictions of data.
The code's MIT license does not replace the applicable terms of source datasets.
-## Decorative branding
+## Branding
-`header.png` and `logo.png` are AI-generated decorative artwork created for this
-project. They are separate assets, both designed on white backgrounds; the header
-does not embed the standalone logo. The abstract contour shapes are not observations, model outputs, or
-geographic boundaries. They are deliberately separate from the computed gallery.
+The logo and banner live in [`docs/brand/`](../brand/). They are hand-authored
+SVG line drawings (white background) of the library's signature layout — a map
+with isolines and stippling, a lower-left histogram inset, a right-hand latitude
+profile and a discrete colourbar. They are decorative and contain no data.
+PNG exports (`logo.png`, `banner.png`) were rendered from the SVGs with CairoSVG;
+the banner's text uses Inter and Noto Sans CJK SC when rendered.
+
+## Journal gallery (0.4)
+
+`docs/gallery/*.png` are generated by `examples/journal_figures.py`:
+
+- `fig1_regimes.png`, `fig2_sites.png`, `fig4_regional_significance.png`:
+ **synthetic example data** (seeded random fields, pseudo-ensemble members,
+ random sites and z-test p-values). They demonstrate layout only and carry no
+ scientific meaning; each figure says so in its footnote.
+- `fig3_era5_profile.png`: the bundled ERA5 1978 file, July minus the
+ 12-month mean, with inter-month IQR/range bands on the zonal profile.
## Reproduction
diff --git a/docs/assets/east-asia.png b/docs/assets/east-asia.png
index 41abdf0..10ff306 100644
Binary files a/docs/assets/east-asia.png and b/docs/assets/east-asia.png differ
diff --git a/docs/assets/global-temperature.png b/docs/assets/global-temperature.png
index 5a98d4c..e5b5b73 100644
Binary files a/docs/assets/global-temperature.png and b/docs/assets/global-temperature.png differ
diff --git a/docs/assets/header.png b/docs/assets/header.png
deleted file mode 100644
index a3a5332..0000000
Binary files a/docs/assets/header.png and /dev/null differ
diff --git a/docs/assets/journal-global.png b/docs/assets/journal-global.png
index 465f5e8..d39a2e2 100644
Binary files a/docs/assets/journal-global.png and b/docs/assets/journal-global.png differ
diff --git a/docs/assets/journal-regional.png b/docs/assets/journal-regional.png
index bfc3e33..817163c 100644
Binary files a/docs/assets/journal-regional.png and b/docs/assets/journal-regional.png differ
diff --git a/docs/assets/journal-significance.png b/docs/assets/journal-significance.png
index 6dec415..b150f38 100644
Binary files a/docs/assets/journal-significance.png and b/docs/assets/journal-significance.png differ
diff --git a/docs/assets/logo.png b/docs/assets/logo.png
deleted file mode 100644
index 8d3cb60..0000000
Binary files a/docs/assets/logo.png and /dev/null differ
diff --git a/docs/assets/seasonal-contrast.png b/docs/assets/seasonal-contrast.png
index 8a09f54..0f09600 100644
Binary files a/docs/assets/seasonal-contrast.png and b/docs/assets/seasonal-contrast.png differ
diff --git a/docs/assets/seasons.png b/docs/assets/seasons.png
index 03bbb69..0b2226f 100644
Binary files a/docs/assets/seasons.png and b/docs/assets/seasons.png differ
diff --git a/docs/brand/banner.png b/docs/brand/banner.png
new file mode 100644
index 0000000..57e72f4
Binary files /dev/null and b/docs/brand/banner.png differ
diff --git a/docs/brand/banner.svg b/docs/brand/banner.svg
new file mode 100644
index 0000000..6b1bf2f
--- /dev/null
+++ b/docs/brand/banner.svg
@@ -0,0 +1,56 @@
+
diff --git a/docs/brand/logo.png b/docs/brand/logo.png
new file mode 100644
index 0000000..abbddd7
Binary files /dev/null and b/docs/brand/logo.png differ
diff --git a/docs/brand/logo.svg b/docs/brand/logo.svg
new file mode 100644
index 0000000..d95e083
--- /dev/null
+++ b/docs/brand/logo.svg
@@ -0,0 +1,43 @@
+
diff --git a/docs/gallery/fig1_regimes.png b/docs/gallery/fig1_regimes.png
new file mode 100644
index 0000000..1251584
Binary files /dev/null and b/docs/gallery/fig1_regimes.png differ
diff --git a/docs/gallery/fig2_sites.png b/docs/gallery/fig2_sites.png
new file mode 100644
index 0000000..d1c2e3c
Binary files /dev/null and b/docs/gallery/fig2_sites.png differ
diff --git a/docs/gallery/fig3_era5_profile.png b/docs/gallery/fig3_era5_profile.png
new file mode 100644
index 0000000..fa9b971
Binary files /dev/null and b/docs/gallery/fig3_era5_profile.png differ
diff --git a/docs/gallery/fig4_regional_significance.png b/docs/gallery/fig4_regional_significance.png
new file mode 100644
index 0000000..8470df6
Binary files /dev/null and b/docs/gallery/fig4_regional_significance.png differ
diff --git a/examples/journal_figures.py b/examples/journal_figures.py
new file mode 100644
index 0000000..ebfb081
--- /dev/null
+++ b/examples/journal_figures.py
@@ -0,0 +1,389 @@
+"""Journal-grade figure recipes built with ``plot_function.journal``.
+
+Every figure below is reproducible from this script alone. Figures 1, 2 and 4
+use **synthetic example data** (seeded random fields / sites, clearly labelled in
+the figures); figure 3 uses the bundled ERA5 1978 monthly 2 m temperature file.
+
+ python examples/journal_figures.py # writes docs/gallery/*.png
+ python examples/journal_figures.py --output out/ # elsewhere
+ python examples/journal_figures.py --dpi 600 --pdf # print-ready
+"""
+
+import argparse
+import warnings
+from pathlib import Path
+
+import matplotlib
+
+matplotlib.use("Agg")
+import cartopy.crs as ccrs # noqa: E402
+import matplotlib.pyplot as plt # noqa: E402
+import numpy as np # noqa: E402
+import xarray as xr # noqa: E402
+from matplotlib.lines import Line2D # noqa: E402
+from matplotlib.patches import Patch, Rectangle # noqa: E402
+
+import plot_function.journal as pj # noqa: E402
+
+ROOT = Path(__file__).resolve().parents[1]
+PC = ccrs.PlateCarree()
+
+
+# ----------------------------------------------------------------------------- data helpers
+
+
+def land_mask(lon, lat):
+ """Boolean land mask on a lon/lat grid from Natural Earth 110m land (cached by Cartopy)."""
+ import cartopy.io.shapereader as shpreader
+ import shapely
+ from shapely.ops import unary_union
+
+ geoms = list(shpreader.Reader(shpreader.natural_earth("110m", "physical", "land")).geometries())
+ land = unary_union(geoms)
+ shapely.prepare(land)
+ xx, yy = np.meshgrid(lon, lat)
+ return shapely.contains_xy(land, xx, yy)
+
+
+def smooth_field(rng, lon, lat, n=60, scale=(8, 22)):
+ """Sum of random Gaussian blobs: a smooth synthetic field in roughly [-1, 1]."""
+ xx, yy = np.meshgrid(lon, lat)
+ out = np.zeros_like(xx, dtype=float)
+ for _ in range(n):
+ x0 = rng.uniform(np.min(lon) - 10, np.max(lon) + 10)
+ y0 = rng.uniform(max(np.min(lat), -60) - 5, min(np.max(lat), 80) + 5)
+ sx, sy = rng.uniform(*scale), rng.uniform(*scale) * 0.7
+ dx = (xx - x0 + 180) % 360 - 180
+ out += rng.normal() * np.exp(-((dx / sx) ** 2) - ((yy - y0) / sy) ** 2)
+ return out / np.nanmax(np.abs(out))
+
+
+def footnote(fig, text, y=-0.01):
+ fig.text(0.01, y, text, fontsize=5.5, color=pj.MUTED, ha="left", va="top")
+
+
+# ----------------------------------------------------------------------------- figure 1
+
+
+def fig_regimes(out, dpi, pdf):
+ """Two-regime diverging map + stippling + stacked inset bars + aligned zonal profile."""
+ rng = np.random.default_rng(42)
+ lon = np.arange(-179.5, 180, 1.0)
+ lat = np.arange(-59.5, 90, 1.0)
+ land = land_mask(lon, lat)
+ land &= lat[:, None] > -58
+ n_mem = 12
+ signal = 12 * smooth_field(rng, lon, lat, n=160, scale=(5, 16))
+ shared = 5 * smooth_field(rng, lon, lat, n=200, scale=(2, 6)) # small-scale texture
+ members = np.stack(
+ [signal + shared + 9 * smooth_field(rng, lon, lat, n=90, scale=(6, 20))
+ for _ in range(n_mem)]
+ )
+ members[:, ~land] = np.nan
+ delta = np.nanmean(members, axis=0)
+ agree = (np.sign(members) == np.sign(delta)).sum(0)
+ low_agree = land & (agree < 8)
+ regime = smooth_field(rng, lon, lat, n=50, scale=(15, 35)) > 0 # True: blue, False: green
+ blue = np.where(land & regime, delta, np.nan)
+ green = np.where(land & ~regime, -delta, np.nan)
+
+ levels = np.arange(-10, 10.1, 2.5)
+ brown_blue = ["#8c6d46", "#a68a62", "#cdb994", "#ece2c8", "#fbf6e6",
+ "#e9f3fb", "#b6d5f0", "#7cb2e6", "#3f86d9", "#2a6fbf"]
+ green_purple = ["#1f7a1f", "#3d9a3d", "#77b96f", "#b9dbb0", "#eef7ea",
+ "#e7e0ef", "#c8b5da", "#a68cc4", "#8566ad", "#6a4b95"]
+ cmap_b, norm = pj.discrete_cmap(brown_blue, levels, extend="both")
+ cmap_g, _ = pj.discrete_cmap(green_purple, levels, extend="both")
+
+ with pj.journal_style():
+ fig = plt.figure(figsize=pj.figsize("double", 0.49))
+ ax = fig.add_axes([0.065, 0.255, 0.79, 0.72], projection=PC)
+ ax.set_extent([-180, 180, -60, 90], crs=PC)
+ ax.coastlines("110m", lw=0.35, color=pj.INK, zorder=3)
+ ax.pcolormesh(lon, lat, blue, cmap=cmap_b, norm=norm, transform=PC, zorder=1)
+ mg = ax.pcolormesh(lon, lat, green, cmap=cmap_g, norm=norm, transform=PC, zorder=1)
+ mb = ax.collections[-2]
+ pj.geo_ticks(ax, xticks=range(-180, 181, 60), yticks=range(-40, 81, 20))
+ pj.add_stippling(ax, lon, lat, low_agree, stride=3, size=0.9)
+ pj.add_panel_label(ax, "a")
+
+ # lower-left inset: land area per regime, hatched part = low agreement
+ area = np.cos(np.deg2rad(lat))[:, None] * np.ones_like(delta)
+ tot = np.nansum(area[land])
+ cats = {
+ "Greener": land & ~regime & (green < 0),
+ "Less\ngreen": land & ~regime & (green >= 0),
+ "Less\nblue": land & regime & (blue < 0),
+ "Bluer": land & regime & (blue >= 0),
+ }
+ vals = [100 * area[m].sum() / tot for m in cats.values()]
+ low = [100 * area[m & low_agree].sum() / tot for m in cats.values()]
+ pj.add_inset_bars(
+ ax,
+ list(cats),
+ vals,
+ hatched=low,
+ colors=["#4f9a4f", "#a68cc4", "#c9b58f", "#5ea3e6"],
+ bounds=(0.075, 0.13, 0.22, 0.36),
+ ylabel="Land area [%]",
+ )
+
+ # right marginal: zonal mean with ensemble IQR and min-max envelope
+ with warnings.catch_warnings():
+ warnings.simplefilter("ignore", RuntimeWarning)
+ zonal = np.nanmean(members, axis=2) # (member, lat)
+ with warnings.catch_warnings(), np.errstate(all="ignore"):
+ warnings.simplefilter("ignore", RuntimeWarning)
+ q25, q50, q75 = np.nanpercentile(zonal, [25, 50, 75], axis=0)
+ mn, mx = np.nanmin(zonal, 0), np.nanmax(zonal, 0)
+ ok = np.isfinite(q50)
+ pj.add_latitude_profile(
+ ax, lat[ok], q50[ok], lower=q25[ok], upper=q75[ok], outer=(mn[ok], mx[ok]),
+ width=0.115, pad=0.03, color="#2b5d8c", xlabel="Zonal mean\nΔBGWS [ppts]",
+ yticks=range(-40, 81, 20), xlim=(-8, 8),
+ )
+
+ # legend for stippling (lower right, inside map)
+ ax.legend(
+ [pj.hatch_patch("....", facecolor="white")],
+ ["Low ensemble\nsign agreement (<8/12)"],
+ loc="lower right", bbox_to_anchor=(0.995, 0.01), frameon=True, fancybox=False,
+ edgecolor="none", framealpha=0.9, handlelength=2.0, handleheight=1.4,
+ prop={"weight": "bold", "size": 6.5},
+ )
+ pj.add_colorbar(mb, ax, title="Historical blue water regime", label="ΔBGWS [ppts]",
+ bounds=[0.10, 0.11, 0.30, 0.03], tick_every=2)
+ pj.add_colorbar(mg, ax, title="Historical green water regime", label="ΔBGWS [ppts]",
+ bounds=[0.52, 0.11, 0.30, 0.03], tick_every=2)
+ footnote(fig, "Synthetic example data (seeded random fields, 12 pseudo-members); "
+ "layout after Nature-style blue/green-water regime maps.", y=0.0)
+ path = pj.save_figure(fig, out / "fig1_regimes.png", dpi=dpi)
+ if pdf:
+ pj.save_figure(fig, out / "fig1_regimes.pdf")
+ plt.close(fig)
+ return path
+
+
+# ----------------------------------------------------------------------------- figure 2
+
+
+def fig_sites(out, dpi, pdf):
+ """Cream-land site map, sized markers, log histogram + cumulative inset, legends below."""
+ rng = np.random.default_rng(7)
+ hubs = np.array([[-95, 37], [-75, 41], [10, 47], [32, 30], [77, 24], [105, 30],
+ [-47, -18], [28, -26], [140, -33], [3, 10], [120, 35], [55, 30],
+ [-99, 19], [78, 28], [110, -7], [-58, -34], [36, 0], [-2, 52],
+ [67, 33], [-77, 4], [44, 40], [126, 37], [100, 14]])
+
+ def sites(n, spread):
+ idx = rng.integers(0, len(hubs), n)
+ pts = hubs[idx] + rng.normal(0, spread, (n, 2))
+ return pts[:, 0], np.clip(pts[:, 1], -50, 70)
+
+ groups = {
+ "managed": (*sites(108, 6), rng.lognormal(0.3, 1.0, 108)),
+ "dumping": (*sites(43, 4), rng.lognormal(0.9, 0.8, 43)),
+ "urban": (*sites(46, 5), rng.lognormal(1.2, 0.6, 46)),
+ }
+ scale = lambda q: 9 + 22 * np.log10(1 + q) # marker area in pt^2 # noqa: E731
+ orange, purple = "#f28e1c", "#9b4fc2"
+
+ with pj.journal_style():
+ fig = plt.figure(figsize=pj.figsize("double", 0.5))
+ ax = fig.add_axes([0.0, 0.22, 1.0, 0.77], projection=ccrs.Robinson())
+ ax.set_extent([-180, 180, -57, 84], crs=PC)
+ pj.add_land(ax, borders=True, coast_width=0.3)
+ ax.spines["geo"].set_visible(False)
+ x, y, q = groups["urban"]
+ ax.scatter(x, y, s=scale(q) * 1.3, marker="D", facecolor="#a9a9a9", edgecolor="#6f6f6f",
+ lw=0.6, transform=PC, zorder=5)
+ x, y, q = groups["managed"]
+ ax.scatter(x, y, s=scale(q), marker="^", facecolor="#ffd9a8", edgecolor=orange, lw=0.8,
+ transform=PC, zorder=6)
+ x, y, q = groups["dumping"]
+ ax.scatter(x, y, s=scale(q), marker="o", facecolor="#e2c4f0", edgecolor=purple, lw=0.8,
+ transform=PC, zorder=6)
+
+ rates = np.r_[groups["managed"][2], groups["dumping"][2]]
+ pj.add_inset_histogram(
+ ax,
+ {"Managed landfills": groups["managed"][2], "Dumping sites": groups["dumping"][2]},
+ bins=18, log=True, colors=["#fbb260", "#b37fd6"], alpha=0.85,
+ total_label="All waste disposal sites",
+ cumulative=rates, cumulative_label="Cumulative emissions\n(t h$^{-1}$)",
+ bounds=(0.10, 0.11, 0.165, 0.38),
+ xlabel="Methane emission rate (t h$^{-1}$)", ylabel="Number of sites",
+ )
+ # legends below the map (three tidy blocks)
+ lax = fig.add_axes([0.0, 0.0, 1.0, 0.2])
+ lax.axis("off")
+ h1 = [Patch(facecolor="none", edgecolor=pj.INK, lw=0.8),
+ Patch(facecolor="#fbb260", alpha=0.85), Patch(facecolor="#b37fd6", alpha=0.85)]
+ leg1 = lax.legend(h1, ["All waste disposal sites", "Managed landfills", "Dumping sites"],
+ loc="upper left", bbox_to_anchor=(0.06, 1.0))
+ lax.add_artist(leg1)
+ mk = dict(ls="none", markersize=6.5, markeredgewidth=0.8)
+ h2 = [Line2D([], [], marker="D", markerfacecolor="#a9a9a9", markeredgecolor="#6f6f6f", **mk),
+ Line2D([], [], marker="^", markerfacecolor="#ffd9a8", markeredgecolor=orange, **mk),
+ Line2D([], [], marker="o", markerfacecolor="#e2c4f0", markeredgecolor=purple, **mk)]
+ lax.legend(h2, [f"Urban plume clusters (N = {len(groups['urban'][2])})",
+ f"Managed landfills (N = {len(groups['managed'][2])})",
+ f"Dumping sites (N = {len(groups['dumping'][2])})"],
+ loc="upper left", bbox_to_anchor=(0.38, 1.0))
+ sax = fig.add_axes([0.74, 0.0, 0.25, 0.2])
+ sax.axis("off")
+ pj.add_size_legend(sax, [1, 10], ["1 t h$^{-1}$", "10 t h$^{-1}$"], scale=scale,
+ marker="o", edgecolor=purple, title="Site emissions",
+ loc="upper left", ncol=2, columnspacing=1.2)
+ pj.add_panel_label(ax, "b", x=0.0, y=1.0)
+ footnote(fig, "Synthetic example data (seeded random sites and emission rates); "
+ "layout after Nature-style point-source maps.", y=0.05)
+ path = pj.save_figure(fig, out / "fig2_sites.png", dpi=dpi)
+ if pdf:
+ pj.save_figure(fig, out / "fig2_sites.pdf")
+ plt.close(fig)
+ return path
+
+
+# ----------------------------------------------------------------------------- figure 3
+
+
+def load_era5():
+ with xr.open_dataset(ROOT / "data/ERA5temp_1978_monthly.nc") as ds:
+ t = ds.t2m.isel(latitude=slice(None, None, 4), longitude=slice(None, None, 4)).load()
+ t = (t - 273.15).rename(latitude="lat", longitude="lon")
+ t = t.assign_coords(lon=((t.lon + 180) % 360) - 180).sortby("lon").sortby("lat")
+ return t
+
+
+def fig_era5(out, dpi, pdf):
+ """Real ERA5 1978: July anomaly from the annual mean, Robinson, aligned zonal profile."""
+ t = load_era5()
+ anom = t.isel(time=6) - t.mean("time")
+ monthly_anom = t - t.mean("time")
+ levels = np.arange(-20, 20.1, 4)
+ cmap, norm = pj.discrete_cmap("RdBu_r", levels, extend="both")
+ with pj.journal_style():
+ fig = plt.figure(figsize=pj.figsize("double", 0.53))
+ ax = fig.add_axes([0.04, 0.21, 0.75, 0.76], projection=ccrs.Robinson())
+ ax.set_global()
+ m = ax.pcolormesh(anom.lon, anom.lat, anom, cmap=cmap, norm=norm, transform=PC,
+ rasterized=True)
+ ax.coastlines("110m", lw=0.35, color=pj.INK)
+ # latitude labels live on the aligned profile panel, so the map only labels longitude
+ pj.geo_ticks(ax, xticks=range(-180, 181, 60), yticks=range(-60, 61, 30), lat_labels=False)
+ pj.add_panel_label(ax, "c", x=0.0, y=1.0)
+ w = np.cos(np.deg2rad(anom.lat))
+ # zonal mean of the July anomaly and the spread of all 12 months' zonal means
+ zm = anom.mean("lon")
+ q = monthly_anom.mean("lon").quantile([0.25, 0.75], "time")
+ mm = monthly_anom.mean("lon")
+ pj.add_latitude_profile(
+ ax, anom.lat, zm, lower=q.sel(quantile=0.25), upper=q.sel(quantile=0.75),
+ outer=(mm.min("time"), mm.max("time")), color="#b2182b", width=0.16, pad=0.02,
+ xlabel="Zonal mean (°C)", yticks=range(-60, 61, 30), xlim=(-22, 22),
+ title="Jul (line) · monthly IQR / range",
+ )
+ vals = anom.values.ravel()
+ okv = np.isfinite(vals)
+ child, _ = pj.add_inset_histogram(
+ ax, vals[okv], bins=np.arange(-24, 24.1, 2), colors=["#9e9e9e"], total=False,
+ stacked=False, bounds=(0.035, 0.06, 0.17, 0.25), xlabel="Anomaly (°C)",
+ ylabel="Cells", alpha=1.0,
+ )
+ # colour the histogram bars like the map
+ for patch in child.patches:
+ if not isinstance(patch, Rectangle):
+ continue
+ xc = patch.get_x() + patch.get_width() / 2
+ patch.set_facecolor(cmap(norm(xc)))
+ area_warm = float((w * (anom > 0)).sum() / (w * anom.notnull()).sum())
+ child.text(1.0, 1.02, f"{100 * area_warm:.0f}% of area warmer", transform=child.transAxes,
+ ha="right", va="bottom", fontsize=6)
+ pj.add_colorbar(m, ax, label="July minus annual-mean 2 m temperature (°C)",
+ bounds=[0.205, 0.095, 0.42, 0.026])
+ footnote(fig, "Data: ERA5 monthly 2 m temperature, 1978 (bundled file, 1° display grid). "
+ "Profile band: inter-month IQR; light band: inter-month range.", y=0.01)
+ path = pj.save_figure(fig, out / "fig3_era5_profile.png", dpi=dpi)
+ if pdf:
+ pj.save_figure(fig, out / "fig3_era5_profile.pdf")
+ plt.close(fig)
+ return path
+
+
+# ----------------------------------------------------------------------------- figure 4
+
+
+def fig_regional(out, dpi, pdf):
+ """Two regional panels sharing one discrete scale: stippling vs hatching (synthetic p)."""
+ rng = np.random.default_rng(2026)
+ lon = np.arange(70.25, 140, 0.5)
+ lat = np.arange(15.25, 55, 0.5)
+ xx, yy = np.meshgrid(lon, lat)
+ trend = (1.2 * np.exp(-((xx - 110) / 14) ** 2 - ((yy - 33) / 8) ** 2)
+ - 0.9 * np.exp(-((xx - 88) / 9) ** 2 - ((yy - 42) / 6) ** 2)
+ + 0.45 * smooth_field(rng, lon, lat, n=160, scale=(2, 6))
+ + 0.35 * smooth_field(rng, lon, lat, n=60, scale=(6, 14)))
+ se = 0.35 + 0.1 * rng.random(trend.shape)
+ from math import erfc
+ p = np.vectorize(lambda z: erfc(abs(z) / np.sqrt(2)))(trend / se)
+ land = land_mask(lon, lat)
+ trend[~land] = np.nan
+ levels = np.arange(-1.5, 1.51, 0.3)
+ cmap, norm = pj.discrete_cmap("BrBG", levels, extend="both")
+ with pj.journal_style():
+ fig = plt.figure(figsize=pj.figsize("double", 0.36))
+ axes = [fig.add_axes([0.045 + 0.505 * i, 0.25, 0.355, 0.72], projection=PC)
+ for i in range(2)]
+ for i, (ax, style) in enumerate(zip(axes, ["stipple", "hatch"])):
+ ax.set_extent([70, 140, 15, 55], crs=PC)
+ m = ax.pcolormesh(lon, lat, trend, cmap=cmap, norm=norm, transform=PC)
+ pj.add_land(ax, color="none", borders=True, coast_width=0.35, resolution="50m",
+ zorder=2)
+ pj.geo_ticks(ax, xticks=range(80, 131, 20), yticks=range(20, 51, 10),
+ lat_labels=(i == 0)) # panel e shares latitude with d's profile
+ sig = (p < 0.05) & land
+ if style == "stipple":
+ pj.add_stippling(ax, lon, lat, sig, stride=3, size=0.7)
+ handle = Line2D([], [], ls="none", marker="o", ms=1.6, color=pj.INK)
+ label = "p < 0.05 (dots)"
+ else:
+ pj.add_hatching(ax, lon, lat, sig, hatch="////", linewidth=0.35)
+ handle = pj.hatch_patch("////", facecolor="white")
+ label = "p < 0.05 (hatch)"
+ ax.legend([handle], [label], loc="lower right", frameon=True, edgecolor="none",
+ framealpha=0.85)
+ zm = np.nanmean(trend, axis=1)
+ sd = np.nanstd(trend, axis=1)
+ prof = pj.add_latitude_profile(ax, lat, zm, lower=zm - sd, upper=zm + sd, width=0.17,
+ pad=0.045, color="#01665e", xlabel="Mean ± SD",
+ yticks=range(20, 51, 10), xlim=(-1.6, 1.6))
+ prof.set_xticks([-1, 0, 1])
+ pj.add_panel_label(ax, "de"[i], x=0.01, y=0.98)
+ pj.add_colorbar(m, axes[0], label="Trend (mm d$^{-1}$ decade$^{-1}$)",
+ bounds=[0.3, 0.10, 0.4, 0.04])
+ footnote(fig, "Synthetic example data (seeded trend field and z-test p-values).", y=0.0)
+ path = pj.save_figure(fig, out / "fig4_regional_significance.png", dpi=dpi)
+ if pdf:
+ pj.save_figure(fig, out / "fig4_regional_significance.pdf")
+ plt.close(fig)
+ return path
+
+
+FIGURES = {"regimes": fig_regimes, "sites": fig_sites, "era5": fig_era5, "regional": fig_regional}
+
+
+def main():
+ parser = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawTextHelpFormatter)
+ parser.add_argument("--output", type=Path, default=ROOT / "docs/gallery")
+ parser.add_argument("--dpi", type=int, default=300)
+ parser.add_argument("--pdf", action="store_true", help="also write vector PDFs")
+ parser.add_argument("--only", choices=list(FIGURES), nargs="*")
+ args = parser.parse_args()
+ args.output.mkdir(parents=True, exist_ok=True)
+ for name in args.only or FIGURES:
+ print("wrote", FIGURES[name](args.output, args.dpi, args.pdf))
+
+
+if __name__ == "__main__":
+ main()
diff --git a/examples/journal_gallery.py b/examples/journal_gallery.py
index 8aaa34f..1f031b8 100644
--- a/examples/journal_gallery.py
+++ b/examples/journal_gallery.py
@@ -19,12 +19,14 @@
from plot_function import Distribution, MapFeatures, Profile, Significance, plot_map
ROOT = Path(__file__).resolve().parents[1]
-INK, MUTED, TEAL = "#163438", "#647b80", "#267f83"
+INK, MUTED, TEAL = "#1a1a1a", "#555555", "#2b5d8c"
def footer(fig, text):
- fig.text(0.055, 0.025, text, fontsize=8, color=MUTED)
- fig.text(0.95, 0.025, "plot-function / research atlas", ha="right", fontsize=8, color=INK)
+ from plot_function.journal import resolved_font
+
+ fig.text(0.01, 0.005, text, fontsize=5.5, color=MUTED, ha="left", va="bottom",
+ fontfamily=resolved_font())
def synthetic_experiment(path):
@@ -72,26 +74,24 @@ def gallery(output, *, coastlines=True, map_details=True):
reduce="time",
cmap="RdYlBu_r",
levels=np.arange(-40, 41, 5),
- title="Temperature, in space and in distribution",
- subtitle="ERA5 / 1978 · Mean of 12 monthly means · 1° display grid",
- figsize=(13, 7.5),
+ title="Annual-mean 2 m air temperature",
+ subtitle="ERA5 1978 · mean of 12 monthly means · 1° display grid",
+ label="Temperature (°C)",
+ figsize=(7.2, 4.0),
coastlines=coastlines,
- profiles=Profile(style="band", spread="std", reference=0, width=0.23),
+ profiles=Profile(style="band", spread="std", reference=0, width=0.15,
+ label="Zonal mean ± 1 s.d."),
distribution=Distribution(
style="bars",
weights="coslat",
color="map",
- bounds=(0.035, 0.07, 0.26, 0.25),
- label="Cos-lat weighted density",
+ bounds=(0.03, 0.06, 0.17, 0.22),
),
panel_label="a",
)
- global_map.figure.set_facecolor("white")
- footer(
- global_map.figure,
- "SPATIAL SUMMARY / Ribbon: ±1 spatial SD, not a confidence interval",
- )
- global_map.save(output / "journal-global.png", dpi=180)
+ footer(global_map.figure, "Band: ±1 spatial s.d. across longitudes (not a confidence "
+ "interval). Inset: cos-latitude-weighted density.")
+ global_map.save(output / "journal-global.png", dpi=300)
plt.close(global_map.figure)
regional = plot_map(
@@ -100,59 +100,46 @@ def gallery(output, *, coastlines=True, map_details=True):
extent=[92, 142, 8, 53],
projection=ccrs.PlateCarree(),
cmap="RdYlBu_r",
- levels=np.arange(-10, 37, 2),
- title="East Asia / climate in context",
- subtitle="ERA5 / July 1978 · Regional summaries within the displayed bounds",
- figsize=(10, 9),
+ levels=np.arange(-10, 37, 3),
+ title="East Asia, July 1978",
+ subtitle="ERA5 2 m air temperature · summaries within the displayed bounds",
+ label="Temperature (°C)",
+ figsize=(5.2, 5.0),
coastlines=coastlines,
profiles=[
- Profile(style="line", color=TEAL),
+ Profile(style="band", spread="std", color=TEAL, label="Zonal mean ± 1 s.d."),
Profile(
position="top",
style="band",
statistic="median",
spread="iqr",
- color="#ad7843",
- width=0.18,
+ color="#a6611a",
+ width=0.16,
+ label="Meridional median (IQR)",
),
],
distribution=Distribution(
- style="step", bins=20, weights="coslat", bounds=(0.045, 0.09, 0.30, 0.24)
+ style="step", bins=20, weights="coslat", bounds=(0.06, 0.07, 0.3, 0.2)
),
features=MapFeatures(
- resolution="50m", rivers=map_details, lakes=map_details, borders=map_details
+ resolution="50m", rivers=map_details, lakes=map_details, borders=map_details,
+ river_color="#4a90b8", border_color="#8c8c8c",
),
panel_label="b",
)
- regional.figure.set_facecolor("white")
- footer(
- regional.figure,
- "REGIONAL ATLAS / Top: median + spatial IQR · Right: zonal mean",
- )
- regional.save(output / "journal-regional.png", dpi=180)
+ regional.save(output / "journal-regional.png", dpi=300)
plt.close(regional.figure)
synthetic = ROOT / "examples/output/synthetic-significance.nc"
synthetic_experiment(synthetic)
- fig, axes = plt.subplots(
- 1, 3, figsize=(16, 6.0), subplot_kw={"projection": ccrs.PlateCarree()}, facecolor="white"
- )
- fig.subplots_adjust(left=0.04, right=0.94, bottom=0.24, top=0.77, wspace=0.54)
- fig.text(
- 0.04,
- 0.94,
- "One statistical decision. Three visual languages.",
- fontsize=23,
- fontweight="bold",
- color=INK,
- )
- fig.text(
- 0.04,
- 0.865,
- "SYNTHETIC EXPERIMENT / 40 independent realizations · Two-sided z-test, known σ · BH FDR q = 0.05",
- fontsize=10,
- color=MUTED,
- )
+ import plot_function.journal as pj
+
+ with pj.journal_style():
+ fig, axes = plt.subplots(
+ 1, 3, figsize=(7.2, 2.7), subplot_kw={"projection": ccrs.PlateCarree()},
+ facecolor="white",
+ )
+ fig.subplots_adjust(left=0.06, right=0.9, bottom=0.3, top=0.9, wspace=0.55)
for ax, style, profile_style, dist_style, letter in zip(
axes,
["stipple", "hatch", "contour"],
@@ -167,43 +154,40 @@ def gallery(output, *, coastlines=True, map_details=True):
extent=[80, 150, 0, 60],
cmap="RdBu_r",
levels=np.linspace(-1.8, 1.8, 13),
- title=style.capitalize(),
+ title={"stipple": "Stippling", "hatch": "Hatching", "contour": "Contour"}[style],
coastlines=False,
colorbar=False,
- panel_label=letter,
significance=Significance(
synthetic,
variable="p_value",
style=style,
correction="fdr_bh",
stride=2,
- size=2.8,
- hatch="....",
+ size=2.2,
+ hatch="///",
legend=False,
),
profiles=Profile(
style=profile_style,
- width=0.19,
+ width=0.2,
reference=0,
- label="Mean" if profile_style != "band" else "Mean ± SD",
+ label="Mean" if profile_style != "band" else "Mean ± s.d.",
),
distribution=Distribution(
- style=dist_style, bins=14, bounds=(0.075, 0.065, 0.41, 0.25), show_mean=False
+ style=dist_style, bins=14, bounds=(0.1, 0.16, 0.36, 0.22), show_mean=False
),
)
- bar = fig.colorbar(
- result.artist, cax=fig.add_axes([0.32, 0.14, 0.36, 0.022]), orientation="horizontal"
- )
- bar.set_label("Synthetic response / arbitrary units", fontsize=9, color=INK)
- bar.ax.tick_params(labelsize=8, colors=MUTED, length=2)
- bar.outline.set_visible(False)
- footer(
- fig,
- "METHOD STUDY / Identical p-values and FDR family in all panels · No observed climate significance is implied",
- )
- fig.savefig(
- output / "journal-significance.png", dpi=180, bbox_inches="tight", facecolor="white"
- )
+ ax.annotate(letter, xy=(0, 0), xycoords=ax._left_title, xytext=(-4, 0),
+ textcoords="offset points", fontsize=10, fontweight="bold", ha="right",
+ fontfamily=pj.resolved_font())
+ if letter != "a": # latitude is labelled once, on panel a and on each profile
+ ax.tick_params(labelleft=False)
+ with pj.journal_style():
+ pj.add_colorbar(result.artist, axes[1], label="Synthetic response (arbitrary units)",
+ bounds=[0.3, 0.13, 0.4, 0.03], tick_every=2)
+ footer(fig, "Synthetic experiment: 40 realizations, two-sided z-test (known σ), "
+ "BH FDR q = 0.05. Identical decisions in all panels.")
+ pj.save_figure(fig, output / "journal-significance.png", dpi=300)
plt.close(fig)
print(f"Rendered three journal compositions to {output}")
diff --git a/examples/journal_quickstart.py b/examples/journal_quickstart.py
new file mode 100644
index 0000000..8df1039
--- /dev/null
+++ b/examples/journal_quickstart.py
@@ -0,0 +1,48 @@
+"""Ten-line journal figure: discrete map, stippling, inset bars, aligned zonal profile.
+
+Uses an analytic **example field** (no downloads except Natural Earth coastlines):
+ python examples/journal_quickstart.py
+"""
+
+from pathlib import Path
+
+import matplotlib
+
+matplotlib.use("Agg")
+import cartopy.crs as ccrs # noqa: E402
+import matplotlib.pyplot as plt # noqa: E402
+import numpy as np # noqa: E402
+
+import plot_function.journal as pj # noqa: E402
+
+
+def main(coastlines=True):
+ lon, lat = np.arange(-179.5, 180, 1.0), np.arange(-59.5, 90, 1.0)
+ field = (8 * np.sin(np.deg2rad(2 * lon))[None, :] * np.cos(np.deg2rad(lat))[:, None]
+ + 4 * np.sin(np.deg2rad(3 * lat))[:, None]) # example data
+ low_agreement = np.abs(field) < 1.5
+
+ with pj.journal_style():
+ fig = plt.figure(figsize=pj.figsize("double", 0.45))
+ ax = fig.add_axes([0.06, 0.24, 0.78, 0.72], projection=ccrs.PlateCarree())
+ ax.set_extent([-180, 180, -60, 90], crs=ccrs.PlateCarree())
+ cmap, norm = pj.discrete_cmap("BrBG", np.arange(-10, 10.1, 2.5), extend="both")
+ mesh = ax.pcolormesh(lon, lat, field, cmap=cmap, norm=norm,
+ transform=ccrs.PlateCarree())
+ if coastlines:
+ ax.coastlines(lw=0.3)
+ pj.geo_ticks(ax, xticks=range(-180, 181, 60), yticks=range(-40, 81, 20))
+ pj.add_stippling(ax, lon, lat, low_agreement, stride=3)
+ pj.add_inset_bars(ax, ["Drier", "Wetter"], [42, 58], hatched=[9, 12],
+ colors=["#a6611a", "#018571"], ylabel="Area [%]")
+ pj.add_latitude_profile(ax, lat, field.mean(1), lower=np.percentile(field, 25, 1),
+ upper=np.percentile(field, 75, 1), xlabel="Zonal mean")
+ pj.add_colorbar(mesh, ax, title="Example regime", label="Δ [units]")
+ pj.add_panel_label(ax, "a")
+ out = Path(__file__).resolve().parent / "output" / "journal_quickstart.png"
+ pj.save_figure(fig, out, dpi=300)
+ print(f"Created {out}")
+
+
+if __name__ == "__main__":
+ main()
diff --git a/plot_function/__init__.py b/plot_function/__init__.py
index 920ab23..3ec1f66 100644
--- a/plot_function/__init__.py
+++ b/plot_function/__init__.py
@@ -1,11 +1,13 @@
"""NetCDF-first plotting, built on the original plot-function map utilities."""
+from . import journal
from .data import open_field
from .maps import MapResult, plot_map
from .options import Distribution, MapFeatures, Profile, Significance
from .statistics import histogram, significance_mask, spatial_profile
__all__ = [
+ "journal",
"MapResult",
"open_field",
"plot_map",
@@ -17,4 +19,4 @@
"histogram",
"significance_mask",
]
-__version__ = "0.3.0"
+__version__ = "0.4.0"
diff --git a/plot_function/journal.py b/plot_function/journal.py
new file mode 100644
index 0000000..9864469
--- /dev/null
+++ b/plot_function/journal.py
@@ -0,0 +1,988 @@
+"""Publication ("journal") styling helpers that work with any Matplotlib/Cartopy axes.
+
+Every helper is functional: it receives the axes it should decorate and returns
+the artists or axes it created, so it composes with ``plot_map`` *and* with
+hand-built figures. Nothing here mutates global ``rcParams``; wrap your figure
+code in :func:`journal_style` to apply the typography only while drawing.
+
+Typical use::
+
+ import plot_function.journal as pj
+
+ with pj.journal_style():
+ fig = plt.figure(figsize=pj.figsize("double", 0.5))
+ ax = fig.add_subplot(projection=ccrs.PlateCarree())
+ cmap, norm = pj.discrete_cmap("BrBG", np.arange(-10, 11, 2.5), extend="both")
+ mesh = ax.pcolormesh(lon, lat, field, cmap=cmap, norm=norm, transform=ccrs.PlateCarree())
+ pj.geo_ticks(ax, xticks=range(-180, 181, 60), yticks=range(-40, 81, 20))
+ pj.add_stippling(ax, lon, lat, low_agreement)
+ pj.add_colorbar(mesh, ax, label="ΔBGWS [ppts]", title="Blue water regime")
+ pj.add_latitude_profile(ax, lat, zonal_mean, lower=q25, upper=q75)
+ pj.add_panel_label(ax, "a")
+ pj.save_figure(fig, "figure1.png", dpi=600)
+"""
+
+from __future__ import annotations
+
+import contextlib
+from pathlib import Path
+
+import cartopy.crs as ccrs
+import matplotlib as mpl
+import matplotlib.pyplot as plt
+import numpy as np
+from matplotlib.colors import BoundaryNorm, ListedColormap
+from matplotlib.lines import Line2D
+from matplotlib.patches import Patch
+from matplotlib.ticker import FixedLocator, FuncFormatter, MaxNLocator, NullLocator
+
+__all__ = [
+ "INK",
+ "MUTED",
+ "GRID",
+ "LAND",
+ "JOURNAL_RC",
+ "journal_style",
+ "figsize",
+ "discrete_cmap",
+ "clean_formatter",
+ "add_colorbar",
+ "geo_ticks",
+ "format_lat",
+ "format_lon",
+ "add_land",
+ "add_stippling",
+ "add_hatching",
+ "hatch_patch",
+ "add_inset_bars",
+ "add_inset_histogram",
+ "add_latitude_profile",
+ "add_longitude_profile",
+ "add_panel_label",
+ "add_size_legend",
+ "save_figure",
+]
+
+INK = "#1a1a1a"
+MUTED = "#555555"
+GRID = "#b4b4b4"
+LAND = "#fbf6e4" # cream land (Nature-style site maps)
+SANS = [
+ "Arial",
+ "Helvetica",
+ "Helvetica Neue",
+ "Liberation Sans",
+ "Nimbus Sans",
+ "TeX Gyre Heros",
+ "DejaVu Sans",
+]
+
+#: rcParams for 7–9 pt figures sized for a 89 mm (single) or 183 mm (double) column.
+JOURNAL_RC = {
+ # An explicit family list (not the generic "sans-serif") is stored on every Text
+ # object, so lazily created tick labels keep Arial/Helvetica after the context ends.
+ "font.family": SANS,
+ "font.sans-serif": SANS,
+ "font.size": 7,
+ "axes.titlesize": 8,
+ "axes.titleweight": "bold",
+ "axes.labelsize": 7,
+ "axes.linewidth": 0.5,
+ "axes.edgecolor": INK,
+ "axes.labelcolor": INK,
+ "axes.titlecolor": INK,
+ "axes.labelpad": 2.0,
+ "xtick.labelsize": 6.5,
+ "ytick.labelsize": 6.5,
+ "xtick.color": INK,
+ "ytick.color": INK,
+ "xtick.major.size": 2.5,
+ "ytick.major.size": 2.5,
+ "xtick.minor.size": 1.5,
+ "ytick.minor.size": 1.5,
+ "xtick.major.width": 0.5,
+ "ytick.major.width": 0.5,
+ "xtick.minor.width": 0.4,
+ "ytick.minor.width": 0.4,
+ "xtick.major.pad": 1.8,
+ "ytick.major.pad": 1.8,
+ "xtick.direction": "out",
+ "ytick.direction": "out",
+ "legend.fontsize": 6.5,
+ "legend.frameon": False,
+ "legend.handlelength": 1.4,
+ "legend.handletextpad": 0.5,
+ "legend.borderaxespad": 0.3,
+ "lines.linewidth": 0.9,
+ "patch.linewidth": 0.5,
+ "hatch.linewidth": 0.4,
+ "hatch.color": INK,
+ "figure.dpi": 150,
+ "savefig.dpi": 600,
+ "savefig.bbox": "tight",
+ "savefig.pad_inches": 0.02,
+ "savefig.facecolor": "white",
+ "figure.facecolor": "white",
+ "pdf.fonttype": 42, # editable text in Illustrator / Inkscape
+ "ps.fonttype": 42,
+ "svg.fonttype": "none",
+ "axes.unicode_minus": True,
+ "mathtext.default": "regular",
+}
+
+_MM = 1 / 25.4
+_WIDTHS = {"single": 89 * _MM, "onehalf": 120 * _MM, "double": 183 * _MM}
+
+
+def installed_fonts():
+ """The installed subset of the Arial/Helvetica fallback chain (DejaVu Sans last)."""
+ from matplotlib import font_manager
+
+ names = {f.name for f in font_manager.fontManager.ttflist}
+ found = [n for n in SANS if n in names]
+ return found if found and found[-1] == "DejaVu Sans" else found + ["DejaVu Sans"]
+
+
+def resolved_font():
+ """Name of the first installed font in the Arial/Helvetica fallback chain."""
+ return installed_fonts()[0]
+
+
+def journal_style(**overrides):
+ """Context manager applying :data:`JOURNAL_RC` (plus overrides) temporarily."""
+ rc = dict(JOURNAL_RC)
+ rc["font.family"] = rc["font.sans-serif"] = installed_fonts()
+ rc.update(overrides)
+ return mpl.rc_context(rc)
+
+
+def figsize(width="double", ratio=0.55):
+ """Figure size in inches for a journal column: ``'single'``, ``'onehalf'``, ``'double'``
+ or a width in millimetres; ``ratio`` is height / width."""
+ w = _WIDTHS[width] if isinstance(width, str) else float(width) * _MM
+ return (w, w * ratio)
+
+
+# --------------------------------------------------------------------------- colour
+
+
+def discrete_cmap(cmap, levels, *, extend="both", under=None, over=None, bad="none"):
+ """Return ``(cmap, norm)`` with one colour per interval and distinct extend colours.
+
+ ``levels`` are the interval boundaries. With ``extend='both'`` the colormap is
+ sampled at ``len(levels) + 1`` evenly spaced points so the triangles at the
+ ends of the colourbar get the darkest colours, exactly as in printed atlases.
+ A list of colours can be given instead of a colormap name.
+ """
+ levels = np.asarray(levels, dtype=float)
+ if levels.ndim != 1 or len(levels) < 2 or np.any(np.diff(levels) <= 0):
+ raise ValueError("levels must be at least two strictly increasing values.")
+ if extend not in ("neither", "min", "max", "both"):
+ raise ValueError("extend must be 'neither', 'min', 'max', or 'both'.")
+ n_int = len(levels) - 1
+ n_extra = {"neither": 0, "min": 1, "max": 1, "both": 2}[extend]
+ if isinstance(cmap, (list, tuple)):
+ colors = list(cmap)
+ if len(colors) != n_int + n_extra:
+ raise ValueError(f"Expected {n_int + n_extra} colours for {n_int} intervals.")
+ else:
+ base = mpl.colormaps[cmap] if isinstance(cmap, str) else cmap
+ colors = list(base(np.linspace(0, 1, n_int + n_extra)))
+ lo = 1 if extend in ("min", "both") else 0
+ inner = colors[lo : lo + n_int]
+ extremes = {"bad": bad}
+ if extend in ("min", "both"):
+ extremes["under"] = under if under is not None else colors[0]
+ if extend in ("max", "both"):
+ extremes["over"] = over if over is not None else colors[-1]
+ out = ListedColormap(inner, name=f"{getattr(cmap, 'name', cmap)}_discrete")
+ out = out.with_extremes(**extremes)
+ norm = BoundaryNorm(levels, n_int, extend="neither")
+ out.colorbar_extend = extend
+ return out, norm
+
+
+def add_colorbar(
+ mappable,
+ ax,
+ *,
+ label=None,
+ title=None,
+ orientation="horizontal",
+ extend=None,
+ ticks=None,
+ bounds=None,
+ size=0.035,
+ pad=0.12,
+ length=0.6,
+ anchor=0.5,
+ labelsize=None,
+ outline=True,
+ tick_every=1,
+):
+ """Slim colourbar with triangle ends, an optional bold title and a label.
+
+ The colourbar is placed relative to ``ax`` (axes fraction) unless explicit
+ figure-fraction ``bounds=[x, y, w, h]`` are given. ``title`` is drawn above a
+ horizontal bar (bold, like "Historical blue water regime") and ``label`` below.
+ """
+ fig = ax.figure
+ horizontal = orientation == "horizontal"
+ if bounds is not None:
+ cax = fig.add_axes(bounds)
+ elif horizontal:
+ pad = pad + (0.07 if title else 0.0)
+ cax = ax.inset_axes([anchor - length / 2, -pad - size, length, size], transform=ax.transAxes)
+ else:
+ cax = ax.inset_axes([1 + pad, 0.5 - length / 2, size, length], transform=ax.transAxes)
+ if extend is None:
+ extend = getattr(getattr(mappable, "cmap", None), "colorbar_extend", "neither")
+ extend = extend if extend in ("neither", "min", "max", "both") else "neither"
+ cbar = fig.colorbar(
+ mappable,
+ cax=cax,
+ orientation=orientation,
+ extend=extend,
+ extendfrac=0.05 if horizontal else 0.04,
+ drawedges=False,
+ )
+ norm = getattr(mappable, "norm", None)
+ if ticks is None and isinstance(norm, BoundaryNorm):
+ ticks = norm.boundaries[::tick_every]
+ if ticks is not None:
+ cbar.set_ticks(ticks)
+ cbar.ax.minorticks_off()
+ fmt = _tick_formatter(np.asarray(ticks, dtype=float))
+ (cbar.ax.xaxis if horizontal else cbar.ax.yaxis).set_major_formatter(fmt)
+ cbar.outline.set_visible(outline)
+ cbar.outline.set_linewidth(0.5)
+ cbar.outline.set_edgecolor(INK)
+ size_pt = labelsize or mpl.rcParams["xtick.labelsize"]
+ cbar.ax.tick_params(labelsize=size_pt, length=2.0, width=0.5, pad=1.5, color=INK)
+ for patch in getattr(cbar, "_extend_patches", []):
+ patch.set_edgecolor(INK)
+ patch.set_linewidth(0.5)
+ if label:
+ cbar.set_label(label, fontsize=mpl.rcParams["axes.labelsize"] + 0.5, labelpad=1.5)
+ if title:
+ if horizontal:
+ cax.set_title(title, fontsize=mpl.rcParams["axes.titlesize"], fontweight="bold", pad=3)
+ else:
+ cax.set_title(title, fontsize=mpl.rcParams["axes.labelsize"], pad=4, loc="left")
+ return cbar
+
+
+def clean_formatter():
+ """``%g`` tick labels with a true minus sign and no ``-0``."""
+ return FuncFormatter(lambda v, _: "0" if abs(v) < 1e-12 else f"{v:g}".replace("-", "\u2212"))
+
+
+def _tick_formatter(ticks):
+ finite = ticks[np.isfinite(ticks)]
+ decimals = 0
+ for d in range(0, 4):
+ if np.allclose(finite, np.round(finite, d)):
+ decimals = d
+ break
+ else:
+ decimals = 3
+ return FuncFormatter(
+ lambda v, _: "0" if abs(v) < 1e-12 else f"{v:.{decimals}f}".replace("-", "\u2212")
+ )
+
+
+# --------------------------------------------------------------------------- geography
+
+
+def format_lat(value, _=None):
+ if np.isclose(value, 0):
+ return "0°"
+ return f"{abs(value):g}°{'N' if value > 0 else 'S'}"
+
+
+def format_lon(value, _=None):
+ value = (value + 180) % 360 - 180 if abs(value) > 180 else value
+ if np.isclose(abs(value), 180):
+ return "180°"
+ if np.isclose(value, 0):
+ return "0°"
+ return f"{abs(value):g}°{'E' if value > 0 else 'W'}"
+
+
+def geo_ticks(
+ ax,
+ *,
+ xticks=None,
+ yticks=None,
+ gridlines=True,
+ grid_color=GRID,
+ grid_style=(0, (3, 3)),
+ grid_width=0.35,
+ labels=True,
+ lat_labels=True,
+ lon_labels=True,
+ frame=True,
+):
+ """Clean ``60°E`` / ``40°N`` ticks and light dashed gridlines on a GeoAxes.
+
+ Rectangular projections (PlateCarree, Mercator) get real outward ticks on
+ the left/bottom frame; other projections fall back to Cartopy gridline labels.
+ """
+ proj = ax.projection
+ xticks = list(xticks) if xticks is not None else None
+ yticks = list(yticks) if yticks is not None else None
+ rectangular = isinstance(proj, (ccrs.PlateCarree, ccrs.Mercator))
+ artists = []
+ if rectangular and labels:
+ if xticks is not None and lon_labels:
+ ax.set_xticks(xticks, crs=ccrs.PlateCarree())
+ ax.xaxis.set_major_formatter(FuncFormatter(format_lon))
+ if yticks is not None and lat_labels:
+ ax.set_yticks(yticks, crs=ccrs.PlateCarree())
+ ax.yaxis.set_major_formatter(FuncFormatter(format_lat))
+ ax.tick_params(
+ which="major", length=2.5, width=0.5, labelsize=mpl.rcParams["xtick.labelsize"] + 0.5
+ )
+ if gridlines:
+ grid = ax.gridlines(
+ crs=ccrs.PlateCarree(),
+ draw_labels=(labels and not rectangular),
+ xlocs=xticks,
+ ylocs=yticks,
+ linewidth=grid_width,
+ color=grid_color,
+ linestyle=grid_style,
+ alpha=0.9,
+ zorder=0.6,
+ )
+ if labels and not rectangular:
+ grid.top_labels = grid.right_labels = False
+ grid.left_labels = lat_labels
+ grid.bottom_labels = lon_labels
+ grid.rotate_labels = False
+ grid.x_inline = grid.y_inline = False
+ grid.xpadding = grid.ypadding = 3
+ style = {"size": mpl.rcParams["xtick.labelsize"], "color": INK}
+ grid.xlabel_style = grid.ylabel_style = style
+ artists.append(grid)
+ spine = ax.spines.get("geo")
+ if spine is not None:
+ spine.set_visible(frame)
+ spine.set_linewidth(0.5)
+ spine.set_edgecolor(INK if rectangular else "#7a7a7a")
+ return artists
+
+
+def add_land(
+ ax,
+ *,
+ color=LAND,
+ edgecolor="#8c8c8c",
+ linewidth=0.25,
+ borders=True,
+ coast_color=INK,
+ coast_width=0.3,
+ resolution="110m",
+ zorder=0.5,
+):
+ """Cream land, thin grey country borders and a dark coastline (Natural Earth)."""
+ import cartopy.feature as cfeature
+
+ artists = [
+ ax.add_feature(
+ cfeature.LAND.with_scale(resolution), facecolor=color, edgecolor="none", zorder=zorder
+ )
+ ]
+ if borders:
+ artists.append(
+ ax.add_feature(
+ cfeature.BORDERS.with_scale(resolution),
+ facecolor="none",
+ edgecolor=edgecolor,
+ linewidth=linewidth,
+ zorder=zorder + 0.1,
+ )
+ )
+ if coast_width:
+ artists.append(
+ ax.add_feature(
+ cfeature.COASTLINE.with_scale(resolution),
+ facecolor="none",
+ edgecolor=coast_color,
+ linewidth=coast_width,
+ zorder=zorder + 0.2,
+ )
+ )
+ return artists
+
+
+# --------------------------------------------------------------------------- significance
+
+
+def _grid(lon, lat, mask):
+ lon, lat, mask = np.asarray(lon), np.asarray(lat), np.asarray(mask)
+ if lon.ndim == 1 and lat.ndim == 1:
+ lon, lat = np.meshgrid(lon, lat)
+ if mask.shape != lon.shape:
+ raise ValueError(f"mask shape {mask.shape} must match the lat/lon grid {lon.shape}.")
+ return lon, lat, mask.astype(bool) & np.isfinite(lon) & np.isfinite(lat)
+
+
+def add_stippling(
+ ax,
+ lon,
+ lat,
+ mask,
+ *,
+ stride=2,
+ size=1.2,
+ color=INK,
+ marker="o",
+ alpha=1.0,
+ offset=True,
+ transform=None,
+ zorder=4,
+):
+ """Black dots where ``mask`` is True (e.g. low ensemble agreement or p < 0.05).
+
+ ``stride`` thins the dot lattice for readability; ``offset`` staggers every
+ other row (a hexagonal pattern looks much calmer than a square lattice).
+ Thinning only changes the display, never the underlying decision.
+ """
+ lon2, lat2, m = _grid(lon, lat, mask)
+ keep = np.zeros_like(m)
+ for k, row in enumerate(range(0, m.shape[0], stride)):
+ start = stride // 2 if (offset and stride > 1 and k % 2) else 0
+ keep[row, start::stride] = True
+ show = m & keep
+ if not show.any():
+ return None
+ return ax.scatter(
+ lon2[show],
+ lat2[show],
+ s=size,
+ c=color,
+ marker=marker,
+ linewidths=0,
+ alpha=alpha,
+ transform=transform or ccrs.PlateCarree(),
+ zorder=zorder,
+ rasterized=False,
+ )
+
+
+def add_hatching(
+ ax, lon, lat, mask, *, hatch="....", color=INK, linewidth=0.4, transform=None, zorder=4
+):
+ """Hatch (``'....'``, ``'///'``, ``'xxx'``) the region where ``mask`` is True."""
+ lon2, lat2, m = _grid(lon, lat, mask)
+ if not m.any():
+ return None
+ with mpl.rc_context({"hatch.color": color, "hatch.linewidth": linewidth}):
+ cs = ax.contourf(
+ lon2,
+ lat2,
+ m.astype(float),
+ levels=[0.5, 1.5],
+ colors="none",
+ hatches=[hatch],
+ transform=transform or ccrs.PlateCarree(),
+ zorder=zorder,
+ )
+ for coll in [cs] if hasattr(cs, "set_edgecolor") else cs.collections:
+ coll.set_edgecolor(color)
+ coll.set_linewidth(0)
+ if hasattr(coll, "set_hatch_linewidth"):
+ coll.set_hatch_linewidth(linewidth)
+ return cs
+
+
+def hatch_patch(hatch="....", *, color=INK, facecolor="white", edgecolor=INK, linewidth=0.5):
+ """A legend handle for a hatched region (use with ``ax.legend``)."""
+ patch = Patch(facecolor=facecolor, edgecolor=edgecolor, hatch=hatch, linewidth=linewidth)
+ patch._hatch_color = mpl.colors.to_rgba(color)
+ return patch
+
+
+# --------------------------------------------------------------------------- insets
+
+
+def _halo(ax, width=2.2):
+ """White halo behind inset tick labels and axis labels, so they stay legible on maps."""
+ from matplotlib import patheffects
+
+ effect = [patheffects.withStroke(linewidth=width, foreground="white")]
+ original = ax.draw
+
+ def draw(renderer):
+ texts = [ax.xaxis.label, ax.yaxis.label]
+ texts += ax.get_xticklabels() + ax.get_yticklabels()
+ for text in texts:
+ text.set_path_effects(effect)
+ return original(renderer)
+
+ ax.draw = draw
+ return ax
+
+
+def _clean_inset(ax, frame):
+ if frame == "box":
+ for spine in ax.spines.values():
+ spine.set_visible(True)
+ spine.set_linewidth(0.5)
+ spine.set_color(INK)
+ ax.patch.set_facecolor("white")
+ ax.patch.set_alpha(0.92)
+ else: # "open": transparent background, left + bottom spines only
+ for name, spine in ax.spines.items():
+ spine.set_visible(name in ("left", "bottom"))
+ spine.set_linewidth(0.5)
+ spine.set_color(INK)
+ ax.patch.set_alpha(0)
+ ax.tick_params(length=2.0, width=0.5, pad=1.5)
+ _halo(ax)
+ return ax
+
+
+def _inset(ax, bounds, zorder):
+ child = ax.inset_axes(bounds, transform=ax.transAxes, zorder=zorder)
+ # GeoAxes children inherit nothing from the projection; they are plain axes.
+ return child
+
+
+def add_inset_bars(
+ ax,
+ labels,
+ values,
+ *,
+ colors=None,
+ hatched=None,
+ hatch="....",
+ hatch_label="Low agreement",
+ bounds=(0.07, 0.10, 0.2, 0.36),
+ ylabel=None,
+ frame="box",
+ edgecolor=INK,
+ width=0.72,
+ ylim=None,
+ rotation=0,
+ zorder=8,
+ backdrop=True,
+):
+ """Category bar chart inset (e.g. land area per regime) for the lower-left corner.
+
+ ``hatched`` gives, per bar, the part of ``values`` that is uncertain (e.g. low
+ ensemble agreement). It is stacked *on top* of the confident part with dot
+ hatching and separated by a black line, as in Nature/Science regime maps.
+ """
+ labels = list(labels)
+ values = np.asarray(values, dtype=float)
+ child = _inset(ax, list(bounds), zorder)
+ _clean_inset(child, frame)
+ x = np.arange(len(values))
+ colors = colors if colors is not None else ["#7f7f7f"] * len(values)
+ artists = []
+ if hatched is None:
+ artists.append(
+ child.bar(x, values, width, color=colors, edgecolor=edgecolor, linewidth=0.5)
+ )
+ else:
+ hatched = np.asarray(hatched, dtype=float)
+ solid = values - hatched
+ artists.append(child.bar(x, solid, width, color=colors, edgecolor="none"))
+ with mpl.rc_context({"hatch.color": INK, "hatch.linewidth": 0.6}):
+ artists.append(
+ child.bar(
+ x,
+ hatched,
+ width,
+ bottom=solid,
+ color=colors,
+ edgecolor=INK,
+ linewidth=0,
+ hatch=hatch,
+ )
+ )
+ for xi, s, v in zip(x, solid, values):
+ child.plot([xi - width / 2, xi + width / 2], [s, s], color=INK, lw=0.6, zorder=3)
+ child.add_patch(
+ mpl.patches.Rectangle(
+ (xi - width / 2, 0), width, v, fill=False, edgecolor=edgecolor, lw=0.5, zorder=3
+ )
+ )
+ if hatch_label:
+ child.legend(
+ [hatch_patch(hatch, facecolor="white")],
+ [hatch_label],
+ loc="upper right",
+ bbox_to_anchor=(1.02, 1.04),
+ handlelength=1.6,
+ handleheight=0.9,
+ borderpad=0.2,
+ fontsize=mpl.rcParams["legend.fontsize"],
+ )
+ child.set_xticks(x)
+ child.set_xticklabels(labels, rotation=rotation, ha="right" if rotation else "center")
+ child.tick_params(axis="x", length=0)
+ child.set_xlim(-0.6, len(values) - 0.4)
+ top = ylim[1] if ylim else float(np.nanmax(values)) * (1.38 if hatch_label else 1.12)
+ child.set_ylim(ylim[0] if ylim else 0, top)
+ child.yaxis.set_major_locator(MaxNLocator(4, integer=True))
+ if ylabel:
+ child.set_ylabel(ylabel)
+ if backdrop and frame == "box":
+ child.patch.set_alpha(0.92)
+ return child
+
+
+def add_inset_histogram(
+ ax,
+ groups,
+ *,
+ bins=20,
+ colors=None,
+ labels=None,
+ total=True,
+ total_label="All",
+ stacked=True,
+ log=False,
+ cumulative=None,
+ cumulative_color="#e8161b",
+ cumulative_label=None,
+ cumulative_reverse=True,
+ bounds=(0.09, 0.08, 0.22, 0.3),
+ xlabel=None,
+ ylabel="Count",
+ frame="open",
+ zorder=8,
+ alpha=0.55,
+ mean_line=False,
+ legend=False,
+):
+ """Histogram inset with optional stacked groups, a black total outline and a
+ cumulative curve on a second, outward-offset y axis (the Nature "site map" look).
+
+ ``groups`` is an array or a dict ``{label: values}``. ``cumulative`` may be
+ ``True`` (cumulative count) or an array of weights the same length as the
+ concatenated values (e.g. emissions), plotted as a red line on the outer axis.
+ ``log=True`` uses logarithmic bins and axis.
+ Returns ``(hist_axes, cumulative_axes_or_None)``.
+ """
+ if not isinstance(groups, dict):
+ groups = {total_label: np.asarray(groups, dtype=float)}
+ total = False if len(groups) == 1 and not stacked else total
+ names = list(groups)
+ arrays = [np.asarray(groups[n], dtype=float).ravel() for n in names]
+ arrays = [a[np.isfinite(a)] for a in arrays]
+ allv = np.concatenate(arrays)
+ if allv.size == 0:
+ raise ValueError("The histogram has no finite values.")
+ if np.isscalar(bins):
+ if log:
+ pos = allv[allv > 0]
+ edges = np.logspace(np.log10(pos.min()), np.log10(pos.max()), int(bins) + 1)
+ else:
+ edges = np.histogram_bin_edges(allv, bins=int(bins))
+ else:
+ edges = np.asarray(bins, dtype=float)
+ colors = colors or ["#f7a541", "#b56ccf", "#4f9fd6", "#61a85c"][: len(names)]
+ child = _inset(ax, list(bounds), zorder)
+ _clean_inset(child, frame)
+ single = len(names) == 1
+ child.hist(
+ arrays if not single else arrays[0],
+ bins=edges,
+ stacked=stacked and not single,
+ color=colors if not single else colors[0],
+ alpha=alpha,
+ edgecolor="none",
+ label=labels or names,
+ )
+ if total or single:
+ counts, _ = np.histogram(allv, bins=edges)
+ child.stairs(counts, edges, color=INK, lw=0.8, label=total_label if total else None)
+ if log:
+ child.set_xscale("log")
+ if mean_line:
+ child.axvline(np.mean(allv), color=INK, lw=0.6, ls=(0, (2, 1.5)))
+ child.set_xlim(edges[0], edges[-1])
+ child.set_ylim(bottom=0)
+ child.yaxis.set_major_locator(MaxNLocator(4, integer=True))
+ child.xaxis.set_major_locator(MaxNLocator(4) if not log else mpl.ticker.LogLocator())
+ if xlabel:
+ child.set_xlabel(xlabel)
+ if ylabel:
+ child.set_ylabel(ylabel)
+ twin = None
+ if cumulative is not None and cumulative is not False:
+ w = np.ones_like(allv) if cumulative is True else np.asarray(cumulative, dtype=float)
+ order = np.argsort(allv)
+ xs, ws = allv[order], w[order]
+ if cumulative_reverse: # amount from sources >= x (Nature style)
+ ys = np.cumsum(ws[::-1])[::-1]
+ else:
+ ys = np.cumsum(ws)
+ twin = child.twinx()
+ _halo(twin)
+ twin.patch.set_alpha(0)
+ twin.plot(xs, ys, color=cumulative_color, lw=0.9, zorder=5)
+ for name, spine in twin.spines.items():
+ spine.set_visible(name == "left")
+ twin.spines["left"].set_position(("outward", 26))
+ twin.spines["left"].set_linewidth(0.5)
+ twin.yaxis.set_ticks_position("left")
+ twin.yaxis.set_label_position("left")
+ twin.tick_params(length=2.0, width=0.5, pad=1.5)
+ twin.set_ylim(bottom=0)
+ twin.yaxis.set_major_locator(MaxNLocator(4))
+ if cumulative_label:
+ twin.set_ylabel(cumulative_label)
+ twin.set_zorder(child.get_zorder() + 0.1)
+ if legend:
+ child.legend(loc="upper right", fontsize=mpl.rcParams["legend.fontsize"] - 0.5)
+ return child, twin
+
+
+# --------------------------------------------------------------------------- marginal profiles
+
+
+def _proj_coord(ax, lat=None, lon=None, ref=None):
+ """Map latitudes (or longitudes) onto the projected y (or x) axis of ``ax``."""
+ proj = ax.projection
+ geo = ccrs.PlateCarree()
+ if lat is not None:
+ lat = np.asarray(lat, dtype=float)
+ ref = ref if ref is not None else getattr(proj, "proj4_params", {}).get("lon_0", 0.0)
+ pts = proj.transform_points(geo, np.full_like(lat, ref), lat)
+ return pts[:, 1]
+ lon = np.asarray(lon, dtype=float)
+ ref = 0.0 if ref is None else ref
+ pts = proj.transform_points(geo, lon, np.full_like(lon, ref))
+ return pts[:, 0]
+
+
+def add_latitude_profile(
+ ax,
+ lat,
+ center,
+ *,
+ lower=None,
+ upper=None,
+ members=None,
+ outer=None,
+ width=0.13,
+ pad=0.03,
+ color="#2a6f97",
+ band_alpha=0.22,
+ outer_alpha=0.10,
+ linewidth=1.0,
+ xlabel=None,
+ xlim=None,
+ reference=0.0,
+ yticks=None,
+ ticklabels="right",
+ gridlines=True,
+ title=None,
+ ref_lon=None,
+ zorder=6,
+ extra=None,
+):
+ """Right-hand marginal panel: zonal-mean line with a shaded spread, *aligned* with
+ the map's latitudes for any projection.
+
+ The y axis of the panel lives in the map's projected coordinates (latitudes are
+ pushed through ``ax.projection``), and its limits are linked to the map
+ ``ylim`` at draw time, so 40°N on the profile sits exactly level with 40°N on
+ the map — on Robinson and Equal Earth too.
+
+ ``lower``/``upper`` shade an inner band (e.g. IQR or ±1 SD); ``outer`` =
+ ``(low, high)`` adds a lighter outer band (e.g. ensemble min–max).
+ ``members`` (n_members × n_lat) draws thin individual lines instead.
+ ``extra`` is a list of ``(center, dict(**line_kwargs))`` for further lines.
+ """
+ lat = np.asarray(lat, dtype=float)
+ center = np.asarray(center, dtype=float)
+ child = ax.inset_axes([1 + pad, 0, width, 1], transform=ax.transAxes, zorder=zorder)
+ y = _proj_coord(ax, lat=lat, ref=ref_lon)
+ if outer is not None:
+ child.fill_betweenx(
+ y, outer[0], outer[1], color=color, alpha=outer_alpha, lw=0, zorder=1
+ )
+ if members is not None:
+ for m in np.asarray(members, dtype=float):
+ child.plot(m, y, color=color, lw=0.3, alpha=0.35, zorder=2)
+ if lower is not None and upper is not None:
+ child.fill_betweenx(y, lower, upper, color=color, alpha=band_alpha, lw=0, zorder=2)
+ for line, kw in extra or []:
+ child.plot(np.asarray(line, dtype=float), y, **{"lw": 0.8, "zorder": 3, **kw})
+ child.plot(center, y, color=color, lw=linewidth, zorder=4, solid_capstyle="round")
+ if reference is not None:
+ child.axvline(reference, color=INK, lw=0.5, ls=(0, (2.5, 2)), zorder=1.5)
+ for name, spine in child.spines.items():
+ spine.set_visible(name in ("bottom", ticklabels if ticklabels in ("left", "right") else ""))
+ spine.set_linewidth(0.5)
+ spine.set_color(INK)
+ child.patch.set_alpha(0)
+ if yticks is None:
+ span = np.nanmax(lat) - np.nanmin(lat)
+ step = 30 if span > 120 else (20 if span > 60 else (10 if span > 25 else 5))
+ yticks = np.arange(np.ceil(np.nanmin(lat) / step) * step, np.nanmax(lat) + 1e-9, step)
+ yticks = np.asarray(yticks, dtype=float)
+ ypos = _proj_coord(ax, lat=yticks, ref=ref_lon)
+ child.yaxis.set_major_locator(FixedLocator(ypos))
+
+ def _lat_label(value, _pos=None):
+ i = int(np.argmin(np.abs(ypos - value)))
+ return format_lat(yticks[i])
+
+ child.yaxis.set_major_formatter(FuncFormatter(_lat_label))
+ if ticklabels == "right":
+ child.yaxis.tick_right()
+ elif ticklabels in (None, False, "none"):
+ child.tick_params(axis="y", left=False, labelleft=False)
+ child.yaxis.set_minor_locator(NullLocator())
+ if gridlines:
+ for t in child.get_yticks():
+ child.axhline(t, color=GRID, lw=0.35, ls=(0, (3, 3)), zorder=0)
+ child.xaxis.set_major_locator(MaxNLocator(3, symmetric=reference == 0))
+ child.tick_params(length=2.0, width=0.5, pad=1.5)
+ if xlim is not None:
+ child.set_xlim(xlim)
+ if xlabel:
+ child.set_xlabel(xlabel)
+ if title:
+ child.set_title(title, fontsize=mpl.rcParams["axes.labelsize"], fontweight="normal", pad=3)
+
+ _orig_draw = child.draw
+
+ def draw(renderer):
+ # Keep the profile's latitude axis locked to the map, even after set_extent.
+ child.set_ylim(ax.get_ylim())
+ return _orig_draw(renderer)
+
+ child.set_ylim(ax.get_ylim())
+ child.draw = draw
+ return child
+
+
+def add_longitude_profile(
+ ax,
+ lon,
+ center,
+ *,
+ lower=None,
+ upper=None,
+ height=0.16,
+ pad=0.02,
+ color="#2a6f97",
+ band_alpha=0.22,
+ linewidth=1.0,
+ ylabel=None,
+ reference=None,
+ ref_lat=0.0,
+ zorder=6,
+):
+ """Top marginal panel (meridional mean vs longitude), aligned with the map x axis."""
+ lon = np.asarray(lon, dtype=float)
+ center = np.asarray(center, dtype=float)
+ child = ax.inset_axes([0, 1 + pad, 1, height], transform=ax.transAxes, zorder=zorder)
+ x = _proj_coord(ax, lon=lon, ref=ref_lat)
+ if lower is not None and upper is not None:
+ child.fill_between(x, lower, upper, color=color, alpha=band_alpha, lw=0)
+ child.plot(x, center, color=color, lw=linewidth)
+ if reference is not None:
+ child.axhline(reference, color=INK, lw=0.5, ls=(0, (2.5, 2)))
+ for name, spine in child.spines.items():
+ spine.set_visible(name in ("left", "bottom"))
+ spine.set_linewidth(0.5)
+ child.patch.set_alpha(0)
+ child.tick_params(axis="x", bottom=False, labelbottom=False)
+ child.yaxis.set_major_locator(MaxNLocator(3))
+ child.tick_params(length=2.0, width=0.5, pad=1.5)
+ if ylabel:
+ child.set_ylabel(ylabel)
+ _orig_draw = child.draw
+
+ def draw(renderer):
+ child.set_xlim(ax.get_xlim())
+ return _orig_draw(renderer)
+
+ child.set_xlim(ax.get_xlim())
+ child.draw = draw
+ return child
+
+
+# --------------------------------------------------------------------------- labels & export
+
+
+def add_panel_label(ax, label, *, style="({})", x=0.008, y=0.985, size=None, inside=True, **kw):
+ """Bold panel label such as ``(a)`` (``style='{}'`` gives a bare ``a``)."""
+ text = style.format(label)
+ opts = dict(
+ transform=ax.transAxes,
+ fontsize=size or mpl.rcParams["axes.titlesize"] + 1,
+ fontweight="bold",
+ va="top" if inside else "bottom",
+ ha="left",
+ color=INK,
+ zorder=20,
+ )
+ opts.update(kw)
+ return ax.text(x, y if inside else 1.01, text, **opts)
+
+
+def add_size_legend(
+ ax,
+ values,
+ labels,
+ *,
+ scale,
+ marker="o",
+ facecolor="none",
+ edgecolor=INK,
+ linewidth=0.8,
+ title=None,
+ **legend_kw,
+):
+ """Legend for markers sized by value; ``scale(value) -> marker area (pt²)``."""
+ handles = [
+ Line2D(
+ [],
+ [],
+ ls="none",
+ marker=marker,
+ markersize=np.sqrt(scale(v)),
+ markerfacecolor=facecolor,
+ markeredgecolor=edgecolor,
+ markeredgewidth=linewidth,
+ )
+ for v in values
+ ]
+ opts = dict(frameon=False, title=title, labelspacing=0.9, borderpad=0.2, handletextpad=0.4)
+ opts.update(legend_kw)
+ leg = ax.legend(handles, labels, **opts)
+ if title:
+ leg.get_title().set_fontsize(mpl.rcParams["legend.fontsize"] + 0.5)
+ return leg
+
+
+def save_figure(fig, path, *, dpi=600, formats=None, **kwargs):
+ """Save at print resolution; ``formats=('png', 'pdf')`` writes several files."""
+ path = Path(path)
+ path.parent.mkdir(parents=True, exist_ok=True)
+ opts = {"dpi": dpi, "bbox_inches": "tight", "pad_inches": 0.02, "facecolor": "white"}
+ opts.update(kwargs)
+ paths = []
+ for fmt in formats or [path.suffix.lstrip(".") or "png"]:
+ target = path.with_suffix("." + fmt)
+ fig.savefig(target, **opts)
+ paths.append(target)
+ return paths[0] if len(paths) == 1 else paths
+
+
+@contextlib.contextmanager
+def figure(width="double", ratio=0.55, **subplot_kw):
+ """Convenience: ``with journal.figure('double', .5) as fig:`` in journal style."""
+ with journal_style():
+ fig = plt.figure(figsize=figsize(width, ratio), **subplot_kw)
+ yield fig
diff --git a/plot_function/layers.py b/plot_function/layers.py
index 55c5939..ed2edf0 100644
--- a/plot_function/layers.py
+++ b/plot_function/layers.py
@@ -9,8 +9,9 @@
import xarray as xr
from matplotlib.lines import Line2D
from matplotlib.patches import Patch
-from matplotlib.ticker import MaxNLocator
+from matplotlib.ticker import FixedLocator, FuncFormatter, MaxNLocator
+from . import journal as pj
from .data import _select_variable, open_field
from .legacy import _cyclic_if_global
from .options import Distribution, MapFeatures, Profile, Significance
@@ -45,19 +46,32 @@ def _scope(result, weights=None):
def _style_axes(ax, result):
- ink = "#163438" if result.theme == "light" else "#edf5ef"
- muted = "#647b80" if result.theme == "light" else "#a3bec2"
+ """Journal inset styling: transparent background, thin left/bottom spines."""
+ light = result.theme == "light"
+ ink = pj.INK if light else "#edf5ef"
+ muted = pj.MUTED if light else "#a3bec2"
for name, spine in ax.spines.items():
spine.set_visible(name in ("left", "bottom"))
- spine.set_color(muted)
- spine.set_linewidth(0.55)
- ax.set_facecolor("white" if result.theme == "light" else "#10252e")
- ax.tick_params(labelsize=7, colors=muted, length=2, width=0.5)
+ spine.set_color(ink)
+ spine.set_linewidth(0.5)
+ ax.patch.set_alpha(0)
+ ax.tick_params(labelsize=6, colors=ink, length=2, width=0.5, pad=1.5)
ax.xaxis.label.set_color(ink)
ax.yaxis.label.set_color(ink)
return ink, muted
+def _attach_draw_hook(child, sync):
+ original = child.draw
+
+ def draw(renderer):
+ sync()
+ return original(renderer)
+
+ child.draw = draw
+ sync()
+
+
def add_profile(result, options=None, **kwargs):
spec = options if options is not None else Profile(**kwargs)
if not isinstance(spec, Profile):
@@ -79,51 +93,104 @@ def add_profile(result, options=None, **kwargs):
)
if not bool(np.isfinite(stats.center).any()):
raise ValueError("The profile has no finite positive-weight samples.")
+ map_ax = result.axes
bounds = [1 + spec.pad, 0, spec.width, 1] if right else [0, 1 + spec.pad, 1, spec.width]
- ax = result.axes.inset_axes(bounds, transform=result.axes.transAxes, zorder=6)
+ ax = map_ax.inset_axes(bounds, transform=map_ax.transAxes, zorder=6)
ink, muted = _style_axes(ax, result)
coord, center = stats[coordinate].values, stats.center.values
+ # Push geographic coordinates through the map projection so that, e.g., 40°N on
+ # the profile is level with 40°N on the map, for Robinson/Equal Earth as well.
+ if result.extent is not None:
+ west, east, south, north = result.extent
+ ref = (west + east) / 2 if right else (south + north) / 2
+ else:
+ ref = None if right else 0.0
+ try:
+ pos = (
+ pj._proj_coord(map_ax, lat=coord, ref=ref)
+ if right
+ else pj._proj_coord(map_ax, lon=coord, ref=ref)
+ )
+ aligned = bool(np.isfinite(pos).all())
+ except Exception: # pragma: no cover - exotic projections
+ aligned = False
+ if not aligned:
+ pos = coord
artists = []
if spec.style == "band" and spread is not None:
fill = ax.fill_betweenx if right else ax.fill_between
artists.append(
- fill(coord, stats.lower, stats.upper, color=spec.color, alpha=0.16, linewidth=0)
+ fill(pos, stats.lower, stats.upper, color=spec.color, alpha=0.22, linewidth=0)
)
if spec.style == "bars":
- spacing = float(np.min(np.diff(coord))) * 0.75 if len(coord) > 1 else 0.6
+ spacing = float(np.min(np.abs(np.diff(pos)))) * 0.8 if len(pos) > 1 else 0.6
if right:
- artists.extend(ax.barh(coord, center, height=spacing, color=spec.color, alpha=0.75))
+ artists.extend(ax.barh(pos, center, height=spacing, color=spec.color, alpha=0.8))
else:
- artists.extend(ax.bar(coord, center, width=spacing, color=spec.color, alpha=0.75))
+ artists.extend(ax.bar(pos, center, width=spacing, color=spec.color, alpha=0.8))
else:
artists.extend(
- ax.plot(center, coord, color=spec.color, lw=1.35)
+ ax.plot(center, pos, color=spec.color, lw=1.0)
if right
- else ax.plot(coord, center, color=spec.color, lw=1.35)
+ else ax.plot(pos, center, color=spec.color, lw=1.0)
)
if spec.reference is not None:
reference = ax.axvline if right else ax.axhline
- artists.append(reference(spec.reference, color=muted, lw=0.65, ls="--"))
- if result.extent is not None:
- limits = result.extent[2:] if right else result.extent[:2]
- else:
- limits = [float(coord.min()), float(coord.max())]
- if limits[0] == limits[1]:
- limits = [limits[0] - 0.5, limits[1] + 0.5]
+ artists.append(reference(spec.reference, color=ink, lw=0.5, ls=(0, (2.5, 2))))
unit = data.attrs.get("units", "")
+ if aligned:
+ if result.extent is not None:
+ lo, hi = result.extent[2:] if right else result.extent[:2]
+ else:
+ lo, hi = (-90, 90) if right else (-180, 180)
+ span = hi - lo
+ step = next((v for v in (5, 10, 15, 20, 30, 60) if span / v <= 6), 60)
+ ticks = np.arange(np.ceil(lo / step) * step, hi + 1e-9, step)
+ if right:
+ ticks = ticks[np.abs(ticks) < 90]
+ tick_pos = (
+ pj._proj_coord(map_ax, lat=ticks, ref=ref)
+ if right
+ else pj._proj_coord(map_ax, lon=ticks, ref=ref)
+ )
+ fmt = pj.format_lat if right else pj.format_lon
+ axis = ax.yaxis if right else ax.xaxis
+ axis.set_major_locator(FixedLocator(tick_pos))
+ axis.set_major_formatter(
+ FuncFormatter(lambda v, _p: fmt(ticks[int(np.argmin(np.abs(tick_pos - v)))]))
+ )
+ for t in tick_pos:
+ (ax.axhline if right else ax.axvline)(t, color=pj.GRID, lw=0.35, ls=(0, (3, 3)), zorder=0)
+
+ def sync():
+ if right:
+ ax.set_ylim(map_ax.get_ylim())
+ else:
+ ax.set_xlim(map_ax.get_xlim())
+
+ _attach_draw_hook(ax, sync)
+ else:
+ limits = (
+ (result.extent[2:] if right else result.extent[:2])
+ if result.extent is not None
+ else [float(coord.min()), float(coord.max())]
+ )
+ if limits[0] == limits[1]:
+ limits = [limits[0] - 0.5, limits[1] + 0.5]
+ (ax.set_ylim if right else ax.set_xlim)(limits)
+ band_label = {"std": " ± 1 s.d.", "iqr": " (IQR)", None: ""}[spread]
+ title = spec.label if spec.label is not None else spec.statistic.capitalize() + band_label
if right:
- ax.set_ylim(limits)
ax.yaxis.tick_right()
ax.spines["right"].set_visible(True)
ax.spines["left"].set_visible(False)
- ax.set_xlabel(unit, fontsize=8)
- ax.yaxis.set_label_position("right")
- ax.set_ylabel("Latitude / °", fontsize=8, labelpad=5)
+ # The description goes under the axis, where it cannot collide with titles.
+ ax.set_xlabel(f"{title}\n({unit})" if unit else title, fontsize=6.5)
ax.xaxis.set_major_locator(MaxNLocator(3))
+ ax.xaxis.set_major_formatter(pj.clean_formatter())
else:
- ax.set_xlim(limits)
- ax.set_ylabel(unit, fontsize=8)
- ax.tick_params(axis="x", labelbottom=False)
+ ax.set_ylabel(unit, fontsize=7)
+ ax.tick_params(axis="x", labelbottom=False, bottom=False)
ax.yaxis.set_major_locator(MaxNLocator(3))
if result._title is not None:
result._title = result.axes.set_title(
@@ -131,14 +198,12 @@ def add_profile(result, options=None, **kwargs):
loc="left",
fontproperties=result._title.get_fontproperties(),
color=result._title.get_color(),
- y=1 + spec.pad + spec.width + 0.07,
- pad=28 if result._subtitle is not None else 14,
+ y=1 + spec.pad + spec.width + 0.08,
+ pad=14 if result._subtitle is not None else 5,
)
if result._subtitle is not None:
- result._subtitle.set_y(1 + spec.pad + spec.width + 0.07)
- band_label = {"std": " ±1 spatial SD", "iqr": " · spatial IQR", None: ""}[spread]
- title = spec.label if spec.label is not None else spec.statistic.capitalize() + band_label
- ax.set_title(title, loc="left", fontsize=8, color=ink, pad=7)
+ result._subtitle.set_y(1 + spec.pad + spec.width + 0.08)
+ ax.set_title(title, loc="left", fontsize=6.5, color=ink, pad=2, fontweight="normal")
layer = LayerResult(ax, stats, artists)
result.profiles[spec.position] = layer
return layer
@@ -173,10 +238,11 @@ def add_distribution(result, options=None, **kwargs):
stats = histogram(data, bins=spec.bins, weights=weights, density=spec.density)
ax = result.axes.inset_axes(spec.bounds, transform=result.axes.transAxes, zorder=8)
ink, muted = _style_axes(ax, result)
- # The solid white panel separates the data summary from the map beneath it.
- ax.patch.set_alpha(0.96)
+ if spec.background is not None:
+ ax.patch.set_facecolor(spec.background)
+ ax.patch.set_alpha(spec.background_alpha)
artists = []
- color = spec.color if spec.color != "map" else "#267f83"
+ color = spec.color if spec.color != "map" else "#3b6f8f"
if spec.style == "bars":
colors = (
result.artist.cmap(result.artist.norm(stats.bin.values))
@@ -187,52 +253,50 @@ def add_distribution(result, options=None, **kwargs):
ax.bar(
stats.bin,
stats.height,
- width=(stats.right - stats.left) * 0.88,
+ width=(stats.right - stats.left),
color=colors,
- edgecolor="none",
- alpha=0.9,
+ edgecolor="white" if result.theme == "light" else "none",
+ linewidth=0.3,
)
)
+ edges = np.r_[stats.left.values, stats.right.values[-1]]
+ artists.append(ax.stairs(stats.height, edges, color=ink, lw=0.6, fill=False))
elif spec.style == "step":
edges = np.r_[stats.left.values, stats.right.values[-1]]
- artists.append(ax.stairs(stats.height, edges, color=color, lw=1.2, fill=False))
+ artists.append(ax.stairs(stats.height, edges, color=color, lw=0.9, fill=False))
elif spec.style == "line":
- artists.extend(ax.plot(stats.bin, stats.height, color=color, lw=1.35))
- artists.append(ax.fill_between(stats.bin, stats.height, color=color, alpha=0.10))
+ artists.extend(ax.plot(stats.bin, stats.height, color=color, lw=1.0))
+ artists.append(ax.fill_between(stats.bin, stats.height, color=color, alpha=0.15, lw=0))
else:
artists.extend(
- ax.step(stats.value, stats["cumulative"], where="post", color=color, lw=1.35)
+ ax.step(stats.value, stats["cumulative"], where="post", color=color, lw=1.0)
)
if spec.show_mean:
- artists.append(ax.axvline(stats.attrs["mean"], color=ink, lw=0.75, ls=(0, (3, 2))))
+ mean = stats.attrs["mean"]
+ artists.append(
+ ax.plot(
+ [mean], [0], marker="^", ms=4, color=ink, clip_on=False, zorder=5,
+ transform=ax.get_xaxis_transform(),
+ )[0]
+ )
ax.text(
- 0.98,
- 0.90,
- f"Mean {stats.attrs['mean']:.2g}",
- ha="right",
- va="top",
- transform=ax.transAxes,
- fontsize=6.5,
- color=ink,
+ mean, 1.0, f"mean {mean:.3g}", transform=ax.get_xaxis_transform(), ha="center",
+ va="bottom", fontsize=5.5, color=ink,
)
weighted = spec.weights is not None
label = (
- "Cumulative fraction"
+ "CDF"
if spec.style == "ecdf"
- else ("Density" if spec.density else ("Weight" if weighted else "Cell count"))
- )
- ax.set_title(
- spec.label or ("Weighted distribution" if weighted else "Cell distribution"),
- loc="left",
- fontsize=8,
- color=ink,
- pad=5,
+ else ("Density" if spec.density else ("Weight" if weighted else "Cells"))
)
- ax.set_xlabel(data.attrs.get("units", ""), fontsize=7, labelpad=2)
- ax.set_ylabel(label, fontsize=7, labelpad=3)
+ if spec.label:
+ ax.set_title(spec.label, loc="left", fontsize=6.5, color=ink, pad=8, fontweight="normal")
+ ax.set_xlabel(data.attrs.get("units", ""), fontsize=6.5, labelpad=1)
+ ax.set_ylabel(label, fontsize=6.5, labelpad=2)
ax.xaxis.set_major_locator(MaxNLocator(3))
ax.yaxis.set_major_locator(MaxNLocator(3))
- ax.tick_params(labelsize=6.5)
+ ax.xaxis.set_major_formatter(pj.clean_formatter())
+ ax.yaxis.set_major_formatter(pj.clean_formatter())
ax.set_ylim(bottom=0)
if spec.style == "ecdf":
ax.set_ylim(0, 1.04)
@@ -351,12 +415,12 @@ def load_overlay(source):
legend = ax.legend(
[handle],
[spec.label or label],
- loc="upper right",
- fontsize=7,
+ loc="lower right",
+ fontsize=6.5,
frameon=True,
facecolor="white",
- edgecolor="0.85",
- framealpha=0.96,
+ edgecolor="none",
+ framealpha=0.85,
)
legend.set_zorder(10)
artists.append(legend)
diff --git a/plot_function/maps.py b/plot_function/maps.py
index ad9da56..fc3a7d6 100644
--- a/plot_function/maps.py
+++ b/plot_function/maps.py
@@ -9,12 +9,14 @@
import numpy as np
import xarray as xr
+from . import journal as pj
from .data import Source, open_field
from .legacy import one_map_flat
from .options import MapFeatures, Profile
+# "light" is the publication theme: white paper, near-black ink, light dashed grid.
_THEMES = {
- "light": {"paper": "#f8faf9", "ink": "#163438", "muted": "#647b80", "water": "#edf3f2"},
+ "light": {"paper": "#ffffff", "ink": pj.INK, "muted": pj.MUTED, "water": "#ffffff"},
"dark": {"paper": "#10252e", "ink": "#edf5ef", "muted": "#a3bec2", "water": "#16333e"},
}
_PROJECTIONS = {
@@ -67,13 +69,19 @@ def add_features(self, options=None, **kwargs):
return add_features(self, options, **kwargs)
- def save(self, path: str | Path, *, dpi: int = 180, **kwargs) -> Path:
+ def save(self, path: str | Path, *, dpi: int = 300, **kwargs) -> Path:
"""Save PNG, PDF, SVG, or another Matplotlib format; create parent directories."""
path = Path(path)
path.parent.mkdir(parents=True, exist_ok=True)
- options = {"dpi": dpi, "bbox_inches": "tight", "facecolor": self.figure.get_facecolor()}
+ options = {
+ "dpi": dpi,
+ "bbox_inches": "tight",
+ "pad_inches": 0.03,
+ "facecolor": self.figure.get_facecolor(),
+ }
options.update(kwargs)
- self.figure.savefig(path, **options)
+ with pj.journal_style():
+ self.figure.savefig(path, **options)
return path
@@ -110,9 +118,9 @@ def plot_map(
features=None,
panel_label=None,
ax=None,
- figsize=(10, 5.8),
+ figsize=(7.2, 4.2),
output=None,
- dpi=180,
+ dpi=300,
) -> MapResult:
"""Plot a NetCDF field with one function call; return editable Matplotlib objects.
@@ -121,6 +129,10 @@ def plot_map(
Choose ``coastlines=False`` for a completely offline, NetCDF-only plot.
Coastlines otherwise use Cartopy's cached/downloaded Natural Earth 110m data.
A supplied Cartopy ``ax`` owns the projection; ``projection`` is then ignored.
+
+ Figures use the journal typography of :mod:`plot_function.journal` (Arial /
+ Helvetica-like sans serif, 7–9 pt, thin 0.5 pt lines) and are exported at
+ 300 dpi by default; pass ``dpi=600`` for print.
"""
if theme not in _THEMES:
raise ValueError(f"Unknown theme {theme!r}; choose 'light' or 'dark'.")
@@ -163,15 +175,15 @@ def plot_map(
)
colors = _THEMES[theme]
owns_figure = ax is None
- with mpl.rc_context({"font.size": 10}):
+ with pj.journal_style():
if owns_figure:
figure, ax = plt.subplots(figsize=figsize, subplot_kw={"projection": projection})
figure.patch.set_facecolor(colors["paper"])
figure.subplots_adjust(
- left=0.055,
- right=0.78 if any(p.position == "right" for p in profiles) else 0.945,
- top=0.67 if any(p.position == "top" for p in profiles) else 0.83,
- bottom=0.14,
+ left=0.06,
+ right=0.80 if any(p.position == "right" for p in profiles) else 0.97,
+ top=0.72 if any(p.position == "top" for p in profiles) else 0.88,
+ bottom=0.13,
)
else:
figure = ax.figure
@@ -194,41 +206,29 @@ def plot_map(
ax.coastlines(
resolution=features.resolution if features else "110m",
color=colors["ink"],
- linewidth=0.45,
- alpha=0.75,
+ linewidth=0.35,
)
- ax.spines["geo"].set_edgecolor(colors["muted"])
- ax.spines["geo"].set_linewidth(0.6)
+ ax.spines["geo"].set_edgecolor(colors["ink"] if theme == "light" else colors["muted"])
+ ax.spines["geo"].set_linewidth(0.5)
if gridlines:
- grid = ax.gridlines(
- draw_labels=extent is not None,
- x_inline=False,
- y_inline=False,
- linewidth=0.4,
- color=colors["muted"],
- alpha=0.45,
- linestyle=(0, (3, 5)),
- )
- grid.top_labels = grid.right_labels = False
- grid.rotate_labels = False
- grid.xlabel_style = grid.ylabel_style = {"color": colors["muted"], "size": 8}
+ _journal_grid(ax, extent, colors, theme)
heading = title if title is not None else data.attrs.get("long_name", data.name or "")
title_artist = ax.set_title(
heading,
loc="left",
- fontsize=17 if owns_figure else 12,
+ fontsize=9 if owns_figure else 8,
fontweight="bold",
color=colors["ink"],
- pad=28 if subtitle else 14,
+ pad=14 if subtitle else 5,
)
subtitle_artist = None
if subtitle:
subtitle_artist = ax.text(
0,
- 1.015,
+ 1.012,
subtitle,
transform=ax.transAxes,
- fontsize=9,
+ fontsize=7,
color=colors["muted"],
ha="left",
va="bottom",
@@ -239,18 +239,25 @@ def plot_map(
artist,
ax=ax,
orientation="horizontal",
- pad=0.09,
- fraction=0.045,
- shrink=0.72,
- aspect=36,
+ pad=0.06 if extent is None else 0.09,
+ fraction=0.04,
+ shrink=0.6,
+ aspect=38,
+ extendfrac=0.04,
)
cbar.set_label(
label if label is not None else data.attrs.get("units", ""),
color=colors["ink"],
- labelpad=7,
+ fontsize=7.5,
+ labelpad=2,
+ )
+ cbar.ax.tick_params(
+ colors=colors["ink"], labelsize=6.5, length=2, width=0.5, pad=1.5
)
- cbar.ax.tick_params(colors=colors["muted"], labelsize=8, length=3)
- cbar.outline.set_visible(False)
+ cbar.outline.set_linewidth(0.5)
+ cbar.outline.set_edgecolor(colors["ink"])
+ cbar.ax.minorticks_off()
+ cbar.ax.xaxis.set_major_formatter(pj.clean_formatter())
result = MapResult(
figure,
ax,
@@ -271,18 +278,20 @@ def plot_map(
if distribution is not None:
result.add_distribution(distribution)
if panel_label is not None:
- ax.text(
- 0.015,
- 0.975,
- panel_label,
- transform=ax.transAxes,
- va="top",
- fontsize=12,
- fontweight="bold",
- color=colors["ink"],
- zorder=11,
- bbox={"facecolor": colors["paper"], "edgecolor": "none", "pad": 3},
+ label_style = dict(
+ fontsize=10, fontweight="bold", color=colors["ink"], zorder=11, ha="right"
)
+ if heading:
+ # Nature style: bold lowercase letter immediately left of the title.
+ ax.annotate(
+ panel_label, xy=(0, 0), xycoords=result._title, xytext=(-5, 0),
+ textcoords="offset points", va="bottom", **label_style,
+ )
+ else:
+ ax.text(
+ -0.01, 1.01, panel_label, transform=ax.transAxes, va="bottom",
+ **label_style,
+ )
if output is not None:
result.save(output, dpi=dpi)
return result
@@ -290,3 +299,34 @@ def plot_map(
if owns_figure:
plt.close(figure)
raise
+
+
+def _journal_grid(ax, extent, colors, theme):
+ """Light dashed graticule; outward degree ticks on rectangular regional maps."""
+ rectangular = isinstance(ax.projection, (ccrs.PlateCarree, ccrs.Mercator))
+ if extent is not None:
+ west, east, south, north = extent
+ xstep = _nice_step(east - west)
+ ystep = _nice_step(north - south)
+ xticks = np.arange(np.ceil(west / xstep) * xstep, east + 1e-9, xstep)
+ yticks = np.arange(np.ceil(south / ystep) * ystep, north + 1e-9, ystep)
+ else:
+ xticks, yticks = np.arange(-180, 181, 60), np.arange(-60, 61, 30)
+ pj.geo_ticks(
+ ax,
+ xticks=xticks,
+ yticks=yticks,
+ labels=extent is not None or rectangular,
+ grid_color=pj.GRID if theme == "light" else colors["muted"],
+ )
+ # Axes created outside journal_style() have DejaVu tick labels; restyle them.
+ ax.tick_params(labelfontfamily=pj.resolved_font())
+ if theme != "light":
+ ax.tick_params(colors=colors["muted"], labelcolor=colors["muted"])
+
+
+def _nice_step(span):
+ for step in (1, 2, 5, 10, 20, 30, 60):
+ if span / step <= 5:
+ return step
+ return 60
diff --git a/plot_function/options.py b/plot_function/options.py
index b186e61..cff86a3 100644
--- a/plot_function/options.py
+++ b/plot_function/options.py
@@ -11,9 +11,9 @@ class Profile:
statistic: str = "mean"
spread: str | None = "std"
weights: Any = None
- color: str = "#267f83"
- width: float = 0.20
- pad: float = 0.045
+ color: str = "#2b5d8c"
+ width: float = 0.16
+ pad: float = 0.03
label: str | None = None
reference: float | None = None
@@ -24,10 +24,12 @@ class Distribution:
bins: Any = 18
weights: Any = None
density: bool = True
- bounds: tuple = (0.045, 0.075, 0.27, 0.25)
- color: str = "#267f83"
+ bounds: tuple = (0.06, 0.1, 0.2, 0.22)
+ color: str = "#3b6f8f"
label: str | None = None
show_mean: bool = True
+ background: str | None = "white" # translucent backdrop over busy maps; None = clear
+ background_alpha: float = 0.8
@dataclass(frozen=True)
@@ -39,8 +41,8 @@ class Significance:
alpha: float = 0.05
correction: str = "none"
stride: int = 3
- color: str = "#243e44"
- size: float = 3.0
+ color: str = "#1a1a1a"
+ size: float = 1.2
hatch: str = "...."
label: str | None = None
legend: bool = True
diff --git a/pyproject.toml b/pyproject.toml
index 6772f49..f1739b8 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -4,16 +4,25 @@ build-backend = "setuptools.build_meta"
[project]
name = "plot-function"
-version = "0.3.0"
-description = "From NetCDF to beautiful maps, with xarray, Cartopy, and Matplotlib."
+version = "0.4.0"
+description = "Journal-grade maps from NetCDF: xarray + Cartopy + Matplotlib, with insets, marginal profiles and significance layers."
readme = "README.md"
requires-python = ">=3.10"
license = "MIT"
authors = [{name = "Longhao Wang"}]
-keywords = ["netcdf", "climate", "cartopy", "xarray", "visualization"]
+keywords = ["netcdf", "climate", "hydrology", "cartopy", "xarray", "matplotlib", "visualization", "scientific-figures", "geoscience"]
+classifiers = [
+ "Programming Language :: Python :: 3",
+ "Programming Language :: Python :: 3 :: Only",
+ "Intended Audience :: Science/Research",
+ "Topic :: Scientific/Engineering :: Atmospheric Science",
+ "Topic :: Scientific/Engineering :: GIS",
+ "Topic :: Scientific/Engineering :: Visualization",
+ "Framework :: Matplotlib",
+]
dependencies = [
"numpy>=1.24",
- "matplotlib>=3.7",
+ "matplotlib>=3.8",
"cartopy>=0.22",
"xarray>=2023.1",
"netCDF4>=1.6",
@@ -27,7 +36,10 @@ dev = ["pytest>=8", "ruff>=0.9", "build>=1"]
[project.urls]
Homepage = "https://github.com/GISWLH/plot-function"
+Gallery = "https://github.com/GISWLH/plot-function#gallery"
Issues = "https://github.com/GISWLH/plot-function/issues"
+Changelog = "https://github.com/GISWLH/plot-function/blob/main/CHANGELOG.md"
+Author = "https://giswlh.github.io/"
[project.scripts]
plot-function = "plot_function.cli:main"
diff --git a/tests/test_journal.py b/tests/test_journal.py
new file mode 100644
index 0000000..2898a8f
--- /dev/null
+++ b/tests/test_journal.py
@@ -0,0 +1,116 @@
+import cartopy.crs as ccrs
+import matplotlib as mpl
+import matplotlib.pyplot as plt
+import numpy as np
+import pytest
+from matplotlib.colors import BoundaryNorm
+
+import plot_function.journal as pj
+
+
+def test_style_is_temporary_and_uses_sans_list():
+ before = mpl.rcParams.copy()
+ with pj.journal_style():
+ assert mpl.rcParams["font.size"] == 7
+ assert mpl.rcParams["pdf.fonttype"] == 42
+ assert mpl.rcParams["font.family"][-1] == "DejaVu Sans"
+ assert mpl.rcParams == before
+
+
+def test_figsize_columns():
+ w, h = pj.figsize("single", 0.5)
+ assert w == pytest.approx(89 / 25.4)
+ assert h == pytest.approx(w / 2)
+ assert pj.figsize(120)[0] == pytest.approx(120 / 25.4)
+
+
+@pytest.mark.parametrize("extend,n_colors", [("both", 4), ("min", 4), ("neither", 4)])
+def test_discrete_cmap_extends(extend, n_colors):
+ cmap, norm = pj.discrete_cmap("BrBG", [-2, -1, 0, 1, 2], extend=extend)
+ assert cmap.N == n_colors
+ assert isinstance(norm, BoundaryNorm)
+ assert cmap.colorbar_extend == extend
+ with pytest.raises(ValueError):
+ pj.discrete_cmap("BrBG", [0, 0, 1])
+
+
+def test_discrete_cmap_from_colour_list():
+ colors = ["#000000", "#111111", "#222222", "#333333"]
+ cmap, _ = pj.discrete_cmap(colors, [0, 1, 2], extend="both")
+ assert mpl.colors.to_hex(cmap.get_under()) == "#000000"
+ assert mpl.colors.to_hex(cmap.get_over()) == "#333333"
+
+
+def test_formatters():
+ assert pj.format_lat(40) == "40°N"
+ assert pj.format_lat(-20) == "20°S"
+ assert pj.format_lat(0) == "0°"
+ assert pj.format_lon(-120) == "120°W"
+ assert pj.format_lon(180) == "180°"
+
+
+@pytest.fixture
+def geo():
+ lon = np.arange(-179, 180, 2.0)
+ lat = np.arange(-59, 90, 2.0)
+ field = np.cos(np.deg2rad(lat))[:, None] * np.sin(np.deg2rad(lon))[None, :] * 10
+ return lon, lat, field
+
+
+@pytest.mark.parametrize("projection", [ccrs.PlateCarree(), ccrs.Robinson()])
+def test_full_composition_offline(geo, projection, tmp_path):
+ lon, lat, field = geo
+ with pj.journal_style():
+ fig = plt.figure(figsize=pj.figsize("double", 0.5))
+ ax = fig.add_axes([0.05, 0.25, 0.8, 0.7], projection=projection)
+ ax.set_global()
+ cmap, norm = pj.discrete_cmap("BrBG", np.arange(-10, 11, 2.5))
+ mesh = ax.pcolormesh(lon, lat, field, cmap=cmap, norm=norm, transform=ccrs.PlateCarree())
+ pj.geo_ticks(ax, xticks=range(-180, 181, 60), yticks=range(-60, 61, 30))
+ dots = pj.add_stippling(ax, lon, lat, np.abs(field) < 2, stride=2)
+ assert dots is not None and len(dots.get_offsets()) > 0
+ assert pj.add_hatching(ax, lon, lat, field > 8) is not None
+ assert pj.add_stippling(ax, lon, lat, np.zeros_like(field, bool)) is None
+ bars = pj.add_inset_bars(ax, ["A", "B"], [30, 20], hatched=[10, 5])
+ hist, twin = pj.add_inset_histogram(
+ ax, {"x": np.abs(field).ravel() + 0.1, "y": np.ones(5)}, log=True, cumulative=True
+ )
+ assert twin is not None
+ profile = pj.add_latitude_profile(
+ ax, lat, field.mean(1), lower=field.min(1), upper=field.max(1)
+ )
+ top = pj.add_longitude_profile(ax, lon, field.mean(0))
+ cbar = pj.add_colorbar(mesh, ax, label="Δ", title="Regime")
+ assert cbar.extend == "both"
+ pj.add_panel_label(ax, "a")
+ out = pj.save_figure(fig, tmp_path / "fig.png", dpi=120)
+ assert out.stat().st_size > 10_000
+ # The profile's latitude axis is locked to the map's projected y-range.
+ assert profile.get_ylim() == pytest.approx(ax.get_ylim())
+ assert top.get_xlim() == pytest.approx(ax.get_xlim())
+ assert bars in ax.child_axes and hist in ax.child_axes
+ plt.close(fig)
+
+
+def test_profile_alignment_matches_projection(geo):
+ lon, lat, field = geo
+ fig = plt.figure()
+ ax = fig.add_subplot(projection=ccrs.Robinson())
+ ax.set_global()
+ child = pj.add_latitude_profile(ax, lat, field.mean(1), yticks=[-30, 0, 30, 60])
+ expected = ccrs.Robinson().transform_point(0, 60, ccrs.PlateCarree())[1]
+ assert np.min(np.abs(np.asarray(child.get_yticks()) - expected)) < 1e-6
+ fig.canvas.draw()
+ labels = [t.get_text() for t in child.get_yticklabels()]
+ assert "60°N" in labels and "30°S" in labels
+
+
+def test_size_legend_and_bad_mask(geo):
+ lon, lat, field = geo
+ fig, ax = plt.subplots()
+ leg = pj.add_size_legend(ax, [1, 10], ["1", "10"], scale=lambda v: 10 * v)
+ assert len(leg.legend_handles) == 2
+ gfig = plt.figure()
+ gax = gfig.add_subplot(projection=ccrs.PlateCarree())
+ with pytest.raises(ValueError, match="shape"):
+ pj.add_stippling(gax, lon, lat, np.ones((3, 3), bool))