# Author: LKouadio <etanoyau@gmail.com>
# License: LGPL-3.0
"""2-D station-map API for pyCSAMT surveys."""
from __future__ import annotations
from dataclasses import replace
from typing import Any
import numpy as np
from ._backends import (
require_matplotlib_pyplot,
require_plotly,
)
from ._core import (
MapData,
ensure_map_data,
skin_depth_at_frequency,
station_dataframe,
value_at_frequency,
)
from .config import StationMapOptions
from .overlays import (
build_basemap_layout,
build_profile_line_overlay,
)
from .styles import (
is_dark_theme,
theme_colors,
to_plotly_cmap,
)
[docs]
class StationMap:
"""Builder object for station-map figures."""
def __init__(
self,
data: Any,
*,
options: StationMapOptions | None = None,
**ensure_kwargs: Any,
) -> None:
self.data: MapData = ensure_map_data(
data,
**ensure_kwargs,
)
self.options = options or StationMapOptions()
[docs]
def with_overlay(
self,
overlay: str,
**kwargs: Any,
) -> StationMap:
"""Return a copy with a different overlay."""
opts = replace(
self.options,
overlay=overlay,
**kwargs,
)
new = object.__new__(StationMap)
new.data = self.data
new.options = opts
return new
[docs]
def with_options(self, **kwargs: Any) -> StationMap:
"""Return a copy with updated options."""
new = object.__new__(StationMap)
new.data = self.data
new.options = replace(self.options, **kwargs)
return new
[docs]
def plot_station_map(
data: Any,
*,
options: StationMapOptions | None = None,
**kwargs: Any,
) -> Any:
"""Build a 2-D station map."""
return StationMap(
data,
options=options,
**kwargs,
).figure()
[docs]
def build_station_map(
data: MapData,
options: StationMapOptions | None = None,
) -> Any:
"""Build a Plotly 2-D station map."""
opts = options or StationMapOptions()
df = _station_values(data, opts)
colors = theme_colors(opts.theme)
if opts.backend == "matplotlib":
return _matplotlib_station_map(df, opts, colors)
if opts.backend != "plotly":
msg = f"Unknown station map backend: {opts.backend}"
raise ValueError(msg)
dark = is_dark_theme(opts.theme)
if df.empty:
return _empty_station_figure(colors)
if _has_geo(df):
return _geo_station_map(df, opts, colors, dark)
return _profile_station_map(df, opts, colors)
def _station_values(
data: MapData,
opts: StationMapOptions,
):
df = station_dataframe(data)
if df.empty:
return df
overlay = opts.overlay.lower()
if overlay in {"index", "station"}:
df["_value"] = df["Index"]
df["_label"] = "Station index"
elif overlay == "elevation":
df["_value"] = df["Elevation"]
df["_label"] = "Elevation (m)"
elif overlay in {"rho", "resistivity"}:
freq = getattr(opts, "frequency", None) or 1.0
vals = value_at_frequency(
data,
frequency=float(freq),
quantity="rho",
component=opts.component,
tolerance=opts.frequency_tolerance,
)
df["_value"] = df["ID"].map(vals)
df["_label"] = "ρ (Ω·m)"
elif overlay == "phase":
freq = getattr(opts, "frequency", None) or 1.0
vals = value_at_frequency(
data,
frequency=float(freq),
quantity="phase",
component=opts.component,
tolerance=opts.frequency_tolerance,
)
df["_value"] = df["ID"].map(vals)
df["_label"] = "φ (°)"
elif overlay in {"depth", "skin_depth"}:
freq = getattr(opts, "frequency", None) or 1.0
vals = skin_depth_at_frequency(
data,
frequency=float(freq),
component=opts.component,
tolerance=opts.frequency_tolerance,
)
df["_value"] = df["ID"].map(vals)
df["_label"] = "δ skin depth (m)"
else:
if opts.overlay in df:
df["_value"] = df[opts.overlay]
else:
df["_value"] = df["Index"]
df["_label"] = opts.overlay
df["_plot_value"] = _plot_values(df["_value"], opts)
return df
def _has_geo(df) -> bool:
lat = np.asarray(df["Latitude"], dtype=float)
lon = np.asarray(df["Longitude"], dtype=float)
return bool(np.isfinite(lat).any() and np.isfinite(lon).any())
def _geo_station_map(df, opts, colors, dark: bool):
go = require_plotly()
fig = go.Figure()
for line, group in df.groupby("Line", dropna=False):
line_name = str(line or "line")
if opts.line_filter and line_name not in opts.line_filter:
continue
lat = group["Latitude"].astype(float)
lon = group["Longitude"].astype(float)
marker_size = _marker_sizes(group["ID"], opts)
if opts.show_contours:
_add_density_layer(fig, group, opts)
fig.add_trace(
_scatter_map_trace(
go,
lat=lat,
lon=lon,
mode=_marker_mode(opts.show_labels),
text=group["ID"],
customdata=group["ID"],
name=line_name,
marker=dict(
size=marker_size,
color=group["_plot_value"],
colorscale=to_plotly_cmap(opts.cmap),
opacity=opts.opacity,
showscale=True,
cmin=_cmin(opts),
cmax=_cmax(opts),
colorbar=dict(
title=dict(
text=_color_title(group, opts),
side="right",
),
),
),
textposition="top right",
hovertemplate=(
"<b>%{text}</b><br>"
"lat=%{lat:.5f}<br>lon=%{lon:.5f}"
"<extra></extra>"
),
)
)
if opts.show_profiles and len(group) > 1:
fig.add_trace(
build_profile_line_overlay(
lon,
lat,
geo=True,
name=f"{line_name} profile",
color=colors["accent"],
width=2,
)
)
basemap = build_basemap_layout(
df["Longitude"],
df["Latitude"],
dark=dark,
style=opts.basemap,
bearing=opts.bearing,
)
layers = list(basemap.layers)
contour_layer = _contour_image_layer(df, opts)
if contour_layer is not None:
layers.append(contour_layer)
map_layout = dict(
style=basemap.style,
center=basemap.center,
zoom=basemap.zoom,
bearing=basemap.bearing,
)
if layers:
map_layout["layers"] = layers
fig.update_layout(
**{_map_layout_key(fig): map_layout},
paper_bgcolor=colors["paper"],
font=dict(color=colors["text"]),
margin=dict(l=0, r=0, t=30, b=0),
legend=dict(orientation="h"),
title=opts.title or None,
)
return fig
def _contour_image_layer(df, opts):
"""Build a Surfer-style filled-contour image map layer, or None."""
if not opts.contour_image:
return None
from .overlays import build_geo_contour_image
overlay = build_geo_contour_image(
df["Longitude"],
df["Latitude"],
df["_value"],
cmap=opts.cmap,
n_levels=opts.contour_levels,
opacity=float(opts.contour_opacity),
mode=opts.contour_mode,
log_scale=bool(opts.log_color),
grid_res=int(opts.contour_grid_res),
interp=opts.contour_interp,
smooth_sigma=float(opts.contour_smooth),
)
if overlay is None:
return None
return {
"below": "traces",
"sourcetype": "image",
"source": overlay["image"],
"coordinates": overlay["coordinates"],
}
def _scatter_map_trace(go, **kwargs):
if hasattr(go, "Scattermap"):
return go.Scattermap(**kwargs)
return go.Scattermapbox(**kwargs)
def _density_map_trace(go, **kwargs):
if hasattr(go, "Densitymap"):
return go.Densitymap(**kwargs)
return go.Densitymapbox(**kwargs)
def _map_layout_key(fig) -> str:
modern = {"scattermap", "densitymap"}
if any(getattr(trace, "type", "") in modern for trace in fig.data):
return "map"
return "mapbox"
def _matplotlib_station_map(df, opts, colors):
plt = require_matplotlib_pyplot()
fig, ax = plt.subplots(figsize=(8, 5))
fig.patch.set_facecolor(colors["paper"])
ax.set_facecolor(colors["plot"])
ax.tick_params(colors=colors["text"])
if df.empty:
ax.text(
0.5,
0.5,
"No stations available",
ha="center",
va="center",
color=colors["text"],
)
return fig
x_col = "Longitude" if _has_geo(df) else "Index"
y_col = "Latitude" if _has_geo(df) else "_value"
sizes = (
np.asarray(
_marker_sizes(df["ID"], opts),
dtype=float,
)
** 2
)
sc = ax.scatter(
df[x_col],
df[y_col],
c=df["_plot_value"],
s=sizes,
cmap=opts.cmap,
alpha=float(opts.opacity),
vmin=_cmin(opts),
vmax=_cmax(opts),
)
if opts.show_labels:
for _, row in df.iterrows():
ax.text(
row[x_col],
row[y_col],
str(row["ID"]),
color=colors["text"],
fontsize=8,
)
if opts.show_profiles and _has_geo(df):
for _, group in df.groupby("Line", dropna=False):
ax.plot(
group["Longitude"],
group["Latitude"],
color=colors["accent"],
linewidth=1.2,
)
cbar = fig.colorbar(sc, ax=ax)
cbar.set_label(
_color_title(df, opts),
color=colors["text"],
)
ax.set_xlabel(x_col, color=colors["text"])
ax.set_ylabel(y_col, color=colors["text"])
if opts.title:
ax.set_title(opts.title, color=colors["text"])
return fig
def _profile_station_map(df, opts, colors):
go = require_plotly()
fig = go.Figure()
fig.add_trace(
go.Scatter(
x=df["Index"],
y=df["_value"],
mode=_marker_mode(opts.show_labels),
text=df["ID"],
marker=dict(
size=_marker_sizes(df["ID"], opts),
color=df["_plot_value"],
colorscale=to_plotly_cmap(opts.cmap),
showscale=True,
cmin=_cmin(opts),
cmax=_cmax(opts),
colorbar=dict(
title=dict(text=_color_title(df, opts), side="right")
),
),
)
)
fig.update_layout(
xaxis_title="Station index",
yaxis_title=str(df["_label"].iloc[0]),
paper_bgcolor=colors["paper"],
plot_bgcolor=colors["plot"],
font=dict(color=colors["text"]),
title=opts.title or None,
)
return fig
def _plot_values(values, opts):
arr = np.asarray(values, dtype=float)
if opts.log_color:
arr = np.where(arr > 0, np.log10(arr), np.nan)
return arr
def _color_title(df, opts) -> str:
label = str(df["_label"].iloc[0])
return f"log₁₀ {label}" if opts.log_color else label
def _cmin(opts) -> float | None:
if not opts.value_range:
return None
lo, _ = opts.value_range
if opts.log_color and lo > 0:
return float(np.log10(lo))
return float(lo)
def _cmax(opts) -> float | None:
if not opts.value_range:
return None
_, hi = opts.value_range
if opts.log_color and hi > 0:
return float(np.log10(hi))
return float(hi)
def _add_density_layer(fig, group, opts) -> None:
go = require_plotly()
values = np.asarray(group["_plot_value"], dtype=float)
lat = np.asarray(group["Latitude"], dtype=float)
lon = np.asarray(group["Longitude"], dtype=float)
good = np.isfinite(values) & np.isfinite(lat) & np.isfinite(lon)
if good.sum() < 3:
return
fig.add_trace(
_density_map_trace(
go,
lat=lat[good],
lon=lon[good],
z=values[good],
radius=max(10, int(opts.marker_size * 2)),
colorscale=to_plotly_cmap(opts.cmap),
opacity=float(opts.contour_opacity),
showscale=False,
hoverinfo="skip",
name="overlay density",
)
)
def _marker_sizes(ids, opts) -> list[int]:
selected = str(opts.selected_id) if opts.selected_id else None
return [
int(opts.marker_size * 1.8)
if selected == str(sid)
else int(opts.marker_size)
for sid in ids
]
def _marker_mode(show_labels: bool) -> str:
return "markers+text" if show_labels else "markers"
def _empty_station_figure(colors):
go = require_plotly()
fig = go.Figure()
fig.add_annotation(
text="No stations available",
x=0.5,
y=0.5,
xref="paper",
yref="paper",
showarrow=False,
font=dict(color=colors["text"]),
)
fig.update_layout(
paper_bgcolor=colors["paper"],
plot_bgcolor=colors["plot"],
)
return fig