diff --git a/apps/geolibre-desktop/src/components/layout/toolbar/ProcessingMenu.tsx b/apps/geolibre-desktop/src/components/layout/toolbar/ProcessingMenu.tsx index 89773aad6..e6d570dff 100644 --- a/apps/geolibre-desktop/src/components/layout/toolbar/ProcessingMenu.tsx +++ b/apps/geolibre-desktop/src/components/layout/toolbar/ProcessingMenu.tsx @@ -309,13 +309,16 @@ export function ProcessingMenu({ - {t("toolbar.item.subGroupH3")} + {t("toolbar.item.subGroupDggs")} - setVectorToolOpen("h3-grid")}> - {t("toolbar.vectorTool.h3Grid")} + setVectorToolOpen("dggs-grid")}> + {t("toolbar.vectorTool.dggsGenerator")} - setVectorToolOpen("h3-bin-points")}> - {t("toolbar.vectorTool.h3BinPoints")} + setVectorToolOpen("dggs-bin")}> + {t("toolbar.vectorTool.dggsBinning")} + + setVectorToolOpen("dggs-compact")}> + {t("toolbar.vectorTool.dggsCompact")} diff --git a/apps/geolibre-desktop/src/components/layout/toolbar/constants.ts b/apps/geolibre-desktop/src/components/layout/toolbar/constants.ts index d7bb7206a..a2fcc8a2d 100644 --- a/apps/geolibre-desktop/src/components/layout/toolbar/constants.ts +++ b/apps/geolibre-desktop/src/components/layout/toolbar/constants.ts @@ -171,8 +171,9 @@ export const VECTOR_TOOL_COMMANDS: Array<{ { kind: "smooth", titleKey: "toolbar.vectorTool.smooth" }, { kind: "grid", titleKey: "toolbar.vectorTool.grid" }, { kind: "voronoi", titleKey: "toolbar.vectorTool.voronoi" }, - { kind: "h3-grid", titleKey: "toolbar.vectorTool.h3Grid" }, - { kind: "h3-bin-points", titleKey: "toolbar.vectorTool.h3BinPoints" }, + { kind: "dggs-grid", titleKey: "toolbar.vectorTool.dggsGenerator" }, + { kind: "dggs-bin", titleKey: "toolbar.vectorTool.dggsBinning" }, + { kind: "dggs-compact", titleKey: "toolbar.vectorTool.dggsCompact" }, ]; export const RASTER_TOOL_COMMANDS: Array<{ diff --git a/apps/geolibre-desktop/src/components/processing/ParameterField.tsx b/apps/geolibre-desktop/src/components/processing/ParameterField.tsx index 85a49775c..6d8b2f89d 100644 --- a/apps/geolibre-desktop/src/components/processing/ParameterField.tsx +++ b/apps/geolibre-desktop/src/components/processing/ParameterField.tsx @@ -106,16 +106,21 @@ export function ParameterField({ if (param.type === "boolean") { return ( - +
+ + {param.description ? ( +

{param.description}

+ ) : null} +
); } diff --git a/apps/geolibre-desktop/src/components/processing/ProcessingHistoryDialog.tsx b/apps/geolibre-desktop/src/components/processing/ProcessingHistoryDialog.tsx index 6e70a9796..809673c38 100644 --- a/apps/geolibre-desktop/src/components/processing/ProcessingHistoryDialog.tsx +++ b/apps/geolibre-desktop/src/components/processing/ProcessingHistoryDialog.tsx @@ -6,6 +6,7 @@ import { type StatisticsToolKind, type VectorToolKind, } from "@geolibre/core"; +import { resolveVectorRerun } from "@geolibre/processing"; import { allAlgorithms } from "../../lib/scripting/scriptingApi"; import { Button, @@ -127,16 +128,22 @@ export function ProcessingHistoryDialog(): ReactElement { // the reopened dialog is visible. const openForRun = useCallback( (run: ProcessingRun, autoRun: boolean) => { + // Pre-DGGS history may still store h3-grid / h3-bin-points — map those + // onto dggs-* before queuing the re-run and opening the dialog. + const vectorResolved = + run.kind === "vector" ? resolveVectorRerun(run.toolId, run.parameters) : null; + const toolId = vectorResolved?.toolId ?? run.toolId; + const parameters = vectorResolved?.parameters ?? run.parameters; setProcessingRerun({ kind: run.kind, - toolId: run.toolId, - parameters: run.parameters, + toolId, + parameters, engine: run.engine, autoRun, }); switch (run.kind) { case "vector": - setVectorToolOpen(run.toolId as VectorToolKind); + setVectorToolOpen(toolId as VectorToolKind); break; case "statistics": setStatisticsToolOpen(run.toolId as StatisticsToolKind); diff --git a/apps/geolibre-desktop/src/components/processing/VectorToolsDialog.tsx b/apps/geolibre-desktop/src/components/processing/VectorToolsDialog.tsx index 03a868915..9dd52b083 100644 --- a/apps/geolibre-desktop/src/components/processing/VectorToolsDialog.tsx +++ b/apps/geolibre-desktop/src/components/processing/VectorToolsDialog.tsx @@ -3,8 +3,10 @@ import { detectGeometryProfile, type MapController } from "@geolibre/map"; import { VECTOR_TOOLS, getVectorTool, + resolveVectorRerun, runVectorTool, fetchVectorStatus, + maxResolutionForDggs, type AlgorithmParameter, type GeometryFamily, type ProcessingAlgorithm, @@ -109,9 +111,10 @@ export function VectorToolsDialog({ mapControllerRef }: VectorToolsDialogProps): // win over the defaults when both effects fire in the same commit. useEffect(() => { if (!open || !rerun || rerun.kind !== "vector") return; - // A saved-project history entry can reference a tool that was renamed or - // removed since; drop the request instead of leaving it pending forever. - if (!getVectorTool(rerun.toolId)) { + // History may still carry pre-DGGS H3 tool ids — map them before lookup so + // we do not clear the request as "unavailable" and leave it pending. + const resolved = resolveVectorRerun(rerun.toolId, rerun.parameters); + if (!getVectorTool(resolved.toolId)) { setLog((prev) => [ ...prev, `Error: ${t("processing.history.toolUnavailable", { toolId: rerun.toolId })}`, @@ -119,8 +122,13 @@ export function VectorToolsDialog({ mapControllerRef }: VectorToolsDialogProps): setProcessingRerun(null); return; } - if (rerun.toolId !== tool.id) return; - setParams({ ...rerun.parameters }); + if (resolved.toolId !== tool.id) { + // Alias mapped to a different tool than the one currently selected + // (e.g. History opened with a stale VectorToolKind). Switch over. + if (selectedId !== resolved.toolId) setSelectedId(resolved.toolId); + return; + } + setParams({ ...resolved.parameters }); if (rerun.engine === "client" || rerun.engine === "sidecar" || rerun.engine === "pyodide") { // A history entry recorded on a sidecar run re-runs on Pyodide in the // Mac App Store build, where the sidecar engine does not exist. @@ -128,14 +136,14 @@ export function VectorToolsDialog({ mapControllerRef }: VectorToolsDialogProps): } setProcessingRerun(null); if (rerun.autoRun) setAutoRunPending(true); - }, [open, rerun, tool, setProcessingRerun, t]); + }, [open, rerun, tool, selectedId, setProcessingRerun, t]); - // Prefill the H3 grid's manual bounding-box fields from the current map + // Prefill the DGGS Generator's manual bounding-box fields from the current map // viewport when the user first switches to that source, so they can tweak the // box rather than type it from scratch. Only fills empty fields, so it never // clobbers manual edits. Keyed on the source value, not every keystroke. useEffect(() => { - if (selectedId !== "h3-grid" || params.source !== "bbox") return; + if (selectedId !== "dggs-grid" || params.source !== "bbox") return; if ( params.west !== undefined || params.south !== undefined || @@ -475,16 +483,43 @@ export function VectorToolsDialog({ mapControllerRef }: VectorToolsDialogProps):

{tool.description}

- {tool.parameters.filter(isParamVisible).map((param) => ( - handleParamChange(param.id, value)} - /> - ))} + {tool.parameters.filter(isParamVisible).map((param) => { + // Narrow the resolution spinner to the selected DGGS type's range. + const fieldParam = + (tool.id === "dggs-grid" || + tool.id === "dggs-bin" || + tool.id === "dggs-compact") && + param.id === "resolution" + ? (() => { + const dggsType = + params.dggsType === "s2" || + params.dggsType === "a5" || + params.dggsType === "dggrid" || + params.dggsType === "dggal" + ? params.dggsType + : "h3"; + const rawSubtype = + dggsType === "dggal" ? params.dggalType : params.dggridType; + const subtype = typeof rawSubtype === "string" ? rawSubtype : undefined; + const max = maxResolutionForDggs(dggsType, subtype); + return { + ...param, + max, + label: t("processing.vectorTools.resolutionRange", { max }), + }; + })() + : param; + return ( + handleParamChange(param.id, value)} + /> + ); + })}
{tool.supportsSidecar || tool.requiresSidecar ? ( diff --git a/apps/geolibre-desktop/src/i18n/locales/ar.json b/apps/geolibre-desktop/src/i18n/locales/ar.json index 40fa6ace7..59d27a09c 100644 --- a/apps/geolibre-desktop/src/i18n/locales/ar.json +++ b/apps/geolibre-desktop/src/i18n/locales/ar.json @@ -2405,8 +2405,6 @@ "grid": "شبكة منتظمة", "voronoi": "Voronoi / Delaunay", "cellSectors": "تغطية المواقع الخلوية", - "h3Grid": "إنشاء شبكة H3", - "h3BinPoints": "تجميع النقاط في خلايا H3", "trajectorySpeed": "سرعة مسار الحركة", "detectStops": "اكتشاف التوقفات", "spaceTimeProximity": "التقارب الزماني المكاني", diff --git a/apps/geolibre-desktop/src/i18n/locales/de.json b/apps/geolibre-desktop/src/i18n/locales/de.json index 9859150e5..91b78b75b 100644 --- a/apps/geolibre-desktop/src/i18n/locales/de.json +++ b/apps/geolibre-desktop/src/i18n/locales/de.json @@ -2234,8 +2234,6 @@ "grid": "Regelmäßiges Raster", "voronoi": "Voronoi / Delaunay", "cellSectors": "Mobilfunkzellen-Abdeckung", - "h3Grid": "H3-Raster erstellen", - "h3BinPoints": "Punkte zu H3 binnen", "trajectorySpeed": "Trajektoriengeschwindigkeit", "detectStops": "Stopps erkennen", "spaceTimeProximity": "Raum-Zeit-Nähe", diff --git a/apps/geolibre-desktop/src/i18n/locales/en.json b/apps/geolibre-desktop/src/i18n/locales/en.json index 14ccbf6a4..d9545feb5 100644 --- a/apps/geolibre-desktop/src/i18n/locales/en.json +++ b/apps/geolibre-desktop/src/i18n/locales/en.json @@ -2235,8 +2235,9 @@ "grid": "Regular grid", "voronoi": "Voronoi / Delaunay", "cellSectors": "Cell-site coverage", - "h3Grid": "Create H3 grid", - "h3BinPoints": "Bin points to H3", + "dggsGenerator": "DGGS Generator", + "dggsBinning": "DGGS Binning", + "dggsCompact": "DGGS Compact", "trajectorySpeed": "Trajectory speed", "detectStops": "Detect stops", "spaceTimeProximity": "Space-time proximity", @@ -2421,7 +2422,8 @@ "subGroupOverlay": "Overlay", "subGroupJoin": "Join", "subGroupSelect": "Select", - "subGroupH3": "H3", + "subGroupDggs": "DGGS", + "subGroupH3": "DGGS", "subGroupMovement": "Movement & time", "subGroupDataQuality": "Data quality", "subGroupTerrain": "Terrain", @@ -3002,9 +3004,9 @@ "title": "OLC", "controlTitle": "OLC settings", "autoResolution": "Automatic resolution", - "resolution": "Code length", + "resolution": "Resolution", "cellCount": "{{count}} cells in view", - "tooManyCells": "This view exceeds the {{limit}} cell limit. Zoom in or lower the code length.", + "tooManyCells": "This view exceeds the {{limit}} cell limit. Zoom in or lower the resolution.", "fillColor": "Fill color", "fillOpacity": "Fill opacity", "lineColor": "Outline color", @@ -3029,9 +3031,9 @@ "title": "Geohash", "controlTitle": "Geohash settings", "autoResolution": "Automatic resolution", - "resolution": "Precision", + "resolution": "Resolution", "cellCount": "{{count}} cells in view", - "tooManyCells": "This view exceeds the {{limit}} cell limit. Zoom in or lower the precision.", + "tooManyCells": "This view exceeds the {{limit}} cell limit. Zoom in or lower the resolution.", "fillColor": "Fill color", "fillOpacity": "Fill opacity", "lineColor": "Outline color", @@ -3056,9 +3058,9 @@ "title": "Tilecode", "controlTitle": "Tilecode settings", "autoResolution": "Automatic resolution", - "resolution": "Zoom level", + "resolution": "Resolution", "cellCount": "{{count}} tiles in view", - "tooManyCells": "This view exceeds the {{limit}} tile limit. Zoom in or lower the zoom level.", + "tooManyCells": "This view exceeds the {{limit}} tile limit. Zoom in or lower the resolution.", "fillColor": "Fill color", "fillOpacity": "Fill opacity", "lineColor": "Outline color", @@ -3993,7 +3995,8 @@ "sidecarUnavailablePyodide": "The GeoPandas sidecar is not available. Start the sidecar with the vector extra, or switch to Python (Pyodide).", "pyodideNote": "Runs GeoPandas in your browser. The first run downloads the Python runtime (one-time, needs an internet connection).", "run": "Run", - "outputPlaceholder": "Output will appear here." + "outputPlaceholder": "Output will appear here.", + "resolutionRange": "Resolution (0-{{max}})" } }, "segmentation": { diff --git a/apps/geolibre-desktop/src/i18n/locales/es.json b/apps/geolibre-desktop/src/i18n/locales/es.json index 8578c1bab..08832516b 100644 --- a/apps/geolibre-desktop/src/i18n/locales/es.json +++ b/apps/geolibre-desktop/src/i18n/locales/es.json @@ -2234,8 +2234,6 @@ "grid": "Cuadrícula regular", "voronoi": "Voronoi / Delaunay", "cellSectors": "Cobertura de sitios celulares", - "h3Grid": "Crear cuadrícula H3", - "h3BinPoints": "Agrupar puntos en H3", "trajectorySpeed": "Velocidad de trayectoria", "detectStops": "Detectar paradas", "spaceTimeProximity": "Proximidad espacio-temporal", diff --git a/apps/geolibre-desktop/src/i18n/locales/fr.json b/apps/geolibre-desktop/src/i18n/locales/fr.json index e34ae0dec..0d44e9aad 100644 --- a/apps/geolibre-desktop/src/i18n/locales/fr.json +++ b/apps/geolibre-desktop/src/i18n/locales/fr.json @@ -2234,8 +2234,6 @@ "grid": "Grille régulière", "voronoi": "Voronoï / Delaunay", "cellSectors": "Couverture de site cellulaire", - "h3Grid": "Créer une grille H3", - "h3BinPoints": "Regrouper les points en H3", "trajectorySpeed": "Vitesse de trajectoire", "detectStops": "Détecter les arrêts", "spaceTimeProximity": "Proximité spatio-temporelle", diff --git a/apps/geolibre-desktop/src/i18n/locales/hi.json b/apps/geolibre-desktop/src/i18n/locales/hi.json index c6c2b4981..ec789bfe3 100644 --- a/apps/geolibre-desktop/src/i18n/locales/hi.json +++ b/apps/geolibre-desktop/src/i18n/locales/hi.json @@ -2234,8 +2234,6 @@ "grid": "नियमित ग्रिड", "voronoi": "Voronoi / Delaunay", "cellSectors": "सेल-साइट कवरेज", - "h3Grid": "H3 ग्रिड बनाएं", - "h3BinPoints": "बिंदुओं को H3 में बिन करें", "trajectorySpeed": "प्रक्षेपवक्र गति", "detectStops": "रुकावटें पहचानें", "spaceTimeProximity": "स्थान-समय निकटता", diff --git a/apps/geolibre-desktop/src/i18n/locales/id.json b/apps/geolibre-desktop/src/i18n/locales/id.json index 910af0db2..c3588d319 100644 --- a/apps/geolibre-desktop/src/i18n/locales/id.json +++ b/apps/geolibre-desktop/src/i18n/locales/id.json @@ -2191,8 +2191,6 @@ "grid": "Grid teratur", "voronoi": "Voronoi / Delaunay", "cellSectors": "Cakupan sel BTS", - "h3Grid": "Buat grid H3", - "h3BinPoints": "Bin titik ke H3", "trajectorySpeed": "Kecepatan lintasan", "detectStops": "Deteksi pemberhentian", "spaceTimeProximity": "Kedekatan ruang-waktu", diff --git a/apps/geolibre-desktop/src/i18n/locales/it.json b/apps/geolibre-desktop/src/i18n/locales/it.json index e351cf80f..487371769 100644 --- a/apps/geolibre-desktop/src/i18n/locales/it.json +++ b/apps/geolibre-desktop/src/i18n/locales/it.json @@ -2234,8 +2234,6 @@ "grid": "Griglia regolare", "voronoi": "Voronoi / Delaunay", "cellSectors": "Copertura delle celle radio", - "h3Grid": "Crea griglia H3", - "h3BinPoints": "Assegna punti a celle H3", "trajectorySpeed": "Velocità di traiettoria", "detectStops": "Rileva soste", "spaceTimeProximity": "Prossimità spazio-temporale", diff --git a/apps/geolibre-desktop/src/i18n/locales/ja.json b/apps/geolibre-desktop/src/i18n/locales/ja.json index 13fd0dba0..6a97f39cd 100644 --- a/apps/geolibre-desktop/src/i18n/locales/ja.json +++ b/apps/geolibre-desktop/src/i18n/locales/ja.json @@ -2191,8 +2191,6 @@ "grid": "正方格子", "voronoi": "ボロノイ / ドロネー", "cellSectors": "セルサイトのカバレッジ", - "h3Grid": "H3グリッドを作成", - "h3BinPoints": "ポイントをH3にビン化", "trajectorySpeed": "軌跡速度", "detectStops": "停止を検出", "spaceTimeProximity": "時空間近接", diff --git a/apps/geolibre-desktop/src/i18n/locales/ka.json b/apps/geolibre-desktop/src/i18n/locales/ka.json index 6a03c9f57..ba00f87ad 100644 --- a/apps/geolibre-desktop/src/i18n/locales/ka.json +++ b/apps/geolibre-desktop/src/i18n/locales/ka.json @@ -2234,8 +2234,6 @@ "grid": "რეგულარული ბადე", "voronoi": "ვორონოი / დელონე", "cellSectors": "ფიჭური საიტის დაფარვა", - "h3Grid": "H3 ბადის შექმნა", - "h3BinPoints": "წერტილების დაჯგუფება H3-ში", "trajectorySpeed": "ტრაექტორიის სიჩქარე", "detectStops": "გაჩერებების ამოცნობა", "spaceTimeProximity": "სივრცე-დროითი სიახლოვე", diff --git a/apps/geolibre-desktop/src/i18n/locales/ko.json b/apps/geolibre-desktop/src/i18n/locales/ko.json index 3e088f0ba..bad3789e0 100644 --- a/apps/geolibre-desktop/src/i18n/locales/ko.json +++ b/apps/geolibre-desktop/src/i18n/locales/ko.json @@ -2191,8 +2191,6 @@ "grid": "정규 격자", "voronoi": "Voronoi / Delaunay", "cellSectors": "기지국 커버리지", - "h3Grid": "H3 격자 생성", - "h3BinPoints": "포인트를 H3 셀로 구간화", "trajectorySpeed": "궤적 속도", "detectStops": "정지 지점 감지", "spaceTimeProximity": "시공간 근접성", diff --git a/apps/geolibre-desktop/src/i18n/locales/nl.json b/apps/geolibre-desktop/src/i18n/locales/nl.json index ff7080f3c..05e9531f8 100644 --- a/apps/geolibre-desktop/src/i18n/locales/nl.json +++ b/apps/geolibre-desktop/src/i18n/locales/nl.json @@ -2234,8 +2234,6 @@ "grid": "Regelmatig raster", "voronoi": "Voronoi / Delaunay", "cellSectors": "Dekkingsgebied zendmast", - "h3Grid": "H3-raster maken", - "h3BinPoints": "Punten binnen in H3", "trajectorySpeed": "Trajectsnelheid", "detectStops": "Stops detecteren", "spaceTimeProximity": "Ruimte-tijdnabijheid", diff --git a/apps/geolibre-desktop/src/i18n/locales/pt.json b/apps/geolibre-desktop/src/i18n/locales/pt.json index e6812924a..074d5d3e8 100644 --- a/apps/geolibre-desktop/src/i18n/locales/pt.json +++ b/apps/geolibre-desktop/src/i18n/locales/pt.json @@ -2234,8 +2234,6 @@ "grid": "Grade regular", "voronoi": "Voronoi / Delaunay", "cellSectors": "Cobertura de estações de célula", - "h3Grid": "Criar grade H3", - "h3BinPoints": "Agrupar pontos em H3", "trajectorySpeed": "Velocidade da trajetória", "detectStops": "Detectar paradas", "spaceTimeProximity": "Proximidade espaço-temporal", diff --git a/apps/geolibre-desktop/src/i18n/locales/ru.json b/apps/geolibre-desktop/src/i18n/locales/ru.json index 2680e5414..70adaf5c7 100644 --- a/apps/geolibre-desktop/src/i18n/locales/ru.json +++ b/apps/geolibre-desktop/src/i18n/locales/ru.json @@ -2320,8 +2320,6 @@ "grid": "Регулярная сетка", "voronoi": "Вороного / Делоне", "cellSectors": "Зоны покрытия базовых станций", - "h3Grid": "Создать сетку H3", - "h3BinPoints": "Группировать точки в H3", "trajectorySpeed": "Скорость траектории", "detectStops": "Определить остановки", "spaceTimeProximity": "Пространственно-временная близость", diff --git a/apps/geolibre-desktop/src/i18n/locales/th.json b/apps/geolibre-desktop/src/i18n/locales/th.json index 2947f7ec8..554468e2e 100644 --- a/apps/geolibre-desktop/src/i18n/locales/th.json +++ b/apps/geolibre-desktop/src/i18n/locales/th.json @@ -2192,8 +2192,6 @@ "grid": "กริดสม่ำเสมอ", "voronoi": "โวโรนอย / เดอโลเน", "cellSectors": "พื้นที่ครอบคลุมของสถานีฐาน", - "h3Grid": "สร้างกริด H3", - "h3BinPoints": "จัดกลุ่มจุดลงกริด H3", "trajectorySpeed": "ความเร็วตามเส้นทางการเคลื่อนที่", "detectStops": "ตรวจหาจุดหยุด", "spaceTimeProximity": "ความใกล้เคียงเชิงพื้นที่และเวลา", diff --git a/apps/geolibre-desktop/src/i18n/locales/tr.json b/apps/geolibre-desktop/src/i18n/locales/tr.json index 3ba0ed3cc..25ce3e231 100644 --- a/apps/geolibre-desktop/src/i18n/locales/tr.json +++ b/apps/geolibre-desktop/src/i18n/locales/tr.json @@ -2234,8 +2234,6 @@ "grid": "Düzenli ızgara", "voronoi": "Voronoi / Delaunay", "cellSectors": "Baz istasyonu kapsama alanı", - "h3Grid": "H3 ızgarası oluştur", - "h3BinPoints": "Noktaları H3'e grupla", "trajectorySpeed": "Güzergah hızı", "detectStops": "Durakları algıla", "spaceTimeProximity": "Uzay-zaman yakınlığı", diff --git a/apps/geolibre-desktop/src/i18n/locales/zh.json b/apps/geolibre-desktop/src/i18n/locales/zh.json index f23660a6f..c2bd84a3d 100644 --- a/apps/geolibre-desktop/src/i18n/locales/zh.json +++ b/apps/geolibre-desktop/src/i18n/locales/zh.json @@ -2191,8 +2191,6 @@ "grid": "规则格网", "voronoi": "泰森多边形 / 德劳内三角网", "cellSectors": "基站覆盖扇区", - "h3Grid": "创建 H3 网格", - "h3BinPoints": "将点聚合到 H3", "trajectorySpeed": "轨迹速度", "detectStops": "检测停留点", "spaceTimeProximity": "时空邻近性", diff --git a/apps/geolibre-desktop/src/lib/duckdb-processing.ts b/apps/geolibre-desktop/src/lib/duckdb-processing.ts index a028d5180..5409a3488 100644 --- a/apps/geolibre-desktop/src/lib/duckdb-processing.ts +++ b/apps/geolibre-desktop/src/lib/duckdb-processing.ts @@ -2,6 +2,8 @@ import type { DuckDbCapability, DuckDbGeoJsonSource } from "@geolibre/processing import type { FeatureCollection } from "geojson"; import { stripAutoFidColumn } from "./duckdb-geometry"; import { + ensureA5Extension, + ensureDuckDggsExtension, ensureH3Extension, ensureSpatialExtension, getDatabase, @@ -25,6 +27,8 @@ export function createDuckDbCapability(): DuckDbCapability { try { if (names.includes("spatial")) await ensureSpatialExtension(db, connection); if (names.includes("h3")) await ensureH3Extension(connection); + if (names.includes("a5")) await ensureA5Extension(connection); + if (names.includes("duck_dggs")) await ensureDuckDggsExtension(connection); } finally { await connection.close(); } diff --git a/apps/geolibre-desktop/src/lib/duckdb-vector-loader.ts b/apps/geolibre-desktop/src/lib/duckdb-vector-loader.ts index c20bdeb1d..45c20403d 100644 --- a/apps/geolibre-desktop/src/lib/duckdb-vector-loader.ts +++ b/apps/geolibre-desktop/src/lib/duckdb-vector-loader.ts @@ -238,30 +238,47 @@ export async function ensureSpatialExtension( } } -let h3ExtensionPromise: Promise | null = null; +/** + * Memoized INSTALL/LOAD for a DuckDB community extension. `name` must stay a + * module-internal literal — never accept caller-supplied SQL identifiers. + */ +function createCommunityExtensionLoader( + name: "h3" | "a5" | "duck_dggs", +): (connection: duckdb.AsyncDuckDBConnection) => Promise { + let promise: Promise | null = null; + return async (connection) => { + promise ??= (async () => { + await connection.query(`INSTALL ${name} FROM community`); + await connection.query(`LOAD ${name}`); + })(); + try { + await promise; + } catch (error) { + promise = null; + throw error; + } + }; +} /** * Install and load the DuckDB `h3` community extension once per database - * instance. Mirrors {@link ensureSpatialExtension}: memoized as a promise so - * concurrent callers share one INSTALL/LOAD, and cleared on failure so a later - * call can retry. `h3` is published for the bundled DuckDB version (v1.5.1) on - * all WASM platforms. + * instance. Concurrent callers share one INSTALL/LOAD; failures clear the + * memo so a later call can retry. `h3` is published for the bundled DuckDB + * version (v1.5.1) on all WASM platforms. */ -export async function ensureH3Extension(connection: duckdb.AsyncDuckDBConnection): Promise { - h3ExtensionPromise ??= (async () => { - // Unlike `ensureSpatialExtension`, no `beforeLoad` warm-up is needed here: - // the duckdb-wasm v1.33.1-dev45 remote-read bug only affects `spatial`. If a - // similar issue ever surfaces for `h3`, add a `beforeLoad` hook to match. - await connection.query("INSTALL h3 FROM community"); - await connection.query("LOAD h3"); - })(); - try { - await h3ExtensionPromise; - } catch (error) { - h3ExtensionPromise = null; - throw error; - } -} +export const ensureH3Extension = createCommunityExtensionLoader("h3"); + +/** + * Install and load the DuckDB `a5` community extension once per database + * instance. Mirrors {@link ensureH3Extension}. + */ +export const ensureA5Extension = createCommunityExtensionLoader("a5"); + +/** + * Install and load the DuckDB `duck_dggs` community extension (DGGRID v8) + * once per database instance. Mirrors {@link ensureH3Extension}. + */ +export const ensureDuckDggsExtension = createCommunityExtensionLoader("duck_dggs"); async function createDatabase(): Promise { const bundle = await selectDuckDbBundle(); diff --git a/package-lock.json b/package-lock.json index d9f97f0ba..5e5019724 100644 --- a/package-lock.json +++ b/package-lock.json @@ -22601,10 +22601,12 @@ "@turf/tin": "^7.3.5", "@turf/union": "^7.3.5", "@turf/voronoi": "^7.3.5", + "dggal": "^0.0.6", "fflate": "^0.8.3", "geolibre-wasm": "^1.4.2", "geotiff": "^3.0.5", - "onnxruntime-web": "1.27.0" + "onnxruntime-web": "1.27.0", + "s2js": "^1.44.0" }, "devDependencies": { "@types/geojson": "^7946.0.16", diff --git a/packages/core/src/project.ts b/packages/core/src/project.ts index 77c6f94e6..bbccfd75a 100644 --- a/packages/core/src/project.ts +++ b/packages/core/src/project.ts @@ -552,6 +552,15 @@ const PROCESSING_RUN_KINDS = new Set([ "algorithm", ]); +/** + * Old H3 vector-tool IDs from projects saved before the DGGS rename. + * Mapped onto current tool ids during project load. + */ +const LEGACY_H3_PROCESSING_TOOL_IDS: Readonly> = { + "h3-grid": "dggs-grid", + "h3-bin-points": "dggs-bin", +}; + /** * Coerce an untrusted (possibly hand-edited) `processingHistory` array into * valid {@link ProcessingRun} records. Drops entries without a usable id, tool @@ -576,7 +585,7 @@ export function normalizeProcessingHistory(value: unknown): ProcessingRun[] | nu if (!entry || typeof entry !== "object") continue; const candidate = entry as Partial; const id = normalizeString(candidate.id).trim(); - const toolId = normalizeString(candidate.toolId).trim(); + let toolId = normalizeString(candidate.toolId).trim(); const kind = candidate.kind; if (!id || !toolId || seen.has(id)) continue; if (!kind || !PROCESSING_RUN_KINDS.has(kind)) continue; @@ -592,16 +601,22 @@ export function normalizeProcessingHistory(value: unknown): ProcessingRun[] | nu const outputLayerNames = Array.isArray(candidate.outputLayerNames) ? candidate.outputLayerNames.filter((name): name is string => typeof name === "string") : undefined; + let parameters: Record = + candidate.parameters && typeof candidate.parameters === "object" + ? { ...(candidate.parameters as Record) } + : {}; + const migrated = LEGACY_H3_PROCESSING_TOOL_IDS[toolId]; + if (migrated) { + toolId = migrated; + if (parameters.dggsType == null) parameters = { ...parameters, dggsType: "h3" }; + } runs.push({ id, kind, toolId, toolName: normalizeString(candidate.toolName) || toolId, engine: normalizeString(candidate.engine), - parameters: - candidate.parameters && typeof candidate.parameters === "object" - ? (candidate.parameters as Record) - : {}, + parameters, ...(inputLayerNames && Object.keys(inputLayerNames).length > 0 ? { inputLayerNames } : {}), ...(outputLayerNames?.length ? { outputLayerNames } : {}), ...(normalizeString(candidate.inputPath) diff --git a/packages/core/src/store.ts b/packages/core/src/store.ts index 6bc52f1ca..4d843cfed 100644 --- a/packages/core/src/store.ts +++ b/packages/core/src/store.ts @@ -123,8 +123,9 @@ export type VectorToolKind = | "grid" | "voronoi" | "cell-sectors" - | "h3-grid" - | "h3-bin-points" + | "dggs-grid" + | "dggs-bin" + | "dggs-compact" | "trajectory-speed" | "detect-stops" | "space-time-proximity" @@ -1265,7 +1266,19 @@ export const useAppStore = create()( setProcessingInitialTool: (toolId) => set((s) => ({ ui: { ...s.ui, processingInitialTool: toolId } })), setConversionOpen: (kind) => set((s) => ({ ui: { ...s.ui, conversionOpen: kind } })), - setVectorToolOpen: (kind) => set((s) => ({ ui: { ...s.ui, vectorToolOpen: kind } })), + setVectorToolOpen: (kind) => + set((s) => ({ + ui: { + ...s.ui, + // Pre-DGGS projects / callers may still pass h3-grid / h3-bin-points. + vectorToolOpen: + (kind as string | null) === "h3-grid" + ? "dggs-grid" + : (kind as string | null) === "h3-bin-points" + ? "dggs-bin" + : kind, + }, + })), setNetworkToolOpen: (kind) => set((s) => ({ ui: { ...s.ui, networkToolOpen: kind } })), setStatisticsToolOpen: (kind) => set((s) => ({ ui: { ...s.ui, statisticsToolOpen: kind } })), setRasterToolOpen: (kind) => set((s) => ({ ui: { ...s.ui, rasterToolOpen: kind } })), diff --git a/packages/plugins/src/plugins/maplibre-geohash.ts b/packages/plugins/src/plugins/maplibre-geohash.ts index ce2f09fd7..15a2a7562 100644 --- a/packages/plugins/src/plugins/maplibre-geohash.ts +++ b/packages/plugins/src/plugins/maplibre-geohash.ts @@ -89,10 +89,10 @@ export const DEFAULT_GEOHASH_LABELS: GeohashLabels = { title: "Geohash", controlTitle: "Geohash settings", autoResolution: "Automatic resolution", - resolution: "Precision", + resolution: "Resolution", cellCount: (count) => `${count.toLocaleString()} cells in view`, tooManyCells: (limit) => - `This view exceeds the ${limit.toLocaleString()} cell limit. Zoom in or lower the precision.`, + `This view exceeds the ${limit.toLocaleString()} cell limit. Zoom in or lower the resolution.`, fillColor: "Fill color", fillOpacity: "Fill opacity", lineColor: "Outline color", @@ -739,7 +739,7 @@ function renderPanel(container: HTMLElement): void { labels.addAsLayer, () => { if (currentGrid.features.length) { - appRef?.addGeoJsonLayer(`Geohash (precision ${effectiveResolution()})`, currentGrid); + appRef?.addGeoJsonLayer(`Geohash (res ${effectiveResolution()})`, currentGrid); } }, currentGrid.features.length === 0, diff --git a/packages/plugins/src/plugins/maplibre-olc.ts b/packages/plugins/src/plugins/maplibre-olc.ts index ad549db3b..7c530bf8b 100644 --- a/packages/plugins/src/plugins/maplibre-olc.ts +++ b/packages/plugins/src/plugins/maplibre-olc.ts @@ -98,10 +98,10 @@ export const DEFAULT_OLC_LABELS: OlcLabels = { title: "OLC", controlTitle: "OLC settings", autoResolution: "Automatic resolution", - resolution: "Code length", + resolution: "Resolution", cellCount: (count) => `${count.toLocaleString()} cells in view`, tooManyCells: (limit) => - `This view exceeds the ${limit.toLocaleString()} cell limit. Zoom in or lower the code length.`, + `This view exceeds the ${limit.toLocaleString()} cell limit. Zoom in or lower the resolution.`, fillColor: "Fill color", fillOpacity: "Fill opacity", lineColor: "Outline color", @@ -785,7 +785,7 @@ function renderPanel(container: HTMLElement): void { labels.addAsLayer, () => { if (currentGrid.features.length) { - appRef?.addGeoJsonLayer(`OLC (code length ${effectiveResolution()})`, currentGrid); + appRef?.addGeoJsonLayer(`OLC (res ${effectiveResolution()})`, currentGrid); } }, currentGrid.features.length === 0, diff --git a/packages/plugins/src/plugins/maplibre-tilecode.ts b/packages/plugins/src/plugins/maplibre-tilecode.ts index 3f0aa6d66..62803f66f 100644 --- a/packages/plugins/src/plugins/maplibre-tilecode.ts +++ b/packages/plugins/src/plugins/maplibre-tilecode.ts @@ -95,10 +95,10 @@ export const DEFAULT_TILECODE_LABELS: TilecodeLabels = { title: "Tilecode", controlTitle: "Tilecode settings", autoResolution: "Automatic resolution", - resolution: "Zoom level", + resolution: "Resolution", cellCount: (count) => `${count.toLocaleString()} tiles in view`, tooManyCells: (limit) => - `This view exceeds the ${limit.toLocaleString()} tile limit. Zoom in or lower the zoom level.`, + `This view exceeds the ${limit.toLocaleString()} tile limit. Zoom in or lower the resolution.`, fillColor: "Fill color", fillOpacity: "Fill opacity", lineColor: "Outline color", @@ -810,7 +810,7 @@ function renderPanel(container: HTMLElement): void { labels.addAsLayer, () => { if (currentGrid.features.length) { - appRef?.addGeoJsonLayer(`Tilecode (zoom ${effectiveResolution()})`, currentGrid); + appRef?.addGeoJsonLayer(`Tilecode (res ${effectiveResolution()})`, currentGrid); } }, currentGrid.features.length === 0, diff --git a/packages/processing/package.json b/packages/processing/package.json index f363c805a..8d44cb4de 100644 --- a/packages/processing/package.json +++ b/packages/processing/package.json @@ -32,10 +32,12 @@ "@turf/tin": "^7.3.5", "@turf/union": "^7.3.5", "@turf/voronoi": "^7.3.5", + "dggal": "^0.0.6", "fflate": "^0.8.3", "geolibre-wasm": "^1.4.2", "geotiff": "^3.0.5", - "onnxruntime-web": "1.27.0" + "onnxruntime-web": "1.27.0", + "s2js": "^1.44.0" }, "devDependencies": { "@types/geojson": "^7946.0.16", diff --git a/packages/processing/src/a5-tools.ts b/packages/processing/src/a5-tools.ts new file mode 100644 index 000000000..27b51f995 --- /dev/null +++ b/packages/processing/src/a5-tools.ts @@ -0,0 +1,260 @@ +import type { FeatureCollection, Geometry } from "geojson"; +import { bboxToWktPolygon, normalizeLonLatBbox, sqlIdent, sqlStr } from "./h3-tools"; + +/** + * Approximate average A5 cell area (km²) at resolutions 0..30. + * A5 cells are equal-area; counts grow by 4× per level from 12 root cells + * (earth surface ≈ 5.101×10⁸ km²). Used only for auto-suggest / hard-cap + * estimates before the DuckDB query runs. + */ +export const A5_AVG_AREA_KM2: number[] = Array.from({ length: 31 }, (_, res) => { + return 510_065_621.724 / (12 * 4 ** res); +}); + +/** Soft target used when auto-suggesting a resolution. */ +export const A5_TARGET_CELLS = 10_000; +/** Finest resolution the auto-suggester will pick. */ +export const A5_MAX_SUGGESTED_RES = 12; +/** Hard ceiling: a grid larger than this aborts rather than running away. */ +export const A5_HARD_CAP = 200_000; +/** + * Max resolution offered in the processing dialog. Matches A5's + * `MAX_RESOLUTION` (0–30); level 31 does not exist in the encoding. + */ +export const A5_MAX_TOOL_RES = 30; +/** + * Max longitude span (degrees) per `a5_geometry_to_cells` call. Wider rings + * (especially ±180) return empty or dateline-only cells in DuckDB A5. + */ +export const A5_MAX_POLYFILL_LON_SPAN = 90; + +/** Estimated number of A5 cells covering `areaKm2` at `res`. */ +export function estimateA5CellCount(areaKm2: number, res: number): number { + const cellArea = A5_AVG_AREA_KM2[res]; + if (cellArea === undefined) return Number.POSITIVE_INFINITY; + return areaKm2 / cellArea; +} + +/** Finest resolution whose estimated cell count stays <= the target. */ +export function suggestA5Resolution( + areaKm2: number, + targetCells = A5_TARGET_CELLS, + maxRes = A5_MAX_SUGGESTED_RES, +): number { + const capped = Math.min(maxRes, A5_MAX_TOOL_RES); + for (let res = capped; res >= 0; res -= 1) { + if (estimateA5CellCount(areaKm2, res) <= targetCells) return res; + } + return 0; +} + +const GRID_SELECT = + "SELECT a5_u64_to_hex(cell) AS a5, " + + "CAST(ST_AsGeoJSON(a5_cell_to_geometry(cell)) AS VARCHAR) AS geojson FROM cells"; + +/** + * Expand a compacted covering from `a5_geometry_to_cells` to a uniform + * resolution. Without this, fully covered parents stay at coarser levels and + * the output mixes cell sizes. + */ +function cellsFromGeomExpr(geomSql: string, res: number): string { + return `unnest(a5_uncompact(a5_geometry_to_cells(${geomSql}, ${res}), ${res}))`; +} + +/** + * All cells at `res` from the 12 resolution-0 roots. Prefer this for full + * longitude (±180): `a5_geometry_to_cells` returns an empty list for a single + * world ring, and longitude strips also under-cover (gaps at poles / seams). + */ +function cellsFromRes0Expr(res: number): string { + return `unnest(a5_uncompact(a5_get_res0_cells(), ${res}))`; +} + +/** Wrap a `SELECT … AS cell` query; optionally `a5_compact` the result. */ +function finalizeA5Cells(rawSelect: string, compact: boolean): string { + if (!compact) { + return `WITH cells AS (${rawSelect}) ` + GRID_SELECT; + } + return ( + `WITH raw AS (${rawSelect}), ` + + `arr AS (SELECT list(cell) AS cells FROM raw), ` + + `cells AS (SELECT unnest(a5_compact(cells)) AS cell FROM arr) ` + + GRID_SELECT + ); +} + +/** Grid SQL from a polygon WKT literal (used for bbox / viewport sources). */ +export function buildA5GridFromWktSql(wkt: string, res: number, compact = false): string { + return finalizeA5Cells( + `SELECT ${cellsFromGeomExpr(`ST_GeomFromText(${sqlStr(wkt)})`, res)} AS cell`, + compact, + ); +} + +/** + * Grid SQL from a lon/lat bbox. + * + * - Full longitude (`[-180, 180]` after normalize): enumerate via + * {@link cellsFromRes0Expr}, optionally filtering by latitude. + * - Spans wider than {@link A5_MAX_POLYFILL_LON_SPAN}: slice into strips + * (a single wide ring returns empty or incomplete cells from DuckDB A5). + * - Narrower spans: one WKT polyfill. + */ +export function buildA5GridFromBboxSql( + bbox: [number, number, number, number], + res: number, + compact = false, +): string { + const [w, s, e, n] = normalizeLonLatBbox(bbox); + if (w === -180 && e === 180) { + // Full globe in longitude. Filter by cell centroid latitude when the view + // is not essentially ±90 (e.g. Web Mercator max ~±85). + const fullLat = s <= -89.999 && n >= 89.999; + if (fullLat) { + return finalizeA5Cells(`SELECT ${cellsFromRes0Expr(res)} AS cell`, compact); + } + return finalizeA5Cells( + `SELECT cell FROM (SELECT ${cellsFromRes0Expr(res)} AS cell) ` + + `WHERE list_extract(a5_cell_to_lonlat(cell), 2) BETWEEN ${s} AND ${n}`, + compact, + ); + } + const lonSpan = e - w; + if (lonSpan <= A5_MAX_POLYFILL_LON_SPAN) { + return buildA5GridFromWktSql(bboxToWktPolygon([w, s, e, n]), res, compact); + } + const parts = Math.ceil(lonSpan / A5_MAX_POLYFILL_LON_SPAN); + const step = lonSpan / parts; + const selects: string[] = []; + for (let i = 0; i < parts; i += 1) { + const left = w + i * step; + const right = w + (i + 1) * step; + const wkt = bboxToWktPolygon([left, s, right, n]); + selects.push(`SELECT ${cellsFromGeomExpr(`ST_GeomFromText(${sqlStr(wkt)})`, res)} AS cell`); + } + return finalizeA5Cells(`SELECT DISTINCT cell FROM (${selects.join(" UNION ALL ")})`, compact); +} + +/** + * Grid SQL that unions all geometry from a registered source into one + * (multi)polygon and fills it (used for the polyfill source). `sourceSql` is a + * FROM-able expression whose geometry column is `geom` (DuckDB `ST_Read`). + */ +export function buildA5GridFromSourceSql(sourceSql: string, res: number, compact = false): string { + // Union only polygonal geometries: a mixed layer would otherwise aggregate to + // a GEOMETRYCOLLECTION that a5_geometry_to_cells rejects. The outer select + // filters a NULL union result so a NULL geometry never reaches the a5 function. + return finalizeA5Cells( + `SELECT ${cellsFromGeomExpr("g", res)} AS cell FROM (` + + `SELECT ST_Union_Agg(geom) AS g FROM ${sourceSql} ` + + `WHERE geom IS NOT NULL AND ST_GeometryType(geom) IN ('POLYGON', 'MULTIPOLYGON')` + + `) merged WHERE g IS NOT NULL`, + compact, + ); +} + +/** Supported point-binning aggregate operations (same set as H3 binning). */ +export type A5AggOp = "count" | "sum" | "mean" | "min" | "max"; + +const AGG_FN: Record, string> = { + sum: "sum", + mean: "avg", + min: "min", + max: "max", +}; + +/** + * Aggregate point geometry from `sourceSql` into A5 cells. Mirrors + * {@link buildBinSql} but uses `a5_lonlat_to_cell` (lon, lat order) and + * `a5_cell_to_geometry` for boundaries. + */ +export function buildA5BinSql(sourceSql: string, res: number, op: A5AggOp, field?: string): string { + const fn = op === "count" ? undefined : AGG_FN[op]; + const aggSelect = fn && field ? `, ${fn}(CAST(${sqlIdent(field)} AS DOUBLE)) AS value` : ""; + const aggOut = fn && field ? ", value" : ""; + return ( + `WITH pts AS (SELECT ST_Centroid(geom) AS pt` + + (field ? `, ${sqlIdent(field)}` : "") + + ` FROM ${sourceSql} ` + + `WHERE geom IS NOT NULL AND ST_GeometryType(geom) IN ('POINT', 'MULTIPOINT')), ` + + `binned AS (SELECT a5_lonlat_to_cell(ST_X(pt), ST_Y(pt), ${res}) AS cell, ` + + `count(*) AS count${aggSelect} FROM pts GROUP BY cell) ` + + `SELECT a5_u64_to_hex(cell) AS a5, count${aggOut}, ` + + `CAST(ST_AsGeoJSON(a5_cell_to_geometry(cell)) AS VARCHAR) AS geojson FROM binned` + ); +} + +/** + * Collect A5 cell IDs (hex strings) from `cellField` on `sourceSql` into an array. + */ +function a5CellArrayCte(sourceSql: string, cellField: string): string { + const f = sqlIdent(cellField); + return ( + `input AS (SELECT DISTINCT a5_hex_to_u64(CAST(${f} AS VARCHAR)) AS cell FROM ${sourceSql} ` + + `WHERE ${f} IS NOT NULL AND CAST(${f} AS VARCHAR) <> ''), ` + + `arr AS (SELECT list(cell) AS cells FROM input)` + ); +} + +/** Compact A5 cells from a polygon cell layer (IDs in `cellField`, default `a5`). */ +export function buildA5CompactSql(sourceSql: string, cellField = "a5"): string { + return ( + `WITH ${a5CellArrayCte(sourceSql, cellField)}, ` + + `cells AS (SELECT unnest(a5_compact(cells)) AS cell FROM arr) ` + + GRID_SELECT + ); +} + +/** Expand (uncompact) A5 cells to a uniform `res`. */ +export function buildA5ExpandSql(sourceSql: string, res: number, cellField = "a5"): string { + return ( + `WITH ${a5CellArrayCte(sourceSql, cellField)}, ` + + `cells AS (SELECT unnest(a5_uncompact(cells, ${res})) AS cell FROM arr) ` + + GRID_SELECT + ); +} + +/** Count of cells that {@link buildA5ExpandSql} would emit (for the hard-cap guard). */ +export function buildA5ExpandCountSql(sourceSql: string, res: number, cellField = "a5"): string { + return ( + `WITH ${a5CellArrayCte(sourceSql, cellField)} ` + + `SELECT coalesce(len(a5_uncompact(cells, ${res})), 0) AS n FROM arr` + ); +} + +/** Parse a DuckDB `ST_AsGeoJSON` cell (VARCHAR or already-decoded JSON object). */ +function geometryFromGeoJsonCell(raw: unknown): Geometry | null { + if (typeof raw === "string") { + try { + return JSON.parse(raw) as Geometry; + } catch { + return null; + } + } + if (raw && typeof raw === "object" && "type" in raw) { + return raw as Geometry; + } + return null; +} + +/** Build a FeatureCollection from rows carrying `a5`, optional `count`/`value`, and `geojson`. */ +export function a5RowsToFeatureCollection(rows: Record[]): FeatureCollection { + const features = []; + for (const row of rows) { + const geometry = geometryFromGeoJsonCell(row.geojson); + if (!geometry) continue; + const properties: Record = { a5: String(row.a5) }; + if (row.count !== undefined && row.count !== null) { + properties.count = Number(row.count); + } + if (row.value !== undefined && row.value !== null) { + properties.value = Number(row.value); + } + features.push({ + type: "Feature" as const, + geometry, + properties, + }); + } + return { type: "FeatureCollection", features }; +} diff --git a/packages/processing/src/antimeridian.ts b/packages/processing/src/antimeridian.ts new file mode 100644 index 000000000..6051de37f --- /dev/null +++ b/packages/processing/src/antimeridian.ts @@ -0,0 +1,46 @@ +import type { Geometry } from "geojson"; + +/** + * Unwrap successive ring vertices so consecutive longitudes stay within 180°. + * DuckDB H3 / duck_dggs emit longitudes in [-180, 180]; dateline-straddling + * cells then draw the long way around MapLibre unless vertices are shifted + * into an adjacent world copy. + */ +export function unwrapAntimeridianRing(ring: number[][]): number[][] { + if (ring.length === 0) return ring; + const firstLon = ring[0]![0]!; + const firstLat = ring[0]![1]!; + const last = ring[ring.length - 1]!; + const closed = ring.length > 1 && last[0] === firstLon && last[1] === firstLat; + const limit = closed ? ring.length - 1 : ring.length; + + // Preserve elevation / M and any further components; only longitude shifts. + const out: number[][] = [[...ring[0]!]]; + for (let i = 1; i < limit; i += 1) { + let lon = ring[i]![0]!; + const rest = ring[i]!.slice(1); + const prev = out[i - 1]![0]!; + while (lon - prev > 180) lon -= 360; + while (lon - prev < -180) lon += 360; + out.push([lon, ...rest]); + } + if (closed) out.push([...out[0]!]); + return out; +} + +/** Unwrap Polygon / MultiPolygon rings across ±180°; other geometry types pass through. */ +export function unwrapAntimeridianGeometry(geometry: Geometry): Geometry { + if (geometry.type === "Polygon") { + return { + type: "Polygon", + coordinates: geometry.coordinates.map(unwrapAntimeridianRing), + }; + } + if (geometry.type === "MultiPolygon") { + return { + type: "MultiPolygon", + coordinates: geometry.coordinates.map((poly) => poly.map(unwrapAntimeridianRing)), + }; + } + return geometry; +} diff --git a/packages/processing/src/dggal-tools.ts b/packages/processing/src/dggal-tools.ts new file mode 100644 index 000000000..808214a7b --- /dev/null +++ b/packages/processing/src/dggal-tools.ts @@ -0,0 +1,538 @@ +import type { Feature, FeatureCollection, Geometry, Polygon, Position } from "geojson"; +import booleanPointInPolygon from "@turf/boolean-point-in-polygon"; +import { point as turfPoint } from "@turf/helpers"; + +/** + * Named DGGAL DGGRS types for DGGS Generator / Binning. Keys are stable UI + * values; `className` is passed to `dggal.createDGGRS(...)`. + */ +export const DGGAL_TYPES = { + gnosis: { + minRes: 0, + maxRes: 28, + defaultRes: 16, + className: "GNOSISGlobalGrid", + }, + isea4r: { minRes: 0, maxRes: 25, defaultRes: 12, className: "ISEA4R" }, + isea9r: { minRes: 0, maxRes: 16, defaultRes: 10, className: "ISEA9R" }, + isea3h: { minRes: 0, maxRes: 33, defaultRes: 21, className: "ISEA3H" }, + isea7h: { minRes: 0, maxRes: 19, defaultRes: 11, className: "ISEA7H" }, + isea7h_z7: { minRes: 0, maxRes: 19, defaultRes: 11, className: "ISEA7H_Z7" }, + ivea4r: { minRes: 0, maxRes: 25, defaultRes: 12, className: "IVEA4R" }, + ivea9r: { minRes: 0, maxRes: 16, defaultRes: 10, className: "IVEA9R" }, + ivea3h: { minRes: 0, maxRes: 33, defaultRes: 21, className: "IVEA3H" }, + ivea7h: { minRes: 0, maxRes: 19, defaultRes: 11, className: "IVEA7H" }, + ivea7h_z7: { minRes: 0, maxRes: 19, defaultRes: 11, className: "IVEA7H_Z7" }, + rtea4r: { minRes: 0, maxRes: 25, defaultRes: 12, className: "RTEA4R" }, + rtea9r: { minRes: 0, maxRes: 16, defaultRes: 10, className: "RTEA9R" }, + rtea3h: { minRes: 0, maxRes: 33, defaultRes: 21, className: "RTEA3H" }, + rtea7h: { minRes: 0, maxRes: 19, defaultRes: 11, className: "RTEA7H" }, + rtea7h_z7: { minRes: 0, maxRes: 19, defaultRes: 11, className: "RTEA7H_Z7" }, + healpix: { minRes: 0, maxRes: 26, defaultRes: 18, className: "HEALPix" }, + rhealpix: { minRes: 0, maxRes: 16, defaultRes: 10, className: "rHEALPix" }, +} as const; + +export type DggalGridType = keyof typeof DGGAL_TYPES; + +export const DGGAL_GRID_TYPES = Object.keys(DGGAL_TYPES) as DggalGridType[]; + +export const DEFAULT_DGGAL_GRID_TYPE: DggalGridType = "isea3h"; + +export type DggalGridSpec = (typeof DGGAL_TYPES)[DggalGridType]; + +/** Labels shown in the processing dialog (engine class names). */ +export const DGGAL_GRID_TYPE_OPTIONS: { value: DggalGridType; label: string }[] = + DGGAL_GRID_TYPES.map((value) => ({ value, label: DGGAL_TYPES[value].className })); + +export function resolveDggalGridType(raw: unknown): DggalGridType { + if (typeof raw === "string" && Object.hasOwn(DGGAL_TYPES, raw)) { + return raw as DggalGridType; + } + return DEFAULT_DGGAL_GRID_TYPE; +} + +export function maxResolutionForDggal(gridType: DggalGridType = DEFAULT_DGGAL_GRID_TYPE): number { + return DGGAL_TYPES[gridType].maxRes; +} + +/** Soft target used when auto-suggesting a resolution. */ +export const DGGAL_TARGET_CELLS = 10_000; +/** Finest resolution the auto-suggester will pick. */ +export const DGGAL_MAX_SUGGESTED_RES = 12; +/** Hard ceiling: a grid larger than this aborts rather than running away. */ +export const DGGAL_HARD_CAP = 200_000; +/** Absolute finest resolution across exposed DGGAL types (ISEA3H / IVEA3H / RTEA3H). */ +export const DGGAL_MAX_TOOL_RES = Math.max(...DGGAL_GRID_TYPES.map((t) => DGGAL_TYPES[t].maxRes)); + +const EARTH_AREA_KM2 = 510_065_621.724; +const DEG_PER_RAD = 180 / Math.PI; +const RAD_PER_DEG = Math.PI / 180; + +/** Approximate global zone count at `res` (for suggest / pre-cap without WASM). */ +export function dggalApproxGlobalCount(res: number, gridType: DggalGridType): number { + const name = DGGAL_TYPES[gridType].className; + if (name === "HEALPix" || name === "rHEALPix") return 12 * 4 ** res; + if (name.endsWith("3H")) return 10 * 3 ** res + 2; + if (name.endsWith("7H") || name.endsWith("7H_Z7")) return 10 * 7 ** res + 2; + if (name.endsWith("9R")) return 10 * 9 ** res + 2; + if (name.endsWith("4R")) return 10 * 4 ** res + 2; + // GNOSISGlobalGrid and fallbacks: treat like aperture-4. + return 10 * 4 ** res + 2; +} + +export function estimateDggalCellCount( + areaKm2: number, + res: number, + gridType: DggalGridType = DEFAULT_DGGAL_GRID_TYPE, +): number { + const global = dggalApproxGlobalCount(res, gridType); + if (!Number.isFinite(global) || global <= 0) return Number.POSITIVE_INFINITY; + return (areaKm2 / EARTH_AREA_KM2) * global; +} + +export function suggestDggalResolution( + areaKm2: number, + targetCells = DGGAL_TARGET_CELLS, + maxRes = DGGAL_MAX_SUGGESTED_RES, + gridType: DggalGridType = DEFAULT_DGGAL_GRID_TYPE, +): number { + const capped = Math.min(maxRes, DGGAL_TYPES[gridType].maxRes); + for (let res = capped; res >= 0; res -= 1) { + if (estimateDggalCellCount(areaKm2, res, gridType) <= targetCells) return res; + } + return 0; +} + +/** Geographic point in radians (DGGAL native). */ +interface GeoPoint { + lat: number; + lon: number; +} + +/** Subset of a DGGAL DGGRS instance used by the processing tools. */ +export interface DggalDggrs { + getZoneFromTextID(zoneId: string): bigint; + getZoneTextID(zone: bigint): string; + getZoneLevel(zone: bigint): number; + getZoneWGS84Centroid(zone: bigint): GeoPoint; + getZoneRefinedWGS84Vertices(zone: bigint, edgeRefinement: number): GeoPoint[]; + listZones(level: number, bbox: { ll: GeoPoint; ur: GeoPoint }): bigint[]; + getZoneFromWGS84Centroid(level: number, geoPoint: GeoPoint): bigint; + countZones(level: number): bigint; + /** Recursively replace full child sets with parents (mutates conceptually; returns new list). */ + compactZones(zones: bigint[]): bigint[]; + /** Sub-zones of `zone` at relative `depth` (1 = immediate children). */ + getSubZones(zone: bigint, depth: number): bigint[]; + /** Number of sub-zones at relative `depth` (bigint from WASM). */ + countSubZones(zone: bigint, depth: number): bigint | number; + delete(): void; +} + +export interface DggalEngine { + createDGGRS(name: string): DggalDggrs; +} + +let dggalPromise: Promise | null = null; + +/** Load the DGGAL WASM module once (dynamic import keeps it out of cold paths). */ +export function loadDggal(): Promise { + dggalPromise ??= import("dggal") + .then((module) => module.DGGAL.init() as unknown as Promise) + .catch((error) => { + dggalPromise = null; + throw error; + }); + return dggalPromise; +} + +/** Run `fn` with a short-lived DGGRS for `gridType`, always deleting the instance. */ +export async function withDggalDggrs( + gridType: DggalGridType, + fn: (engine: DggalDggrs) => T | Promise, +): Promise { + const dggal = await loadDggal(); + const engine = dggal.createDGGRS(DGGAL_TYPES[gridType].className); + try { + return await fn(engine); + } finally { + engine.delete(); + } +} + +function normalizeLon(lon: number): number { + let x = lon; + while (x > 180) x -= 360; + while (x < -180) x += 360; + return x; +} + +function zoneRing(engine: DggalDggrs, zone: bigint): number[][] { + const ring = engine + .getZoneRefinedWGS84Vertices(zone, 0) + .map(({ lat, lon }): number[] => [lon * DEG_PER_RAD, lat * DEG_PER_RAD]); + if (ring.length > 0) { + const [firstLng, firstLat] = ring[0]!; + const [lastLng, lastLat] = ring[ring.length - 1]!; + if (firstLng !== lastLng || firstLat !== lastLat) ring.push([firstLng!, firstLat!]); + } + return ring; +} + +/** Convert a DGGAL zone text ID to a GeoJSON polygon feature. */ +export function dggalZoneFeature(engine: DggalDggrs, cell: string): Feature { + const zone = engine.getZoneFromTextID(cell); + const centroid = engine.getZoneWGS84Centroid(zone); + return { + type: "Feature", + id: cell, + properties: { + dggal: cell, + resolution: engine.getZoneLevel(zone), + center_lat: centroid.lat * DEG_PER_RAD, + center_lng: centroid.lon * DEG_PER_RAD, + }, + geometry: { type: "Polygon", coordinates: [zoneRing(engine, zone)] }, + }; +} + +/** Build polygon features from DGGAL zone text IDs. */ +export function dggalTokensToFeatureCollection( + engine: DggalDggrs, + tokens: Iterable, +): FeatureCollection { + return { + type: "FeatureCollection", + features: [...tokens].map((token) => dggalZoneFeature(engine, token)), + }; +} + +/** Drop null / unreadable zone handles from WASM array paddings. */ +function validZones(engine: DggalDggrs, zones: Iterable): bigint[] { + const out: bigint[] = []; + for (const zone of zones) { + if (zone === 0n) continue; + try { + engine.getZoneLevel(zone); + out.push(zone); + } catch { + // Padding / invalid handle. + } + } + return out; +} + +/** + * Fill a WGS84 bounding box with DGGAL zones via `listZones` (same approach as + * the maplibre-dggal plugin). + */ +export function dggalGridFromBbox( + engine: DggalDggrs, + bounds: [number, number, number, number], + resolution: number, + limit = DGGAL_HARD_CAP, + options: { compact?: boolean } = {}, +): FeatureCollection { + let [west, south, east, north] = bounds; + south = Math.max(-90, Math.min(90, south)); + north = Math.max(-90, Math.min(90, north)); + if (south > north) [south, north] = [north, south]; + if (east - west >= 360) { + west = -180; + east = 180; + } else { + west = normalizeLon(west); + east = normalizeLon(east); + } + if (east < west) { + // Cover each side without compacting, then compact once over the union so + // sibling sets that straddle the antimeridian can still merge. + const left = dggalGridFromBbox(engine, [west, south, 180, north], resolution, limit); + const right = dggalGridFromBbox(engine, [-180, south, east, north], resolution, limit); + const seen = new Set(); + const features: Feature[] = []; + for (const feature of [...left.features, ...right.features]) { + const id = String(feature.properties?.dggal ?? feature.id); + if (seen.has(id)) continue; + seen.add(id); + features.push(feature); + if (features.length > limit) { + throw new RangeError(`DGGAL zone limit exceeded: ${limit}`); + } + } + if (options.compact) { + const tokens = compactDggalTokens( + engine, + features.map((f) => String(f.properties?.dggal ?? f.id)), + ); + return dggalTokensToFeatureCollection(engine, tokens); + } + return { type: "FeatureCollection", features }; + } + + let zones = validZones( + engine, + engine.listZones(resolution, { + ll: { lat: south * RAD_PER_DEG, lon: west * RAD_PER_DEG }, + ur: { lat: north * RAD_PER_DEG, lon: east * RAD_PER_DEG }, + }), + ); + if (zones.length > limit) { + throw new RangeError(`DGGAL zone limit exceeded: ${limit}`); + } + if (options.compact) { + zones = validZones(engine, engine.compactZones(zones)); + } + return { + type: "FeatureCollection", + features: zones.map((zone) => dggalZoneFeature(engine, engine.getZoneTextID(zone))), + }; +} + +/** + * Polyfill polygon geometry: list zones over each feature's bbox, keep those + * whose centroid falls inside the polygon (DGGAL has no native polyfill). + */ +export function dggalGridFromFeatureCollection( + engine: DggalDggrs, + fc: FeatureCollection, + resolution: number, + limit = DGGAL_HARD_CAP, + options: { compact?: boolean } = {}, +): FeatureCollection { + const seen = new Set(); + const features: Feature[] = []; + + const consider = (poly: Polygon) => { + let minX = Infinity; + let minY = Infinity; + let maxX = -Infinity; + let maxY = -Infinity; + for (const ring of poly.coordinates) { + for (const pos of ring) { + const x = pos[0]; + const y = pos[1]; + if (x === undefined || y === undefined) continue; + if (x < minX) minX = x; + if (y < minY) minY = y; + if (x > maxX) maxX = x; + if (y > maxY) maxY = y; + } + } + if (!Number.isFinite(minX) || minX >= maxX || minY >= maxY) return; + const candidates = dggalGridFromBbox(engine, [minX, minY, maxX, maxY], resolution, limit); + for (const feature of candidates.features) { + const id = String(feature.properties?.dggal ?? feature.id); + if (seen.has(id)) continue; + const lng = Number(feature.properties?.center_lng); + const lat = Number(feature.properties?.center_lat); + if (!Number.isFinite(lng) || !Number.isFinite(lat)) continue; + if (!booleanPointInPolygon(turfPoint([lng, lat]), poly)) continue; + seen.add(id); + features.push(feature); + if (features.length > limit) { + throw new RangeError(`DGGAL zone limit exceeded: ${limit}`); + } + } + }; + + for (const feature of fc.features) { + const g = feature.geometry; + if (!g) continue; + if (g.type === "Polygon") consider(g); + else if (g.type === "MultiPolygon") { + for (const coords of g.coordinates) { + consider({ type: "Polygon", coordinates: coords }); + } + } + } + if (options.compact) { + const tokens = compactDggalTokens( + engine, + features.map((f) => String(f.properties?.dggal ?? f.id)), + ); + return dggalTokensToFeatureCollection(engine, tokens); + } + return { type: "FeatureCollection", features }; +} + +/** Collect DGGAL zone text IDs from a feature property (default `dggal`). */ +export function tokensFromDggalLayer(fc: FeatureCollection, cellField = "dggal"): string[] { + const out: string[] = []; + for (const feature of fc.features) { + const raw = feature.properties?.[cellField]; + if (raw === undefined || raw === null) continue; + const token = String(raw).trim(); + if (token) out.push(token); + } + return out; +} + +/** + * Compact zone text IDs with {@link DggalDggrs.compactZones}: complete child + * sets become their parents (mixed resolutions; behaviour is DGGRS-specific). + */ +export function compactDggalTokens(engine: DggalDggrs, tokens: Iterable): string[] { + const zones: bigint[] = []; + for (const token of tokens) { + zones.push(engine.getZoneFromTextID(token)); + } + if (zones.length === 0) return []; + return validZones(engine, engine.compactZones(zones)).map((z) => engine.getZoneTextID(z)); +} + +/** + * Expand zone text IDs to a uniform `level` via {@link DggalDggrs.getSubZones}. + * Cells finer than `level` throw. + */ +export function expandDggalTokens( + engine: DggalDggrs, + tokens: Iterable, + level: number, +): string[] { + const seen = new Set(); + const out: string[] = []; + for (const token of tokens) { + const zone = engine.getZoneFromTextID(token); + const L = engine.getZoneLevel(zone); + if (L > level) { + throw new RangeError( + `DGGAL cell ${token} is finer than target level ${level}; choose a finer target or compact first`, + ); + } + if (L === level) { + const id = engine.getZoneTextID(zone); + if (!seen.has(id)) { + seen.add(id); + out.push(id); + } + continue; + } + for (const sub of validZones(engine, engine.getSubZones(zone, level - L))) { + const id = engine.getZoneTextID(sub); + if (seen.has(id)) continue; + seen.add(id); + out.push(id); + } + } + return out; +} + +/** Exact cell count after expanding `tokens` to `level`. */ +export function estimateDggalExpandCount( + engine: DggalDggrs, + tokens: Iterable, + level: number, +): number { + let n = 0; + for (const token of tokens) { + const zone = engine.getZoneFromTextID(token); + const L = engine.getZoneLevel(zone); + if (L > level) return 0; + if (L === level) { + n += 1; + continue; + } + n += Number(engine.countSubZones(zone, level - L)); + } + return n; +} + +/** Compact a polygon cell layer's `dggal` (or other) ID field. */ +export function compactDggalFeatureCollection( + engine: DggalDggrs, + fc: FeatureCollection, + options: { cellField?: string } = {}, +): FeatureCollection { + const tokens = compactDggalTokens(engine, tokensFromDggalLayer(fc, options.cellField ?? "dggal")); + return dggalTokensToFeatureCollection(engine, tokens); +} + +/** Expand a polygon cell layer to a uniform DGGAL level. */ +export function expandDggalFeatureCollection( + engine: DggalDggrs, + fc: FeatureCollection, + level: number, + options: { cellField?: string } = {}, +): FeatureCollection { + const tokens = expandDggalTokens( + engine, + tokensFromDggalLayer(fc, options.cellField ?? "dggal"), + level, + ); + return dggalTokensToFeatureCollection(engine, tokens); +} + +/** Supported point-binning aggregate operations (same set as H3). */ +export type DggalAggOp = "count" | "sum" | "mean" | "min" | "max"; + +function eachPointCoord( + geometry: Geometry | null | undefined, + visit: (pos: Position) => void, +): void { + if (!geometry) return; + if (geometry.type === "Point") { + visit(geometry.coordinates); + return; + } + if (geometry.type === "MultiPoint") { + for (const c of geometry.coordinates) visit(c); + } +} + +/** Aggregate point geometry into DGGAL zones (client-side). */ +export function binPointsToDggal( + engine: DggalDggrs, + fc: FeatureCollection, + resolution: number, + op: DggalAggOp, + field?: string, +): FeatureCollection { + type Acc = { count: number; sum: number; min: number; max: number }; + const byCell = new Map(); + + for (const feature of fc.features) { + eachPointCoord(feature.geometry, (pos) => { + const lng = pos[0]; + const lat = pos[1]; + if ( + lng === undefined || + lat === undefined || + !Number.isFinite(lng) || + !Number.isFinite(lat) + ) { + return; + } + const zone = engine.getZoneFromWGS84Centroid(resolution, { + lat: lat * RAD_PER_DEG, + lon: lng * RAD_PER_DEG, + }); + const token = engine.getZoneTextID(zone); + let acc = byCell.get(token); + if (!acc) { + acc = { count: 0, sum: 0, min: Number.POSITIVE_INFINITY, max: Number.NEGATIVE_INFINITY }; + byCell.set(token, acc); + } + acc.count += 1; + if (op !== "count" && field) { + const raw = feature.properties?.[field]; + const n = typeof raw === "number" ? raw : Number(raw); + if (!Number.isFinite(n)) return; + acc.sum += n; + if (n < acc.min) acc.min = n; + if (n > acc.max) acc.max = n; + } + }); + } + + const features: Feature[] = []; + for (const [token, acc] of byCell) { + const feature = dggalZoneFeature(engine, token); + const properties: Record = { ...feature.properties, count: acc.count }; + if (op === "sum") properties.value = acc.sum; + else if (op === "mean") properties.value = acc.count > 0 ? acc.sum / acc.count : 0; + else if (op === "min") properties.value = Number.isFinite(acc.min) ? acc.min : null; + else if (op === "max") properties.value = Number.isFinite(acc.max) ? acc.max : null; + features.push({ ...feature, properties }); + } + return { type: "FeatureCollection", features }; +} diff --git a/packages/processing/src/dggrid-tools.ts b/packages/processing/src/dggrid-tools.ts new file mode 100644 index 000000000..6ed3ce2f2 --- /dev/null +++ b/packages/processing/src/dggrid-tools.ts @@ -0,0 +1,464 @@ +import type { FeatureCollection, Geometry } from "geojson"; +import { unwrapAntimeridianGeometry } from "./antimeridian"; +import { sqlIdent, sqlStr } from "./h3-tools"; + +/** + * Default DGGRID orientation (ISEA/FULLER pole) used by duck_dggs examples: + * https://duckdb.org/community_extensions/extensions/duck_dggs + */ +const DGGRID_DEFAULT_ORIENT = "0.0, 58.3971459, 11.2" as const; + +/** + * Named DGGRID types exposed in DGGS Generator / Binning. Ranges match the + * standard DGGRID_TYPES table; `dggs_params` maps each name onto duck_dggs. + */ +export const DGGRID_GRID_TYPES = [ + "SUPERFUND", + "PLANETRISK", + "ISEA3H", + "ISEA4H", + "ISEA4T", + "ISEA4D", + "ISEA43H", + "ISEA7H", + "IGEO7", + "FULLER3H", + "FULLER4H", + "FULLER4T", + "FULLER4D", + "FULLER43H", + "FULLER7H", +] as const; +export type DggridGridType = (typeof DGGRID_GRID_TYPES)[number]; + +export const DEFAULT_DGGRID_GRID_TYPE: DggridGridType = "ISEA4H"; + +/** Resolution bounds + duck_dggs configuration for a named DGGRID type. */ +export type DggridGridSpec = { + minRes: number; + maxRes: number; + defaultRes: number; + /** Pure aperture (ignored when `apertureSequence` is set, except as a fallback). */ + aperture: 3 | 4 | 7; + projection: "ISEA" | "FULLER"; + topology: "HEXAGON" | "TRIANGLE" | "DIAMOND"; + /** + * Mixed-aperture sequence (digits 3/4/7). When set, SQL uses the 8-arg + * `dggs_params(..., true, sequence)` overload. + */ + apertureSequence?: string; +}; + +/** + * Specs for every named type. SUPERFUND = two aperture-4 then fifteen + * aperture-3 (Appendix E); PLANETRISK = `433347…7` (Appendix F); ISEA43H / + * FULLER43H default to MIXED43 with zero leading aperture-4 levels (pure + * aperture 3 geometrically). IGEO7 uses the canonical ISEA aperture-7 hex grid. + */ +export const DGGRID_GRID_SPECS: Record = { + SUPERFUND: { + minRes: 0, + maxRes: 17, + defaultRes: 9, + aperture: 3, + projection: "FULLER", + topology: "HEXAGON", + // Two aperture-4 + fifteen aperture-3 (EPA Superfund_500m). + apertureSequence: "44333333333333333", + }, + PLANETRISK: { + minRes: 0, + maxRes: 22, + defaultRes: 13, + aperture: 7, + projection: "ISEA", + topology: "HEXAGON", + apertureSequence: "43334777777777777777777", + }, + ISEA3H: { + minRes: 0, + maxRes: 35, + defaultRes: 20, + aperture: 3, + projection: "ISEA", + topology: "HEXAGON", + }, + ISEA4H: { + minRes: 0, + maxRes: 29, + defaultRes: 16, + aperture: 4, + projection: "ISEA", + topology: "HEXAGON", + }, + ISEA4T: { + minRes: 0, + maxRes: 29, + defaultRes: 15, + aperture: 4, + projection: "ISEA", + topology: "TRIANGLE", + }, + ISEA4D: { + minRes: 0, + maxRes: 29, + defaultRes: 16, + aperture: 4, + projection: "ISEA", + topology: "DIAMOND", + }, + ISEA43H: { + minRes: 0, + maxRes: 18, + defaultRes: 10, + // MIXED43 with dggs_num_aperture_4_res = 0 → all aperture 3. + aperture: 3, + projection: "ISEA", + topology: "HEXAGON", + }, + ISEA7H: { + minRes: 0, + maxRes: 21, + defaultRes: 11, + aperture: 7, + projection: "ISEA", + topology: "HEXAGON", + }, + IGEO7: { + minRes: 0, + maxRes: 20, + defaultRes: 12, + aperture: 7, + projection: "ISEA", + topology: "HEXAGON", + }, + FULLER3H: { + minRes: 0, + maxRes: 35, + defaultRes: 20, + aperture: 3, + projection: "FULLER", + topology: "HEXAGON", + }, + FULLER4H: { + minRes: 0, + maxRes: 30, + defaultRes: 16, + aperture: 4, + projection: "FULLER", + topology: "HEXAGON", + }, + FULLER4T: { + minRes: 0, + maxRes: 29, + defaultRes: 15, + aperture: 4, + projection: "FULLER", + topology: "TRIANGLE", + }, + FULLER4D: { + minRes: 0, + maxRes: 30, + defaultRes: 16, + aperture: 4, + projection: "FULLER", + topology: "DIAMOND", + }, + FULLER43H: { + minRes: 0, + maxRes: 18, + defaultRes: 10, + aperture: 3, + projection: "FULLER", + topology: "HEXAGON", + }, + FULLER7H: { + minRes: 0, + maxRes: 21, + defaultRes: 11, + aperture: 7, + projection: "FULLER", + topology: "HEXAGON", + }, +}; + +/** Aperture used by each preset (drives pure-grid cell-count estimates). */ +export const DGGRID_GRID_APERTURE: Record = Object.fromEntries( + DGGRID_GRID_TYPES.map((t) => [t, DGGRID_GRID_SPECS[t].aperture]), +) as Record; + +function dggsParamsSql(spec: DggridGridSpec): string { + const head = `'${spec.projection}', ${spec.aperture}, '${spec.topology}', ${DGGRID_DEFAULT_ORIENT}`; + if (spec.apertureSequence) { + return `dggs_params(${head}, true, '${spec.apertureSequence}')`; + } + return `dggs_params(${head})`; +} + +/** SQL `dggs_params(...)` fragment for each named type. */ +export const DGGRID_GRID_PARAMS_SQL: Record = Object.fromEntries( + DGGRID_GRID_TYPES.map((t) => [t, dggsParamsSql(DGGRID_GRID_SPECS[t])]), +) as Record; + +const DGGRID_GRID_TYPE_LABEL: Record = { + SUPERFUND: "SUPERFUND", + PLANETRISK: "PLANETRISK", + ISEA3H: "ISEA3H", + ISEA4H: "ISEA4H", + ISEA4T: "ISEA4T", + ISEA4D: "ISEA4D", + ISEA43H: "ISEA43H", + ISEA7H: "ISEA7H", + IGEO7: "IGEO7", + FULLER3H: "FULLER3H", + FULLER4H: "FULLER4H", + FULLER4T: "FULLER4T", + FULLER4D: "FULLER4D", + FULLER43H: "FULLER43H", + FULLER7H: "FULLER7H", +}; + +/** Labels shown in the processing dialog. */ +export const DGGRID_GRID_TYPE_OPTIONS: { value: DggridGridType; label: string }[] = + DGGRID_GRID_TYPES.map((value) => ({ value, label: DGGRID_GRID_TYPE_LABEL[value] })); + +export function resolveDggridGridType(raw: unknown): DggridGridType { + if (typeof raw === "string" && (DGGRID_GRID_TYPES as readonly string[]).includes(raw)) { + return raw as DggridGridType; + } + return DEFAULT_DGGRID_GRID_TYPE; +} + +/** Max resolution allowed for a DGGRID named type. */ +export function maxResolutionForDggrid( + gridType: DggridGridType = DEFAULT_DGGRID_GRID_TYPE, +): number { + return DGGRID_GRID_SPECS[gridType].maxRes; +} + +/** Soft target used when auto-suggesting a resolution. */ +export const DGGRID_TARGET_CELLS = 10_000; +/** Finest resolution the auto-suggester will pick. */ +export const DGGRID_MAX_SUGGESTED_RES = 12; +/** Hard ceiling: a grid larger than this aborts rather than running away. */ +export const DGGRID_HARD_CAP = 200_000; +/** + * Absolute finest resolution across exposed DGGRID types (ISEA3H / FULLER3H). + * The dialog narrows this per selected {@link DggridGridType}. + */ +export const DGGRID_MAX_TOOL_RES = Math.max( + ...DGGRID_GRID_TYPES.map((t) => DGGRID_GRID_SPECS[t].maxRes), +); +/** + * Cap on sample-grid axes when polyfilling. duck_dggs has no polygon cover + * function, so the generator densifies the envelope and maps points → cells. + */ +export const DGGRID_SAMPLE_AXIS_CAP = 500; + +const EARTH_AREA_KM2 = 510_065_621.724; + +/** Approximate global cell count at `res` for a named type. */ +export function dggridCellCountAtRes(res: number, gridType: DggridGridType): number { + const spec = DGGRID_GRID_SPECS[gridType]; + if (spec.apertureSequence) { + let prod = 1; + const seq = spec.apertureSequence; + for (let i = 0; i < res; i += 1) { + const digit = Number(seq[i] ?? seq[seq.length - 1] ?? spec.aperture); + prod *= digit; + } + return 10 * prod + 2; + } + return 10 * spec.aperture ** res + 2; +} + +/** Approximate average cell area (km²) at `res` for an aperture-N icosahedral grid. */ +export function dggridAvgCellAreaKm2(res: number, aperture: 3 | 4 | 7): number { + return EARTH_AREA_KM2 / (10 * aperture ** res + 2); +} + +/** Approximate average cell area (km²) at res for aperture-4 icosahedral grids. */ +export const DGGRID_AVG_AREA_KM2_A4: number[] = Array.from( + { length: DGGRID_MAX_TOOL_RES + 1 }, + (_, res) => dggridAvgCellAreaKm2(res, 4), +); + +/** Approximate average cell area (km²) at res for aperture-3 hexagon grids. */ +export const DGGRID_AVG_AREA_KM2_A3: number[] = Array.from( + { length: DGGRID_MAX_TOOL_RES + 1 }, + (_, res) => dggridAvgCellAreaKm2(res, 3), +); + +/** @deprecated Prefer {@link DGGRID_AVG_AREA_KM2_A4}; kept for existing imports. */ +export const DGGRID_AVG_AREA_KM2 = DGGRID_AVG_AREA_KM2_A4; + +/** Estimated number of DGGRID cells covering `areaKm2` at `res`. */ +export function estimateDggridCellCount( + areaKm2: number, + res: number, + gridType: DggridGridType = DEFAULT_DGGRID_GRID_TYPE, +): number { + const spec = DGGRID_GRID_SPECS[gridType]; + if (!Number.isInteger(res) || res < spec.minRes || res > spec.maxRes) { + return Number.POSITIVE_INFINITY; + } + return areaKm2 / (EARTH_AREA_KM2 / dggridCellCountAtRes(res, gridType)); +} + +/** Finest resolution whose estimated cell count stays <= the target. */ +export function suggestDggridResolution( + areaKm2: number, + targetCells = DGGRID_TARGET_CELLS, + maxRes = DGGRID_MAX_SUGGESTED_RES, + gridType: DggridGridType = DEFAULT_DGGRID_GRID_TYPE, +): number { + const capped = Math.min(maxRes, maxResolutionForDggrid(gridType)); + for (let res = capped; res >= 0; res -= 1) { + if (estimateDggridCellCount(areaKm2, res, gridType) <= targetCells) return res; + } + return 0; +} + +function paramsSql(gridType: DggridGridType): string { + return DGGRID_GRID_PARAMS_SQL[gridType]; +} + +/** + * Sample-based covering: densify the geometry envelope, keep points that + * intersect it, map each to a seqnum via `geo_to_seqnum`, then emit boundaries. + * duck_dggs ([docs](https://duckdb.org/community_extensions/extensions/duck_dggs)) + * only converts POINT → cell; there is no H3-style polygon polyfill. + * + * When the envelope would need more than {@link DGGRID_SAMPLE_AXIS_CAP} samples + * on an axis, the step is scaled up so the series still spans the full bbox + * (large areas may be under-sampled rather than clipped mid-extent). + * + * @param areaSelectSql SQL for the `_dggs_area` CTE body, e.g. + * `SELECT ST_GeomFromText(...) AS g` or `SELECT g FROM (...) WHERE g IS NOT NULL`. + */ +function cellsCteFromGeom(areaSelectSql: string, res: number, gridType: DggridGridType): string { + const cap = DGGRID_SAMPLE_AXIS_CAP; + const p = paramsSql(gridType); + return ( + `WITH _dggs_area AS (${areaSelectSql}), ` + + `_dggs_meta AS (` + + `SELECT ST_XMin(g) AS w, ST_YMin(g) AS s, ST_XMax(g) AS e, ST_YMax(g) AS n, g, ` + + // Half the characteristic length scale (km) → degrees (~111.32 km/deg), + // raised enough that ≤ cap samples cover each axis end-to-end. + `GREATEST(dggs_cls_km(${res}, ${p}) / 222.64, (e - w) / ${cap}, (n - s) / ${cap}, 1e-5) AS step ` + + `FROM _dggs_area), ` + + `_dggs_grid AS (` + + `SELECT ST_Point(w + i * step, s + j * step) AS pt, g FROM _dggs_meta, ` + + `generate_series(0, CAST(CEIL((e - w) / step) AS BIGINT)) AS t(i), ` + + `generate_series(0, CAST(CEIL((n - s) / step) AS BIGINT)) AS u(j)), ` + + `_dggs_pts AS (` + + `SELECT pt FROM _dggs_grid WHERE ST_Intersects(pt, g) ` + + `UNION ALL SELECT ST_Centroid(g) FROM _dggs_meta WHERE g IS NOT NULL), ` + + `cells AS (SELECT DISTINCT geo_to_seqnum(pt, ${res}, ${p}) AS cell FROM _dggs_pts WHERE pt IS NOT NULL)` + ); +} + +function gridSelect(res: number, gridType: DggridGridType): string { + const p = paramsSql(gridType); + return ( + `SELECT CAST(cell AS VARCHAR) AS dggrid, ` + + `CAST(ST_AsGeoJSON(seqnum_to_boundary(cell, ${res}, ${p})) AS VARCHAR) AS geojson FROM cells` + ); +} + +/** Grid SQL from a polygon WKT literal (bbox / viewport sources). */ +export function buildDggridGridFromWktSql( + wkt: string, + res: number, + gridType: DggridGridType = DEFAULT_DGGRID_GRID_TYPE, +): string { + return ( + `${cellsCteFromGeom(`SELECT ST_GeomFromText(${sqlStr(wkt)}) AS g`, res, gridType)} ` + + gridSelect(res, gridType) + ); +} + +/** + * Grid SQL that unions polygonal geometry from a registered source, then + * sample-covers it (polyfill source). + */ +export function buildDggridGridFromSourceSql( + sourceSql: string, + res: number, + gridType: DggridGridType = DEFAULT_DGGRID_GRID_TYPE, +): string { + const merged = + `SELECT g FROM (` + + `SELECT ST_Union_Agg(geom) AS g FROM ${sourceSql} ` + + `WHERE geom IS NOT NULL AND ST_GeometryType(geom) IN ('POLYGON', 'MULTIPOLYGON')` + + `) WHERE g IS NOT NULL`; + return `${cellsCteFromGeom(merged, res, gridType)} ` + gridSelect(res, gridType); +} + +/** Supported point-binning aggregate operations. */ +export type DggridAggOp = "count" | "sum" | "mean" | "min" | "max"; + +const AGG_FN: Record, string> = { + sum: "sum", + mean: "avg", + min: "min", + max: "max", +}; + +/** + * Aggregate point geometry into DGGRID cells via `geo_to_seqnum`. Boundaries + * come from `seqnum_to_boundary`. + */ +export function buildDggridBinSql( + sourceSql: string, + res: number, + op: DggridAggOp, + field?: string, + gridType: DggridGridType = DEFAULT_DGGRID_GRID_TYPE, +): string { + const fn = op === "count" ? undefined : AGG_FN[op]; + const aggSelect = fn && field ? `, ${fn}(CAST(${sqlIdent(field)} AS DOUBLE)) AS value` : ""; + const aggOut = fn && field ? ", value" : ""; + const p = paramsSql(gridType); + return ( + `WITH pts AS (SELECT ST_Centroid(geom) AS pt` + + (field ? `, ${sqlIdent(field)}` : "") + + ` FROM ${sourceSql} ` + + `WHERE geom IS NOT NULL AND ST_GeometryType(geom) IN ('POINT', 'MULTIPOINT')), ` + + `binned AS (SELECT geo_to_seqnum(pt, ${res}, ${p}) AS cell, ` + + `count(*) AS count${aggSelect} FROM pts GROUP BY cell) ` + + `SELECT CAST(cell AS VARCHAR) AS dggrid, count${aggOut}, ` + + `CAST(ST_AsGeoJSON(seqnum_to_boundary(cell, ${res}, ${p})) AS VARCHAR) AS geojson FROM binned` + ); +} + +/** Build a FeatureCollection from rows carrying `dggrid`, optional aggregates, and `geojson`. */ +export function dggridRowsToFeatureCollection( + rows: Record[], + fixAntimeridian = true, +): FeatureCollection { + const features = []; + for (const row of rows) { + const raw = row.geojson; + if (typeof raw !== "string") continue; + let geometry: Geometry; + try { + geometry = JSON.parse(raw) as Geometry; + } catch { + continue; + } + const properties: Record = { dggrid: String(row.dggrid) }; + if (row.count !== undefined && row.count !== null) { + properties.count = Number(row.count); + } + if (row.value !== undefined && row.value !== null) { + properties.value = Number(row.value); + } + features.push({ + type: "Feature" as const, + geometry: fixAntimeridian ? unwrapAntimeridianGeometry(geometry) : geometry, + properties, + }); + } + return { type: "FeatureCollection", features }; +} diff --git a/packages/processing/src/dggs-tools.ts b/packages/processing/src/dggs-tools.ts new file mode 100644 index 000000000..ba95702df --- /dev/null +++ b/packages/processing/src/dggs-tools.ts @@ -0,0 +1,1030 @@ +import type { FeatureCollection } from "geojson"; +import bbox from "@turf/bbox"; +import type { GeoLibreLayer } from "@geolibre/core"; +import type { + DuckDbCapability, + DuckDbGeoJsonSource, + ProcessingAlgorithm, + ProcessingContext, +} from "./types"; +import { + a5RowsToFeatureCollection, + buildA5BinSql, + buildA5CompactSql, + buildA5ExpandCountSql, + buildA5ExpandSql, + buildA5GridFromBboxSql, + buildA5GridFromSourceSql, + A5_HARD_CAP, + A5_MAX_TOOL_RES, + estimateA5CellCount, + suggestA5Resolution, + type A5AggOp, +} from "./a5-tools"; +import { + buildDggridBinSql, + buildDggridGridFromSourceSql, + buildDggridGridFromWktSql, + DEFAULT_DGGRID_GRID_TYPE, + DGGRID_GRID_TYPE_OPTIONS, + DGGRID_HARD_CAP, + DGGRID_MAX_TOOL_RES, + dggridRowsToFeatureCollection, + estimateDggridCellCount, + maxResolutionForDggrid, + resolveDggridGridType, + suggestDggridResolution, + type DggridAggOp, + type DggridGridType, +} from "./dggrid-tools"; +import { + binPointsToDggal, + compactDggalFeatureCollection, + dggalGridFromBbox, + dggalGridFromFeatureCollection, + DEFAULT_DGGAL_GRID_TYPE, + DGGAL_GRID_TYPE_OPTIONS, + DGGAL_HARD_CAP, + DGGAL_MAX_TOOL_RES, + DGGAL_TYPES, + estimateDggalCellCount, + estimateDggalExpandCount, + expandDggalFeatureCollection, + maxResolutionForDggal, + resolveDggalGridType, + suggestDggalResolution, + tokensFromDggalLayer, + withDggalDggrs, + type DggalAggOp, + type DggalGridType, +} from "./dggal-tools"; +import { + binPointsToS2, + compactS2FeatureCollection, + estimateS2CellCount, + estimateS2ExpandCount, + expandS2FeatureCollection, + s2GridFromBbox, + s2GridFromFeatureCollection, + suggestS2Resolution, + tokensFromS2Layer, + S2_HARD_CAP, + S2_MAX_TOOL_RES, + type S2AggOp, +} from "./s2-tools"; +import { + bboxAreaKm2, + bboxToWktPolygon, + buildBinSql, + buildGridFromBboxSql, + buildGridFromSourceSql, + buildH3CompactSql, + buildH3ExpandCountSql, + buildH3ExpandSql, + estimateCellCount, + H3_AGG_OPS, + H3_HARD_CAP, + normalizeLonLatBbox, + rowsToFeatureCollection, + suggestResolution, + type H3AggOp, +} from "./h3-tools"; + +/** Supported DGGS backends for the DGGS Generator / Binning tools. */ +export type DggsType = "h3" | "s2" | "a5" | "dggrid" | "dggal"; + +export const DGGS_TYPES: readonly DggsType[] = ["h3", "s2", "a5", "dggrid", "dggal"]; + +const DGGS_TYPE_LABEL: Record = { + h3: "H3", + s2: "S2", + a5: "A5", + dggrid: "DGGRID", + dggal: "DGGAL", +}; + +/** + * Max resolution for the selected DGGS type. `subtype` is the DGGRID or DGGAL + * named type when applicable. + */ +export function maxResolutionForDggs( + type: DggsType, + subtype?: DggridGridType | DggalGridType | string, +): number { + if (type === "s2") return S2_MAX_TOOL_RES; + if (type === "a5") return A5_MAX_TOOL_RES; + if (type === "dggrid") return maxResolutionForDggrid(resolveDggridGridType(subtype)); + if (type === "dggal") return maxResolutionForDggal(resolveDggalGridType(subtype)); + return 15; +} + +/** DuckDB community extension name required for `type` (S2/DGGAL are client-side). */ +export function extensionForDggs(type: DggsType): string | null { + if (type === "s2" || type === "dggal") return null; + if (type === "a5") return "a5"; + if (type === "dggrid") return "duck_dggs"; + return "h3"; +} + +const NO_DUCKDB = "This tool requires DuckDB-WASM, which is unavailable in this environment."; + +function requireDuckDb(ctx: ProcessingContext): DuckDbCapability { + if (!ctx.duckdb) throw new Error(NO_DUCKDB); + return ctx.duckdb; +} + +function getLayer(ctx: ProcessingContext, paramId = "layer"): GeoLibreLayer | undefined { + const id = ctx.parameters[paramId] as string | undefined; + return ctx.layers.find((l) => l.id === id); +} + +function numberParam(ctx: ProcessingContext, id: string): number { + const raw = ctx.parameters[id]; + if (raw === undefined || raw === null || raw === "") return NaN; + return typeof raw === "string" ? Number(raw) : (raw as number); +} + +function bboxFromParams(ctx: ProcessingContext): [number, number, number, number] | null { + const west = numberParam(ctx, "west"); + const south = numberParam(ctx, "south"); + const east = numberParam(ctx, "east"); + const north = numberParam(ctx, "north"); + if ([west, south, east, north].some((n) => !Number.isFinite(n))) { + ctx.log("Error: enter numeric west, south, east, and north values"); + return null; + } + if (west >= east || south >= north) { + ctx.log("Error: bounding box must have west < east and south < north"); + return null; + } + return [west, south, east, north]; +} + +function resolveDggsType(ctx: ProcessingContext): DggsType | null { + const raw = (ctx.parameters.dggsType as string) || "h3"; + if (raw === "h3" || raw === "s2" || raw === "a5" || raw === "dggrid" || raw === "dggal") { + return raw; + } + ctx.log(`Error: unknown DGGS type "${raw}"`); + return null; +} + +function resolveResolution( + ctx: ProcessingContext, + type: DggsType, + areaKm2: number, + dggridType: DggridGridType = DEFAULT_DGGRID_GRID_TYPE, + dggalType: DggalGridType = DEFAULT_DGGAL_GRID_TYPE, +): number | null { + const maxRes = maxResolutionForDggs(type, type === "dggal" ? dggalType : dggridType); + const raw = ctx.parameters.resolution; + if (raw === undefined || raw === null || raw === "") { + const suggested = + type === "s2" + ? suggestS2Resolution(areaKm2) + : type === "a5" + ? suggestA5Resolution(areaKm2) + : type === "dggrid" + ? suggestDggridResolution(areaKm2, undefined, undefined, dggridType) + : type === "dggal" + ? suggestDggalResolution(areaKm2, undefined, undefined, dggalType) + : suggestResolution(areaKm2); + ctx.log(`Using suggested resolution ${suggested}`); + return suggested; + } + const res = typeof raw === "string" ? Number(raw) : (raw as number); + if (!Number.isInteger(res) || res < 0 || res > maxRes) { + ctx.log( + `Error: resolution must be an integer from 0 to ${maxRes} for ${dggsLabel(type, dggridType, dggalType)}`, + ); + return null; + } + return res; +} + +function estimateFor( + type: DggsType, + areaKm2: number, + res: number, + dggridType: DggridGridType = DEFAULT_DGGRID_GRID_TYPE, + dggalType: DggalGridType = DEFAULT_DGGAL_GRID_TYPE, +): number { + if (type === "s2") return estimateS2CellCount(areaKm2, res); + if (type === "a5") return estimateA5CellCount(areaKm2, res); + if (type === "dggrid") return estimateDggridCellCount(areaKm2, res, dggridType); + if (type === "dggal") return estimateDggalCellCount(areaKm2, res, dggalType); + return estimateCellCount(areaKm2, res); +} + +function hardCapFor(type: DggsType): number { + if (type === "s2") return S2_HARD_CAP; + if (type === "a5") return A5_HARD_CAP; + if (type === "dggrid") return DGGRID_HARD_CAP; + if (type === "dggal") return DGGAL_HARD_CAP; + return H3_HARD_CAP; +} + +const DGGS_TYPE_PARAM = { + id: "dggsType", + label: "DGGS type", + type: "select" as const, + default: "h3", + options: [ + { value: "h3", label: "H3" }, + { value: "s2", label: "S2" }, + { value: "a5", label: "A5" }, + { value: "dggrid", label: "DGGRID" }, + { value: "dggal", label: "DGGAL" }, + ], +}; + +/** Shown only when DGGS type is DGGRID — duck_dggs presets. */ +const DGGRID_TYPE_PARAM = { + id: "dggridType", + label: "DGGRID type", + type: "select" as const, + default: DEFAULT_DGGRID_GRID_TYPE, + options: DGGRID_GRID_TYPE_OPTIONS, + visibleWhen: { param: "dggsType", in: ["dggrid"] }, + description: + "Grid configuration passed to duck_dggs as dggs_params (ISEA/FULLER × aperture × topology).", +}; + +/** Shown only when DGGS type is DGGAL — DGGRS class presets. */ +const DGGAL_TYPE_PARAM = { + id: "dggalType", + label: "DGGAL type", + type: "select" as const, + default: DEFAULT_DGGAL_GRID_TYPE, + options: DGGAL_GRID_TYPE_OPTIONS, + visibleWhen: { param: "dggsType", in: ["dggal"] }, + description: "DGGRS passed to dggal.createDGGRS (ISEA/IVEA/RTEA/HEALPix/GNOSIS).", +}; + +/** + * Unwrap dateline-straddling cell rings for MapLibre. H3, S2, and DGGRID only; + * A5 and DGGAL emit dateline-safe geometry natively. + */ +const FIX_ANTIMERIDIAN_PARAM = { + id: "fixAntimeridian", + label: "Fix antimeridian", + type: "boolean" as const, + default: true, + visibleWhen: { param: "dggsType", in: ["h3", "s2", "dggrid"] }, + description: "Unwrap cell rings that cross ±180° longitude.", +}; + +/** + * After filling at the requested resolution, merge complete child sets into + * parents (mixed resolutions). H3, A5, S2, and DGGAL; default off (uniform cells). + */ +const COMPACT_CELLS_PARAM = { + id: "compactCells", + label: "Compact cells", + type: "boolean" as const, + default: false, + visibleWhen: { param: "dggsType", in: ["h3", "a5", "s2", "dggal"] }, + description: + "Merge complete sets of sibling cells into coarser parents (fewer features, mixed resolutions).", +}; + +function resolveFixAntimeridian(ctx: ProcessingContext, type: DggsType): boolean { + if (type !== "h3" && type !== "s2" && type !== "dggrid") { + return false; + } + const raw = ctx.parameters.fixAntimeridian; + // Default checked when the param was never set (e.g. scripted runs). + if (raw === undefined || raw === null || raw === "") return true; + return Boolean(raw); +} + +function dggsLabel(type: DggsType, dggridType: DggridGridType, dggalType: DggalGridType): string { + if (type === "dggrid") return dggridType; + if (type === "dggal") return DGGAL_TYPES[dggalType].className; + return DGGS_TYPE_LABEL[type]; +} + +function resolveCompactCells(ctx: ProcessingContext, type: DggsType): boolean { + if (type !== "h3" && type !== "a5" && type !== "s2" && type !== "dggal") return false; + return Boolean(ctx.parameters.compactCells); +} + +/** + * Fill an area with DGGS cells. H3/A5/DGGRID use DuckDB-WASM community + * extensions; S2 (s2js) and DGGAL (dggal WASM) run client-side. + */ +export const createDggsGridTool: ProcessingAlgorithm = { + id: "dggs-grid", + name: "DGGS Generator", + description: + "Fill an area with DGGS cells (H3, S2, A5, DGGRID, DGGAL). Source: a layer's geometry, a layer's extent, the current map view, or a manual bounding box.", + group: "DGGS", + parameters: [ + DGGS_TYPE_PARAM, + DGGRID_TYPE_PARAM, + DGGAL_TYPE_PARAM, + { + id: "source", + label: "Area source", + type: "select", + default: "polyfill", + options: [ + { value: "polyfill", label: "Layer geometry (polyfill)" }, + { value: "extent", label: "Layer extent (bbox)" }, + { value: "viewport", label: "Map viewport" }, + { value: "bbox", label: "Manual bounding box" }, + ], + }, + { + id: "layer", + label: "Input layer", + type: "layer", + required: true, + visibleWhen: { param: "source", in: ["polyfill", "extent"] }, + }, + { + id: "west", + label: "West (min lon)", + type: "number", + required: true, + min: -180, + max: 180, + visibleWhen: { param: "source", in: ["bbox"] }, + }, + { + id: "south", + label: "South (min lat)", + type: "number", + required: true, + min: -90, + max: 90, + visibleWhen: { param: "source", in: ["bbox"] }, + }, + { + id: "east", + label: "East (max lon)", + type: "number", + required: true, + min: -180, + max: 180, + visibleWhen: { param: "source", in: ["bbox"] }, + }, + { + id: "north", + label: "North (max lat)", + type: "number", + required: true, + min: -90, + max: 90, + visibleWhen: { param: "source", in: ["bbox"] }, + }, + { + // Absolute ceiling is the finest supported DGGS (currently ISEA3H); the + // dialog narrows the input max for the selected type / subtype. + id: "resolution", + label: "Resolution", + type: "number", + min: 0, + max: Math.max(S2_MAX_TOOL_RES, A5_MAX_TOOL_RES, DGGRID_MAX_TOOL_RES, DGGAL_MAX_TOOL_RES), + step: 1, + description: "Range depends on DGGS type. Leave blank to auto-pick from the area.", + }, + COMPACT_CELLS_PARAM, + FIX_ANTIMERIDIAN_PARAM, + ], + run: async (ctx) => { + const type = resolveDggsType(ctx); + if (!type) return; + const source = (ctx.parameters.source as string) || "polyfill"; + + let areaKm2: number; + let areaBbox: [number, number, number, number] | null = null; + let inputGeojson: FeatureCollection | null = null; + if (source === "viewport") { + const bounds = ctx.viewportBounds?.(); + if (!bounds) { + ctx.log("Error: map viewport is unavailable"); + return; + } + if (bounds[0] >= bounds[2]) { + ctx.log( + "Error: the map view crosses the antimeridian; pan so it doesn't wrap +/-180, or use a manual bounding box", + ); + return; + } + areaBbox = normalizeLonLatBbox(bounds); + areaKm2 = bboxAreaKm2(areaBbox); + } else if (source === "bbox") { + const bounds = bboxFromParams(ctx); + if (!bounds) return; + areaBbox = normalizeLonLatBbox(bounds); + areaKm2 = bboxAreaKm2(areaBbox); + } else { + const layer = getLayer(ctx, "layer"); + if (!layer?.geojson?.features?.length) { + ctx.log('Error: parameter "layer" has no GeoJSON features'); + return; + } + if (source === "polyfill") { + const hasPolygon = layer.geojson.features.some( + (f) => f.geometry?.type === "Polygon" || f.geometry?.type === "MultiPolygon", + ); + if (!hasPolygon) { + ctx.log( + 'Error: polyfill needs a polygon layer; use the "Layer extent" source for point or line layers', + ); + return; + } + } + inputGeojson = layer.geojson; + const bb = normalizeLonLatBbox(bbox(layer.geojson) as [number, number, number, number]); + areaKm2 = bboxAreaKm2(bb); + if (source === "extent") areaBbox = bb; + // A5 `geometry_to_cells` returns [] for a ±180° world ring; use the bbox + // path (res0 enumeration) when the layer extent is full longitude. + if (source === "polyfill" && type === "a5" && bb[0] === -180 && bb[2] === 180) { + areaBbox = bb; + } + } + + const dggridType = resolveDggridGridType(ctx.parameters.dggridType); + const dggalType = resolveDggalGridType(ctx.parameters.dggalType); + const res = resolveResolution(ctx, type, areaKm2, dggridType, dggalType); + if (res === null) return; + + const estimate = estimateFor(type, areaKm2, res, dggridType, dggalType); + const hardCap = hardCapFor(type); + if (estimate > hardCap) { + ctx.log( + `Error: resolution ${res} would generate about ${Math.round( + estimate, + ).toLocaleString()} cells (cap ${hardCap.toLocaleString()}). Choose a coarser resolution.`, + ); + return; + } + + const fixAntimeridian = resolveFixAntimeridian(ctx, type); + const compactCells = resolveCompactCells(ctx, type); + const label = dggsLabel(type, dggridType, dggalType); + + // S2 is entirely client-side (s2js); no DuckDB extension. + if (type === "s2") { + try { + const gridOpts = { + limit: hardCap, + unwrap: fixAntimeridian, + compact: compactCells, + }; + const fc = + areaBbox != null + ? s2GridFromBbox(areaBbox, res, gridOpts) + : s2GridFromFeatureCollection(inputGeojson!, res, gridOpts); + if (fc.features.length === 0) { + ctx.log( + `No S2 cells were produced at resolution ${res}. Try a finer resolution or a larger area.`, + ); + return; + } + ctx.log( + `Created ${fc.features.length} S2 cell(s) at resolution ${res}` + + (compactCells ? " (compacted)" : ""), + ); + ctx.addResultLayer?.( + compactCells ? `S2 grid (res ${res}, compact)` : `S2 grid (res ${res})`, + fc, + ); + } catch (error) { + const message = error instanceof Error ? error.message : String(error); + ctx.log(`Error: ${message}`); + } + return; + } + + // DGGAL is client-side WASM (dggal); no DuckDB. + if (type === "dggal") { + try { + const fc = await withDggalDggrs(dggalType, (engine) => + areaBbox != null + ? dggalGridFromBbox(engine, areaBbox, res, hardCap, { compact: compactCells }) + : dggalGridFromFeatureCollection(engine, inputGeojson!, res, hardCap, { + compact: compactCells, + }), + ); + if (fc.features.length === 0) { + ctx.log( + `No ${label} cells were produced at resolution ${res}. Try a finer resolution or a larger area.`, + ); + return; + } + ctx.log( + `Created ${fc.features.length} ${label} cell(s) at resolution ${res}` + + (compactCells ? " (compacted)" : ""), + ); + ctx.addResultLayer?.( + compactCells ? `${label} grid (res ${res}, compact)` : `${label} grid (res ${res})`, + fc, + ); + } catch (error) { + const message = error instanceof Error ? error.message : String(error); + ctx.log(`Error: ${message}`); + } + return; + } + + const duckdb = requireDuckDb(ctx); + const extension = extensionForDggs(type); + let registered: DuckDbGeoJsonSource | null = null; + try { + await duckdb.ensureExtensions(["spatial", extension!]); + let sql: string; + if (areaBbox) { + sql = + type === "a5" + ? buildA5GridFromBboxSql(areaBbox, res, compactCells) + : type === "dggrid" + ? buildDggridGridFromWktSql(bboxToWktPolygon(areaBbox), res, dggridType) + : buildGridFromBboxSql(areaBbox, res, compactCells); + } else { + registered = await duckdb.registerGeoJson(inputGeojson!); + sql = + type === "a5" + ? buildA5GridFromSourceSql(registered.sql, res, compactCells) + : type === "dggrid" + ? buildDggridGridFromSourceSql(registered.sql, res, dggridType) + : buildGridFromSourceSql(registered.sql, res, compactCells); + } + const rows = await duckdb.query(sql); + const fc = + type === "a5" + ? a5RowsToFeatureCollection(rows) + : type === "dggrid" + ? dggridRowsToFeatureCollection(rows, fixAntimeridian) + : rowsToFeatureCollection(rows, fixAntimeridian); + if (fc.features.length === 0) { + ctx.log( + `No ${label} cells were produced at resolution ${res}. Try a finer resolution or a larger area.`, + ); + return; + } + ctx.log( + `Created ${fc.features.length} ${label} cell(s) at resolution ${res}` + + (compactCells ? " (compacted)" : ""), + ); + ctx.addResultLayer?.( + compactCells ? `${label} grid (res ${res}, compact)` : `${label} grid (res ${res})`, + fc, + ); + } catch (error) { + const message = error instanceof Error ? error.message : String(error); + ctx.log(`Error: ${message}`); + } finally { + await registered?.release(); + } + }, +}; + +/** + * Aggregate a point layer into DGGS cells (count, or sum/mean/min/max of a + * numeric field). Same dialog shape as {@link createDggsGridTool}'s type picker. + */ +export const dggsBinPointsTool: ProcessingAlgorithm = { + id: "dggs-bin", + name: "DGGS Binning", + description: + "Aggregate a point layer into DGGS cells (count, or sum/mean/min/max of a numeric field).", + group: "DGGS", + parameters: [ + DGGS_TYPE_PARAM, + DGGRID_TYPE_PARAM, + DGGAL_TYPE_PARAM, + { + id: "layer", + label: "Input point layer", + type: "layer", + required: true, + geometryFilter: ["point"], + }, + { + id: "aggOp", + label: "Aggregate", + type: "select", + default: "count", + options: [ + { value: "count", label: "Count" }, + { value: "sum", label: "Sum" }, + { value: "mean", label: "Mean" }, + { value: "min", label: "Min" }, + { value: "max", label: "Max" }, + ], + }, + { + id: "field", + label: "Field", + type: "field", + fieldSource: "layer", + required: true, + visibleWhen: { param: "aggOp", notIn: ["count"] }, + description: "Numeric field to aggregate.", + }, + { + id: "resolution", + label: "Resolution", + type: "number", + min: 0, + max: Math.max(S2_MAX_TOOL_RES, A5_MAX_TOOL_RES, DGGRID_MAX_TOOL_RES, DGGAL_MAX_TOOL_RES), + step: 1, + description: "Range depends on DGGS type. Leave blank to auto-pick from the area.", + }, + FIX_ANTIMERIDIAN_PARAM, + ], + run: async (ctx) => { + const type = resolveDggsType(ctx); + if (!type) return; + const layer = getLayer(ctx, "layer"); + if (!layer?.geojson?.features?.length) { + ctx.log('Error: parameter "layer" has no GeoJSON features'); + return; + } + const op = (ctx.parameters.aggOp as string) || "count"; + if (!H3_AGG_OPS.includes(op as H3AggOp)) { + ctx.log(`Error: unknown aggregate "${op}"`); + return; + } + const field = ctx.parameters.field as string | undefined; + if (op !== "count" && !field) { + ctx.log(`Error: select a numeric field to ${op}`); + return; + } + + const bb = bbox(layer.geojson) as [number, number, number, number]; + const dggridType = resolveDggridGridType(ctx.parameters.dggridType); + const dggalType = resolveDggalGridType(ctx.parameters.dggalType); + const res = resolveResolution(ctx, type, bboxAreaKm2(bb), dggridType, dggalType); + if (res === null) return; + + const fixAntimeridian = resolveFixAntimeridian(ctx, type); + const label = dggsLabel(type, dggridType, dggalType); + + if (type === "s2") { + const fc = binPointsToS2(layer.geojson, res, op as S2AggOp, field, { + unwrap: fixAntimeridian, + }); + if (fc.features.length === 0) { + ctx.log( + `No points fell into S2 cells at resolution ${res}. Check the layer has point geometries.`, + ); + return; + } + ctx.log(`Binned points into ${fc.features.length} S2 cell(s) at resolution ${res}`); + ctx.addResultLayer?.(`S2 bins (res ${res})`, fc); + return; + } + + if (type === "dggal") { + try { + const fc = await withDggalDggrs(dggalType, (engine) => + binPointsToDggal(engine, layer.geojson!, res, op as DggalAggOp, field), + ); + if (fc.features.length === 0) { + ctx.log( + `No points fell into ${label} cells at resolution ${res}. Check the layer has point geometries.`, + ); + return; + } + ctx.log(`Binned points into ${fc.features.length} ${label} cell(s) at resolution ${res}`); + ctx.addResultLayer?.(`${label} bins (res ${res})`, fc); + } catch (error) { + const message = error instanceof Error ? error.message : String(error); + ctx.log(`Error: ${message}`); + } + return; + } + + const duckdb = requireDuckDb(ctx); + const extension = extensionForDggs(type); + let registered: DuckDbGeoJsonSource | null = null; + try { + await duckdb.ensureExtensions(["spatial", extension!]); + registered = await duckdb.registerGeoJson(layer.geojson); + const sql = + type === "a5" + ? buildA5BinSql(registered.sql, res, op as A5AggOp, field) + : type === "dggrid" + ? buildDggridBinSql(registered.sql, res, op as DggridAggOp, field, dggridType) + : buildBinSql(registered.sql, res, op as H3AggOp, field); + const rows = await duckdb.query(sql); + const fc = + type === "a5" + ? a5RowsToFeatureCollection(rows) + : type === "dggrid" + ? dggridRowsToFeatureCollection(rows, fixAntimeridian) + : rowsToFeatureCollection(rows, fixAntimeridian); + if (fc.features.length === 0) { + ctx.log( + `No points fell into ${label} cells at resolution ${res}. Check the layer has point geometries.`, + ); + return; + } + ctx.log(`Binned points into ${fc.features.length} ${label} cell(s) at resolution ${res}`); + ctx.addResultLayer?.(`${label} bins (res ${res})`, fc); + } catch (error) { + const message = error instanceof Error ? error.message : String(error); + ctx.log(`Error: ${message}`); + } finally { + await registered?.release(); + } + }, +}; + +/** Compact / expand modes for {@link dggsCompactTool}. */ +export type DggsCompactMode = "compact" | "expand"; + +type DggsCompactType = "h3" | "a5" | "s2" | "dggal"; + +/** + * Compact a polygon DGGS cell layer, or expand (uncompact) it to a target + * resolution. Reads cell IDs from a property (`h3` / `a5` / `s2` / `dggal` by + * default). H3 and A5 use DuckDB; S2 and DGGAL run client-side. + */ +export const dggsCompactTool: ProcessingAlgorithm = { + id: "dggs-compact", + name: "DGGS Compact", + description: + "Compact DGGS polygon cells, or expand them to a uniform resolution (H3, S2, A5, DGGAL). Input must be a cell layer with an ID property.", + group: "DGGS", + parameters: [ + { + id: "dggsType", + label: "DGGS type", + type: "select", + default: "h3", + options: [ + { value: "h3", label: "H3" }, + { value: "s2", label: "S2" }, + { value: "a5", label: "A5" }, + { value: "dggal", label: "DGGAL" }, + ], + }, + DGGAL_TYPE_PARAM, + { + id: "mode", + label: "Mode", + type: "select", + default: "compact", + options: [ + { value: "compact", label: "Compact" }, + { value: "expand", label: "Expand" }, + ], + }, + { + id: "layer", + label: "Input DGGS layer", + type: "layer", + required: true, + geometryFilter: ["polygon"], + description: "Polygon cell layer from DGGS Generator / Binning (or equivalent).", + }, + { + id: "cellField", + label: "Cell ID field", + type: "field", + fieldSource: "layer", + description: "Defaults to h3, a5, s2, or dggal for the selected type when left blank.", + }, + { + id: "resolution", + label: "Target resolution", + type: "number", + min: 0, + max: Math.max(15, A5_MAX_TOOL_RES, S2_MAX_TOOL_RES, DGGAL_MAX_TOOL_RES), + step: 1, + required: true, + visibleWhen: { param: "mode", in: ["expand"] }, + description: "Resolution to expand to. Must be at least as fine as the input cells.", + }, + { + id: "fixAntimeridian", + label: "Fix antimeridian", + type: "boolean", + default: true, + visibleWhen: { param: "dggsType", in: ["h3", "s2"] }, + description: "Unwrap cell rings that cross ±180° longitude.", + }, + ], + run: async (ctx) => { + const typeRaw = (ctx.parameters.dggsType as string) || "h3"; + if (typeRaw !== "h3" && typeRaw !== "a5" && typeRaw !== "s2" && typeRaw !== "dggal") { + ctx.log( + `Error: DGGS Compact currently supports H3, A5, S2, and DGGAL only (got "${typeRaw}")`, + ); + return; + } + const type: DggsCompactType = typeRaw; + const dggalType = resolveDggalGridType(ctx.parameters.dggalType); + const mode = ((ctx.parameters.mode as string) || "compact") as DggsCompactMode; + if (mode !== "compact" && mode !== "expand") { + ctx.log(`Error: unknown mode "${mode}"`); + return; + } + + const layer = getLayer(ctx, "layer"); + if (!layer?.geojson?.features?.length) { + ctx.log('Error: parameter "layer" has no GeoJSON features'); + return; + } + const hasPolygon = layer.geojson.features.some( + (f) => f.geometry?.type === "Polygon" || f.geometry?.type === "MultiPolygon", + ); + if (!hasPolygon) { + ctx.log("Error: input must be a polygon DGGS cell layer"); + return; + } + + const defaultField = + type === "a5" ? "a5" : type === "s2" ? "s2" : type === "dggal" ? "dggal" : "h3"; + const cellField = + typeof ctx.parameters.cellField === "string" && ctx.parameters.cellField.trim() + ? ctx.parameters.cellField.trim() + : defaultField; + + const sample = layer.geojson.features.find((f) => f.properties?.[cellField] != null); + if (!sample) { + ctx.log( + `Error: no features have a "${cellField}" property. Pick the cell ID field, or run DGGS Generator first.`, + ); + return; + } + + let res = 0; + if (mode === "expand") { + const maxRes = maxResolutionForDggs(type, dggalType); + const raw = ctx.parameters.resolution; + if (raw === undefined || raw === null || raw === "") { + ctx.log("Error: enter a target resolution to expand to"); + return; + } + res = typeof raw === "string" ? Number(raw) : (raw as number); + if (!Number.isInteger(res) || res < 0 || res > maxRes) { + ctx.log( + `Error: resolution must be an integer from 0 to ${maxRes} for ${dggsLabel( + type, + DEFAULT_DGGRID_GRID_TYPE, + dggalType, + )}`, + ); + return; + } + } + + const label = dggsLabel(type, DEFAULT_DGGRID_GRID_TYPE, dggalType); + const emitResult = (fc: FeatureCollection) => { + if (fc.features.length === 0) { + ctx.log( + mode === "compact" + ? `No ${label} cells were produced by compact. Check the cell ID field.` + : `No ${label} cells were produced by expand. Check the cell ID field and target resolution.`, + ); + return; + } + const verb = mode === "compact" ? "Compacted" : "Expanded"; + const suffix = mode === "expand" ? ` to res ${res}` : ""; + ctx.log(`${verb} to ${fc.features.length} ${label} cell(s)${suffix}`); + ctx.addResultLayer?.( + mode === "compact" ? `${label} compact` : `${label} expand (res ${res})`, + fc, + ); + }; + + // S2 is client-side (s2js); no DuckDB. + if (type === "s2") { + try { + const fixAntimeridian = resolveFixAntimeridian(ctx, "s2"); + if (mode === "expand") { + const tokens = tokensFromS2Layer(layer.geojson, cellField); + const n = estimateS2ExpandCount(tokens, res); + const hardCap = hardCapFor("s2"); + if (!Number.isFinite(n) || n <= 0) { + ctx.log( + `No S2 cells to expand. Check the cell ID field and that cells are coarser than resolution ${res}.`, + ); + return; + } + if (n > hardCap) { + ctx.log( + `Error: expanding to resolution ${res} would generate ${Math.round(n).toLocaleString()} cells (cap ${hardCap.toLocaleString()}). Choose a coarser target.`, + ); + return; + } + emitResult( + expandS2FeatureCollection(layer.geojson, res, { + cellField, + unwrap: fixAntimeridian, + }), + ); + } else { + emitResult( + compactS2FeatureCollection(layer.geojson, { + cellField, + unwrap: fixAntimeridian, + }), + ); + } + } catch (error) { + const message = error instanceof Error ? error.message : String(error); + ctx.log(`Error: ${message}`); + } + return; + } + + // DGGAL is client-side WASM; no DuckDB. + if (type === "dggal") { + try { + await withDggalDggrs(dggalType, (engine) => { + if (mode === "expand") { + const tokens = tokensFromDggalLayer(layer.geojson!, cellField); + const n = estimateDggalExpandCount(engine, tokens, res); + const hardCap = hardCapFor("dggal"); + if (!Number.isFinite(n) || n <= 0) { + ctx.log( + `No ${label} cells to expand. Check the cell ID field, DGGAL type, and that cells are coarser than resolution ${res}.`, + ); + return; + } + if (n > hardCap) { + ctx.log( + `Error: expanding to resolution ${res} would generate ${Math.round(n).toLocaleString()} cells (cap ${hardCap.toLocaleString()}). Choose a coarser target.`, + ); + return; + } + emitResult(expandDggalFeatureCollection(engine, layer.geojson!, res, { cellField })); + } else { + emitResult(compactDggalFeatureCollection(engine, layer.geojson!, { cellField })); + } + }); + } catch (error) { + const message = error instanceof Error ? error.message : String(error); + ctx.log(`Error: ${message}`); + } + return; + } + + const duckdb = requireDuckDb(ctx); + let registered: DuckDbGeoJsonSource | null = null; + try { + await duckdb.ensureExtensions(["spatial", extensionForDggs(type)!]); + registered = await duckdb.registerGeoJson(layer.geojson); + if (mode === "expand") { + const countSql = + type === "a5" + ? buildA5ExpandCountSql(registered.sql, res, cellField) + : buildH3ExpandCountSql(registered.sql, res, cellField); + const countRows = await duckdb.query(countSql); + const n = Number(countRows[0]?.n ?? 0); + const hardCap = hardCapFor(type); + if (!Number.isFinite(n) || n <= 0) { + ctx.log( + `No ${label} cells to expand. Check the cell ID field and that cells are coarser than resolution ${res}.`, + ); + return; + } + if (n > hardCap) { + ctx.log( + `Error: expanding to resolution ${res} would generate ${Math.round(n).toLocaleString()} cells (cap ${hardCap.toLocaleString()}). Choose a coarser target.`, + ); + return; + } + } + + const sql = + mode === "compact" + ? type === "a5" + ? buildA5CompactSql(registered.sql, cellField) + : buildH3CompactSql(registered.sql, cellField) + : type === "a5" + ? buildA5ExpandSql(registered.sql, res, cellField) + : buildH3ExpandSql(registered.sql, res, cellField); + + const rows = await duckdb.query(sql); + const fc = + type === "a5" + ? a5RowsToFeatureCollection(rows) + : rowsToFeatureCollection(rows, resolveFixAntimeridian(ctx, type)); + emitResult(fc); + } catch (error) { + const message = error instanceof Error ? error.message : String(error); + ctx.log(`Error: ${message}`); + } finally { + await registered?.release(); + } + }, +}; + +export const DGGS_TOOLS: ProcessingAlgorithm[] = [ + createDggsGridTool, + dggsBinPointsTool, + dggsCompactTool, +]; + +export function getDggsTool(id: string): ProcessingAlgorithm | undefined { + return DGGS_TOOLS.find((tool) => tool.id === id); +} diff --git a/packages/processing/src/h3-tools.ts b/packages/processing/src/h3-tools.ts index d3deeb129..d93da2028 100644 --- a/packages/processing/src/h3-tools.ts +++ b/packages/processing/src/h3-tools.ts @@ -1,12 +1,6 @@ import type { FeatureCollection, Geometry } from "geojson"; -import bbox from "@turf/bbox"; -import type { GeoLibreLayer } from "@geolibre/core"; -import type { - DuckDbCapability, - DuckDbGeoJsonSource, - ProcessingAlgorithm, - ProcessingContext, -} from "./types"; +import { unwrapAntimeridianGeometry } from "./antimeridian"; +import type { ProcessingAlgorithm } from "./types"; /** Average area (km^2) of an H3 cell at each resolution 0..15 (official values). */ export const H3_AVG_AREA_KM2: number[] = [ @@ -63,11 +57,11 @@ export function suggestResolution( return 0; } -function sqlStr(value: string): string { +export function sqlStr(value: string): string { return `'${value.replaceAll("'", "''")}'`; } -function sqlIdent(value: string): string { +export function sqlIdent(value: string): string { return `"${value.replaceAll('"', '""')}"`; } @@ -77,33 +71,107 @@ export function bboxToWktPolygon(bbox: [number, number, number, number]): string return `POLYGON((${w} ${s}, ${e} ${s}, ${e} ${n}, ${w} ${n}, ${w} ${s}))`; } +/** + * Clamp a lon/lat bbox into WGS84. Spans wider than 180° of longitude (typical + * zoomed-out MapLibre viewports / world copies) collapse to [-180, 180]. + * Edges outside ±180° also collapse to the full-width path so clipping does + * not silently drop coverage (e.g. west=-190, east=-10). + * Callers that polyfill must still split that full-width ring — DuckDB H3 + * treats `POLYGON((-180 … 180 …))` as a dateline sliver (~tens of cells). + */ +export function normalizeLonLatBbox( + bbox: [number, number, number, number], +): [number, number, number, number] { + let [west, south, east, north] = bbox; + south = Math.max(-90, Math.min(90, south)); + north = Math.max(-90, Math.min(90, north)); + if (south > north) [south, north] = [north, south]; + + let lonSpan = east - west; + if (lonSpan < 0) lonSpan += 360; + if (lonSpan > 180 || west < -180 || east > 180 || west > 180 || east < -180) { + return [-180, south, 180, north]; + } + + west = Math.max(-180, Math.min(180, west)); + east = Math.max(-180, Math.min(180, east)); + return [west, south, east, north]; +} + const GRID_SELECT = "SELECT h3_h3_to_string(cell) AS h3, " + "ST_AsGeoJSON(ST_GeomFromText(h3_cell_to_boundary_wkt(cell))) AS geojson FROM cells"; -/** Grid SQL from a polygon WKT literal (used for bbox / viewport sources). */ -export function buildGridFromWktSql(wkt: string, res: number): string { +/** + * Prefer the experimental polyfill with overlap containment. The legacy + * `h3_polygon_wkt_to_cells` (center containment) returns an empty list when the + * polygon is smaller than a cell, and also for near-global ±180° rings that + * H3's center algorithm mishandles. + */ +function polyfillUnnest(wktSql: string, res: number): string { + return `unnest(h3_polygon_wkt_to_cells_experimental(${wktSql}, ${res}, 'overlap'))`; +} + +/** Wrap a `SELECT … AS cell` query; optionally `h3_compact_cells` the result. */ +function finalizeH3Cells(rawSelect: string, compact: boolean): string { + if (!compact) { + return `WITH cells AS (${rawSelect}) ` + GRID_SELECT; + } return ( - `WITH cells AS (SELECT unnest(h3_polygon_wkt_to_cells(${sqlStr(wkt)}, ${res})) AS cell) ` + + `WITH raw AS (${rawSelect}), ` + + `arr AS (SELECT list(cell) AS cells FROM raw), ` + + `cells AS (SELECT unnest(h3_compact_cells(cells)) AS cell FROM arr) ` + GRID_SELECT ); } +/** Grid SQL from a polygon WKT literal (used for bbox / viewport sources). */ +export function buildGridFromWktSql(wkt: string, res: number, compact = false): string { + return finalizeH3Cells(`SELECT ${polyfillUnnest(sqlStr(wkt), res)} AS cell`, compact); +} + +/** + * Grid SQL from a lon/lat bbox. After {@link normalizeLonLatBbox}, a whole-world + * request is [-180, s, 180, n]; a single ring polyfills to a dateline sliver, so + * that case is split into western/eastern hemispheres and unioned. + */ +export function buildGridFromBboxSql( + bbox: [number, number, number, number], + res: number, + compact = false, +): string { + const [w, s, e, n] = normalizeLonLatBbox(bbox); + if (w === -180 && e === 180) { + const left = bboxToWktPolygon([-180, s, 0, n]); + const right = bboxToWktPolygon([0, s, 180, n]); + return finalizeH3Cells( + `SELECT DISTINCT cell FROM (` + + `SELECT ${polyfillUnnest(sqlStr(left), res)} AS cell ` + + `UNION ALL ` + + `SELECT ${polyfillUnnest(sqlStr(right), res)} AS cell` + + `)`, + compact, + ); + } + return buildGridFromWktSql(bboxToWktPolygon([w, s, e, n]), res, compact); +} + /** * Grid SQL that unions all geometry from a registered source into one * (multi)polygon and fills it (used for the polyfill source). `sourceSql` is a * FROM-able expression whose geometry column is `geom` (DuckDB `ST_Read`). */ -export function buildGridFromSourceSql(sourceSql: string, res: number): string { +export function buildGridFromSourceSql(sourceSql: string, res: number, compact = false): string { // Union only polygonal geometries: a mixed layer would otherwise aggregate to - // a GEOMETRYCOLLECTION that `h3_polygon_wkt_to_cells` rejects. The `cells` CTE + // a GEOMETRYCOLLECTION that the H3 polyfill rejects. The `cells` CTE // filters a NULL union result (no polygons survived) so a NULL WKT never // reaches the h3 function, which can throw on NULL. - return ( - `WITH merged AS (SELECT ST_AsText(ST_Union_Agg(geom)) AS wkt FROM ${sourceSql} ` + - `WHERE geom IS NOT NULL AND ST_GeometryType(geom) IN ('POLYGON', 'MULTIPOLYGON')), ` + - `cells AS (SELECT unnest(h3_polygon_wkt_to_cells(wkt, ${res})) AS cell FROM merged WHERE wkt IS NOT NULL) ` + - GRID_SELECT + return finalizeH3Cells( + `SELECT ${polyfillUnnest("wkt", res)} AS cell FROM (` + + `SELECT ST_AsText(ST_Union_Agg(geom)) AS wkt FROM ${sourceSql} ` + + `WHERE geom IS NOT NULL AND ST_GeometryType(geom) IN ('POLYGON', 'MULTIPOLYGON')` + + `) merged WHERE wkt IS NOT NULL`, + compact, ); } @@ -144,8 +212,53 @@ export function buildBinSql(sourceSql: string, res: number, op: H3AggOp, field?: ); } +/** + * Collect H3 cell IDs from `cellField` on `sourceSql` into an array. + * `sourceSql` is a FROM-able expression (DuckDB `ST_Read`) with that column. + */ +function h3CellArrayCte(sourceSql: string, cellField: string): string { + const f = sqlIdent(cellField); + return ( + `input AS (SELECT DISTINCT h3_string_to_h3(CAST(${f} AS VARCHAR)) AS cell FROM ${sourceSql} ` + + `WHERE ${f} IS NOT NULL AND CAST(${f} AS VARCHAR) <> ''), ` + + `arr AS (SELECT list(cell) AS cells FROM input)` + ); +} + +/** Compact H3 cells from a polygon cell layer (IDs in `cellField`, default `h3`). */ +export function buildH3CompactSql(sourceSql: string, cellField = "h3"): string { + return ( + `WITH ${h3CellArrayCte(sourceSql, cellField)}, ` + + `cells AS (SELECT unnest(h3_compact_cells(cells)) AS cell FROM arr) ` + + GRID_SELECT + ); +} + +/** + * Expand (uncompact) H3 cells to a uniform `res`. Cells already finer than + * `res` are rejected by the H3 extension. + */ +export function buildH3ExpandSql(sourceSql: string, res: number, cellField = "h3"): string { + return ( + `WITH ${h3CellArrayCte(sourceSql, cellField)}, ` + + `cells AS (SELECT unnest(h3_uncompact_cells(cells, ${res})) AS cell FROM arr) ` + + GRID_SELECT + ); +} + +/** Count of cells that {@link buildH3ExpandSql} would emit (for the hard-cap guard). */ +export function buildH3ExpandCountSql(sourceSql: string, res: number, cellField = "h3"): string { + return ( + `WITH ${h3CellArrayCte(sourceSql, cellField)} ` + + `SELECT coalesce(len(h3_uncompact_cells(cells, ${res})), 0) AS n FROM arr` + ); +} + /** Build a FeatureCollection from rows carrying `h3`, optional `count`/`value`, and `geojson`. */ -export function rowsToFeatureCollection(rows: Record[]): FeatureCollection { +export function rowsToFeatureCollection( + rows: Record[], + fixAntimeridian = true, +): FeatureCollection { const features = []; for (const row of rows) { const raw = row.geojson; @@ -165,326 +278,16 @@ export function rowsToFeatureCollection(rows: Record[]): Featur if (row.value !== undefined && row.value !== null) { properties.value = Number(row.value); } - features.push({ type: "Feature" as const, geometry, properties }); + features.push({ + type: "Feature" as const, + geometry: fixAntimeridian ? unwrapAntimeridianGeometry(geometry) : geometry, + properties, + }); } return { type: "FeatureCollection", features }; } -const NO_DUCKDB = "This tool requires DuckDB-WASM, which is unavailable in this environment."; - -function requireDuckDb(ctx: ProcessingContext): DuckDbCapability { - if (!ctx.duckdb) throw new Error(NO_DUCKDB); - return ctx.duckdb; -} - -// Mirrors the same helper in vector-tools.ts and registry.ts; intentionally -// duplicated because vector-tools.ts imports from this file, so importing the -// other direction would create a cycle. Keep the three copies in sync. -function getLayer(ctx: ProcessingContext, paramId = "layer"): GeoLibreLayer | undefined { - const id = ctx.parameters[paramId] as string | undefined; - return ctx.layers.find((l) => l.id === id); -} - -/** Read a numeric parameter, returning NaN when missing or non-numeric. */ -function numberParam(ctx: ProcessingContext, id: string): number { - const raw = ctx.parameters[id]; - if (raw === undefined || raw === null || raw === "") return NaN; - return typeof raw === "string" ? Number(raw) : (raw as number); -} - -/** - * Read and validate the manual [west, south, east, north] bbox parameters. - * Logs a clear error and returns null when any value is missing or the box is - * degenerate (west >= east or south >= north). - */ -function bboxFromParams(ctx: ProcessingContext): [number, number, number, number] | null { - const west = numberParam(ctx, "west"); - const south = numberParam(ctx, "south"); - const east = numberParam(ctx, "east"); - const north = numberParam(ctx, "north"); - if ([west, south, east, north].some((n) => !Number.isFinite(n))) { - ctx.log("Error: enter numeric west, south, east, and north values"); - return null; - } - if (west >= east || south >= north) { - ctx.log("Error: bounding box must have west < east and south < north"); - return null; - } - return [west, south, east, north]; -} - -/** Parse the `resolution` param, or auto-suggest from area. Logs + returns null on bad input. */ -function resolveResolution(ctx: ProcessingContext, areaKm2: number): number | null { - const raw = ctx.parameters.resolution; - if (raw === undefined || raw === null || raw === "") { - const suggested = suggestResolution(areaKm2); - ctx.log(`Using suggested resolution ${suggested}`); - return suggested; - } - const res = typeof raw === "string" ? Number(raw) : (raw as number); - if (!Number.isInteger(res) || res < 0 || res > 15) { - ctx.log("Error: resolution must be an integer from 0 to 15"); - return null; - } - return res; -} - -export const createH3GridTool: ProcessingAlgorithm = { - id: "h3-grid", - name: "Create H3 grid", - description: - "Fill an area with H3 hexagons (DuckDB h3 extension). Source: a layer's geometry, a layer's extent, the current map view, or a manual bounding box.", - group: "H3", - parameters: [ - { - id: "source", - label: "Area source", - type: "select", - default: "polyfill", - options: [ - { value: "polyfill", label: "Layer geometry (polyfill)" }, - { value: "extent", label: "Layer extent (bbox)" }, - { value: "viewport", label: "Map viewport" }, - { value: "bbox", label: "Manual bounding box" }, - ], - }, - { - id: "layer", - label: "Input layer", - type: "layer", - required: true, - // No geometry filter: "extent" fills any layer's bounding box, while - // "polyfill" needs polygons (validated at run time below). The layer is - // only required for the layer-based sources, so it stays hidden (and - // skips required validation) for the viewport and bbox sources. - visibleWhen: { param: "source", in: ["polyfill", "extent"] }, - }, - { - id: "west", - label: "West (min lon)", - type: "number", - required: true, - min: -180, - max: 180, - visibleWhen: { param: "source", in: ["bbox"] }, - }, - { - id: "south", - label: "South (min lat)", - type: "number", - required: true, - min: -90, - max: 90, - visibleWhen: { param: "source", in: ["bbox"] }, - }, - { - id: "east", - label: "East (max lon)", - type: "number", - required: true, - min: -180, - max: 180, - visibleWhen: { param: "source", in: ["bbox"] }, - }, - { - id: "north", - label: "North (max lat)", - type: "number", - required: true, - min: -90, - max: 90, - visibleWhen: { param: "source", in: ["bbox"] }, - }, - { - id: "resolution", - label: "Resolution (0-15)", - type: "number", - min: 0, - max: 15, - step: 1, - description: "Leave blank to auto-pick from the area.", - }, - ], - run: async (ctx) => { - const duckdb = requireDuckDb(ctx); - const source = (ctx.parameters.source as string) || "polyfill"; - - let areaKm2: number; - let wkt: string | null = null; - let inputGeojson: FeatureCollection | null = null; - if (source === "viewport") { - const bounds = ctx.viewportBounds?.(); - if (!bounds) { - ctx.log("Error: map viewport is unavailable"); - return; - } - if (bounds[0] >= bounds[2]) { - // west >= east means the viewport wraps the antimeridian; the rectangle - // WKT would self-cross and fill the wrong (340deg) span. Bail with a - // clear message rather than producing wrong cells. - ctx.log( - "Error: the map view crosses the antimeridian; pan so it doesn't wrap +/-180, or use a manual bounding box", - ); - return; - } - areaKm2 = bboxAreaKm2(bounds); - wkt = bboxToWktPolygon(bounds); - } else if (source === "bbox") { - const bounds = bboxFromParams(ctx); - if (!bounds) return; - areaKm2 = bboxAreaKm2(bounds); - wkt = bboxToWktPolygon(bounds); - } else { - const layer = getLayer(ctx, "layer"); - if (!layer?.geojson?.features?.length) { - ctx.log('Error: parameter "layer" has no GeoJSON features'); - return; - } - if (source === "polyfill") { - const hasPolygon = layer.geojson.features.some( - (f) => f.geometry?.type === "Polygon" || f.geometry?.type === "MultiPolygon", - ); - if (!hasPolygon) { - ctx.log( - 'Error: polyfill needs a polygon layer; use the "Layer extent" source for point or line layers', - ); - return; - } - } - inputGeojson = layer.geojson; - const bb = bbox(layer.geojson) as [number, number, number, number]; - areaKm2 = bboxAreaKm2(bb); - if (source === "extent") wkt = bboxToWktPolygon(bb); - } - - const res = resolveResolution(ctx, areaKm2); - if (res === null) return; - - const estimate = estimateCellCount(areaKm2, res); - if (estimate > H3_HARD_CAP) { - ctx.log( - `Error: resolution ${res} would generate about ${Math.round( - estimate, - ).toLocaleString()} cells (cap ${H3_HARD_CAP.toLocaleString()}). Choose a coarser resolution.`, - ); - return; - } - - await duckdb.ensureExtensions(["spatial", "h3"]); - let registered: DuckDbGeoJsonSource | null = null; - try { - let sql: string; - if (wkt) { - sql = buildGridFromWktSql(wkt, res); - } else { - registered = await duckdb.registerGeoJson(inputGeojson!); // non-null: polyfill path only runs after the layer guard above set inputGeojson - sql = buildGridFromSourceSql(registered.sql, res); - } - const rows = await duckdb.query(sql); - const fc = rowsToFeatureCollection(rows); - if (fc.features.length === 0) { - ctx.log( - `No H3 cells were produced at resolution ${res}. Try a coarser resolution or a larger area.`, - ); - return; - } - ctx.log(`Created ${fc.features.length} H3 cell(s) at resolution ${res}`); - ctx.addResultLayer?.(`H3 grid (res ${res})`, fc); - } finally { - await registered?.release(); - } - }, -}; - -export const binPointsTool: ProcessingAlgorithm = { - id: "h3-bin-points", - name: "Bin points to H3", - description: - "Aggregate a point layer into H3 cells (count, or sum/mean/min/max of a numeric field).", - group: "H3", - parameters: [ - { - id: "layer", - label: "Input point layer", - type: "layer", - required: true, - geometryFilter: ["point"], - }, - { - id: "aggOp", - label: "Aggregate", - type: "select", - default: "count", - options: [ - { value: "count", label: "Count" }, - { value: "sum", label: "Sum" }, - { value: "mean", label: "Mean" }, - { value: "min", label: "Min" }, - { value: "max", label: "Max" }, - ], - }, - { - id: "field", - label: "Field", - type: "field", - fieldSource: "layer", - required: true, - visibleWhen: { param: "aggOp", notIn: ["count"] }, - description: "Numeric field to aggregate.", - }, - { - id: "resolution", - label: "Resolution (0-15)", - type: "number", - min: 0, - max: 15, - step: 1, - description: "Leave blank to auto-pick from the area.", - }, - ], - run: async (ctx) => { - const duckdb = requireDuckDb(ctx); - const layer = getLayer(ctx, "layer"); - if (!layer?.geojson?.features?.length) { - ctx.log('Error: parameter "layer" has no GeoJSON features'); - return; - } - const op = (ctx.parameters.aggOp as string) || "count"; - if (!H3_AGG_OPS.includes(op as H3AggOp)) { - ctx.log(`Error: unknown aggregate "${op}"`); - return; - } - const field = ctx.parameters.field as string | undefined; - if (op !== "count" && !field) { - ctx.log(`Error: select a numeric field to ${op}`); - return; - } - - const bb = bbox(layer.geojson) as [number, number, number, number]; - const res = resolveResolution(ctx, bboxAreaKm2(bb)); - if (res === null) return; - - await duckdb.ensureExtensions(["spatial", "h3"]); - const registered = await duckdb.registerGeoJson(layer.geojson); - try { - const sql = buildBinSql(registered.sql, res, op as H3AggOp, field); - const rows = await duckdb.query(sql); - const fc = rowsToFeatureCollection(rows); - if (fc.features.length === 0) { - ctx.log( - `No points fell into H3 cells at resolution ${res}. Check the layer has point geometries.`, - ); - return; - } - ctx.log(`Binned points into ${fc.features.length} H3 cell(s) at resolution ${res}`); - ctx.addResultLayer?.(`H3 bins (res ${res})`, fc); - } finally { - await registered.release(); - } - }, -}; - -export const H3_TOOLS: ProcessingAlgorithm[] = [createH3GridTool, binPointsTool]; +export const H3_TOOLS: ProcessingAlgorithm[] = []; export function getH3Tool(id: string): ProcessingAlgorithm | undefined { return H3_TOOLS.find((tool) => tool.id === id); diff --git a/packages/processing/src/index.ts b/packages/processing/src/index.ts index d06afd599..c2b17e27a 100644 --- a/packages/processing/src/index.ts +++ b/packages/processing/src/index.ts @@ -8,6 +8,7 @@ export { export { VECTOR_TOOLS, getVectorTool, + resolveVectorRerun, matchFeaturesByLocation, MAX_CLIENT_PAIRS, SELECT_LOCATION_PREDICATES, @@ -45,7 +46,89 @@ export { emergingHotSpotTool, emergingPattern, } from "./statistics-tools"; -export { H3_TOOLS, getH3Tool, createH3GridTool, binPointsTool } from "./h3-tools"; +export { + H3_TOOLS, + getH3Tool, + buildBinSql, + buildGridFromBboxSql, + buildH3CompactSql, + buildH3ExpandSql, + buildH3ExpandCountSql, + H3_AGG_OPS, + normalizeLonLatBbox, + type H3AggOp, +} from "./h3-tools"; +export { + buildA5GridFromWktSql, + buildA5GridFromBboxSql, + buildA5GridFromSourceSql, + buildA5BinSql, + buildA5CompactSql, + buildA5ExpandSql, + buildA5ExpandCountSql, + a5RowsToFeatureCollection, + suggestA5Resolution, + estimateA5CellCount, + A5_MAX_TOOL_RES, +} from "./a5-tools"; +export { unwrapAntimeridianGeometry, unwrapAntimeridianRing } from "./antimeridian"; +export { + buildDggridGridFromWktSql, + buildDggridGridFromSourceSql, + buildDggridBinSql, + dggridRowsToFeatureCollection, + suggestDggridResolution, + estimateDggridCellCount, + resolveDggridGridType, + maxResolutionForDggrid, + DGGRID_MAX_TOOL_RES, + DGGRID_GRID_TYPES, + DGGRID_GRID_TYPE_OPTIONS, + DGGRID_GRID_SPECS, + DEFAULT_DGGRID_GRID_TYPE, + type DggridGridType, +} from "./dggrid-tools"; +export { + DGGS_TOOLS, + getDggsTool, + createDggsGridTool, + dggsBinPointsTool, + dggsCompactTool, + maxResolutionForDggs, + extensionForDggs, + type DggsType, +} from "./dggs-tools"; +export { + s2GridFromBbox, + s2GridFromFeatureCollection, + binPointsToS2, + compactS2Tokens, + expandS2Tokens, + compactS2FeatureCollection, + expandS2FeatureCollection, + suggestS2Resolution, + estimateS2CellCount, + S2_MAX_TOOL_RES, +} from "./s2-tools"; +export { + DGGAL_TYPES, + DGGAL_GRID_TYPES, + DGGAL_GRID_TYPE_OPTIONS, + DEFAULT_DGGAL_GRID_TYPE, + DGGAL_MAX_TOOL_RES, + resolveDggalGridType, + maxResolutionForDggal, + suggestDggalResolution, + estimateDggalCellCount, + dggalGridFromBbox, + dggalGridFromFeatureCollection, + binPointsToDggal, + compactDggalTokens, + expandDggalTokens, + compactDggalFeatureCollection, + expandDggalFeatureCollection, + type DggalGridType, +} from "./dggal-tools"; export { RASTER_TOOLS, getRasterTool, diff --git a/packages/processing/src/s2-tools.ts b/packages/processing/src/s2-tools.ts new file mode 100644 index 000000000..713520558 --- /dev/null +++ b/packages/processing/src/s2-tools.ts @@ -0,0 +1,371 @@ +import type { Feature, FeatureCollection, Geometry, Polygon, Position } from "geojson"; +import { geojson as s2geojson, s1, s2 } from "s2js"; +import { unwrapAntimeridianGeometry } from "./antimeridian"; + +/** + * Approximate average S2 cell area (km²) at levels 0..30. + * Six level-0 faces, each subdividing 4× per level (earth ≈ 5.101×10⁸ km²). + */ +export const S2_AVG_AREA_KM2: number[] = Array.from({ length: 31 }, (_, level) => { + return 510_065_621.724 / (6 * 4 ** level); +}); + +/** Soft target used when auto-suggesting a resolution. */ +export const S2_TARGET_CELLS = 10_000; +/** Finest resolution the auto-suggester will pick. */ +export const S2_MAX_SUGGESTED_RES = 12; +/** Hard ceiling: a grid larger than this aborts rather than running away. */ +export const S2_HARD_CAP = 200_000; +/** Max S2 level offered in the processing dialog (0–30). */ +export const S2_MAX_TOOL_RES = 30; + +/** + * Max longitude span (degrees) per RegionCoverer call. Wider rings are + * ambiguous in GeoJSON / s2js, so bounds are chunked and cells deduplicated. + */ +const MAX_COVER_SPAN_DEGREES = 120; + +/** Estimated number of S2 cells covering `areaKm2` at `res`. */ +export function estimateS2CellCount(areaKm2: number, res: number): number { + const cellArea = S2_AVG_AREA_KM2[res]; + if (cellArea === undefined) return Number.POSITIVE_INFINITY; + return areaKm2 / cellArea; +} + +/** Finest resolution whose estimated cell count stays <= the target. */ +export function suggestS2Resolution( + areaKm2: number, + targetCells = S2_TARGET_CELLS, + maxRes = S2_MAX_SUGGESTED_RES, +): number { + const capped = Math.min(maxRes, S2_MAX_TOOL_RES); + for (let res = capped; res >= 0; res -= 1) { + if (estimateS2CellCount(areaKm2, res) <= targetCells) return res; + } + return 0; +} + +function cellIdFromToken(token: string): bigint { + return s2.cellid.fromToken(token); +} + +function cellCenter(id: bigint): [number, number] { + const latLng = s2.cellid.latLng(id); + return [s1.angle.degrees(latLng.lng), s1.angle.degrees(latLng.lat)]; +} + +/** S2 cell token at `level` containing lon/lat (lon first, matching GeoJSON). */ +export function s2CellAtLonLat(lng: number, lat: number, level: number): string { + const leaf = s2.cellid.fromLatLng(s2.LatLng.fromDegrees(lat, lng)); + return s2.cellid.toToken(s2.cellid.parent(leaf, level)); +} + +/** + * Four corners as a closed lon/lat ring. When `unwrap` is true, vertices are + * shifted relative to the first so dateline-straddling cells stay contiguous + * for MapLibre (same idea as {@link unwrapAntimeridianGeometry}). + */ +function cellRing(id: bigint, unwrap: boolean): [number, number][] { + const cell = s2.Cell.fromCellID(id); + const ring: [number, number][] = []; + for (let i = 0; i <= 4; i += 1) { + const vertex = s2.LatLng.fromPoint(cell.vertex(i % 4)); + let lng = s1.angle.degrees(vertex.lng); + const lat = s1.angle.degrees(vertex.lat); + if (unwrap && ring.length > 0) { + const reference = ring[0]![0]!; + if (lng - reference > 180) lng -= 360; + if (lng - reference < -180) lng += 360; + } + ring.push([lng, lat]); + } + return ring; +} + +/** Convert an S2 token to a GeoJSON polygon feature. */ +export function s2CellFeature(token: string, unwrap = true): Feature { + const id = cellIdFromToken(token); + const [lng, lat] = cellCenter(id); + let geometry: Geometry = { type: "Polygon", coordinates: [cellRing(id, unwrap)] }; + // Belt-and-braces: if the ring was built without per-vertex unwrap, still + // offer the shared antimeridian helper when the caller asked for fix-on. + if (unwrap) geometry = unwrapAntimeridianGeometry(geometry); + return { + type: "Feature", + id: token, + properties: { + s2: token, + resolution: s2.cellid.level(id), + center_lat: lat, + center_lng: lng, + }, + geometry: geometry as Polygon, + }; +} + +/** Longitude chunks in [-180, 180] covering `west`+`span` (handles wrap). */ +function lonChunks(west: number, span: number): Array<[number, number]> { + const chunks: Array<[number, number]> = []; + let cursor = (((west % 360) + 540) % 360) - 180; + let remaining = Math.min(360, Math.max(0, span)); + while (remaining > 1e-9) { + const step = Math.min(remaining, MAX_COVER_SPAN_DEGREES, 180 - cursor); + chunks.push([cursor, cursor + step]); + cursor = cursor + step >= 180 ? -180 : cursor + step; + remaining -= step; + } + return chunks; +} + +function coverPolygonTokens( + polygon: Polygon, + level: number, + limit: number, + into: Set, +): void { + const coverer = new s2geojson.RegionCoverer({ minLevel: level, maxLevel: level }); + for (const id of coverer.covering(polygon)) { + into.add(s2.cellid.toToken(id)); + if (into.size > limit) { + throw new RangeError(`S2 cell limit exceeded: ${limit}`); + } + } +} + +/** Build polygon features from S2 tokens (property field `s2`). */ +export function s2TokensToFeatureCollection( + tokens: Iterable, + unwrap = true, +): FeatureCollection { + return { + type: "FeatureCollection", + features: [...tokens].map((token) => s2CellFeature(token, unwrap)), + }; +} + +/** + * Fill a WGS84 bounding box with S2 cells at one level. Mirrors the maplibre-s2 + * plugin covering (chunked longitude, token dedupe). + */ +export function s2GridFromBbox( + bounds: [number, number, number, number], + level: number, + options: { limit?: number; unwrap?: boolean; compact?: boolean } = {}, +): FeatureCollection { + const limit = options.limit ?? S2_HARD_CAP; + const unwrap = options.unwrap !== false; + const [west, southRaw, east, northRaw] = bounds; + const south = Math.max(-89.999999, Math.min(89.999999, southRaw)); + const north = Math.max(-89.999999, Math.min(89.999999, northRaw)); + const span = Math.min(360, east >= west ? east - west : east + 360 - west); + + const cells = new Set(); + for (const [left, right] of lonChunks(west, span)) { + const polygon: Polygon = { + type: "Polygon", + coordinates: [ + [ + [left, south], + [right, south], + [right, north], + [left, north], + [left, south], + ], + ], + }; + coverPolygonTokens(polygon, level, limit, cells); + } + const tokens = options.compact ? compactS2Tokens(cells) : [...cells]; + return s2TokensToFeatureCollection(tokens, unwrap); +} + +/** + * Polyfill polygon / multipolygon features with S2 cells at `level`. + * Non-polygonal geometries are skipped. + */ +export function s2GridFromFeatureCollection( + fc: FeatureCollection, + level: number, + options: { limit?: number; unwrap?: boolean; compact?: boolean } = {}, +): FeatureCollection { + const limit = options.limit ?? S2_HARD_CAP; + const unwrap = options.unwrap !== false; + const cells = new Set(); + for (const feature of fc.features) { + const g = feature.geometry; + if (!g) continue; + if (g.type === "Polygon") { + coverPolygonTokens(g, level, limit, cells); + } else if (g.type === "MultiPolygon") { + for (const coords of g.coordinates) { + coverPolygonTokens({ type: "Polygon", coordinates: coords }, level, limit, cells); + } + } + } + const tokens = options.compact ? compactS2Tokens(cells) : [...cells]; + return s2TokensToFeatureCollection(tokens, unwrap); +} + +/** Collect S2 cell tokens from a feature property (default `s2`). */ +export function tokensFromS2Layer(fc: FeatureCollection, cellField = "s2"): string[] { + const out: string[] = []; + for (const feature of fc.features) { + const raw = feature.properties?.[cellField]; + if (raw === undefined || raw === null) continue; + const token = String(raw).trim(); + if (token) out.push(token); + } + return out; +} + +/** + * Compact S2 tokens with {@link s2.CellUnion.normalize}: complete sets of four + * siblings become their parent (mixed levels). + */ +export function compactS2Tokens(tokens: Iterable): string[] { + const ids: bigint[] = []; + for (const token of tokens) { + ids.push(cellIdFromToken(token)); + } + if (ids.length === 0) return []; + const union = new s2.CellUnion(...ids); + union.normalize(); + const out: string[] = []; + for (let i = 0; i < union.length; i += 1) { + out.push(s2.cellid.toToken(union[i]!)); + } + return out; +} + +/** + * Expand (denormalize) S2 tokens to a uniform `level`. Cells already finer than + * `level` throw; use {@link estimateS2ExpandCount} for the hard-cap guard. + */ +export function expandS2Tokens(tokens: Iterable, level: number): string[] { + const ids: bigint[] = []; + for (const token of tokens) { + const id = cellIdFromToken(token); + if (s2.cellid.level(id) > level) { + throw new RangeError( + `S2 cell ${token} is finer than target level ${level}; choose a finer target or compact first`, + ); + } + ids.push(id); + } + if (ids.length === 0) return []; + const union = new s2.CellUnion(...ids); + union.normalize(); + union.denormalize(level, 1); + const out: string[] = []; + for (let i = 0; i < union.length; i += 1) { + out.push(s2.cellid.toToken(union[i]!)); + } + return out; +} + +/** Exact cell count after expanding `tokens` to `level` (4× per level). */ +export function estimateS2ExpandCount(tokens: Iterable, level: number): number { + let n = 0; + for (const token of tokens) { + const id = cellIdFromToken(token); + const L = s2.cellid.level(id); + if (L > level) return 0; + n += 4 ** (level - L); + } + return n; +} + +/** Compact a polygon cell layer's `s2` (or other) ID field. */ +export function compactS2FeatureCollection( + fc: FeatureCollection, + options: { cellField?: string; unwrap?: boolean } = {}, +): FeatureCollection { + const unwrap = options.unwrap !== false; + const tokens = compactS2Tokens(tokensFromS2Layer(fc, options.cellField ?? "s2")); + return s2TokensToFeatureCollection(tokens, unwrap); +} + +/** Expand a polygon cell layer to a uniform S2 level. */ +export function expandS2FeatureCollection( + fc: FeatureCollection, + level: number, + options: { cellField?: string; unwrap?: boolean } = {}, +): FeatureCollection { + const unwrap = options.unwrap !== false; + const tokens = expandS2Tokens(tokensFromS2Layer(fc, options.cellField ?? "s2"), level); + return s2TokensToFeatureCollection(tokens, unwrap); +} + +/** Supported point-binning aggregate operations (same set as H3). */ +export type S2AggOp = "count" | "sum" | "mean" | "min" | "max"; + +function eachPointCoord( + geometry: Geometry | null | undefined, + visit: (pos: Position) => void, +): void { + if (!geometry) return; + if (geometry.type === "Point") { + visit(geometry.coordinates); + return; + } + if (geometry.type === "MultiPoint") { + for (const c of geometry.coordinates) visit(c); + } +} + +/** + * Aggregate point geometry into S2 cells (client-side; no DuckDB). + */ +export function binPointsToS2( + fc: FeatureCollection, + level: number, + op: S2AggOp, + field?: string, + options: { unwrap?: boolean } = {}, +): FeatureCollection { + const unwrap = options.unwrap !== false; + type Acc = { count: number; sum: number; min: number; max: number }; + const byCell = new Map(); + + for (const feature of fc.features) { + eachPointCoord(feature.geometry, (pos) => { + const lng = pos[0]; + const lat = pos[1]; + if ( + lng === undefined || + lat === undefined || + !Number.isFinite(lng) || + !Number.isFinite(lat) + ) { + return; + } + const token = s2CellAtLonLat(lng, lat, level); + let acc = byCell.get(token); + if (!acc) { + acc = { count: 0, sum: 0, min: Number.POSITIVE_INFINITY, max: Number.NEGATIVE_INFINITY }; + byCell.set(token, acc); + } + acc.count += 1; + if (op !== "count" && field) { + const raw = feature.properties?.[field]; + const n = typeof raw === "number" ? raw : Number(raw); + if (!Number.isFinite(n)) return; + acc.sum += n; + if (n < acc.min) acc.min = n; + if (n > acc.max) acc.max = n; + } + }); + } + + const features: Feature[] = []; + for (const [token, acc] of byCell) { + const feature = s2CellFeature(token, unwrap); + const properties: Record = { ...feature.properties, count: acc.count }; + if (op === "sum") properties.value = acc.sum; + else if (op === "mean") properties.value = acc.count > 0 ? acc.sum / acc.count : 0; + else if (op === "min") properties.value = Number.isFinite(acc.min) ? acc.min : null; + else if (op === "max") properties.value = Number.isFinite(acc.max) ? acc.max : null; + features.push({ ...feature, properties }); + } + return { type: "FeatureCollection", features }; +} diff --git a/packages/processing/src/vector-tools.ts b/packages/processing/src/vector-tools.ts index e8ea016c9..e09bd7813 100644 --- a/packages/processing/src/vector-tools.ts +++ b/packages/processing/src/vector-tools.ts @@ -32,7 +32,7 @@ import type { } from "geojson"; import { layerJoinKey, type GeoLibreLayer } from "@geolibre/core"; import type { GeometryFamily, ProcessingAlgorithm, ProcessingContext } from "./types"; -import { createH3GridTool, binPointsTool } from "./h3-tools"; +import { createDggsGridTool, dggsBinPointsTool, dggsCompactTool } from "./dggs-tools"; import { TOPOLOGY_TOOLS } from "./topology-tools"; /** Upper bound on input×overlay pairs for the main-thread pairwise loops. */ @@ -2691,9 +2691,10 @@ export const VECTOR_TOOLS: ProcessingAlgorithm[] = [ gridTool, voronoiTool, cellSectorsTool, - createH3GridTool, - binPointsTool, - // Movement & time tools come after H3 so the dialog's group order (derived + createDggsGridTool, + dggsBinPointsTool, + dggsCompactTool, + // Movement & time tools come after DGGS so the dialog's group order (derived // from this array) matches the Processing → Vector menu order. trajectorySpeedTool, detectStopsTool, @@ -2705,3 +2706,32 @@ export const VECTOR_TOOLS: ProcessingAlgorithm[] = [ export function getVectorTool(id: string): ProcessingAlgorithm | undefined { return VECTOR_TOOLS.find((tool) => tool.id === id); } + +/** + * Old H3 processing tool IDs from history entries written before the DGGS + * rename. Map them onto the current tools and default `dggsType` to `"h3"`. + */ +const H3_VECTOR_TOOL_ALIASES: Readonly> = { + "h3-grid": "dggs-grid", + "h3-bin-points": "dggs-bin", +}; + +/** + * Resolve a vector History re-run's tool id (and parameters) for today's + * registry. Unknown ids pass through unchanged so the dialog can still report + * "tool unavailable". + */ +export function resolveVectorRerun( + toolId: string, + parameters: Record = {}, +): { toolId: string; parameters: Record } { + const mapped = H3_VECTOR_TOOL_ALIASES[toolId]; + if (!mapped) return { toolId, parameters }; + return { + toolId: mapped, + parameters: { + ...parameters, + dggsType: parameters.dggsType ?? "h3", + }, + }; +} diff --git a/tests/a5-tools.test.ts b/tests/a5-tools.test.ts new file mode 100644 index 000000000..54404cbea --- /dev/null +++ b/tests/a5-tools.test.ts @@ -0,0 +1,151 @@ +import assert from "node:assert/strict"; +import { describe, it } from "node:test"; +import { + A5_AVG_AREA_KM2, + A5_HARD_CAP, + A5_MAX_TOOL_RES, + a5RowsToFeatureCollection, + buildA5BinSql, + buildA5GridFromBboxSql, + buildA5GridFromSourceSql, + buildA5GridFromWktSql, + estimateA5CellCount, + suggestA5Resolution, +} from "../packages/processing/src/a5-tools"; +import { bboxAreaKm2 } from "../packages/processing/src/h3-tools"; + +describe("a5 resolution math", () => { + it("exposes 31 average-area entries (res 0..30), strictly decreasing", () => { + assert.equal(A5_AVG_AREA_KM2.length, 31); + assert.equal(A5_MAX_TOOL_RES, 30); + for (let r = 1; r < 31; r += 1) { + assert.ok(A5_AVG_AREA_KM2[r] < A5_AVG_AREA_KM2[r - 1]); + } + }); + + it("suggests a coarser resolution for larger areas", () => { + const big = bboxAreaKm2([-10, -10, 10, 10]); + const tiny = bboxAreaKm2([0, 0, 0.001, 0.001]); + const rBig = suggestA5Resolution(big); + const rTiny = suggestA5Resolution(tiny); + assert.ok(rBig < rTiny); + assert.ok(rTiny <= 12); + assert.ok(estimateA5CellCount(big, rBig) <= 10_000); + }); + + it("fails safe (Infinity) for an out-of-range resolution", () => { + const area = bboxAreaKm2([0, 0, 1, 1]); + assert.equal(estimateA5CellCount(area, 31), Number.POSITIVE_INFINITY); + assert.equal(estimateA5CellCount(area, -1), Number.POSITIVE_INFINITY); + assert.ok(estimateA5CellCount(area, 31) > A5_HARD_CAP); + // Res 16 is in range for A5 (0–30). + assert.ok(Number.isFinite(estimateA5CellCount(area, 16))); + }); +}); + +describe("a5 SQL builders", () => { + it("builds grid SQL from a WKT literal, escaping single quotes", () => { + const sql = buildA5GridFromWktSql("POLYGON((0 0, 1 0, 1 1, 0 0))'x", 7); + assert.match( + sql, + /a5_uncompact\(a5_geometry_to_cells\(ST_GeomFromText\('POLYGON\(\(0 0, 1 0, 1 1, 0 0\)\)''x'\), 7\), 7\)/, + ); + assert.match(sql, /a5_u64_to_hex\(cell\) AS a5/); + assert.match(sql, /a5_cell_to_geometry\(cell\)/); + }); + + it("enumerates full-longitude A5 bboxes from res0 cells", () => { + const narrow = buildA5GridFromBboxSql([0, 0, 1, 1], 5); + assert.doesNotMatch(narrow, /UNION ALL/); + assert.doesNotMatch(narrow, /a5_get_res0_cells/); + + const world = buildA5GridFromBboxSql([-180, -90, 180, 90], 4); + assert.match(world, /a5_uncompact\(a5_get_res0_cells\(\), 4\)/); + assert.doesNotMatch(world, /a5_geometry_to_cells/); + assert.doesNotMatch(world, /list_extract/); + + const band = buildA5GridFromBboxSql([-180, -60, 180, 60], 2); + assert.match(band, /a5_get_res0_cells\(\)/); + assert.match(band, /list_extract\(a5_cell_to_lonlat\(cell\), 2\) BETWEEN -60 AND 60/); + + const wide = buildA5GridFromBboxSql([-100, -10, 20, 10], 3); + // 120° span → two 60° strips (still under the full-lon res0 path) + assert.match(wide, /UNION ALL/); + assert.equal((wide.match(/a5_geometry_to_cells/g) ?? []).length, 2); + }); + + it("builds polyfill grid SQL that unions only polygon geometry", () => { + const sql = buildA5GridFromSourceSql("ST_Read('x.geojson')", 8); + assert.match(sql, /ST_Union_Agg\(geom\)/); + assert.match(sql, /'POLYGON', 'MULTIPOLYGON'/); + assert.match(sql, /a5_uncompact\(a5_geometry_to_cells\(g, 8\), 8\)/); + assert.doesNotMatch(sql, /a5_compact\(cells\)/); + }); + + it("optionally compacts A5 grid cells after polyfill", () => { + const sql = buildA5GridFromBboxSql([0, 0, 1, 1], 5, true); + assert.match(sql, /a5_compact\(cells\)/); + assert.match(sql, /a5_uncompact\(a5_geometry_to_cells/); + }); + + it("builds bin SQL with lon/lat cell lookup and optional aggregates", () => { + const countSql = buildA5BinSql("ST_Read('p.geojson')", 5, "count"); + assert.match(countSql, /a5_lonlat_to_cell\(ST_X\(pt\), ST_Y\(pt\), 5\)/); + assert.match(countSql, /a5_u64_to_hex\(cell\) AS a5/); + assert.doesNotMatch(countSql, / AS value/); + + const sumSql = buildA5BinSql("ST_Read('p.geojson')", 5, "sum", 'pop"x'); + assert.match(sumSql, /sum\(CAST\("pop""x" AS DOUBLE\)\) AS value/); + }); + + it("converts result rows to a FeatureCollection with a5 props", () => { + const fc = a5RowsToFeatureCollection([ + { + a5: "abc", + count: 2, + value: 4.5, + geojson: '{"type":"Polygon","coordinates":[[[0,0],[1,0],[1,1],[0,0]]]}', + }, + { + a5: "obj", + geojson: { + type: "Polygon", + coordinates: [ + [ + [0, 0], + [1, 0], + [1, 1], + [0, 0], + ], + ], + }, + }, + { a5: "skip", geojson: 12 }, + ]); + assert.equal(fc.features.length, 2); + assert.equal(fc.features[0].properties?.a5, "abc"); + assert.equal(fc.features[0].properties?.count, 2); + assert.equal(fc.features[0].properties?.value, 4.5); + assert.equal(fc.features[1].properties?.a5, "obj"); + }); + + it("leaves A5 cell geometry unchanged (native dateline handling)", () => { + const raw = [ + [170, 0], + [-170, 0], + [-170, 1], + [170, 1], + [170, 0], + ]; + const fc = a5RowsToFeatureCollection([ + { + a5: "x", + geojson: JSON.stringify({ type: "Polygon", coordinates: [raw] }), + }, + ]); + assert.deepEqual( + (fc.features[0].geometry as { coordinates: number[][][] }).coordinates[0], + raw, + ); + }); +}); diff --git a/tests/dggal-tools.test.ts b/tests/dggal-tools.test.ts new file mode 100644 index 000000000..86c08d881 --- /dev/null +++ b/tests/dggal-tools.test.ts @@ -0,0 +1,124 @@ +import assert from "node:assert/strict"; +import { describe, it } from "node:test"; +import { + DGGAL_GRID_TYPES, + DGGAL_HARD_CAP, + DGGAL_MAX_TOOL_RES, + DGGAL_TYPES, + DEFAULT_DGGAL_GRID_TYPE, + binPointsToDggal, + compactDggalTokens, + dggalApproxGlobalCount, + dggalGridFromBbox, + estimateDggalCellCount, + estimateDggalExpandCount, + expandDggalTokens, + maxResolutionForDggal, + resolveDggalGridType, + suggestDggalResolution, + withDggalDggrs, +} from "../packages/processing/src/dggal-tools"; +import { bboxAreaKm2 } from "../packages/processing/src/h3-tools"; + +describe("dggal type table", () => { + it("exposes the documented DGGRS types with class names and max resolutions", () => { + assert.equal(DEFAULT_DGGAL_GRID_TYPE, "isea3h"); + assert.equal(DGGAL_TYPES.isea3h.className, "ISEA3H"); + assert.equal(DGGAL_TYPES.isea3h.maxRes, 33); + assert.equal(DGGAL_TYPES.gnosis.className, "GNOSISGlobalGrid"); + assert.equal(DGGAL_TYPES.gnosis.maxRes, 28); + assert.equal(DGGAL_TYPES.healpix.className, "HEALPix"); + assert.equal(DGGAL_TYPES.rhealpix.className, "rHEALPix"); + assert.equal(DGGAL_MAX_TOOL_RES, 33); + assert.equal(DGGAL_GRID_TYPES.length, 18); + assert.equal(maxResolutionForDggal("isea4r"), 25); + assert.equal(resolveDggalGridType("healpix"), "healpix"); + assert.equal(resolveDggalGridType("nope"), "isea3h"); + }); + + it("suggests coarser resolutions for larger areas", () => { + const big = bboxAreaKm2([-10, -10, 10, 10]); + const tiny = bboxAreaKm2([0, 0, 0.001, 0.001]); + const rBig = suggestDggalResolution(big, undefined, undefined, "isea3h"); + const rTiny = suggestDggalResolution(tiny, undefined, undefined, "isea3h"); + assert.ok(rBig < rTiny); + assert.ok(estimateDggalCellCount(big, rBig, "isea3h") <= 10_000); + assert.ok(dggalApproxGlobalCount(0, "isea3h") >= 10); + assert.ok(estimateDggalCellCount(510_065_621.724, 25, "isea3h") > DGGAL_HARD_CAP); + }); +}); + +describe("dggal compact / expand (WASM)", () => { + it("compacts four ISEA4R siblings into their parent and expands back", async () => { + await withDggalDggrs("isea4r", (engine) => { + const parent = engine.getZoneFromWGS84Centroid(3, { + lat: (10 * Math.PI) / 180, + lon: (10 * Math.PI) / 180, + }); + const kids = [...engine.getSubZones(parent, 1)].map((z) => engine.getZoneTextID(z)); + assert.equal(kids.length, 4); + const compacted = compactDggalTokens(engine, kids); + assert.equal(compacted.length, 1); + assert.equal(compacted[0], engine.getZoneTextID(parent)); + const expanded = expandDggalTokens(engine, compacted, 4); + assert.equal(expanded.length, 4); + assert.equal(estimateDggalExpandCount(engine, compacted, 4), 4); + assert.deepEqual(new Set(expanded), new Set(kids)); + }); + }); + + it("optionally compacts a bbox grid after listZones", async () => { + await withDggalDggrs("isea4r", (engine) => { + const plain = dggalGridFromBbox(engine, [-10, -10, 10, 10], 5); + const compacted = dggalGridFromBbox(engine, [-10, -10, 10, 10], 5, DGGAL_HARD_CAP, { + compact: true, + }); + assert.ok(compacted.features.length > 0); + assert.ok(compacted.features.length <= plain.features.length); + }); + }); +}); + +describe("dggal grid / bin (WASM)", () => { + it("covers a small bbox and bins a point", async () => { + await withDggalDggrs("isea3h", (engine) => { + const fc = dggalGridFromBbox(engine, [0, 0, 1, 1], 4); + assert.ok(fc.features.length > 0); + assert.ok(fc.features.length < 5_000); + assert.equal(fc.features[0]!.geometry.type, "Polygon"); + assert.ok(typeof fc.features[0]!.properties?.dggal === "string"); + + const bins = binPointsToDggal( + engine, + { + type: "FeatureCollection", + features: [ + { + type: "Feature", + properties: { pop: 3 }, + geometry: { type: "Point", coordinates: [0.5, 0.5] }, + }, + ], + }, + 5, + "sum", + "pop", + ); + assert.equal(bins.features.length, 1); + assert.equal(bins.features[0]!.properties?.count, 1); + assert.equal(bins.features[0]!.properties?.value, 3); + }); + }); + + it("returns native contiguous rings for antimeridian-crossing bboxes", async () => { + await withDggalDggrs("isea3h", (engine) => { + const fc = dggalGridFromBbox(engine, [179.6, -0.3, 180.4, 0.3], 8); + assert.ok(fc.features.length > 0); + for (const feature of fc.features) { + const lons = feature.geometry.coordinates[0]!.map(([lng]) => lng); + assert.ok(lons.length >= 4, "DGGAL cell ring must contain a closed polygon"); + assert.ok(Math.max(...lons) - Math.min(...lons) < 180); + } + }); + }); +}); diff --git a/tests/dggrid-tools.test.ts b/tests/dggrid-tools.test.ts new file mode 100644 index 000000000..fe1dab9d0 --- /dev/null +++ b/tests/dggrid-tools.test.ts @@ -0,0 +1,243 @@ +import assert from "node:assert/strict"; +import { describe, it } from "node:test"; +import { + DEFAULT_DGGRID_GRID_TYPE, + DGGRID_AVG_AREA_KM2, + DGGRID_GRID_PARAMS_SQL, + DGGRID_GRID_TYPES, + DGGRID_HARD_CAP, + DGGRID_MAX_TOOL_RES, + buildDggridBinSql, + buildDggridGridFromSourceSql, + buildDggridGridFromWktSql, + dggridRowsToFeatureCollection, + estimateDggridCellCount, + maxResolutionForDggrid, + resolveDggridGridType, + suggestDggridResolution, +} from "../packages/processing/src/dggrid-tools"; +import { bboxAreaKm2 } from "../packages/processing/src/h3-tools"; + +const ISEA4H_PARAMS = DGGRID_GRID_PARAMS_SQL.ISEA4H; +const ISEA3H_PARAMS = DGGRID_GRID_PARAMS_SQL.ISEA3H; + +describe("dggrid resolution math", () => { + it("exposes average-area entries through the finest DGGRID type (ISEA3H 35)", () => { + assert.equal(DGGRID_AVG_AREA_KM2.length, DGGRID_MAX_TOOL_RES + 1); + assert.equal(DGGRID_MAX_TOOL_RES, 35); + for (let r = 1; r < DGGRID_AVG_AREA_KM2.length; r += 1) { + assert.ok(DGGRID_AVG_AREA_KM2[r]! < DGGRID_AVG_AREA_KM2[r - 1]!); + } + }); + + it("suggests a coarser resolution for larger areas", () => { + const big = bboxAreaKm2([-10, -10, 10, 10]); + const tiny = bboxAreaKm2([0, 0, 0.001, 0.001]); + const rBig = suggestDggridResolution(big); + const rTiny = suggestDggridResolution(tiny); + assert.ok(rBig < rTiny); + assert.ok(rTiny <= 12); + assert.ok(estimateDggridCellCount(big, rBig) <= 10_000); + }); + + it("estimates fewer cells for aperture-3 than aperture-4 at the same res", () => { + const area = bboxAreaKm2([0, 0, 1, 1]); + // 3^res grows slower than 4^res, so ISEA3H cells are larger → lower count. + assert.ok( + estimateDggridCellCount(area, 5, "ISEA3H") < estimateDggridCellCount(area, 5, "ISEA4H"), + ); + }); + + it("fails safe (Infinity) for an out-of-range resolution", () => { + const area = bboxAreaKm2([0, 0, 1, 1]); + assert.equal(estimateDggridCellCount(area, 30, "ISEA4H"), Number.POSITIVE_INFINITY); + assert.equal(estimateDggridCellCount(area, 36, "ISEA3H"), Number.POSITIVE_INFINITY); + assert.ok(estimateDggridCellCount(area, 30, "ISEA4H") > DGGRID_HARD_CAP); + }); +}); + +describe("dggrid grid type presets", () => { + it("lists all DGGRID named types with matching dggs_params SQL", () => { + assert.deepEqual( + [...DGGRID_GRID_TYPES], + [ + "SUPERFUND", + "PLANETRISK", + "ISEA3H", + "ISEA4H", + "ISEA4T", + "ISEA4D", + "ISEA43H", + "ISEA7H", + "IGEO7", + "FULLER3H", + "FULLER4H", + "FULLER4T", + "FULLER4D", + "FULLER43H", + "FULLER7H", + ], + ); + assert.equal(DEFAULT_DGGRID_GRID_TYPE, "ISEA4H"); + assert.match(DGGRID_GRID_PARAMS_SQL.ISEA4H, /'ISEA', 4, 'HEXAGON'/); + assert.match(DGGRID_GRID_PARAMS_SQL.ISEA3H, /'ISEA', 3, 'HEXAGON'/); + assert.match(DGGRID_GRID_PARAMS_SQL.ISEA4T, /'ISEA', 4, 'TRIANGLE'/); + assert.match(DGGRID_GRID_PARAMS_SQL.ISEA4D, /'ISEA', 4, 'DIAMOND'/); + assert.match(DGGRID_GRID_PARAMS_SQL.ISEA7H, /'ISEA', 7, 'HEXAGON'/); + assert.match(DGGRID_GRID_PARAMS_SQL.IGEO7, /'ISEA', 7, 'HEXAGON'/); + assert.match(DGGRID_GRID_PARAMS_SQL.FULLER4H, /'FULLER', 4, 'HEXAGON'/); + assert.match(DGGRID_GRID_PARAMS_SQL.FULLER3H, /'FULLER', 3, 'HEXAGON'/); + assert.match(DGGRID_GRID_PARAMS_SQL.FULLER4T, /'FULLER', 4, 'TRIANGLE'/); + assert.match(DGGRID_GRID_PARAMS_SQL.FULLER4D, /'FULLER', 4, 'DIAMOND'/); + assert.match(DGGRID_GRID_PARAMS_SQL.FULLER7H, /'FULLER', 7, 'HEXAGON'/); + assert.match( + DGGRID_GRID_PARAMS_SQL.SUPERFUND, + /'FULLER', 3, 'HEXAGON'.*true, '44333333333333333'/, + ); + assert.match( + DGGRID_GRID_PARAMS_SQL.PLANETRISK, + /'ISEA', 7, 'HEXAGON'.*true, '43334777777777777777777'/, + ); + }); + + it("exposes per-type max resolutions from the DGGRID named-type table", () => { + assert.equal(maxResolutionForDggrid("SUPERFUND"), 17); + assert.equal(maxResolutionForDggrid("PLANETRISK"), 22); + assert.equal(maxResolutionForDggrid("ISEA3H"), 35); + assert.equal(maxResolutionForDggrid("ISEA4H"), 29); + assert.equal(maxResolutionForDggrid("ISEA4T"), 29); + assert.equal(maxResolutionForDggrid("ISEA4D"), 29); + assert.equal(maxResolutionForDggrid("ISEA43H"), 18); + assert.equal(maxResolutionForDggrid("ISEA7H"), 21); + assert.equal(maxResolutionForDggrid("IGEO7"), 20); + assert.equal(maxResolutionForDggrid("FULLER3H"), 35); + assert.equal(maxResolutionForDggrid("FULLER4H"), 30); + assert.equal(maxResolutionForDggrid("FULLER4T"), 29); + assert.equal(maxResolutionForDggrid("FULLER4D"), 30); + assert.equal(maxResolutionForDggrid("FULLER43H"), 18); + assert.equal(maxResolutionForDggrid("FULLER7H"), 21); + }); + + it("resolves unknown values to ISEA4H", () => { + assert.equal(resolveDggridGridType(undefined), "ISEA4H"); + assert.equal(resolveDggridGridType("nope"), "ISEA4H"); + assert.equal(resolveDggridGridType("ISEA3H"), "ISEA3H"); + assert.equal(resolveDggridGridType("PLANETRISK"), "PLANETRISK"); + }); +}); + +describe("dggrid SQL builders", () => { + it("builds sample-cover grid SQL from a WKT literal with default ISEA4H params", () => { + const sql = buildDggridGridFromWktSql("POLYGON((0 0, 1 0, 1 1, 0 0))'x", 5); + assert.match(sql, new RegExp(`geo_to_seqnum\\(pt, 5, ${escapeRegex(ISEA4H_PARAMS)}\\)`)); + assert.match(sql, new RegExp(`seqnum_to_boundary\\(cell, 5, ${escapeRegex(ISEA4H_PARAMS)}\\)`)); + assert.match(sql, new RegExp(`dggs_cls_km\\(5, ${escapeRegex(ISEA4H_PARAMS)}\\)`)); + assert.match(sql, /POLYGON\(\(0 0, 1 0, 1 1, 0 0\)\)''x/); + assert.match(sql, /CAST\(cell AS VARCHAR\) AS dggrid/); + assert.match(sql, /CAST\(ST_AsGeoJSON\(seqnum_to_boundary/); + assert.doesNotMatch(sql, /LEAST\(/); + }); + + it("passes ISEA3H dggs_params into grid and bin SQL", () => { + const grid = buildDggridGridFromWktSql("POLYGON((0 0, 1 0, 1 1, 0 0))", 4, "ISEA3H"); + assert.match(grid, new RegExp(escapeRegex(ISEA3H_PARAMS))); + assert.doesNotMatch(grid, /TRIANGLE/); + + const bin = buildDggridBinSql("ST_Read('p.geojson')", 5, "count", undefined, "FULLER4H"); + assert.match(bin, /dggs_params\('FULLER', 4, 'HEXAGON'/); + }); + + it("builds polyfill grid SQL that unions polygons then sample-covers", () => { + const sql = buildDggridGridFromSourceSql("ST_Read('x.geojson')", 6); + assert.match(sql, /ST_Union_Agg\(geom\)/); + assert.match(sql, /'POLYGON', 'MULTIPOLYGON'/); + assert.match(sql, new RegExp(`geo_to_seqnum\\(pt, 6, ${escapeRegex(ISEA4H_PARAMS)}\\)`)); + }); + + it("builds bin SQL with geo_to_seqnum and optional aggregates", () => { + const countSql = buildDggridBinSql("ST_Read('p.geojson')", 5, "count"); + assert.match(countSql, new RegExp(`geo_to_seqnum\\(pt, 5, ${escapeRegex(ISEA4H_PARAMS)}\\)`)); + assert.match( + countSql, + new RegExp(`seqnum_to_boundary\\(cell, 5, ${escapeRegex(ISEA4H_PARAMS)}\\)`), + ); + assert.doesNotMatch(countSql, / AS value/); + + const sumSql = buildDggridBinSql("ST_Read('p.geojson')", 5, "sum", 'pop"x'); + assert.match(sumSql, /sum\(CAST\("pop""x" AS DOUBLE\)\) AS value/); + }); + + it("converts result rows to a FeatureCollection with dggrid props", () => { + const fc = dggridRowsToFeatureCollection([ + { + dggrid: "2380", + count: 3, + geojson: '{"type":"Polygon","coordinates":[[[0,0],[1,0],[1,1],[0,0]]]}', + }, + { dggrid: "skip", geojson: 12 }, + ]); + assert.equal(fc.features.length, 1); + assert.equal(fc.features[0].properties?.dggrid, "2380"); + assert.equal(fc.features[0].properties?.count, 3); + }); + + it("unwraps antimeridian-crossing cell rings for MapLibre by default", () => { + // duck_dggs-style ring with longitudes clamped to [-180, 180]. + const fc = dggridRowsToFeatureCollection([ + { + dggrid: "dateline", + geojson: JSON.stringify({ + type: "Polygon", + coordinates: [ + [ + [170, 0], + [-170, 0], + [-170, 1], + [170, 1], + [170, 0], + ], + ], + }), + }, + ]); + const ring = (fc.features[0].geometry as { coordinates: number[][][] }).coordinates[0]; + assert.deepEqual(ring, [ + [170, 0], + [190, 0], + [190, 1], + [170, 1], + [170, 0], + ]); + // Contiguous: no edge jumps more than 180° of longitude. + for (let i = 1; i < ring.length; i += 1) { + assert.ok(Math.abs(ring[i]![0]! - ring[i - 1]![0]!) < 180); + } + }); + + it("leaves wrapped rings alone when fixAntimeridian is false", () => { + const raw = [ + [170, 0], + [-170, 0], + [-170, 1], + [170, 1], + [170, 0], + ]; + const fc = dggridRowsToFeatureCollection( + [ + { + dggrid: "dateline", + geojson: JSON.stringify({ type: "Polygon", coordinates: [raw] }), + }, + ], + false, + ); + assert.deepEqual( + (fc.features[0].geometry as { coordinates: number[][][] }).coordinates[0], + raw, + ); + }); +}); + +function escapeRegex(value: string): string { + return value.replace(/[.*+?^${}()|[\]\\]/g, "\\$&"); +} diff --git a/tests/dggs-compact.test.ts b/tests/dggs-compact.test.ts new file mode 100644 index 000000000..cf584a792 --- /dev/null +++ b/tests/dggs-compact.test.ts @@ -0,0 +1,332 @@ +import assert from "node:assert/strict"; +import { describe, it } from "node:test"; +import { s2 } from "s2js"; +import { + buildA5CompactSql, + buildA5ExpandCountSql, + buildA5ExpandSql, +} from "../packages/processing/src/a5-tools"; +import { + buildH3CompactSql, + buildH3ExpandCountSql, + buildH3ExpandSql, +} from "../packages/processing/src/h3-tools"; +import { dggsCompactTool, getDggsTool } from "../packages/processing/src/dggs-tools"; +import { DEFAULT_LAYER_STYLE, type GeoLibreLayer } from "@geolibre/core"; +import type { DuckDbCapability, ProcessingContext } from "../packages/processing/src/types"; +import type { Feature, Polygon } from "geojson"; + +describe("dggs compact/expand SQL", () => { + it("builds H3 compact and expand SQL from a cell ID field", () => { + const compact = buildH3CompactSql("ST_Read('x.geojson')", 'h3"x'); + assert.match(compact, /h3_string_to_h3\(CAST\("h3""x" AS VARCHAR\)\)/); + assert.match(compact, /h3_compact_cells\(cells\)/); + assert.match(compact, /h3_h3_to_string\(cell\) AS h3/); + + const expand = buildH3ExpandSql("ST_Read('x.geojson')", 7); + assert.match(expand, /h3_uncompact_cells\(cells, 7\)/); + assert.match(buildH3ExpandCountSql("ST_Read('x.geojson')", 7), /len\(h3_uncompact_cells/); + }); + + it("builds A5 compact and expand SQL from a cell ID field", () => { + const compact = buildA5CompactSql("ST_Read('x.geojson')"); + assert.match(compact, /a5_hex_to_u64\(CAST\("a5" AS VARCHAR\)\)/); + assert.match(compact, /a5_compact\(cells\)/); + assert.match(compact, /a5_u64_to_hex\(cell\) AS a5/); + + const expand = buildA5ExpandSql("ST_Read('x.geojson')", 8, "a5"); + assert.match(expand, /a5_uncompact\(cells, 8\)/); + assert.match(buildA5ExpandCountSql("ST_Read('x.geojson')", 8), /len\(a5_uncompact/); + }); +}); + +function h3CellLayer(): GeoLibreLayer { + return { + id: "cells", + name: "Cells", + type: "geojson", + source: { type: "geojson" }, + visible: true, + opacity: 1, + style: { ...DEFAULT_LAYER_STYLE }, + metadata: {}, + geojson: { + type: "FeatureCollection", + features: [ + { + type: "Feature", + properties: { h3: "85283473fffffff" }, + geometry: { + type: "Polygon", + coordinates: [ + [ + [0, 0], + [1, 0], + [1, 1], + [0, 0], + ], + ], + }, + }, + ], + }, + }; +} + +function s2SiblingCellLayer(): GeoLibreLayer { + const leaf = s2.cellid.fromLatLng(s2.LatLng.fromDegrees(10, 10)); + const parent = s2.cellid.parent(leaf, 5); + const features: Feature[] = []; + let id = s2.cellid.childBegin(parent); + for (let i = 0; i < 4; i += 1) { + const token = s2.cellid.toToken(id); + features.push({ + type: "Feature", + properties: { s2: token }, + geometry: { + type: "Polygon", + coordinates: [ + [ + [0, 0], + [1, 0], + [1, 1], + [0, 0], + ], + ], + }, + }); + id = s2.cellid.next(id); + } + return { + id: "s2cells", + name: "S2 Cells", + type: "geojson", + source: { type: "geojson" }, + visible: true, + opacity: 1, + style: { ...DEFAULT_LAYER_STYLE }, + metadata: {}, + geojson: { type: "FeatureCollection", features }, + }; +} + +function mockDuckDb(): DuckDbCapability & { queries: string[] } { + const queries: string[] = []; + return { + queries, + ensureExtensions: async () => {}, + registerGeoJson: async () => ({ + sql: "ST_Read('mock.geojson')", + release: async () => {}, + }), + query: async (sql: string) => { + queries.push(sql); + if (sql.includes(" AS n ")) return [{ n: 3 }]; + const geojson = '{"type":"Polygon","coordinates":[[[0,0],[1,0],[1,1],[0,0]]]}'; + if (sql.includes("a5_")) return [{ a5: "1600000000000000", geojson }]; + return [{ h3: "8428347ffffffff", geojson }]; + }, + }; +} + +describe("dggs compact tool", () => { + it("registers under getDggsTool", () => { + assert.equal(getDggsTool("dggs-compact"), dggsCompactTool); + assert.equal(dggsCompactTool.group, "DGGS"); + }); + + it("compacts an H3 cell layer", async () => { + const duckdb = mockDuckDb(); + const logs: string[] = []; + const added: string[] = []; + const ctx: ProcessingContext = { + layers: [h3CellLayer()], + parameters: { dggsType: "h3", mode: "compact", layer: "cells" }, + log: (m) => logs.push(m), + addResultLayer: (name) => added.push(name), + duckdb, + }; + await dggsCompactTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /H3 compact/); + assert.ok(duckdb.queries.some((q) => /h3_compact_cells/.test(q))); + assert.ok(logs.some((l) => /Compacted to \d+ H3 cell/.test(l))); + }); + + it("expands an H3 cell layer after a count guard", async () => { + const duckdb = mockDuckDb(); + const added: string[] = []; + const ctx: ProcessingContext = { + layers: [h3CellLayer()], + parameters: { dggsType: "h3", mode: "expand", layer: "cells", resolution: 6 }, + log: () => {}, + addResultLayer: (name) => added.push(name), + duckdb, + }; + await dggsCompactTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /H3 expand \(res 6\)/); + assert.ok(duckdb.queries.some((q) => /len\(h3_uncompact_cells/.test(q))); + assert.ok(duckdb.queries.some((q) => /h3_uncompact_cells\(cells, 6\)/.test(q))); + }); + + it("rejects expand when the count exceeds the hard cap", async () => { + const duckdb = mockDuckDb(); + duckdb.query = async (sql: string) => { + duckdb.queries.push(sql); + if (sql.includes(" AS n ")) return [{ n: 500_000 }]; + return []; + }; + const logs: string[] = []; + const added: string[] = []; + const ctx: ProcessingContext = { + layers: [h3CellLayer()], + parameters: { dggsType: "h3", mode: "expand", layer: "cells", resolution: 10 }, + log: (m) => logs.push(m), + addResultLayer: (name) => added.push(name), + duckdb, + }; + await dggsCompactTool.run(ctx); + assert.equal(added.length, 0); + assert.ok(logs.some((l) => /cap/.test(l))); + }); + + it("exposes S2 and DGGAL in the type picker; Fix antimeridian is H3/S2 only", () => { + const typeParam = dggsCompactTool.parameters.find((p) => p.id === "dggsType"); + assert.deepEqual( + typeParam?.options?.map((o) => o.value), + ["h3", "s2", "a5", "dggal"], + ); + assert.ok(typeParam?.options?.some((o) => o.value === "s2")); + assert.ok(typeParam?.options?.some((o) => o.value === "dggal")); + const fix = dggsCompactTool.parameters.find((p) => p.id === "fixAntimeridian"); + assert.deepEqual(fix?.visibleWhen, { param: "dggsType", in: ["h3", "s2"] }); + assert.ok(!fix?.visibleWhen?.in.includes("a5")); + assert.ok(!fix?.visibleWhen?.in.includes("dggal")); + const dggalType = dggsCompactTool.parameters.find((p) => p.id === "dggalType"); + assert.deepEqual(dggalType?.visibleWhen, { param: "dggsType", in: ["dggal"] }); + }); + + it("compacts an S2 cell layer without DuckDB", async () => { + const logs: string[] = []; + const added: string[] = []; + let featureCount = 0; + const ctx: ProcessingContext = { + layers: [s2SiblingCellLayer()], + parameters: { dggsType: "s2", mode: "compact", layer: "s2cells" }, + log: (m) => logs.push(m), + addResultLayer: (name, fc) => { + added.push(name); + featureCount = fc.features.length; + }, + }; + await dggsCompactTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0]!, /S2 compact/); + assert.equal(featureCount, 1); + assert.ok(logs.some((l) => /Compacted to 1 S2 cell/.test(l))); + }); + + it("expands an S2 cell layer without DuckDB", async () => { + const compactLayer = s2SiblingCellLayer(); + // Pre-compact to one parent so expand has work to do. + const parentToken = (() => { + const leaf = s2.cellid.fromLatLng(s2.LatLng.fromDegrees(10, 10)); + return s2.cellid.toToken(s2.cellid.parent(leaf, 5)); + })(); + compactLayer.geojson = { + type: "FeatureCollection", + features: [ + { + type: "Feature", + properties: { s2: parentToken }, + geometry: { + type: "Polygon", + coordinates: [ + [ + [0, 0], + [1, 0], + [1, 1], + [0, 0], + ], + ], + }, + }, + ], + }; + const added: string[] = []; + let featureCount = 0; + const ctx: ProcessingContext = { + layers: [compactLayer], + parameters: { dggsType: "s2", mode: "expand", layer: "s2cells", resolution: 6 }, + log: () => {}, + addResultLayer: (name, fc) => { + added.push(name); + featureCount = fc.features.length; + }, + }; + await dggsCompactTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0]!, /S2 expand \(res 6\)/); + assert.equal(featureCount, 4); + }); + + it("compacts a DGGAL ISEA4R cell layer without DuckDB", async () => { + const { withDggalDggrs } = await import("../packages/processing/src/dggal-tools"); + const layer = await withDggalDggrs("isea4r", (engine) => { + const parent = engine.getZoneFromWGS84Centroid(3, { + lat: (10 * Math.PI) / 180, + lon: (10 * Math.PI) / 180, + }); + const features = [...engine.getSubZones(parent, 1)].map((z) => { + const token = engine.getZoneTextID(z); + return { + type: "Feature" as const, + properties: { dggal: token }, + geometry: { + type: "Polygon" as const, + coordinates: [ + [ + [0, 0], + [1, 0], + [1, 1], + [0, 0], + ], + ], + }, + }; + }); + return { + id: "dggalcells", + name: "DGGAL Cells", + type: "geojson" as const, + source: { type: "geojson" as const }, + visible: true, + opacity: 1, + style: { ...DEFAULT_LAYER_STYLE }, + metadata: {}, + geojson: { type: "FeatureCollection" as const, features }, + }; + }); + const added: string[] = []; + let featureCount = 0; + const ctx: ProcessingContext = { + layers: [layer], + parameters: { + dggsType: "dggal", + dggalType: "isea4r", + mode: "compact", + layer: "dggalcells", + }, + log: () => {}, + addResultLayer: (name, fc) => { + added.push(name); + featureCount = fc.features.length; + }, + }; + await dggsCompactTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0]!, /ISEA4R compact/); + assert.equal(featureCount, 1); + }); +}); diff --git a/tests/dggs-tools.test.ts b/tests/dggs-tools.test.ts new file mode 100644 index 000000000..4f3f209ba --- /dev/null +++ b/tests/dggs-tools.test.ts @@ -0,0 +1,617 @@ +import assert from "node:assert/strict"; +import { describe, it } from "node:test"; +import { + createDggsGridTool, + dggsBinPointsTool, + dggsCompactTool, + extensionForDggs, + getDggsTool, + maxResolutionForDggs, +} from "../packages/processing/src/dggs-tools"; +import { DEFAULT_LAYER_STYLE, type GeoLibreLayer } from "@geolibre/core"; +import type { DuckDbCapability, ProcessingContext } from "../packages/processing/src/types"; + +function polygonLayer(): GeoLibreLayer { + return { + id: "poly", + name: "Poly", + type: "geojson", + source: { type: "geojson" }, + visible: true, + opacity: 1, + style: { ...DEFAULT_LAYER_STYLE }, + metadata: {}, + geojson: { + type: "FeatureCollection", + features: [ + { + type: "Feature", + properties: {}, + geometry: { + type: "Polygon", + coordinates: [ + [ + [0, 0], + [1, 0], + [1, 1], + [0, 1], + [0, 0], + ], + ], + }, + }, + ], + }, + }; +} + +function pointLayer(): GeoLibreLayer { + return { + ...polygonLayer(), + id: "pts", + name: "Pts", + geojson: { + type: "FeatureCollection", + features: [ + { + type: "Feature", + properties: {}, + geometry: { type: "Point", coordinates: [0.5, 0.5] }, + }, + ], + }, + }; +} + +function mockDuckDb(): DuckDbCapability & { + queries: string[]; + released: number[]; + extensions: string[][]; +} { + const queries: string[] = []; + const released: number[] = []; + const extensions: string[][] = []; + return { + queries, + released, + extensions, + ensureExtensions: async (names) => { + extensions.push([...names]); + }, + registerGeoJson: async () => ({ + sql: "ST_Read('mock.geojson')", + release: async () => { + released.push(1); + }, + }), + query: async (sql: string) => { + queries.push(sql); + const geojson = '{"type":"Polygon","coordinates":[[[0,0],[1,0],[1,1],[0,0]]]}'; + if (sql.includes("a5_")) { + return [{ a5: "1600000000000000", geojson, count: 1 }]; + } + if (sql.includes("geo_to_seqnum") || sql.includes("seqnum_to_boundary")) { + return [{ dggrid: "2380", geojson, count: 1 }]; + } + return [{ h3: "8928308280fffff", geojson, count: 1 }]; + }, + }; +} + +function baseCtx( + layers: GeoLibreLayer[], + parameters: Record, +): { + ctx: ProcessingContext; + logs: string[]; + added: string[]; + duckdb: ReturnType; +} { + const logs: string[] = []; + const added: string[] = []; + const duckdb = mockDuckDb(); + const ctx: ProcessingContext = { + layers, + parameters, + log: (m) => logs.push(m), + addResultLayer: (name) => added.push(name), + duckdb, + viewportBounds: () => [0, 0, 1, 1], + }; + return { ctx, logs, added, duckdb }; +} + +describe("dggs generator", () => { + it("registers under getDggsTool with DGGS group", () => { + assert.equal(getDggsTool("dggs-grid"), createDggsGridTool); + assert.equal(getDggsTool("dggs-bin"), dggsBinPointsTool); + assert.equal(getDggsTool("dggs-compact"), dggsCompactTool); + assert.equal(getDggsTool("missing"), undefined); + assert.equal(createDggsGridTool.group, "DGGS"); + assert.equal(createDggsGridTool.name, "DGGS Generator"); + assert.equal(dggsBinPointsTool.name, "DGGS Binning"); + }); + + it("exposes Fix antimeridian for H3, S2, and DGGRID only, default checked", () => { + for (const tool of [createDggsGridTool, dggsBinPointsTool]) { + const param = tool.parameters.find((p) => p.id === "fixAntimeridian"); + assert.ok(param); + assert.equal(param.type, "boolean"); + assert.equal(param.default, true); + assert.deepEqual(param.visibleWhen, { + param: "dggsType", + in: ["h3", "s2", "dggrid"], + }); + assert.ok(!param.visibleWhen!.in.includes("a5")); + assert.ok(!param.visibleWhen!.in.includes("dggal")); + } + const typeParam = createDggsGridTool.parameters.find((p) => p.id === "dggsType"); + assert.ok(typeParam && typeParam.type === "select"); + assert.deepEqual( + typeParam.options.map((o) => o.value), + ["h3", "s2", "a5", "dggrid", "dggal"], + ); + }); + + it("exposes H3 max 15, S2/A5 max 30, and per-DGGRID/DGGAL-type maxima", () => { + assert.equal(maxResolutionForDggs("h3"), 15); + assert.equal(maxResolutionForDggs("s2"), 30); + assert.equal(maxResolutionForDggs("a5"), 30); + assert.equal(maxResolutionForDggs("dggrid"), 29); + assert.equal(maxResolutionForDggs("dggrid", "ISEA4H"), 29); + assert.equal(maxResolutionForDggs("dggrid", "ISEA3H"), 35); + assert.equal(maxResolutionForDggs("dggrid", "FULLER4H"), 30); + assert.equal(maxResolutionForDggs("dggrid", "SUPERFUND"), 17); + assert.equal(maxResolutionForDggs("dggrid", "PLANETRISK"), 22); + assert.equal(maxResolutionForDggs("dggrid", "IGEO7"), 20); + assert.equal(maxResolutionForDggs("dggal"), 33); + assert.equal(maxResolutionForDggs("dggal", "isea3h"), 33); + assert.equal(maxResolutionForDggs("dggal", "isea4r"), 25); + assert.equal(maxResolutionForDggs("dggal", "healpix"), 26); + assert.equal(maxResolutionForDggs("dggal", "gnosis"), 28); + assert.equal(extensionForDggs("dggal"), null); + + const dggalParam = createDggsGridTool.parameters.find((p) => p.id === "dggalType"); + assert.ok(dggalParam); + assert.deepEqual(dggalParam.visibleWhen, { param: "dggsType", in: ["dggal"] }); + assert.equal(dggalParam.default, "isea3h"); + }); + + it("throws a clear error when duckdb is unavailable for H3", async () => { + await assert.rejects( + () => + Promise.resolve( + createDggsGridTool.run({ + layers: [], + parameters: { dggsType: "h3", source: "viewport", resolution: 5 }, + log: () => {}, + viewportBounds: () => [0, 0, 1, 1], + }), + ), + /requires DuckDB/, + ); + }); + + it("creates an S2 grid from the map viewport without DuckDB", async () => { + assert.equal(extensionForDggs("s2"), null); + const { ctx, added, duckdb, logs } = baseCtx([], { + dggsType: "s2", + source: "viewport", + resolution: 4, + }); + // S2 must not require DuckDB — clear it to prove the client path. + ctx.duckdb = undefined; + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /S2 grid \(res 4\)/); + assert.equal(duckdb.queries.length, 0); + assert.ok(logs.some((l) => /Created \d+ S2 cell/.test(l))); + }); + + it("compacts an S2 viewport grid when Compact cells is checked", async () => { + let featureCount = 0; + const { ctx, added, duckdb } = baseCtx([], { + dggsType: "s2", + source: "viewport", + resolution: 6, + compactCells: true, + }); + ctx.duckdb = undefined; + ctx.addResultLayer = (name, fc) => { + added.push(name); + featureCount = fc.features.length; + }; + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /S2 grid \(res 6, compact\)/); + assert.equal(duckdb.queries.length, 0); + assert.ok(featureCount > 0); + }); + + it("creates a DGGAL grid from the map viewport without DuckDB", async () => { + const { ctx, added, duckdb } = baseCtx([], { + dggsType: "dggal", + dggalType: "isea3h", + source: "viewport", + resolution: 3, + }); + ctx.duckdb = undefined; + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /ISEA3H grid \(res 3\)/); + assert.equal(duckdb.queries.length, 0); + }); + + it("compacts a DGGAL viewport grid when Compact cells is checked", async () => { + let featureCount = 0; + const { ctx, added, duckdb } = baseCtx([], { + dggsType: "dggal", + dggalType: "isea4r", + source: "viewport", + resolution: 5, + compactCells: true, + }); + ctx.duckdb = undefined; + ctx.addResultLayer = (name, fc) => { + added.push(name); + featureCount = fc.features.length; + }; + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /ISEA4R grid \(res 5, compact\)/); + assert.equal(duckdb.queries.length, 0); + assert.ok(featureCount > 0); + }); + + it("exposes Compact cells for H3, A5, S2, and DGGAL, default off", () => { + const param = createDggsGridTool.parameters.find((p) => p.id === "compactCells"); + assert.ok(param); + assert.equal(param.type, "boolean"); + assert.equal(param.default, false); + assert.deepEqual(param.visibleWhen, { param: "dggsType", in: ["h3", "a5", "s2", "dggal"] }); + }); + + it("creates an H3 grid from the map viewport", async () => { + const { ctx, added, duckdb } = baseCtx([], { + dggsType: "h3", + source: "viewport", + resolution: 5, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /H3 grid \(res 5\)/); + assert.deepEqual(duckdb.extensions[0], ["spatial", "h3"]); + assert.match(duckdb.queries[0], /h3_polygon_wkt_to_cells_experimental/); + assert.doesNotMatch(duckdb.queries[0], /h3_compact_cells/); + }); + + it("compacts an H3 viewport grid when Compact cells is checked", async () => { + const { ctx, added, duckdb } = baseCtx([], { + dggsType: "h3", + source: "viewport", + resolution: 5, + compactCells: true, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /H3 grid \(res 5, compact\)/); + assert.match(duckdb.queries[0], /h3_compact_cells\(cells\)/); + }); + it("creates an A5 grid from the map viewport", async () => { + const { ctx, added, duckdb } = baseCtx([], { + dggsType: "a5", + source: "viewport", + resolution: 5, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /A5 grid \(res 5\)/); + assert.deepEqual(duckdb.extensions[0], ["spatial", "a5"]); + assert.match(duckdb.queries[0], /a5_uncompact\(a5_geometry_to_cells/); + }); + + it("creates a DGGRID grid via duck_dggs sample cover", async () => { + assert.equal(extensionForDggs("dggrid"), "duck_dggs"); + const { ctx, added, duckdb } = baseCtx([], { + dggsType: "dggrid", + source: "viewport", + resolution: 5, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /ISEA4H grid \(res 5\)/); + assert.deepEqual(duckdb.extensions[0], ["spatial", "duck_dggs"]); + assert.match(duckdb.queries[0], /geo_to_seqnum/); + assert.match(duckdb.queries[0], /seqnum_to_boundary/); + assert.match(duckdb.queries[0], /dggs_params\('ISEA', 4, 'HEXAGON'/); + }); + + it("passes the selected DGGRID type into duck_dggs params", async () => { + const dggridTypeParam = createDggsGridTool.parameters.find((p) => p.id === "dggridType"); + assert.ok(dggridTypeParam); + assert.deepEqual(dggridTypeParam.visibleWhen, { param: "dggsType", in: ["dggrid"] }); + assert.equal(dggridTypeParam.default, "ISEA4H"); + + const { ctx, added, duckdb } = baseCtx([], { + dggsType: "dggrid", + dggridType: "ISEA4T", + source: "viewport", + resolution: 4, + }); + await createDggsGridTool.run(ctx); + assert.match(added[0], /ISEA4T grid \(res 4\)/); + assert.match(duckdb.queries[0], /dggs_params\('ISEA', 4, 'TRIANGLE'/); + }); + + it("rejects DGGRID resolutions above the selected type's max", async () => { + const overIsea4 = baseCtx([], { + dggsType: "dggrid", + dggridType: "ISEA4H", + source: "viewport", + resolution: 30, + }); + await createDggsGridTool.run(overIsea4.ctx); + assert.equal(overIsea4.added.length, 0); + assert.ok(overIsea4.logs.some((l) => /0 to 29 for ISEA4H/.test(l))); + + const okIsea3 = baseCtx([], { + dggsType: "dggrid", + dggridType: "ISEA3H", + source: "viewport", + resolution: 30, + }); + // res 30 is valid for ISEA3H but always exceeds the hard cell-count cap on a 1° viewport. + await createDggsGridTool.run(okIsea3.ctx); + assert.equal(okIsea3.added.length, 0); + assert.ok(okIsea3.logs.some((l) => /cap/i.test(l))); + }); + + it("defaults dggsType to h3", async () => { + const { ctx, duckdb } = baseCtx([], { source: "viewport", resolution: 4 }); + await createDggsGridTool.run(ctx); + assert.deepEqual(duckdb.extensions[0], ["spatial", "h3"]); + }); + + it("rejects an antimeridian-crossing viewport", async () => { + const { ctx, added, logs } = baseCtx([], { + dggsType: "h3", + source: "viewport", + resolution: 4, + }); + ctx.viewportBounds = () => [170, 0, -170, 1]; + await createDggsGridTool.run(ctx); + assert.equal(added.length, 0); + assert.ok(logs.some((l) => /antimeridian/i.test(l))); + }); + + it("creates a grid from a manual bounding box", async () => { + const { ctx, added } = baseCtx([], { + dggsType: "a5", + source: "bbox", + west: 0, + south: 0, + east: 1, + north: 1, + resolution: 5, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + }); + + it("rejects a degenerate manual bounding box", async () => { + const { ctx, added, logs } = baseCtx([], { + dggsType: "h3", + source: "bbox", + west: 2, + south: 0, + east: 1, + north: 1, + resolution: 5, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 0); + assert.ok(logs.some((l) => /west < east/i.test(l))); + }); + + it("rejects H3 resolution above 15", async () => { + const { ctx, added, logs } = baseCtx([], { + dggsType: "h3", + source: "viewport", + resolution: 16, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 0); + assert.ok(logs.some((l) => /0 to 15/i.test(l))); + }); + + it("accepts A5 resolution up to 30", async () => { + const { ctx, added, logs, duckdb } = baseCtx([], { + dggsType: "a5", + // Tiny viewport so the hard-cap estimate does not trip at high res. + source: "bbox", + west: 0, + south: 0, + east: 0.000001, + north: 0.000001, + resolution: 30, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.ok(!logs.some((l) => /resolution must be/i.test(l))); + assert.match(duckdb.queries[0], /a5_geometry_to_cells\(.+, 30\), 30\)/); + }); + + it("auto-suggests a resolution when none is given", async () => { + const { ctx, logs } = baseCtx([], { dggsType: "h3", source: "viewport" }); + await createDggsGridTool.run(ctx); + assert.ok(logs.some((l) => /suggested resolution/i.test(l))); + }); + + it("aborts when the requested resolution exceeds the hard cap", async () => { + const { ctx, added, logs } = baseCtx([polygonLayer()], { + dggsType: "h3", + source: "extent", + layer: "poly", + resolution: 15, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 0); + assert.ok(logs.some((l) => /cap/i.test(l))); + }); + + it("polyfills a selected polygon layer and releases the registered source", async () => { + const { ctx, added, duckdb } = baseCtx([polygonLayer()], { + dggsType: "a5", + source: "polyfill", + layer: "poly", + resolution: 6, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + assert.equal(duckdb.released.length, 1); + assert.match(duckdb.queries[0], /ST_Union_Agg/); + }); + + it("rejects polyfill of a non-polygon layer", async () => { + const { ctx, added, logs } = baseCtx([pointLayer()], { + dggsType: "h3", + source: "polyfill", + layer: "pts", + resolution: 6, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 0); + assert.ok(logs.some((l) => /polygon/i.test(l))); + }); + + it("fills the extent of a non-polygon layer", async () => { + const { ctx, added } = baseCtx([pointLayer()], { + dggsType: "h3", + source: "extent", + layer: "pts", + resolution: 6, + }); + await createDggsGridTool.run(ctx); + assert.equal(added.length, 1); + }); + + it("logs a soft message and adds no layer when zero cells are produced", async () => { + const logs: string[] = []; + const added: string[] = []; + const ctx: ProcessingContext = { + layers: [], + parameters: { dggsType: "h3", source: "viewport", resolution: 5 }, + log: (m) => logs.push(m), + addResultLayer: (name) => added.push(name), + viewportBounds: () => [0, 0, 1, 1], + duckdb: { + ensureExtensions: async () => {}, + registerGeoJson: async () => ({ + sql: "ST_Read('mock.geojson')", + release: async () => {}, + }), + query: async () => [], + }, + }; + await createDggsGridTool.run(ctx); + assert.equal(added.length, 0); + assert.ok(logs.some((l) => /no h3 cells/i.test(l))); + }); +}); + +describe("dggs binning", () => { + it("bins points to H3 and requires a field for non-count aggregates", async () => { + const missing = baseCtx([pointLayer()], { + dggsType: "h3", + layer: "pts", + aggOp: "sum", + resolution: 7, + }); + await dggsBinPointsTool.run(missing.ctx); + assert.equal(missing.added.length, 0); + assert.ok(missing.logs.some((l) => /field/i.test(l))); + + const ok = baseCtx([pointLayer()], { + dggsType: "h3", + layer: "pts", + aggOp: "count", + resolution: 7, + }); + await dggsBinPointsTool.run(ok.ctx); + assert.equal(ok.added.length, 1); + assert.match(ok.added[0], /H3 bins/); + assert.deepEqual(ok.duckdb.extensions[0], ["spatial", "h3"]); + assert.match(ok.duckdb.queries[0], /h3_latlng_to_cell/); + assert.equal(ok.duckdb.released.length, 1); + }); + + it("bins points to A5", async () => { + const { ctx, added, duckdb } = baseCtx([pointLayer()], { + dggsType: "a5", + layer: "pts", + aggOp: "count", + resolution: 7, + }); + await dggsBinPointsTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /A5 bins/); + assert.deepEqual(duckdb.extensions[0], ["spatial", "a5"]); + assert.match(duckdb.queries[0], /a5_lonlat_to_cell/); + }); + + it("bins points to S2 without DuckDB", async () => { + const { ctx, added, duckdb } = baseCtx([pointLayer()], { + dggsType: "s2", + layer: "pts", + aggOp: "count", + resolution: 7, + }); + ctx.duckdb = undefined; + await dggsBinPointsTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /S2 bins/); + assert.equal(duckdb.queries.length, 0); + }); + + it("bins points to DGGAL without DuckDB", async () => { + const { ctx, added, duckdb } = baseCtx([pointLayer()], { + dggsType: "dggal", + dggalType: "isea3h", + layer: "pts", + aggOp: "count", + resolution: 5, + }); + ctx.duckdb = undefined; + await dggsBinPointsTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /ISEA3H bins/); + assert.equal(duckdb.queries.length, 0); + }); + + it("bins points to DGGRID", async () => { + const { ctx, added, duckdb } = baseCtx([pointLayer()], { + dggsType: "dggrid", + layer: "pts", + aggOp: "count", + resolution: 5, + }); + await dggsBinPointsTool.run(ctx); + assert.equal(added.length, 1); + assert.match(added[0], /ISEA4H bins/); + assert.deepEqual(duckdb.extensions[0], ["spatial", "duck_dggs"]); + assert.match(duckdb.queries[0], /geo_to_seqnum/); + assert.match(duckdb.queries[0], /dggs_params\('ISEA', 4, 'HEXAGON'/); + }); + + it("rejects an unknown aggregate operation", async () => { + const { ctx, added, logs } = baseCtx([pointLayer()], { + dggsType: "h3", + layer: "pts", + aggOp: "median", + resolution: 7, + }); + await dggsBinPointsTool.run(ctx); + assert.equal(added.length, 0); + assert.ok(logs.some((l) => /unknown aggregate/i.test(l))); + }); +}); diff --git a/tests/h3-tools.test.ts b/tests/h3-tools.test.ts index e7493cbdb..f359e6c9a 100644 --- a/tests/h3-tools.test.ts +++ b/tests/h3-tools.test.ts @@ -5,18 +5,17 @@ import { H3_HARD_CAP, bboxAreaKm2, bboxToWktPolygon, - binPointsTool, buildBinSql, + buildGridFromBboxSql, buildGridFromSourceSql, buildGridFromWktSql, - createH3GridTool, estimateCellCount, getH3Tool, + normalizeLonLatBbox, rowsToFeatureCollection, suggestResolution, } from "../packages/processing/src/h3-tools"; -import { DEFAULT_LAYER_STYLE, type GeoLibreLayer } from "@geolibre/core"; -import type { DuckDbCapability, ProcessingContext } from "../packages/processing/src/types"; +import { getVectorTool, resolveVectorRerun } from "../packages/processing/src/vector-tools"; describe("h3 resolution math", () => { it("exposes 16 average-area entries (res 0..15), strictly decreasing", () => { @@ -71,294 +70,33 @@ describe("h3 resolution math", () => { }); }); -function polygonLayer(): GeoLibreLayer { - return { - id: "poly", - name: "Poly", - type: "geojson", - source: { type: "geojson" }, - visible: true, - opacity: 1, - style: { ...DEFAULT_LAYER_STYLE }, - metadata: {}, - geojson: { - type: "FeatureCollection", - features: [ - { - type: "Feature", - properties: {}, - geometry: { - type: "Polygon", - coordinates: [ - [ - [0, 0], - [1, 0], - [1, 1], - [0, 1], - [0, 0], - ], - ], - }, - }, - ], - }, - }; -} - -function pointLayer(): GeoLibreLayer { - return { - ...polygonLayer(), - id: "pts", - name: "Pts", - geojson: { - type: "FeatureCollection", - features: [ - { - type: "Feature", - properties: { pop: 5 }, - geometry: { type: "Point", coordinates: [0.5, 0.5] }, - }, - ], - }, - }; -} - -/** Capability stub that records queries/releases and returns one canned hex row. */ -function mockDuckDb(): DuckDbCapability & { - queries: string[]; - released: number[]; -} { - const queries: string[] = []; - const released: number[] = []; - return { - queries, - released, - ensureExtensions: async () => {}, - registerGeoJson: async () => ({ - sql: "ST_Read('mock.geojson')", - release: async () => { - released.push(1); - }, - }), - query: async (sql: string) => { - queries.push(sql); - return [ - { - h3: "8928308280fffff", - count: 1, - geojson: '{"type":"Polygon","coordinates":[[[0,0],[1,0],[1,1],[0,0]]]}', - }, - ]; - }, - }; -} - -function baseCtx( - layers: GeoLibreLayer[], - parameters: Record, -): { - ctx: ProcessingContext; - logs: string[]; - added: string[]; - duckdb: ReturnType; -} { - const logs: string[] = []; - const added: string[] = []; - const duckdb = mockDuckDb(); - const ctx: ProcessingContext = { - layers, - parameters, - log: (m) => logs.push(m), - addResultLayer: (name) => added.push(name), - duckdb, - viewportBounds: () => [0, 0, 1, 1], - }; - return { ctx, logs, added, duckdb }; -} - -describe("h3 tools", () => { - it("registers both tools under getH3Tool", () => { - assert.equal(getH3Tool("h3-grid"), createH3GridTool); - assert.equal(getH3Tool("h3-bin-points"), binPointsTool); +describe("h3 tools registry", () => { + it("no longer registers grid/bin tools under getH3Tool (moved to DGGS)", () => { + assert.equal(getH3Tool("h3-bin-points"), undefined); + assert.equal(getH3Tool("h3-grid"), undefined); assert.equal(getH3Tool("missing"), undefined); }); +}); - it("throws a clear error when duckdb is unavailable", async () => { - await assert.rejects( - () => - Promise.resolve( - createH3GridTool.run({ - layers: [], - parameters: { source: "viewport" }, - log: () => {}, - }), - ), - /requires DuckDB/, - ); - }); - - it("creates a grid from the map viewport", async () => { - const { ctx, added } = baseCtx([], { source: "viewport", resolution: 5 }); - await createH3GridTool.run(ctx); - assert.equal(added.length, 1); - assert.match(added[0], /res 5/); - }); - - it("rejects an antimeridian-crossing viewport", async () => { - const { ctx, added, logs } = baseCtx([], { - source: "viewport", - resolution: 4, - }); - ctx.viewportBounds = () => [170, 0, -170, 1]; // west >= east - await createH3GridTool.run(ctx); - assert.equal(added.length, 0); - assert.ok(logs.some((l) => /antimeridian/i.test(l))); - }); - - it("creates a grid from a manual bounding box without a layer", async () => { - const { ctx, added } = baseCtx([], { - source: "bbox", - west: 0, - south: 0, - east: 1, - north: 1, - resolution: 5, - }); - await createH3GridTool.run(ctx); - assert.equal(added.length, 1); - assert.match(added[0], /res 5/); - }); - - it("rejects a degenerate manual bounding box", async () => { - const { ctx, added, logs } = baseCtx([], { - source: "bbox", - west: 2, - south: 0, - east: 1, - north: 1, - resolution: 5, - }); - await createH3GridTool.run(ctx); - assert.equal(added.length, 0); - assert.ok(logs.some((l) => /west < east/i.test(l))); - }); - - it("rejects a manual bounding box with missing values", async () => { - const { ctx, added, logs } = baseCtx([], { - source: "bbox", - west: 0, - south: 0, - east: 1, - resolution: 5, - }); - await createH3GridTool.run(ctx); - assert.equal(added.length, 0); - assert.ok(logs.some((l) => /numeric/i.test(l))); - }); - - it("auto-suggests a resolution when none is given", async () => { - const { ctx, logs } = baseCtx([], { source: "viewport" }); - await createH3GridTool.run(ctx); - assert.ok(logs.some((l) => /suggested resolution/i.test(l))); - }); - - it("aborts when the requested resolution exceeds the hard cap", async () => { - const { ctx, added, logs } = baseCtx([polygonLayer()], { - source: "extent", - layer: "poly", - resolution: 15, - }); - await createH3GridTool.run(ctx); - assert.equal(added.length, 0); - assert.ok(logs.some((l) => /cap/i.test(l))); - }); - - it("polyfills a selected polygon layer and releases the registered source", async () => { - const { ctx, added, duckdb } = baseCtx([polygonLayer()], { - source: "polyfill", - layer: "poly", - resolution: 6, - }); - await createH3GridTool.run(ctx); - assert.equal(added.length, 1); - // The registered temp GeoJSON source must be released after the run. - assert.equal(duckdb.released.length, 1); - }); - - it("rejects polyfill of a non-polygon layer", async () => { - const { ctx, added, logs } = baseCtx([pointLayer()], { - source: "polyfill", - layer: "pts", - resolution: 6, - }); - await createH3GridTool.run(ctx); - assert.equal(added.length, 0); - assert.ok(logs.some((l) => /polygon/i.test(l))); - }); - - it("fills the extent of a non-polygon layer", async () => { - const { ctx, added } = baseCtx([pointLayer()], { - source: "extent", - layer: "pts", - resolution: 6, - }); - await createH3GridTool.run(ctx); - assert.equal(added.length, 1); - }); - - it("bins points and requires a field for non-count aggregates", async () => { - const missing = baseCtx([pointLayer()], { - layer: "pts", - aggOp: "sum", - resolution: 7, - }); - await binPointsTool.run(missing.ctx); - assert.equal(missing.added.length, 0); - assert.ok(missing.logs.some((l) => /field/i.test(l))); - - const ok = baseCtx([pointLayer()], { - layer: "pts", - aggOp: "count", - resolution: 7, - }); - await binPointsTool.run(ok.ctx); - assert.equal(ok.added.length, 1); - // The registered temp GeoJSON source must be released after the run. - assert.equal(ok.duckdb.released.length, 1); - }); - - it("rejects an unknown aggregate operation", async () => { - const { ctx, added, logs } = baseCtx([pointLayer()], { - layer: "pts", - aggOp: "median", - resolution: 7, - }); - await binPointsTool.run(ctx); - assert.equal(added.length, 0); - assert.ok(logs.some((l) => /unknown aggregate/i.test(l))); - }); - - it("logs a soft message and adds no layer when zero cells are produced", async () => { - const logs: string[] = []; - const added: string[] = []; - const ctx: ProcessingContext = { - layers: [], - parameters: { source: "viewport", resolution: 5 }, - log: (m) => logs.push(m), - addResultLayer: (name) => added.push(name), - viewportBounds: () => [0, 0, 1, 1], - duckdb: { - ensureExtensions: async () => {}, - registerGeoJson: async () => ({ - sql: "ST_Read('mock.geojson')", - release: async () => {}, - }), - query: async () => [], // no cells - }, - }; - await createH3GridTool.run(ctx); - assert.equal(added.length, 0); - assert.ok(logs.some((l) => /no h3 cells/i.test(l))); +describe("resolveVectorRerun H3 aliases", () => { + it("maps old H3 tool ids onto DGGS tools with dggsType h3", () => { + const grid = resolveVectorRerun("h3-grid", { resolution: 5, source: "viewport" }); + assert.equal(grid.toolId, "dggs-grid"); + assert.equal(grid.parameters.dggsType, "h3"); + assert.equal(grid.parameters.resolution, 5); + assert.ok(getVectorTool(grid.toolId)); + + const bin = resolveVectorRerun("h3-bin-points", { aggOp: "count" }); + assert.equal(bin.toolId, "dggs-bin"); + assert.equal(bin.parameters.dggsType, "h3"); + assert.ok(getVectorTool(bin.toolId)); + + // Existing dggsType is preserved; unknown ids pass through. + const kept = resolveVectorRerun("h3-grid", { dggsType: "s2" }); + assert.equal(kept.parameters.dggsType, "s2"); + const passthrough = resolveVectorRerun("buffer", { distance: 1 }); + assert.equal(passthrough.toolId, "buffer"); + assert.deepEqual(passthrough.parameters, { distance: 1 }); }); }); @@ -371,7 +109,10 @@ describe("h3 SQL + geometry builders", () => { // Include a single quote in the input so the test actually exercises the // doubling done by sqlStr (a malformed escape would break this assertion). const sql = buildGridFromWktSql("POLYGON((0 0, 1 0, 1 1, 0 0))'x", 7); - assert.match(sql, /h3_polygon_wkt_to_cells\('POLYGON\(\(0 0, 1 0, 1 1, 0 0\)\)''x', 7\)/); + assert.match( + sql, + /h3_polygon_wkt_to_cells_experimental\('POLYGON\(\(0 0, 1 0, 1 1, 0 0\)\)''x', 7, 'overlap'\)/, + ); assert.match(sql, /h3_h3_to_string\(cell\) AS h3/); assert.match( sql, @@ -379,6 +120,24 @@ describe("h3 SQL + geometry builders", () => { ); }); + it("normalizes bboxes wider than 180° of longitude to ±180", () => { + assert.deepEqual(normalizeLonLatBbox([-200, -60, 200, 60]), [-180, -60, 180, 60]); + assert.deepEqual(normalizeLonLatBbox([-100, -40, 100, 40]), [-180, -40, 180, 40]); + assert.deepEqual(normalizeLonLatBbox([10, 20, 30, 40]), [10, 20, 30, 40]); + assert.deepEqual(normalizeLonLatBbox([-190, -10, -10, 10]), [-180, -10, 180, 10]); + }); + + it("splits a full-world bbox into hemispheres for H3 polyfill", () => { + const narrow = buildGridFromBboxSql([0, 0, 1, 1], 5); + assert.match(narrow, /h3_polygon_wkt_to_cells_experimental/); + assert.doesNotMatch(narrow, /UNION ALL/); + + const world = buildGridFromBboxSql([-180, -60, 180, 60], 2); + assert.match(world, /UNION ALL/); + assert.match(world, /SELECT DISTINCT cell/); + assert.match(world, /POLYGON\(\(-180 -60, 0 -60, 0 60, -180 60, -180 -60\)\)/); + assert.match(world, /POLYGON\(\(0 -60, 180 -60, 180 60, 0 60, 0 -60\)\)/); + }); it("builds polyfill grid SQL that unions only polygon geometry and guards NULL", () => { const sql = buildGridFromSourceSql("ST_Read('a.geojson')", 8); assert.match(sql, /ST_Union_Agg\(geom\)/); @@ -388,43 +147,70 @@ describe("h3 SQL + geometry builders", () => { sql, /WHERE geom IS NOT NULL AND ST_GeometryType\(geom\) IN \('POLYGON', 'MULTIPOLYGON'\)/, ); - // A NULL union result (no polygons) is filtered before reaching the h3 fn. - assert.match(sql, /h3_polygon_wkt_to_cells\(wkt, 8\)/); - assert.match(sql, /FROM merged WHERE wkt IS NOT NULL/); + assert.match(sql, /WHERE wkt IS NOT NULL/); + assert.match(sql, /h3_polygon_wkt_to_cells_experimental\(wkt, 8, 'overlap'\)/); + assert.doesNotMatch(sql, /h3_compact_cells/); }); - it("builds bin SQL for count (no field), binning POINT and MULTIPOINT by centroid", () => { - const sql = buildBinSql("ST_Read('p.geojson')", 9, "count"); - assert.match(sql, /h3_latlng_to_cell\(ST_Y\(pt\), ST_X\(pt\), 9\)/); - assert.match(sql, /ST_Centroid\(geom\) AS pt/); - assert.match(sql, /count\(\*\) AS count/); - assert.doesNotMatch(sql, /AS value/); - assert.match(sql, /ST_GeometryType\(geom\) IN \('POINT', 'MULTIPOINT'\)/); + it("optionally compacts H3 grid cells after polyfill", () => { + const sql = buildGridFromBboxSql([0, 0, 1, 1], 5, true); + assert.match(sql, /h3_compact_cells\(cells\)/); + assert.match(sql, /h3_polygon_wkt_to_cells_experimental/); }); - it("builds bin SQL for an aggregate, mapping mean->avg and quoting the field", () => { - const sql = buildBinSql("ST_Read('p.geojson')", 9, "mean", "pop"); - assert.match(sql, /avg\(CAST\("pop" AS DOUBLE\)\) AS value/); - assert.match(sql, /count, value,/); - // The field must be in the SELECT list, not appended after the WHERE clause. - assert.match(sql, /SELECT ST_Centroid\(geom\) AS pt, "pop" FROM/); - assert.doesNotMatch(sql, /MULTIPOINT'\), "pop"/); + it("builds bin SQL for count and for a named aggregate", () => { + const countSql = buildBinSql("ST_Read('p.geojson')", 5, "count"); + assert.match(countSql, /h3_latlng_to_cell/); + assert.match(countSql, /count\(\*\) AS count/); + assert.doesNotMatch(countSql, / AS value/); + + const sumSql = buildBinSql("ST_Read('p.geojson')", 5, "sum", 'pop"x'); + // Field name is double-quote escaped. + assert.match(sumSql, /sum\(CAST\("pop""x" AS DOUBLE\)\) AS value/); }); - it("converts result rows to a FeatureCollection with h3/count/value props", () => { + it("converts result rows to a FeatureCollection with h3 props", () => { const fc = rowsToFeatureCollection([ { - h3: "8928308280fffff", - count: 3n, - value: 12.5, + h3: "abc", + count: 3, + value: 1.5, geojson: '{"type":"Polygon","coordinates":[[[0,0],[1,0],[1,1],[0,0]]]}', }, - { h3: "x", count: 1, geojson: null }, + { h3: "skip", geojson: 12 }, ]); assert.equal(fc.features.length, 1); - assert.equal(fc.features[0].properties?.h3, "8928308280fffff"); + assert.equal(fc.features[0].properties?.h3, "abc"); assert.equal(fc.features[0].properties?.count, 3); - assert.equal(fc.features[0].properties?.value, 12.5); - assert.equal(fc.features[0].geometry?.type, "Polygon"); + assert.equal(fc.features[0].properties?.value, 1.5); + }); + + it("unwraps antimeridian rings by default and can leave them wrapped", () => { + const raw = [ + [170, 0], + [-170, 0], + [-170, 1], + [170, 1], + [170, 0], + ]; + const fixed = rowsToFeatureCollection([ + { h3: "x", geojson: JSON.stringify({ type: "Polygon", coordinates: [raw] }) }, + ]); + assert.deepEqual((fixed.features[0].geometry as { coordinates: number[][][] }).coordinates[0], [ + [170, 0], + [190, 0], + [190, 1], + [170, 1], + [170, 0], + ]); + + const left = rowsToFeatureCollection( + [{ h3: "x", geojson: JSON.stringify({ type: "Polygon", coordinates: [raw] }) }], + false, + ); + assert.deepEqual( + (left.features[0].geometry as { coordinates: number[][][] }).coordinates[0], + raw, + ); }); }); diff --git a/tests/processing-history.test.ts b/tests/processing-history.test.ts index a7116f4d5..81f190037 100644 --- a/tests/processing-history.test.ts +++ b/tests/processing-history.test.ts @@ -194,4 +194,27 @@ describe("normalizeProcessingHistory", () => { assert.equal(runs?.length, MAX_PROCESSING_HISTORY); assert.equal(runs?.[0].id, "run-10"); }); + + it("migrates legacy H3 tool ids to DGGS tools with dggsType h3", () => { + const runs = normalizeProcessingHistory([ + makeRun({ + id: "legacy-grid", + toolId: "h3-grid", + toolName: "H3 Grid", + parameters: { resolution: 5 }, + }), + makeRun({ + id: "legacy-bin", + toolId: "h3-bin-points", + parameters: { resolution: 4, dggsType: "s2" }, + }), + ]); + assert.equal(runs?.length, 2); + assert.equal(runs?.[0].toolId, "dggs-grid"); + assert.equal(runs?.[0].parameters.dggsType, "h3"); + assert.equal(runs?.[0].parameters.resolution, 5); + assert.equal(runs?.[1].toolId, "dggs-bin"); + // Explicit dggsType from the saved run is preserved. + assert.equal(runs?.[1].parameters.dggsType, "s2"); + }); }); diff --git a/tests/s2-tools.test.ts b/tests/s2-tools.test.ts new file mode 100644 index 000000000..d3ff7be8b --- /dev/null +++ b/tests/s2-tools.test.ts @@ -0,0 +1,143 @@ +import assert from "node:assert/strict"; +import { describe, it } from "node:test"; +import { + S2_AVG_AREA_KM2, + S2_HARD_CAP, + S2_MAX_TOOL_RES, + binPointsToS2, + compactS2Tokens, + estimateS2CellCount, + estimateS2ExpandCount, + expandS2Tokens, + s2CellAtLonLat, + s2GridFromBbox, + s2GridFromFeatureCollection, + suggestS2Resolution, +} from "../packages/processing/src/s2-tools"; +import { bboxAreaKm2 } from "../packages/processing/src/h3-tools"; +import { s2 } from "s2js"; + +describe("s2 resolution math", () => { + it("exposes 31 average-area entries (res 0..30), strictly decreasing", () => { + assert.equal(S2_AVG_AREA_KM2.length, 31); + assert.equal(S2_MAX_TOOL_RES, 30); + for (let r = 1; r < 31; r += 1) { + assert.ok(S2_AVG_AREA_KM2[r]! < S2_AVG_AREA_KM2[r - 1]!); + } + }); + + it("suggests a coarser resolution for larger areas", () => { + const big = bboxAreaKm2([-10, -10, 10, 10]); + const tiny = bboxAreaKm2([0, 0, 0.001, 0.001]); + const rBig = suggestS2Resolution(big); + const rTiny = suggestS2Resolution(tiny); + assert.ok(rBig < rTiny); + assert.ok(estimateS2CellCount(big, rBig) <= 10_000); + assert.ok(estimateS2CellCount(big, 30) > S2_HARD_CAP); + assert.ok(Number.isFinite(estimateS2CellCount(big, 20))); + }); +}); + +describe("s2 compact / expand", () => { + function fourSiblingsAtLevel(level: number): string[] { + const leaf = s2.cellid.fromLatLng(s2.LatLng.fromDegrees(10, 10)); + const parent = s2.cellid.parent(leaf, level - 1); + const kids: string[] = []; + let id = s2.cellid.childBegin(parent); + for (let i = 0; i < 4; i += 1) { + kids.push(s2.cellid.toToken(id)); + id = s2.cellid.next(id); + } + return kids; + } + + it("compacts four siblings into their parent", () => { + const kids = fourSiblingsAtLevel(6); + const compacted = compactS2Tokens(kids); + assert.equal(compacted.length, 1); + assert.equal(s2.cellid.level(s2.cellid.fromToken(compacted[0]!)), 5); + }); + + it("expands a parent back to four children at the target level", () => { + const kids = fourSiblingsAtLevel(6); + const [parent] = compactS2Tokens(kids); + const expanded = expandS2Tokens([parent!], 6); + assert.equal(expanded.length, 4); + assert.equal(estimateS2ExpandCount([parent!], 6), 4); + assert.deepEqual(new Set(expanded), new Set(kids)); + }); + + it("optionally compacts a bbox grid after covering", () => { + const plain = s2GridFromBbox([0, 0, 2, 2], 6); + const compacted = s2GridFromBbox([0, 0, 2, 2], 6, { compact: true }); + assert.ok(compacted.features.length > 0); + assert.ok(compacted.features.length <= plain.features.length); + }); +}); + +describe("s2 grid / bin", () => { + it("covers a small bbox with polygon features carrying s2 tokens", () => { + const fc = s2GridFromBbox([0, 0, 1, 1], 8); + assert.ok(fc.features.length > 0); + assert.ok(fc.features.length < 500); + const expected = s2CellAtLonLat(0.5, 0.5, 8); + assert.ok(fc.features.some((f) => f.properties?.s2 === expected)); + assert.equal(fc.features[0]!.geometry.type, "Polygon"); + }); + + it("polyfills a polygon feature collection", () => { + const fc = s2GridFromFeatureCollection( + { + type: "FeatureCollection", + features: [ + { + type: "Feature", + properties: {}, + geometry: { + type: "Polygon", + coordinates: [ + [ + [0, 0], + [1, 0], + [1, 1], + [0, 1], + [0, 0], + ], + ], + }, + }, + ], + }, + 8, + ); + assert.ok(fc.features.length > 0); + }); + + it("bins points into S2 cells with count", () => { + const fc = binPointsToS2( + { + type: "FeatureCollection", + features: [ + { + type: "Feature", + properties: { pop: 10 }, + geometry: { type: "Point", coordinates: [0.5, 0.5] }, + }, + { + type: "Feature", + properties: { pop: 5 }, + geometry: { type: "Point", coordinates: [0.51, 0.51] }, + }, + ], + }, + 10, + "sum", + "pop", + ); + assert.ok(fc.features.length >= 1); + const total = fc.features.reduce((n, f) => n + Number(f.properties?.count ?? 0), 0); + assert.equal(total, 2); + const sum = fc.features.reduce((n, f) => n + Number(f.properties?.value ?? 0), 0); + assert.equal(sum, 15); + }); +});