From 0ee99b7fc8c5e9d7c2c24ae555e110fbe4063581 Mon Sep 17 00:00:00 2001 From: Andrej Fast Date: Tue, 14 Jul 2026 16:31:54 +0200 Subject: [PATCH 1/4] feat(lib): added web-workers for regular grids --- src/lib/grids/regularGridData.worker.ts | 44 +++++++++ src/lib/grids/regularGridDataWorkerClient.ts | 93 +++++++++++++++++++ .../grids/regularGridDataWorkerProtocol.ts | 36 +++++++ src/ui/grids/Regular.vue | 59 +++++++----- 4 files changed, 209 insertions(+), 23 deletions(-) create mode 100644 src/lib/grids/regularGridData.worker.ts create mode 100644 src/lib/grids/regularGridDataWorkerClient.ts create mode 100644 src/lib/grids/regularGridDataWorkerProtocol.ts diff --git a/src/lib/grids/regularGridData.worker.ts b/src/lib/grids/regularGridData.worker.ts new file mode 100644 index 00000000..9d9b206d --- /dev/null +++ b/src/lib/grids/regularGridData.worker.ts @@ -0,0 +1,44 @@ +/// + +import { ZarrDataManager } from "@/lib/data/ZarrDataManager.ts"; +import { + RegularGridDataWorkerMessageType, + type TRegularGridDataWorkerRequest, + type TRegularGridDataWorkerResponse, +} from "@/lib/grids/regularGridDataWorkerProtocol.ts"; + +const workerScope = self as unknown as DedicatedWorkerGlobalScope; + +workerScope.onmessage = async ( + event: MessageEvent +) => { + const { requestId, source, variable, format, selection } = event.data; + try { + const array = await ZarrDataManager.getVariableInfo( + source, + variable, + format + ); + const chunk = await ZarrDataManager.getVariableDataFromArray( + array, + selection + ); + const response: TRegularGridDataWorkerResponse = { + requestId, + type: RegularGridDataWorkerMessageType.RESULT, + data: chunk.data, + }; + const transfer = + ArrayBuffer.isView(chunk.data) && chunk.data.buffer instanceof ArrayBuffer + ? [chunk.data.buffer] + : []; + workerScope.postMessage(response, transfer); + } catch (error) { + const response: TRegularGridDataWorkerResponse = { + requestId, + type: RegularGridDataWorkerMessageType.ERROR, + message: error instanceof Error ? error.message : String(error), + }; + workerScope.postMessage(response); + } +}; diff --git a/src/lib/grids/regularGridDataWorkerClient.ts b/src/lib/grids/regularGridDataWorkerClient.ts new file mode 100644 index 00000000..68bfa7c6 --- /dev/null +++ b/src/lib/grids/regularGridDataWorkerClient.ts @@ -0,0 +1,93 @@ +import type * as zarr from "zarrita"; + +import { + RegularGridDataWorkerMessageType, + type TRegularGridDataWorkerRequest, + type TRegularGridDataWorkerResponse, +} from "./regularGridDataWorkerProtocol.ts"; + +import type { TDataSource, TZarrFormat } from "@/lib/types/GlobeTypes.ts"; + +type TRegularGridDataRequest = { + source: Pick; + variable: string; + format: TZarrFormat; + selection: (number | null | zarr.Slice)[]; +}; + +type TPendingRequest = { + resolve: (data: zarr.TypedArray) => void; + reject: (reason?: unknown) => void; +}; + +let worker: Worker | null = null; +let nextRequestId = 0; +const pendingRequests = new Map(); + +function rejectPendingRequests(error: Error) { + for (const pending of pendingRequests.values()) { + pending.reject(error); + } + pendingRequests.clear(); +} + +function handleWorkerMessage(message: TRegularGridDataWorkerResponse) { + const pending = pendingRequests.get(message.requestId); + if (!pending) { + return; + } + pendingRequests.delete(message.requestId); + if (message.type === RegularGridDataWorkerMessageType.ERROR) { + pending.reject(new Error(message.message)); + return; + } + pending.resolve(message.data); +} + +function getWorker() { + if (worker) { + return worker; + } + worker = new Worker(new URL("./regularGridData.worker.ts", import.meta.url), { + type: "module", + }); + worker.onmessage = (event: MessageEvent) => { + handleWorkerMessage(event.data); + }; + worker.onerror = (event) => { + rejectPendingRequests(new Error(event.message)); + worker?.terminate(); + worker = null; + }; + return worker; +} + +export function getRegularGridVariableData(request: TRegularGridDataRequest) { + const requestId = ++nextRequestId; + const message: TRegularGridDataWorkerRequest = { + requestId, + type: RegularGridDataWorkerMessageType.GET_DATA, + source: { + store: request.source.store, + dataset: request.source.dataset, + }, + variable: request.variable, + format: request.format, + selection: request.selection.map((selection) => + typeof selection === "object" && selection !== null + ? { ...selection } + : selection + ), + }; + + return new Promise>((resolve, reject) => { + pendingRequests.set(requestId, { resolve, reject }); + getWorker().postMessage(message); + }); +} + +export function terminateRegularGridDataWorker() { + worker?.terminate(); + worker = null; + rejectPendingRequests(new Error("Regular grid data worker terminated.")); +} diff --git a/src/lib/grids/regularGridDataWorkerProtocol.ts b/src/lib/grids/regularGridDataWorkerProtocol.ts new file mode 100644 index 00000000..acf6b7f2 --- /dev/null +++ b/src/lib/grids/regularGridDataWorkerProtocol.ts @@ -0,0 +1,36 @@ +import type * as zarr from "zarrita"; + +import type { TDataSource, TZarrFormat } from "@/lib/types/GlobeTypes.ts"; + +export const RegularGridDataWorkerMessageType = { + GET_DATA: "getData", + RESULT: "result", + ERROR: "error", +} as const; + +type TRegularGridDataWorkerMessageType = + (typeof RegularGridDataWorkerMessageType)[keyof typeof RegularGridDataWorkerMessageType]; + +export type TRegularGridDataWorkerRequest = { + requestId: number; + type: typeof RegularGridDataWorkerMessageType.GET_DATA; + source: Pick; + variable: string; + format: TZarrFormat; + selection: (number | null | zarr.Slice)[]; +}; + +type TRegularGridDataWorkerResponseBase = { + requestId: number; + type: TRegularGridDataWorkerMessageType; +}; + +export type TRegularGridDataWorkerResponse = + | (TRegularGridDataWorkerResponseBase & { + type: typeof RegularGridDataWorkerMessageType.RESULT; + data: zarr.TypedArray; + }) + | (TRegularGridDataWorkerResponseBase & { + type: typeof RegularGridDataWorkerMessageType.ERROR; + message: string; + }); diff --git a/src/ui/grids/Regular.vue b/src/ui/grids/Regular.vue index c3ec8949..d9c08e4e 100644 --- a/src/ui/grids/Regular.vue +++ b/src/ui/grids/Regular.vue @@ -27,6 +27,10 @@ import { decodeVariableDataAndGetBounds, } from "@/lib/data/variableDecoding.ts"; import { ZarrDataManager } from "@/lib/data/ZarrDataManager.ts"; +import { + getRegularGridVariableData, + terminateRegularGridDataWorker, +} from "@/lib/grids/regularGridDataWorkerClient.ts"; import { GridTextureExportUserDataKey, getRegularLatLonGridBounds, @@ -124,20 +128,24 @@ const { datasourceUpdate } = useGridDataLoader({ const isLatOnly = ref(false); +function getDimensionData( + grid: TSources["levels"][0]["grid"], + dimensionName: string +) { + return ZarrDataManager.getVariableData( + grid, + ZarrDataManager.resolveVariablePath(varnameSelector.value, dimensionName) + ); +} + async function fetchProjectedXYDims( grid: TSources["levels"][0]["grid"], xDim: string, yDim: string ) { const [xData, yData] = await Promise.all([ - ZarrDataManager.getVariableData( - grid, - ZarrDataManager.resolveVariablePath(varnameSelector.value, xDim) - ), - ZarrDataManager.getVariableData( - grid, - ZarrDataManager.resolveVariablePath(varnameSelector.value, yDim) - ), + getDimensionData(grid, xDim), + getDimensionData(grid, yDim), ]); const crsWkt = await getCRSWkt(props.datasources!, varnameSelector.value); const converted = projectedAxisCoordinatesToLonLat( @@ -172,24 +180,15 @@ async function getDims() { if (isProjectedXY) { await fetchProjectedXYDims(grid, lastDim, secondLastDim); } else if (latOnlyCheck) { - const latitudesData = await ZarrDataManager.getVariableData( - grid, - ZarrDataManager.resolveVariablePath(varnameSelector.value, lastDim) - ); + const latitudesData = await getDimensionData(grid, lastDim); latitudes.value = latitudesData.data as Float32Array; longitudes.value = Float32Array.from({ length: 360 }, (_, i) => i - 179.5); } else { const latName = secondLastDim; const lonName = lastDim; const [latitudesData, longitudesData] = await Promise.all([ - ZarrDataManager.getVariableData( - grid, - ZarrDataManager.resolveVariablePath(varnameSelector.value, latName) - ), - ZarrDataManager.getVariableData( - grid, - ZarrDataManager.resolveVariablePath(varnameSelector.value, lonName) - ), + getDimensionData(grid, latName), + getDimensionData(grid, lonName), ]); const myLongitudes = longitudesData.data as Float32Array; const myLatitudes = latitudesData.data as Float32Array; @@ -796,6 +795,20 @@ async function buildDimensionConfig( ); } +function fetchRegularGridVariableData( + selection: (number | null | zarr.Slice)[] +) { + return getRegularGridVariableData({ + source: ZarrDataManager.getDatasetSource( + props.datasources!, + varnameSelector.value + ), + variable: varnameSelector.value, + format: props.datasources!.zarr_format, + selection, + }); +} + function setMeshMaterials(material: THREE.ShaderMaterial) { if (meshes.length === 0) { disposeMaterial(material); @@ -824,9 +837,8 @@ async function fetchAndRenderData( ) { const { dimensionRanges, indices } = await buildDimensionConfig(datavar); - const rawData = castDataVarToFloat32( - (await ZarrDataManager.getVariableDataFromArray(datavar, indices)).data - ); + const variableData = await fetchRegularGridVariableData(indices); + const rawData = castDataVarToFloat32(variableData); const { min, max, missingValue, fillValue } = decodeVariableDataAndGetBounds( datavar, @@ -864,6 +876,7 @@ onBeforeMount(async () => { }); onBeforeUnmount(() => { + terminateRegularGridDataWorker(); for (const mesh of meshes) { mesh.geometry.dispose(); getScene()?.remove(mesh); From 424145debd297835dc051cdaf9807eb3eab75c8e Mon Sep 17 00:00:00 2001 From: Andrej Fast Date: Wed, 15 Jul 2026 11:01:15 +0200 Subject: [PATCH 2/4] chore: implemented web workers for curvilinear, gaussian reduced and both irregular grids --- src/lib/data/irregularGridHelpers.ts | 25 +- src/lib/grids/curvilinear.worker.ts | 113 +++ src/lib/grids/curvilinearCalculations.ts | 420 ++++++++++ src/lib/grids/curvilinearWorkerClient.ts | 69 ++ src/lib/grids/curvilinearWorkerProtocol.ts | 29 + src/lib/grids/gaussianReduced.worker.ts | 88 ++ src/lib/grids/gaussianReducedCalculations.ts | 195 +++++ src/lib/grids/gaussianReducedWorkerClient.ts | 74 ++ .../grids/gaussianReducedWorkerProtocol.ts | 24 + ...rGridData.worker.ts => gridData.worker.ts} | 20 +- ...orkerClient.ts => gridDataWorkerClient.ts} | 28 +- src/lib/grids/gridDataWorkerProtocol.ts | 36 + src/lib/grids/gridGeometryWorkerClient.ts | 151 ++++ src/lib/grids/gridGeometryWorkerProtocol.ts | 46 + src/lib/grids/gridGeometryWorkerUtils.ts | 98 +++ src/lib/grids/gridWorkerCalculations.ts | 27 + src/lib/grids/gridWorkerTypes.ts | 26 + src/lib/grids/irregular.worker.ts | 78 ++ src/lib/grids/irregularCalculations.ts | 129 +++ src/lib/grids/irregularDelaunay.worker.ts | 114 +++ .../grids/irregularDelaunayCalculations.ts | 312 +++++++ .../grids/irregularDelaunayWorkerClient.ts | 61 ++ .../grids/irregularDelaunayWorkerProtocol.ts | 29 + src/lib/grids/irregularWorkerClient.ts | 70 ++ src/lib/grids/irregularWorkerProtocol.ts | 30 + .../grids/regularGridDataWorkerProtocol.ts | 36 - src/lib/grids/serializedGeoSampleIndex.ts | 37 + src/ui/grids/Curvilinear.vue | 783 +++--------------- src/ui/grids/GaussianReduced.vue | 306 +++---- src/ui/grids/Irregular.vue | 349 +++----- src/ui/grids/IrregularDelaunay.vue | 629 +++----------- src/ui/grids/Regular.vue | 10 +- .../lib/grids/curvilinearCalculations.test.ts | 92 ++ .../grids/gaussianReducedCalculations.test.ts | 47 ++ .../lib/grids/irregularCalculations.test.ts | 63 ++ .../irregularDelaunayCalculations.test.ts | 69 ++ 36 files changed, 3004 insertions(+), 1709 deletions(-) create mode 100644 src/lib/grids/curvilinear.worker.ts create mode 100644 src/lib/grids/curvilinearCalculations.ts create mode 100644 src/lib/grids/curvilinearWorkerClient.ts create mode 100644 src/lib/grids/curvilinearWorkerProtocol.ts create mode 100644 src/lib/grids/gaussianReduced.worker.ts create mode 100644 src/lib/grids/gaussianReducedCalculations.ts create mode 100644 src/lib/grids/gaussianReducedWorkerClient.ts create mode 100644 src/lib/grids/gaussianReducedWorkerProtocol.ts rename src/lib/grids/{regularGridData.worker.ts => gridData.worker.ts} (64%) rename src/lib/grids/{regularGridDataWorkerClient.ts => gridDataWorkerClient.ts} (69%) create mode 100644 src/lib/grids/gridDataWorkerProtocol.ts create mode 100644 src/lib/grids/gridGeometryWorkerClient.ts create mode 100644 src/lib/grids/gridGeometryWorkerProtocol.ts create mode 100644 src/lib/grids/gridGeometryWorkerUtils.ts create mode 100644 src/lib/grids/gridWorkerCalculations.ts create mode 100644 src/lib/grids/gridWorkerTypes.ts create mode 100644 src/lib/grids/irregular.worker.ts create mode 100644 src/lib/grids/irregularCalculations.ts create mode 100644 src/lib/grids/irregularDelaunay.worker.ts create mode 100644 src/lib/grids/irregularDelaunayCalculations.ts create mode 100644 src/lib/grids/irregularDelaunayWorkerClient.ts create mode 100644 src/lib/grids/irregularDelaunayWorkerProtocol.ts create mode 100644 src/lib/grids/irregularWorkerClient.ts create mode 100644 src/lib/grids/irregularWorkerProtocol.ts delete mode 100644 src/lib/grids/regularGridDataWorkerProtocol.ts create mode 100644 src/lib/grids/serializedGeoSampleIndex.ts mode change 100755 => 100644 src/ui/grids/Curvilinear.vue create mode 100644 tests/unit/lib/grids/curvilinearCalculations.test.ts create mode 100644 tests/unit/lib/grids/gaussianReducedCalculations.test.ts create mode 100644 tests/unit/lib/grids/irregularCalculations.test.ts create mode 100644 tests/unit/lib/grids/irregularDelaunayCalculations.test.ts diff --git a/src/lib/data/irregularGridHelpers.ts b/src/lib/data/irregularGridHelpers.ts index 72485604..b9dd625d 100644 --- a/src/lib/data/irregularGridHelpers.ts +++ b/src/lib/data/irregularGridHelpers.ts @@ -46,9 +46,22 @@ export function reconcileCoordinates( ): { latitudes: Float32Array; longitudes: Float32Array } { const latitudes = latitudesVar.data as Float32Array; const longitudes = longitudesVar.data as Float32Array; - const latShape = latitudesVar.shape; - const lonShape = longitudesVar.shape; + return reconcileCoordinateArrays( + latitudes, + longitudes, + latitudesVar.shape, + longitudesVar.shape, + dataLength + ); +} +export function reconcileCoordinateArrays( + latitudes: Float32Array, + longitudes: Float32Array, + latitudeShape: number[], + longitudeShape: number[], + dataLength: number +): { latitudes: Float32Array; longitudes: Float32Array } { // Case 1: All arrays already have same length (typical irregular grid) if (latitudes.length === dataLength && longitudes.length === dataLength) { return { latitudes, longitudes }; @@ -56,9 +69,9 @@ export function reconcileCoordinates( // Case 2: 2D coordinate arrays - just ensure they're flattened and match data // (zarr.Chunk.data should already be flattened) - if (latShape.length === 2 && lonShape.length === 2) { - const latTotal = latShape[0] * latShape[1]; - const lonTotal = lonShape[0] * lonShape[1]; + if (latitudeShape.length === 2 && longitudeShape.length === 2) { + const latTotal = latitudeShape[0] * latitudeShape[1]; + const lonTotal = longitudeShape[0] * longitudeShape[1]; if (latTotal === dataLength && lonTotal === dataLength) { // Already flattened by zarr, should match return { latitudes, longitudes }; @@ -66,7 +79,7 @@ export function reconcileCoordinates( } // Case 3: 1D separate coordinates - expand via meshgrid - if (latShape.length === 1 && lonShape.length === 1) { + if (latitudeShape.length === 1 && longitudeShape.length === 1) { const product = latitudes.length * longitudes.length; if (product === dataLength) { return expandCoordinatesToMeshgrid(latitudes, longitudes, dataLength); diff --git a/src/lib/grids/curvilinear.worker.ts b/src/lib/grids/curvilinear.worker.ts new file mode 100644 index 00000000..308a8250 --- /dev/null +++ b/src/lib/grids/curvilinear.worker.ts @@ -0,0 +1,113 @@ +/// + +import { + buildCurvilinearBatch, + buildCurvilinearHoverIndexData, + detectCurvilinearColumnPeriodicity, + detectCurvilinearLongitudeFlip, + getCurvilinearBatchCount, + type TCurvilinearGrid, +} from "./curvilinearCalculations.ts"; +import type { + TCurvilinearWorkerRequest, + TCurvilinearWorkerResponse, +} from "./curvilinearWorkerProtocol.ts"; + +import { GridGeometryWorkerMessageType } from "@/lib/grids/gridGeometryWorkerProtocol.ts"; +import { + postGridGeometryBatch, + postGridGeometryHoverIndex, + postGridGeometryResponse, +} from "@/lib/grids/gridGeometryWorkerUtils.ts"; +import { ProjectionHelper } from "@/lib/projection/projectionUtils.ts"; + +const workerScope = self as unknown as DedicatedWorkerGlobalScope; + +function postResponse( + response: TCurvilinearWorkerResponse, + transfer: Transferable[] = [] +) { + postGridGeometryResponse(workerScope, response, transfer); +} + +function createGrid(request: TCurvilinearWorkerRequest): TCurvilinearGrid { + return { + latitudes: request.latitudes, + longitudes: request.longitudes, + data: request.data, + nj: request.nj, + ni: request.ni, + shouldFlipLongitude: detectCurvilinearLongitudeFlip( + request.longitudes, + request.latitudes, + request.missingValue, + request.fillValue, + request.nj, + request.ni + ), + isPeriodicI: detectCurvilinearColumnPeriodicity( + request.longitudes, + request.nj, + request.ni + ), + }; +} + +function postHoverIndex(requestId: number, grid: TCurvilinearGrid) { + postGridGeometryHoverIndex( + workerScope, + requestId, + buildCurvilinearHoverIndexData(grid) + ); +} + +function postBatches( + request: TCurvilinearWorkerRequest, + grid: TCurvilinearGrid, + totalBatches: number +) { + const projection = new ProjectionHelper( + request.projectionType, + request.projectionCenter + ); + for (let batchIndex = 0; batchIndex < totalBatches; batchIndex++) { + const batch = buildCurvilinearBatch( + grid, + batchIndex, + request.batchSize, + projection + ); + postGridGeometryBatch(workerScope, request.requestId, batch); + } +} + +function buildGrid(request: TCurvilinearWorkerRequest) { + const grid = createGrid(request); + const totalBatches = getCurvilinearBatchCount(grid, request.batchSize); + postResponse({ + requestId: request.requestId, + type: GridGeometryWorkerMessageType.METADATA, + metadata: { + totalBatches, + shouldFlipLongitude: grid.shouldFlipLongitude, + }, + }); + postHoverIndex(request.requestId, grid); + postBatches(request, grid, totalBatches); + postResponse({ + requestId: request.requestId, + type: GridGeometryWorkerMessageType.DONE, + }); +} + +workerScope.onmessage = (event: MessageEvent) => { + try { + buildGrid(event.data); + } catch (error) { + postResponse({ + requestId: event.data.requestId, + type: GridGeometryWorkerMessageType.ERROR, + message: error instanceof Error ? error.message : String(error), + }); + } +}; diff --git a/src/lib/grids/curvilinearCalculations.ts b/src/lib/grids/curvilinearCalculations.ts new file mode 100644 index 00000000..69d74b94 --- /dev/null +++ b/src/lib/grids/curvilinearCalculations.ts @@ -0,0 +1,420 @@ +import { buildSerializedGeoSampleIndexData } from "./gridWorkerCalculations.ts"; +import type { + TGridGeometryBatch, + TSerializedGeoSampleIndexData, +} from "./gridWorkerTypes.ts"; + +import { ProjectionHelper } from "@/lib/projection/projectionUtils.ts"; + +type TCoordinate = { lat: number; lon: number }; +type TBoundaryFunction = ( + row: number, + fromColumn: number, + toColumn: number +) => TCoordinate; +type TMidpointFunction = (fromColumn: number, toColumn: number) => TCoordinate; + +export type TCurvilinearGrid = { + latitudes: Float32Array; + longitudes: Float32Array; + data: Float32Array; + nj: number; + ni: number; + shouldFlipLongitude: boolean; + isPeriodicI: boolean; +}; + +function isMissingCoordinate( + value: number, + missingValue: number, + fillValue: number +) { + return Number.isNaN(value) || value === missingValue || value === fillValue; +} + +export function detectCurvilinearLongitudeFlip( + longitudes: Float32Array, + latitudes: Float32Array, + missingValue: number, + fillValue: number, + nj: number, + ni: number +) { + for (let row = 0; row < nj - 1; row++) { + const columnLimit = Math.min(ni - 1, 10); + for (let column = 0; column < columnLimit; column++) { + const topLeft = row * ni + column; + const topRight = topLeft + 1; + const bottomLeft = (row + 1) * ni + column; + const bottomRight = bottomLeft + 1; + if ( + isMissingCoordinate(longitudes[topLeft], missingValue, fillValue) || + isMissingCoordinate(longitudes[topRight], missingValue, fillValue) || + isMissingCoordinate(longitudes[bottomLeft], missingValue, fillValue) || + isMissingCoordinate(longitudes[bottomRight], missingValue, fillValue) + ) { + continue; + } + const longitudeI = longitudes[topRight] - longitudes[topLeft]; + const latitudeI = latitudes[topRight] - latitudes[topLeft]; + const longitudeJ = longitudes[bottomLeft] - longitudes[topLeft]; + const latitudeJ = latitudes[bottomLeft] - latitudes[topLeft]; + return longitudeI * latitudeJ - latitudeI * longitudeJ < 0; + } + } + return false; +} + +export function detectCurvilinearColumnPeriodicity( + longitudes: Float32Array, + nj: number, + ni: number +) { + if (ni < 3) { + return false; + } + const sampleRow = Math.floor(nj / 2); + const sampleCount = Math.min(ni - 1, 10); + let spacingTotal = 0; + let validSpacingCount = 0; + for (let column = 0; column < sampleCount; column++) { + const longitude = longitudes[sampleRow * ni + column]; + const nextLongitude = longitudes[sampleRow * ni + column + 1]; + const spacing = Math.abs(((nextLongitude - longitude + 540) % 360) - 180); + if (spacing > 0) { + spacingTotal += spacing; + validSpacingCount++; + } + } + if (validSpacingCount === 0) { + return false; + } + const firstLongitude = longitudes[sampleRow * ni]; + const lastLongitude = longitudes[sampleRow * ni + ni - 1]; + const wrapGap = Math.abs( + ((firstLongitude - lastLongitude + 540) % 360) - 180 + ); + return wrapGap < (spacingTotal / validSpacingCount) * 4; +} + +function getNextColumn(column: number, ni: number, flipLongitude: boolean) { + if (flipLongitude) { + return column === 0 ? ni - 1 : column - 1; + } + return (column + 1) % ni; +} + +function getPreviousColumn(column: number, ni: number, flipLongitude: boolean) { + if (flipLongitude) { + return (column + 1) % ni; + } + return column === 0 ? ni - 1 : column - 1; +} + +function getCircularMeanLongitude( + first: number, + second: number, + third: number, + fourth: number +) { + const radians = Math.PI / 180; + const sine = + Math.sin(first * radians) + + Math.sin(second * radians) + + Math.sin(third * radians) + + Math.sin(fourth * radians); + const cosine = + Math.cos(first * radians) + + Math.cos(second * radians) + + Math.cos(third * radians) + + Math.cos(fourth * radians); + return Math.atan2(sine / 4, cosine / 4) / radians; +} + +function getCellCenter( + grid: TCurvilinearGrid, + topLeft: number, + topRight: number, + bottomLeft: number, + bottomRight: number +) { + return { + lat: + (grid.latitudes[topLeft] + + grid.latitudes[topRight] + + grid.latitudes[bottomLeft] + + grid.latitudes[bottomRight]) / + 4, + lon: getCircularMeanLongitude( + grid.longitudes[topLeft], + grid.longitudes[topRight], + grid.longitudes[bottomLeft], + grid.longitudes[bottomRight] + ), + }; +} + +function getBoundaryPoint( + grid: TCurvilinearGrid, + row: number, + fromColumn: number, + toColumn: number +) { + return getCellCenter( + grid, + row * grid.ni + fromColumn, + row * grid.ni + toColumn, + (row + 1) * grid.ni + fromColumn, + (row + 1) * grid.ni + toColumn + ); +} + +function getRowMidpoint( + grid: TCurvilinearGrid, + row: number, + fromColumn: number, + toColumn: number +) { + return getCellCenter( + grid, + row * grid.ni + fromColumn, + row * grid.ni + toColumn, + row * grid.ni + fromColumn, + row * grid.ni + toColumn + ); +} + +function mirrorBoundaryPoint(center: TCoordinate, boundary: TCoordinate) { + const longitudeOffset = ((boundary.lon - center.lon + 540) % 360) - 180; + return { + lat: 2 * center.lat - boundary.lat, + lon: center.lon - longitudeOffset, + }; +} + +function getPreviousBoundaryCorners( + row: number, + column: number, + nextColumn: number, + boundary: TBoundaryFunction, + midpoint: TMidpointFunction +) { + const forwardNext = boundary(row, column, nextColumn); + const forwardPrevious = mirrorBoundaryPoint( + midpoint(column, column), + forwardNext + ); + const backwardNext = + row === 0 + ? mirrorBoundaryPoint(midpoint(column, nextColumn), forwardNext) + : boundary(row - 1, column, nextColumn); + return { + forwardPrevious, + forwardNext, + backwardPrevious: mirrorBoundaryPoint( + midpoint(column, column), + backwardNext + ), + backwardNext, + }; +} + +function getNextBoundaryCorners( + row: number, + column: number, + previousColumn: number, + boundary: TBoundaryFunction, + midpoint: TMidpointFunction +) { + const forwardPrevious = boundary(row, previousColumn, column); + const forwardNext = mirrorBoundaryPoint( + midpoint(column, column), + forwardPrevious + ); + const backwardPrevious = + row === 0 + ? mirrorBoundaryPoint(midpoint(previousColumn, column), forwardPrevious) + : boundary(row - 1, previousColumn, column); + return { + forwardPrevious, + forwardNext, + backwardPrevious, + backwardNext: mirrorBoundaryPoint( + midpoint(column, column), + backwardPrevious + ), + }; +} + +function getPeriodicCorners( + row: number, + column: number, + previousColumn: number, + nextColumn: number, + boundary: TBoundaryFunction, + midpoint: TMidpointFunction +) { + const forwardPrevious = boundary(row, previousColumn, column); + const forwardNext = boundary(row, column, nextColumn); + return { + forwardPrevious, + forwardNext, + backwardPrevious: + row === 0 + ? mirrorBoundaryPoint(midpoint(previousColumn, column), forwardPrevious) + : boundary(row - 1, previousColumn, column), + backwardNext: + row === 0 + ? mirrorBoundaryPoint(midpoint(column, nextColumn), forwardNext) + : boundary(row - 1, column, nextColumn), + }; +} + +/* eslint-disable-next-line max-lines-per-function */ +function getCenteredCellCorners( + grid: TCurvilinearGrid, + row: number, + column: number +) { + const previousColumn = getPreviousColumn( + column, + grid.ni, + grid.shouldFlipLongitude + ); + const nextColumn = getNextColumn(column, grid.ni, grid.shouldFlipLongitude); + const boundary: TBoundaryFunction = (boundaryRow, fromColumn, toColumn) => + getBoundaryPoint(grid, boundaryRow, fromColumn, toColumn); + const midpoint: TMidpointFunction = (fromColumn, toColumn) => + getRowMidpoint(grid, row, fromColumn, toColumn); + const isFirst = grid.shouldFlipLongitude + ? column === grid.ni - 1 + : column === 0; + const isLast = grid.shouldFlipLongitude + ? column === 0 + : column === grid.ni - 1; + const corners = + !grid.isPeriodicI && isFirst + ? getPreviousBoundaryCorners(row, column, nextColumn, boundary, midpoint) + : !grid.isPeriodicI && isLast + ? getNextBoundaryCorners( + row, + column, + previousColumn, + boundary, + midpoint + ) + : getPeriodicCorners( + row, + column, + previousColumn, + nextColumn, + boundary, + midpoint + ); + return { + latitudes: [ + corners.backwardPrevious.lat, + corners.backwardNext.lat, + corners.forwardNext.lat, + corners.forwardPrevious.lat, + ], + longitudes: [ + corners.backwardPrevious.lon, + corners.backwardNext.lon, + corners.forwardNext.lon, + corners.forwardPrevious.lon, + ], + }; +} + +function fillCell( + grid: TCurvilinearGrid, + projection: ProjectionHelper, + row: number, + column: number, + cellIndex: number, + batch: TGridGeometryBatch +) { + const corners = getCenteredCellCorners(grid, row, column); + const invalid = + corners.latitudes.some((value) => !Number.isFinite(value)) || + corners.longitudes.some((value) => !Number.isFinite(value)); + if (invalid) { + batch.dataValues.fill(NaN, cellIndex * 4, cellIndex * 4 + 4); + } else { + for (let corner = 0; corner < 4; corner++) { + projection.projectLatLonToArrays( + corners.latitudes[corner], + corners.longitudes[corner], + batch.positionValues, + cellIndex * 12 + corner * 3, + batch.latLonValues, + cellIndex * 8 + corner * 2 + ); + } + const dataIndex = row * grid.ni + column; + batch.dataValues.fill( + grid.data[dataIndex], + cellIndex * 4, + cellIndex * 4 + 4 + ); + } + const vertex = cellIndex * 4; + batch.indices.set( + [vertex, vertex + 1, vertex + 2, vertex, vertex + 2, vertex + 3], + cellIndex * 6 + ); +} + +export function getCurvilinearBatchCount( + grid: TCurvilinearGrid, + batchSize: number +) { + return Math.ceil((grid.nj - 1) / batchSize); +} + +export function buildCurvilinearBatch( + grid: TCurvilinearGrid, + batchIndex: number, + batchSize: number, + projection: ProjectionHelper +): TGridGeometryBatch { + const rowStart = batchIndex * batchSize; + const rowEnd = Math.min(rowStart + batchSize, grid.nj - 1); + const cellCount = (rowEnd - rowStart) * grid.ni; + const batch: TGridGeometryBatch = { + batchIndex, + positionValues: new Float32Array(cellCount * 12), + dataValues: new Float32Array(cellCount * 4), + latLonValues: new Float32Array(cellCount * 8), + indices: new Uint32Array(cellCount * 6), + }; + let cellIndex = 0; + for (let row = rowStart; row < rowEnd; row++) { + for (let column = 0; column < grid.ni; column++) { + fillCell(grid, projection, row, column, cellIndex, batch); + cellIndex++; + } + } + return batch; +} + +export function buildCurvilinearHoverIndexData( + grid: TCurvilinearGrid +): TSerializedGeoSampleIndexData { + const sampleCount = (grid.nj - 1) * grid.ni; + const latitudes = new Float64Array(sampleCount); + const longitudes = new Float64Array(sampleCount); + const values = new Float32Array(sampleCount); + let sampleIndex = 0; + for (let row = 0; row < grid.nj - 1; row++) { + for (let column = 0; column < grid.ni; column++) { + const dataIndex = row * grid.ni + column; + latitudes[sampleIndex] = grid.latitudes[dataIndex]; + longitudes[sampleIndex] = grid.longitudes[dataIndex]; + values[sampleIndex] = grid.data[dataIndex]; + sampleIndex++; + } + } + return buildSerializedGeoSampleIndexData(latitudes, longitudes, values); +} diff --git a/src/lib/grids/curvilinearWorkerClient.ts b/src/lib/grids/curvilinearWorkerClient.ts new file mode 100644 index 00000000..c8d0711f --- /dev/null +++ b/src/lib/grids/curvilinearWorkerClient.ts @@ -0,0 +1,69 @@ +import type { + TCurvilinearWorkerMetadata, + TCurvilinearWorkerRequest, +} from "./curvilinearWorkerProtocol.ts"; +import { + copyGridWorkerArray, + createGridGeometryWorkerClient, +} from "./gridGeometryWorkerClient.ts"; +import { GridGeometryWorkerMessageType } from "./gridGeometryWorkerProtocol.ts"; +import type { TGridGeometryBatch } from "./gridWorkerTypes.ts"; + +import type { + TProjectionCenter, + TProjectionType, +} from "@/lib/projection/projectionUtils.ts"; + +type TCurvilinearBuildRequest = { + latitudes: Float32Array; + longitudes: Float32Array; + data: Float32Array; + nj: number; + ni: number; + batchSize: number; + missingValue: number; + fillValue: number; + projectionType: TProjectionType; + projectionCenter: TProjectionCenter; +}; + +const client = createGridGeometryWorkerClient< + TCurvilinearWorkerRequest, + TCurvilinearWorkerMetadata +>( + () => + new Worker(new URL("./curvilinear.worker.ts", import.meta.url), { + type: "module", + }) +); + +export function buildCurvilinearGrid( + request: TCurvilinearBuildRequest, + callbacks: { + onMetadata: (metadata: TCurvilinearWorkerMetadata) => void; + onBatch: (batch: TGridGeometryBatch) => void; + } +) { + return client.build((requestId) => { + const latitudes = copyGridWorkerArray(request.latitudes); + const longitudes = copyGridWorkerArray(request.longitudes); + const data = copyGridWorkerArray(request.data); + return { + message: { + ...request, + requestId, + type: GridGeometryWorkerMessageType.BUILD, + latitudes, + longitudes, + data, + projectionCenter: { + lat: request.projectionCenter.lat, + lon: request.projectionCenter.lon, + }, + }, + transfer: [latitudes.buffer, longitudes.buffer, data.buffer], + }; + }, callbacks); +} + +export const terminateCurvilinearWorker = client.terminate; diff --git a/src/lib/grids/curvilinearWorkerProtocol.ts b/src/lib/grids/curvilinearWorkerProtocol.ts new file mode 100644 index 00000000..b630f580 --- /dev/null +++ b/src/lib/grids/curvilinearWorkerProtocol.ts @@ -0,0 +1,29 @@ +import type { TGridGeometryWorkerResponse } from "./gridGeometryWorkerProtocol.ts"; + +import type { + TProjectionCenter, + TProjectionType, +} from "@/lib/projection/projectionUtils.ts"; + +export type TCurvilinearWorkerRequest = { + requestId: number; + type: "build"; + latitudes: Float32Array; + longitudes: Float32Array; + data: Float32Array; + nj: number; + ni: number; + batchSize: number; + missingValue: number; + fillValue: number; + projectionType: TProjectionType; + projectionCenter: TProjectionCenter; +}; + +export type TCurvilinearWorkerMetadata = { + totalBatches: number; + shouldFlipLongitude: boolean; +}; + +export type TCurvilinearWorkerResponse = + TGridGeometryWorkerResponse; diff --git a/src/lib/grids/gaussianReduced.worker.ts b/src/lib/grids/gaussianReduced.worker.ts new file mode 100644 index 00000000..d1c0b3d4 --- /dev/null +++ b/src/lib/grids/gaussianReduced.worker.ts @@ -0,0 +1,88 @@ +/// + +import { + buildGaussianReducedBatch, + buildGaussianReducedHoverIndexData, + buildGaussianReducedRows, + getGaussianReducedBatchCount, +} from "./gaussianReducedCalculations.ts"; +import { + type TGaussianReducedWorkerRequest, + type TGaussianReducedWorkerResponse, +} from "./gaussianReducedWorkerProtocol.ts"; + +import { GridGeometryWorkerMessageType } from "@/lib/grids/gridGeometryWorkerProtocol.ts"; +import { + postGridGeometryBatch, + postGridGeometryHoverIndex, + postGridGeometryResponse, +} from "@/lib/grids/gridGeometryWorkerUtils.ts"; +import { ProjectionHelper } from "@/lib/projection/projectionUtils.ts"; + +const workerScope = self as unknown as DedicatedWorkerGlobalScope; + +function postResponse( + response: TGaussianReducedWorkerResponse, + transfer: Transferable[] = [] +) { + postGridGeometryResponse(workerScope, response, transfer); +} + +function postBatches( + request: TGaussianReducedWorkerRequest, + grid: ReturnType, + totalBatches: number +) { + const projection = new ProjectionHelper( + request.projectionType, + request.projectionCenter + ); + for (let batchIndex = 0; batchIndex < totalBatches; batchIndex++) { + const batch = buildGaussianReducedBatch( + grid, + batchIndex, + request.batchSize, + request.epsilon, + projection + ); + postGridGeometryBatch(workerScope, request.requestId, batch); + } +} + +function buildGrid(request: TGaussianReducedWorkerRequest) { + const grid = buildGaussianReducedRows( + request.latitudes, + request.longitudes, + request.data + ); + const totalBatches = getGaussianReducedBatchCount(grid, request.batchSize); + postResponse({ + requestId: request.requestId, + type: GridGeometryWorkerMessageType.METADATA, + metadata: { totalBatches }, + }); + postGridGeometryHoverIndex( + workerScope, + request.requestId, + buildGaussianReducedHoverIndexData(grid) + ); + postBatches(request, grid, totalBatches); + postResponse({ + requestId: request.requestId, + type: GridGeometryWorkerMessageType.DONE, + }); +} + +workerScope.onmessage = ( + event: MessageEvent +) => { + try { + buildGrid(event.data); + } catch (error) { + postResponse({ + requestId: event.data.requestId, + type: GridGeometryWorkerMessageType.ERROR, + message: error instanceof Error ? error.message : String(error), + }); + } +}; diff --git a/src/lib/grids/gaussianReducedCalculations.ts b/src/lib/grids/gaussianReducedCalculations.ts new file mode 100644 index 00000000..3be2ca4f --- /dev/null +++ b/src/lib/grids/gaussianReducedCalculations.ts @@ -0,0 +1,195 @@ +import { buildSerializedGeoSampleIndexData } from "./gridWorkerCalculations.ts"; +import type { + TGridGeometryBatch, + TSerializedGeoSampleIndexData, +} from "./gridWorkerTypes.ts"; + +import { ProjectionHelper } from "@/lib/projection/projectionUtils.ts"; + +type TGaussianReducedCell = { + lon: number; + value: number; +}; + +export type TGaussianReducedRows = { + rows: Record; + uniqueLatitudes: number[]; +}; + +export function buildGaussianReducedRows( + latitudes: Float64Array, + longitudes: Float64Array, + data: Float32Array +): TGaussianReducedRows { + const rows: Record = {}; + for (let index = 0; index < latitudes.length; index++) { + const latitude = latitudes[index]; + if (!rows[latitude]) { + rows[latitude] = []; + } + rows[latitude].push({ lon: longitudes[index], value: data[index] }); + } + + const uniqueLatitudes = Object.keys(rows) + .map(Number) + .sort((first, second) => second - first); + return { rows, uniqueLatitudes }; +} + +function getCellData(row: TGaussianReducedCell[], index: number) { + const cell = row[index]; + const nextCell = row[(index + 1) % row.length]; + const longitudeWidth = (nextCell.lon - cell.lon + 360) % 360; + return { cell, longitudeWidth }; +} + +function countCells( + grid: TGaussianReducedRows, + latitudeStart: number, + latitudeEnd: number +) { + let totalCells = 0; + for (let index = latitudeStart; index < latitudeEnd; index++) { + totalCells += grid.rows[grid.uniqueLatitudes[index]].length; + } + return totalCells; +} + +export function getGaussianReducedBatchCount( + grid: TGaussianReducedRows, + batchSize: number +) { + return Math.ceil((grid.uniqueLatitudes.length - 1) / batchSize); +} + +function projectQuad( + projection: ProjectionHelper, + latitudeTop: number, + longitudeLeft: number, + latitudeBottom: number, + longitudeWidth: number, + epsilon: number, + positionValues: Float32Array, + latLonValues: Float32Array, + positionOffset: number, + latLonOffset: number +) { + const longitudeRight = longitudeLeft - longitudeWidth; + projection.projectLatLonToArrays( + latitudeTop, + longitudeLeft + epsilon, + positionValues, + positionOffset, + latLonValues, + latLonOffset + ); + projection.projectLatLonToArrays( + latitudeTop, + longitudeRight - epsilon, + positionValues, + positionOffset + 3, + latLonValues, + latLonOffset + 2 + ); + projection.projectLatLonToArrays( + latitudeBottom - epsilon, + longitudeRight - epsilon, + positionValues, + positionOffset + 6, + latLonValues, + latLonOffset + 4 + ); + projection.projectLatLonToArrays( + latitudeBottom - epsilon, + longitudeLeft + epsilon, + positionValues, + positionOffset + 9, + latLonValues, + latLonOffset + 6 + ); +} + +/* eslint-disable-next-line max-lines-per-function */ +export function buildGaussianReducedBatch( + grid: TGaussianReducedRows, + batchIndex: number, + batchSize: number, + epsilon: number, + projection: ProjectionHelper +): TGridGeometryBatch { + const latitudeStart = batchIndex * batchSize; + const latitudeEnd = Math.min( + latitudeStart + batchSize, + grid.uniqueLatitudes.length - 1 + ); + const totalCells = countCells(grid, latitudeStart, latitudeEnd); + const latLonValues = new Float32Array(totalCells * 8); + const positionValues = new Float32Array(totalCells * 12); + const dataValues = new Float32Array(totalCells * 4); + const indices = new Uint32Array(totalCells * 6); + let cellIndex = 0; + + for ( + let latitudeIndex = latitudeStart; + latitudeIndex < latitudeEnd; + latitudeIndex++ + ) { + const latitudeTop = grid.uniqueLatitudes[latitudeIndex]; + const latitudeBottom = grid.uniqueLatitudes[latitudeIndex + 1]; + const row = grid.rows[latitudeTop]; + + for (let index = 0; index < row.length; index++) { + const { cell, longitudeWidth } = getCellData(row, index); + projectQuad( + projection, + latitudeTop, + cell.lon, + latitudeBottom, + longitudeWidth, + epsilon, + positionValues, + latLonValues, + cellIndex * 12, + cellIndex * 8 + ); + dataValues.fill(cell.value, cellIndex * 4, cellIndex * 4 + 4); + const vertex = cellIndex * 4; + indices.set( + [vertex, vertex + 1, vertex + 2, vertex, vertex + 2, vertex + 3], + cellIndex * 6 + ); + cellIndex++; + } + } + + return { batchIndex, positionValues, dataValues, latLonValues, indices }; +} + +export function buildGaussianReducedHoverIndexData( + grid: TGaussianReducedRows +): TSerializedGeoSampleIndexData { + const rowCount = Math.max(0, grid.uniqueLatitudes.length - 1); + const sampleCount = countCells(grid, 0, rowCount); + const latitudes = new Float64Array(sampleCount); + const longitudes = new Float64Array(sampleCount); + const values = new Float32Array(sampleCount); + let sampleIndex = 0; + + for (let latitudeIndex = 0; latitudeIndex < rowCount; latitudeIndex++) { + const latitudeTop = grid.uniqueLatitudes[latitudeIndex]; + const latitudeBottom = grid.uniqueLatitudes[latitudeIndex + 1]; + const row = grid.rows[latitudeTop]; + for (let indexInRow = 0; indexInRow < row.length; indexInRow++) { + const { cell, longitudeWidth } = getCellData(row, indexInRow); + const latitude = (latitudeTop + latitudeBottom) / 2; + const longitude = ProjectionHelper.normalizeLongitude( + cell.lon - longitudeWidth / 2 + ); + latitudes[sampleIndex] = latitude; + longitudes[sampleIndex] = longitude; + values[sampleIndex] = cell.value; + sampleIndex++; + } + } + return buildSerializedGeoSampleIndexData(latitudes, longitudes, values); +} diff --git a/src/lib/grids/gaussianReducedWorkerClient.ts b/src/lib/grids/gaussianReducedWorkerClient.ts new file mode 100644 index 00000000..87c5246d --- /dev/null +++ b/src/lib/grids/gaussianReducedWorkerClient.ts @@ -0,0 +1,74 @@ +import type { + TGaussianReducedWorkerMetadata, + TGaussianReducedWorkerRequest, +} from "./gaussianReducedWorkerProtocol.ts"; +import { + copyGridWorkerArray, + createGridGeometryWorkerClient, +} from "./gridGeometryWorkerClient.ts"; +import { GridGeometryWorkerMessageType } from "./gridGeometryWorkerProtocol.ts"; +import type { TGridGeometryBatch } from "./gridWorkerTypes.ts"; + +import type { + TProjectionCenter, + TProjectionType, +} from "@/lib/projection/projectionUtils.ts"; + +type TGaussianReducedBuildRequest = { + latitudes: Float64Array; + longitudes: Float64Array; + data: Float32Array; + batchSize: number; + epsilon: number; + projectionType: TProjectionType; + projectionCenter: TProjectionCenter; +}; + +const client = createGridGeometryWorkerClient< + TGaussianReducedWorkerRequest, + TGaussianReducedWorkerMetadata +>( + () => + new Worker(new URL("./gaussianReduced.worker.ts", import.meta.url), { + type: "module", + }) +); + +export function buildGaussianReducedGrid( + request: TGaussianReducedBuildRequest, + callbacks: { + onMetadata: (totalBatches: number) => void; + onBatch: (batch: TGridGeometryBatch) => void; + } +) { + return client + .build( + (requestId) => { + const latitudes = copyGridWorkerArray(request.latitudes); + const longitudes = copyGridWorkerArray(request.longitudes); + const data = copyGridWorkerArray(request.data); + return { + message: { + ...request, + requestId, + type: GridGeometryWorkerMessageType.BUILD, + latitudes, + longitudes, + data, + projectionCenter: { + lat: request.projectionCenter.lat, + lon: request.projectionCenter.lon, + }, + }, + transfer: [latitudes.buffer, longitudes.buffer, data.buffer], + }; + }, + { + onMetadata: ({ totalBatches }) => callbacks.onMetadata(totalBatches), + onBatch: callbacks.onBatch, + } + ) + .then(({ hoverIndexData }) => hoverIndexData); +} + +export const terminateGaussianReducedWorker = client.terminate; diff --git a/src/lib/grids/gaussianReducedWorkerProtocol.ts b/src/lib/grids/gaussianReducedWorkerProtocol.ts new file mode 100644 index 00000000..418a6887 --- /dev/null +++ b/src/lib/grids/gaussianReducedWorkerProtocol.ts @@ -0,0 +1,24 @@ +import type { TGridGeometryWorkerResponse } from "@/lib/grids/gridGeometryWorkerProtocol.ts"; +import type { + TProjectionCenter, + TProjectionType, +} from "@/lib/projection/projectionUtils.ts"; + +export type TGaussianReducedWorkerRequest = { + requestId: number; + type: "build"; + latitudes: Float64Array; + longitudes: Float64Array; + data: Float32Array; + batchSize: number; + epsilon: number; + projectionType: TProjectionType; + projectionCenter: TProjectionCenter; +}; + +export type TGaussianReducedWorkerMetadata = { + totalBatches: number; +}; + +export type TGaussianReducedWorkerResponse = + TGridGeometryWorkerResponse; diff --git a/src/lib/grids/regularGridData.worker.ts b/src/lib/grids/gridData.worker.ts similarity index 64% rename from src/lib/grids/regularGridData.worker.ts rename to src/lib/grids/gridData.worker.ts index 9d9b206d..af94fb7d 100644 --- a/src/lib/grids/regularGridData.worker.ts +++ b/src/lib/grids/gridData.worker.ts @@ -2,16 +2,14 @@ import { ZarrDataManager } from "@/lib/data/ZarrDataManager.ts"; import { - RegularGridDataWorkerMessageType, - type TRegularGridDataWorkerRequest, - type TRegularGridDataWorkerResponse, -} from "@/lib/grids/regularGridDataWorkerProtocol.ts"; + GridDataWorkerMessageType, + type TGridDataWorkerRequest, + type TGridDataWorkerResponse, +} from "@/lib/grids/gridDataWorkerProtocol.ts"; const workerScope = self as unknown as DedicatedWorkerGlobalScope; -workerScope.onmessage = async ( - event: MessageEvent -) => { +workerScope.onmessage = async (event: MessageEvent) => { const { requestId, source, variable, format, selection } = event.data; try { const array = await ZarrDataManager.getVariableInfo( @@ -23,9 +21,9 @@ workerScope.onmessage = async ( array, selection ); - const response: TRegularGridDataWorkerResponse = { + const response: TGridDataWorkerResponse = { requestId, - type: RegularGridDataWorkerMessageType.RESULT, + type: GridDataWorkerMessageType.RESULT, data: chunk.data, }; const transfer = @@ -34,9 +32,9 @@ workerScope.onmessage = async ( : []; workerScope.postMessage(response, transfer); } catch (error) { - const response: TRegularGridDataWorkerResponse = { + const response: TGridDataWorkerResponse = { requestId, - type: RegularGridDataWorkerMessageType.ERROR, + type: GridDataWorkerMessageType.ERROR, message: error instanceof Error ? error.message : String(error), }; workerScope.postMessage(response); diff --git a/src/lib/grids/regularGridDataWorkerClient.ts b/src/lib/grids/gridDataWorkerClient.ts similarity index 69% rename from src/lib/grids/regularGridDataWorkerClient.ts rename to src/lib/grids/gridDataWorkerClient.ts index 68bfa7c6..79c19f31 100644 --- a/src/lib/grids/regularGridDataWorkerClient.ts +++ b/src/lib/grids/gridDataWorkerClient.ts @@ -1,14 +1,14 @@ import type * as zarr from "zarrita"; import { - RegularGridDataWorkerMessageType, - type TRegularGridDataWorkerRequest, - type TRegularGridDataWorkerResponse, -} from "./regularGridDataWorkerProtocol.ts"; + GridDataWorkerMessageType, + type TGridDataWorkerRequest, + type TGridDataWorkerResponse, +} from "./gridDataWorkerProtocol.ts"; import type { TDataSource, TZarrFormat } from "@/lib/types/GlobeTypes.ts"; -type TRegularGridDataRequest = { +export type TGridDataRequest = { source: Pick; variable: string; format: TZarrFormat; @@ -31,13 +31,13 @@ function rejectPendingRequests(error: Error) { pendingRequests.clear(); } -function handleWorkerMessage(message: TRegularGridDataWorkerResponse) { +function handleWorkerMessage(message: TGridDataWorkerResponse) { const pending = pendingRequests.get(message.requestId); if (!pending) { return; } pendingRequests.delete(message.requestId); - if (message.type === RegularGridDataWorkerMessageType.ERROR) { + if (message.type === GridDataWorkerMessageType.ERROR) { pending.reject(new Error(message.message)); return; } @@ -48,10 +48,10 @@ function getWorker() { if (worker) { return worker; } - worker = new Worker(new URL("./regularGridData.worker.ts", import.meta.url), { + worker = new Worker(new URL("./gridData.worker.ts", import.meta.url), { type: "module", }); - worker.onmessage = (event: MessageEvent) => { + worker.onmessage = (event: MessageEvent) => { handleWorkerMessage(event.data); }; worker.onerror = (event) => { @@ -62,11 +62,11 @@ function getWorker() { return worker; } -export function getRegularGridVariableData(request: TRegularGridDataRequest) { +export function getGridVariableData(request: TGridDataRequest) { const requestId = ++nextRequestId; - const message: TRegularGridDataWorkerRequest = { + const message: TGridDataWorkerRequest = { requestId, - type: RegularGridDataWorkerMessageType.GET_DATA, + type: GridDataWorkerMessageType.GET_DATA, source: { store: request.source.store, dataset: request.source.dataset, @@ -86,8 +86,8 @@ export function getRegularGridVariableData(request: TRegularGridDataRequest) { }); } -export function terminateRegularGridDataWorker() { +export function terminateGridDataWorker() { worker?.terminate(); worker = null; - rejectPendingRequests(new Error("Regular grid data worker terminated.")); + rejectPendingRequests(new Error("Grid data worker terminated.")); } diff --git a/src/lib/grids/gridDataWorkerProtocol.ts b/src/lib/grids/gridDataWorkerProtocol.ts new file mode 100644 index 00000000..b93daa3b --- /dev/null +++ b/src/lib/grids/gridDataWorkerProtocol.ts @@ -0,0 +1,36 @@ +import type * as zarr from "zarrita"; + +import type { TDataSource, TZarrFormat } from "@/lib/types/GlobeTypes.ts"; + +export const GridDataWorkerMessageType = { + GET_DATA: "getData", + RESULT: "result", + ERROR: "error", +} as const; + +type TGridDataWorkerMessageType = + (typeof GridDataWorkerMessageType)[keyof typeof GridDataWorkerMessageType]; + +export type TGridDataWorkerRequest = { + requestId: number; + type: typeof GridDataWorkerMessageType.GET_DATA; + source: Pick; + variable: string; + format: TZarrFormat; + selection: (number | null | zarr.Slice)[]; +}; + +type TGridDataWorkerResponseBase = { + requestId: number; + type: TGridDataWorkerMessageType; +}; + +export type TGridDataWorkerResponse = + | (TGridDataWorkerResponseBase & { + type: typeof GridDataWorkerMessageType.RESULT; + data: zarr.TypedArray; + }) + | (TGridDataWorkerResponseBase & { + type: typeof GridDataWorkerMessageType.ERROR; + message: string; + }); diff --git a/src/lib/grids/gridGeometryWorkerClient.ts b/src/lib/grids/gridGeometryWorkerClient.ts new file mode 100644 index 00000000..19150530 --- /dev/null +++ b/src/lib/grids/gridGeometryWorkerClient.ts @@ -0,0 +1,151 @@ +import { + GridGeometryWorkerMessageType, + type TGridGeometryWorkerResponse, +} from "./gridGeometryWorkerProtocol.ts"; +import type { + TGridGeometryBatch, + TSerializedGeoSampleIndexData, + TGridWorkerBatch, +} from "./gridWorkerTypes.ts"; + +type TGridGeometryBuildCallbacks = { + onMetadata: (metadata: TMetadata) => void; + onBatch: (batch: TBatch) => void; +}; + +type TWorkerRequest = { + message: TRequest; + transfer: Transferable[]; +}; + +type TBuildResult = { + metadata: TMetadata; + hoverIndexData: TSerializedGeoSampleIndexData; +}; + +type TActiveBuild = { + requestId: number; + worker: Worker; + callbacks: TGridGeometryBuildCallbacks; + metadata: TMetadata | null; + hoverIndexData: TSerializedGeoSampleIndexData | null; + resolve: (value: TBuildResult) => void; + reject: (reason?: unknown) => void; +}; + +export function copyGridWorkerArray( + array: TArray +): TArray { + return array.slice() as TArray; +} + +/* eslint-disable-next-line max-lines-per-function */ +export function createGridGeometryWorkerClient< + TRequest, + TMetadata, + TBatch extends TGridWorkerBatch = TGridGeometryBatch, +>(createWorker: () => Worker, keepWorkerAlive = false) { + let activeBuild: TActiveBuild | null = null; + let nextRequestId = 0; + let worker: Worker | null = null; + + function terminate() { + if (activeBuild) { + activeBuild.reject(new Error("Grid geometry worker terminated.")); + } + activeBuild = null; + worker?.terminate(); + worker = null; + } + + function fail(error: Error) { + if (!activeBuild) { + return; + } + const { worker: activeWorker, reject } = activeBuild; + activeBuild = null; + activeWorker.terminate(); + worker = null; + reject(error); + } + + function finish() { + if (!activeBuild?.metadata || !activeBuild.hoverIndexData) { + fail(new Error("Grid geometry worker returned incomplete results.")); + return; + } + const { + worker: activeWorker, + metadata, + hoverIndexData, + resolve, + } = activeBuild; + activeBuild = null; + if (!keepWorkerAlive) { + activeWorker.terminate(); + worker = null; + } + resolve({ metadata, hoverIndexData }); + } + + function handleMessage( + message: TGridGeometryWorkerResponse + ) { + if (activeBuild?.requestId !== message.requestId) { + return; + } + if (message.type === GridGeometryWorkerMessageType.METADATA) { + activeBuild.metadata = message.metadata; + activeBuild.callbacks.onMetadata(message.metadata); + } else if (message.type === GridGeometryWorkerMessageType.BATCH) { + activeBuild.callbacks.onBatch(message.batch); + } else if (message.type === GridGeometryWorkerMessageType.HOVER_INDEX) { + activeBuild.hoverIndexData = message.hoverIndexData; + } else if (message.type === GridGeometryWorkerMessageType.ERROR) { + fail(new Error(message.message)); + } else { + finish(); + } + } + + function build( + createRequest: (requestId: number) => TWorkerRequest, + callbacks: TGridGeometryBuildCallbacks + ) { + if (activeBuild) { + terminate(); + } + const requestId = ++nextRequestId; + worker ??= createWorker(); + const activeWorker = worker; + return new Promise>((resolve, reject) => { + activeBuild = { + requestId, + worker: activeWorker, + callbacks, + metadata: null, + hoverIndexData: null, + resolve, + reject, + }; + activeWorker.onmessage = ( + event: MessageEvent> + ) => { + try { + handleMessage(event.data); + } catch (error) { + fail(error instanceof Error ? error : new Error(String(error))); + } + }; + activeWorker.onerror = (event) => { + if (activeBuild?.requestId === requestId) { + fail(new Error(event.message)); + } + }; + const request = createRequest(requestId); + activeWorker.postMessage(request.message, request.transfer); + }); + } + + return { build, terminate }; +} diff --git a/src/lib/grids/gridGeometryWorkerProtocol.ts b/src/lib/grids/gridGeometryWorkerProtocol.ts new file mode 100644 index 00000000..2185d921 --- /dev/null +++ b/src/lib/grids/gridGeometryWorkerProtocol.ts @@ -0,0 +1,46 @@ +import type { + TGridGeometryBatch, + TSerializedGeoSampleIndexData, + TGridWorkerBatch, +} from "./gridWorkerTypes.ts"; + +export const GridGeometryWorkerMessageType = { + BUILD: "build", + METADATA: "metadata", + BATCH: "batch", + HOVER_INDEX: "hoverIndex", + DONE: "done", + ERROR: "error", +} as const; + +type TGridGeometryWorkerMessageType = + (typeof GridGeometryWorkerMessageType)[keyof typeof GridGeometryWorkerMessageType]; + +type TGridGeometryWorkerResponseBase = { + requestId: number; + type: TGridGeometryWorkerMessageType; +}; + +export type TGridGeometryWorkerResponse< + TMetadata, + TBatch extends TGridWorkerBatch = TGridGeometryBatch, +> = + | (TGridGeometryWorkerResponseBase & { + type: typeof GridGeometryWorkerMessageType.METADATA; + metadata: TMetadata; + }) + | (TGridGeometryWorkerResponseBase & { + type: typeof GridGeometryWorkerMessageType.BATCH; + batch: TBatch; + }) + | (TGridGeometryWorkerResponseBase & { + type: typeof GridGeometryWorkerMessageType.HOVER_INDEX; + hoverIndexData: TSerializedGeoSampleIndexData; + }) + | (TGridGeometryWorkerResponseBase & { + type: typeof GridGeometryWorkerMessageType.DONE; + }) + | (TGridGeometryWorkerResponseBase & { + type: typeof GridGeometryWorkerMessageType.ERROR; + message: string; + }); diff --git a/src/lib/grids/gridGeometryWorkerUtils.ts b/src/lib/grids/gridGeometryWorkerUtils.ts new file mode 100644 index 00000000..40c7b84f --- /dev/null +++ b/src/lib/grids/gridGeometryWorkerUtils.ts @@ -0,0 +1,98 @@ +import type { TGridGeometryWorkerResponse } from "./gridGeometryWorkerProtocol.ts"; +import { GridGeometryWorkerMessageType } from "./gridGeometryWorkerProtocol.ts"; +import type { + TGridDataValueBatch, + TGridGeometryBatch, + TGridPointBatch, + TSerializedGeoSampleIndexData, + TGridWorkerBatch, +} from "./gridWorkerTypes.ts"; + +export function postGridGeometryResponse< + TMetadata, + TBatch extends TGridWorkerBatch = TGridGeometryBatch, +>( + workerScope: DedicatedWorkerGlobalScope, + response: TGridGeometryWorkerResponse, + transfer: Transferable[] = [] +) { + workerScope.postMessage(response, transfer); +} + +export function postGridGeometryHoverIndex( + workerScope: DedicatedWorkerGlobalScope, + requestId: number, + hoverIndexData: TSerializedGeoSampleIndexData +) { + postGridGeometryResponse( + workerScope, + { + requestId, + type: GridGeometryWorkerMessageType.HOVER_INDEX, + hoverIndexData, + }, + [ + hoverIndexData.indexData, + hoverIndexData.latitudes.buffer, + hoverIndexData.longitudes.buffer, + hoverIndexData.values.buffer, + ] + ); +} + +export function postGridGeometryBatch( + workerScope: DedicatedWorkerGlobalScope, + requestId: number, + batch: TGridGeometryBatch +) { + postGridGeometryResponse( + workerScope, + { + requestId, + type: GridGeometryWorkerMessageType.BATCH, + batch, + }, + [ + batch.positionValues.buffer, + batch.dataValues.buffer, + batch.latLonValues.buffer, + batch.indices.buffer, + ] + ); +} + +export function postGridPointBatch( + workerScope: DedicatedWorkerGlobalScope, + requestId: number, + batch: TGridPointBatch +) { + postGridGeometryResponse( + workerScope, + { + requestId, + type: GridGeometryWorkerMessageType.BATCH, + batch, + }, + [ + batch.positionValues.buffer, + batch.dataValues.buffer, + batch.latLonValues.buffer, + ] + ); +} + +export function postGridDataValueBatch( + workerScope: DedicatedWorkerGlobalScope, + requestId: number, + batch: TGridDataValueBatch +) { + postGridGeometryResponse( + workerScope, + { + requestId, + type: GridGeometryWorkerMessageType.BATCH, + batch, + }, + [batch.dataValues.buffer] + ); +} diff --git a/src/lib/grids/gridWorkerCalculations.ts b/src/lib/grids/gridWorkerCalculations.ts new file mode 100644 index 00000000..c775c064 --- /dev/null +++ b/src/lib/grids/gridWorkerCalculations.ts @@ -0,0 +1,27 @@ +import KDBush from "kdbush"; + +import type { TSerializedGeoSampleIndexData } from "./gridWorkerTypes.ts"; + +export function buildSerializedGeoSampleIndexData( + latitudes: Float64Array, + longitudes: Float64Array, + values: Float32Array +): TSerializedGeoSampleIndexData { + const index = new KDBush(latitudes.length); + for (let sampleIndex = 0; sampleIndex < latitudes.length; sampleIndex++) { + latitudes[sampleIndex] = Math.max( + -90, + Math.min(90, latitudes[sampleIndex]) + ); + longitudes[sampleIndex] = + (((longitudes[sampleIndex] % 360) + 540) % 360) - 180; + index.add(longitudes[sampleIndex], latitudes[sampleIndex]); + } + index.finish(); + return { + indexData: index.data as ArrayBuffer, + latitudes, + longitudes, + values, + }; +} diff --git a/src/lib/grids/gridWorkerTypes.ts b/src/lib/grids/gridWorkerTypes.ts new file mode 100644 index 00000000..cc01b6b9 --- /dev/null +++ b/src/lib/grids/gridWorkerTypes.ts @@ -0,0 +1,26 @@ +export type TGridGeometryBatch = { + batchIndex: number; + positionValues: Float32Array; + dataValues: Float32Array; + latLonValues: Float32Array; + indices: Uint32Array; +}; + +export type TGridPointBatch = Omit; + +export type TGridDataValueBatch = { + batchIndex: number; + dataValues: Float32Array; +}; + +export type TGridWorkerBatch = + | TGridGeometryBatch + | TGridPointBatch + | TGridDataValueBatch; + +export type TSerializedGeoSampleIndexData = { + indexData: ArrayBuffer; + latitudes: Float64Array; + longitudes: Float64Array; + values: Float32Array; +}; diff --git a/src/lib/grids/irregular.worker.ts b/src/lib/grids/irregular.worker.ts new file mode 100644 index 00000000..fad93311 --- /dev/null +++ b/src/lib/grids/irregular.worker.ts @@ -0,0 +1,78 @@ +/// + +import { + buildIrregularBatch, + buildIrregularGridData, + buildIrregularHoverIndexData, + getIrregularBatchCount, +} from "./irregularCalculations.ts"; +import type { + TIrregularWorkerRequest, + TIrregularWorkerResponse, +} from "./irregularWorkerProtocol.ts"; + +import { GridGeometryWorkerMessageType } from "@/lib/grids/gridGeometryWorkerProtocol.ts"; +import { + postGridGeometryHoverIndex, + postGridGeometryResponse, + postGridPointBatch, +} from "@/lib/grids/gridGeometryWorkerUtils.ts"; +import { ProjectionHelper } from "@/lib/projection/projectionUtils.ts"; + +const workerScope = self as unknown as DedicatedWorkerGlobalScope; + +function postResponse( + response: TIrregularWorkerResponse, + transfer: Transferable[] = [] +) { + postGridGeometryResponse(workerScope, response, transfer); +} + +function buildGrid(request: TIrregularWorkerRequest) { + const projection = new ProjectionHelper( + request.projectionType, + request.projectionCenter + ); + const grid = buildIrregularGridData( + request.latitudes, + request.longitudes, + request.latitudeShape, + request.longitudeShape, + request.data, + projection + ); + const totalBatches = getIrregularBatchCount(grid, request.batchSize); + postResponse({ + requestId: request.requestId, + type: GridGeometryWorkerMessageType.METADATA, + metadata: { totalBatches, estimatedSpacing: grid.estimatedSpacing }, + }); + postGridGeometryHoverIndex( + workerScope, + request.requestId, + buildIrregularHoverIndexData(grid) + ); + for (let batchIndex = 0; batchIndex < totalBatches; batchIndex++) { + postGridPointBatch( + workerScope, + request.requestId, + buildIrregularBatch(grid, batchIndex, request.batchSize) + ); + } + postResponse({ + requestId: request.requestId, + type: GridGeometryWorkerMessageType.DONE, + }); +} + +workerScope.onmessage = (event: MessageEvent) => { + try { + buildGrid(event.data); + } catch (error) { + postResponse({ + requestId: event.data.requestId, + type: GridGeometryWorkerMessageType.ERROR, + message: error instanceof Error ? error.message : String(error), + }); + } +}; diff --git a/src/lib/grids/irregularCalculations.ts b/src/lib/grids/irregularCalculations.ts new file mode 100644 index 00000000..9f2dbf22 --- /dev/null +++ b/src/lib/grids/irregularCalculations.ts @@ -0,0 +1,129 @@ +import { buildSerializedGeoSampleIndexData } from "./gridWorkerCalculations.ts"; +import type { + TGridPointBatch, + TSerializedGeoSampleIndexData, +} from "./gridWorkerTypes.ts"; + +import { reconcileCoordinateArrays } from "@/lib/data/irregularGridHelpers.ts"; +import { ProjectionHelper } from "@/lib/projection/projectionUtils.ts"; + +export type TIrregularGridData = { + positions: Float32Array; + latLonValues: Float32Array; + data: Float32Array; + latitudes: Float32Array; + longitudes: Float32Array; + estimatedSpacing: number; +}; + +export function estimateIrregularAverageSpacing( + positions: Float32Array, + sampleSize = 5000 +) { + const pointCount = positions.length / 3; + const stride = Math.max(1, Math.floor(pointCount / sampleSize)); + let totalDistance = 0; + let totalWeight = 0; + + for (let pointIndex = 0; pointIndex < pointCount; pointIndex += stride) { + const offset = pointIndex * 3; + const x = positions[offset]; + const y = positions[offset + 1]; + const z = positions[offset + 2]; + const weight = Math.cos(Math.asin(z)); + if (weight < 0.1) { + continue; + } + let nearestDistance = Number.POSITIVE_INFINITY; + const start = Math.max(0, pointIndex - 50); + const end = Math.min(pointCount, pointIndex + 50); + for (let neighborIndex = start; neighborIndex < end; neighborIndex++) { + if (pointIndex === neighborIndex) { + continue; + } + const neighborOffset = neighborIndex * 3; + const deltaX = positions[neighborOffset] - x; + const deltaY = positions[neighborOffset + 1] - y; + const deltaZ = positions[neighborOffset + 2] - z; + const distance = Math.sqrt( + deltaX * deltaX + deltaY * deltaY + deltaZ * deltaZ + ); + if (distance > 0 && distance < nearestDistance) { + nearestDistance = distance; + } + } + if (nearestDistance !== Number.POSITIVE_INFINITY) { + totalDistance += nearestDistance * weight; + totalWeight += weight; + } + } + return totalWeight > 0 ? totalDistance / totalWeight : 0.01; +} + +export function buildIrregularGridData( + latitudes: Float32Array, + longitudes: Float32Array, + latitudeShape: number[], + longitudeShape: number[], + data: Float32Array, + projection: ProjectionHelper +): TIrregularGridData { + const coordinates = reconcileCoordinateArrays( + latitudes, + longitudes, + latitudeShape, + longitudeShape, + data.length + ); + const positions = new Float32Array(data.length * 3); + const latLonValues = new Float32Array(data.length * 2); + for (let index = 0; index < data.length; index++) { + projection.projectLatLonToArrays( + coordinates.latitudes[index], + coordinates.longitudes[index], + positions, + index * 3, + latLonValues, + index * 2 + ); + } + return { + positions, + latLonValues, + data, + ...coordinates, + estimatedSpacing: estimateIrregularAverageSpacing(positions), + }; +} + +export function getIrregularBatchCount( + grid: TIrregularGridData, + batchSize: number +) { + return Math.ceil(grid.data.length / batchSize); +} + +export function buildIrregularBatch( + grid: TIrregularGridData, + batchIndex: number, + batchSize: number +): TGridPointBatch { + const start = batchIndex * batchSize; + const end = Math.min(start + batchSize, grid.data.length); + return { + batchIndex, + positionValues: grid.positions.slice(start * 3, end * 3), + latLonValues: grid.latLonValues.slice(start * 2, end * 2), + dataValues: grid.data.slice(start, end), + }; +} + +export function buildIrregularHoverIndexData( + grid: TIrregularGridData +): TSerializedGeoSampleIndexData { + return buildSerializedGeoSampleIndexData( + Float64Array.from(grid.latitudes), + Float64Array.from(grid.longitudes), + grid.data.slice() + ); +} diff --git a/src/lib/grids/irregularDelaunay.worker.ts b/src/lib/grids/irregularDelaunay.worker.ts new file mode 100644 index 00000000..f4d11af9 --- /dev/null +++ b/src/lib/grids/irregularDelaunay.worker.ts @@ -0,0 +1,114 @@ +/// + +import { + buildIrregularDelaunayDataBatch, + buildIrregularDelaunayGeometryBatch, + buildIrregularDelaunayGrid, + buildIrregularDelaunayHoverIndexData, + getIrregularDelaunayBatchCount, + type TIrregularDelaunayGrid, +} from "./irregularDelaunayCalculations.ts"; +import type { + TIrregularDelaunayWorkerRequest, + TIrregularDelaunayWorkerResponse, +} from "./irregularDelaunayWorkerProtocol.ts"; + +import { GridGeometryWorkerMessageType } from "@/lib/grids/gridGeometryWorkerProtocol.ts"; +import { + postGridDataValueBatch, + postGridGeometryHoverIndex, + postGridGeometryResponse, + postGridPointBatch, +} from "@/lib/grids/gridGeometryWorkerUtils.ts"; + +const workerScope = self as unknown as DedicatedWorkerGlobalScope; +let cachedGrid: TIrregularDelaunayGrid | null = null; +let cachedPointCount = 0; + +function postResponse( + response: TIrregularDelaunayWorkerResponse, + transfer: Transferable[] = [] +) { + postGridGeometryResponse(workerScope, response, transfer); +} + +function postBatches( + request: TIrregularDelaunayWorkerRequest, + grid: TIrregularDelaunayGrid, + totalBatches: number, + rebuildGeometry: boolean +) { + for (let batchIndex = 0; batchIndex < totalBatches; batchIndex++) { + if (rebuildGeometry) { + postGridPointBatch( + workerScope, + request.requestId, + buildIrregularDelaunayGeometryBatch( + grid, + request.data, + batchIndex, + request.batchSize + ) + ); + } else { + postGridDataValueBatch( + workerScope, + request.requestId, + buildIrregularDelaunayDataBatch( + grid, + request.data, + batchIndex, + request.batchSize + ) + ); + } + } +} + +function buildGrid(request: TIrregularDelaunayWorkerRequest) { + const rebuildGeometry = + request.forceGeometryRebuild || + !cachedGrid || + cachedPointCount !== request.data.length; + if (rebuildGeometry) { + cachedGrid = buildIrregularDelaunayGrid( + request.latitudes, + request.longitudes, + request.latitudeShape, + request.longitudeShape, + request.data.length + ); + cachedPointCount = request.data.length; + } + const grid = cachedGrid!; + const totalBatches = getIrregularDelaunayBatchCount(grid, request.batchSize); + postResponse({ + requestId: request.requestId, + type: GridGeometryWorkerMessageType.METADATA, + metadata: { totalBatches, rebuildGeometry }, + }); + postGridGeometryHoverIndex( + workerScope, + request.requestId, + buildIrregularDelaunayHoverIndexData(grid, request.data) + ); + postBatches(request, grid, totalBatches, rebuildGeometry); + postResponse({ + requestId: request.requestId, + type: GridGeometryWorkerMessageType.DONE, + }); +} + +workerScope.onmessage = ( + event: MessageEvent +) => { + try { + buildGrid(event.data); + } catch (error) { + postResponse({ + requestId: event.data.requestId, + type: GridGeometryWorkerMessageType.ERROR, + message: error instanceof Error ? error.message : String(error), + }); + } +}; diff --git a/src/lib/grids/irregularDelaunayCalculations.ts b/src/lib/grids/irregularDelaunayCalculations.ts new file mode 100644 index 00000000..0fbe2d8e --- /dev/null +++ b/src/lib/grids/irregularDelaunayCalculations.ts @@ -0,0 +1,312 @@ +import { Delaunay } from "d3-delaunay"; + +import { buildSerializedGeoSampleIndexData } from "./gridWorkerCalculations.ts"; +import type { + TGridDataValueBatch, + TGridPointBatch, + TSerializedGeoSampleIndexData, +} from "./gridWorkerTypes.ts"; + +import { reconcileCoordinateArrays } from "@/lib/data/irregularGridHelpers.ts"; +import { + PROJECTION_TYPES, + ProjectionHelper, +} from "@/lib/projection/projectionUtils.ts"; + +type TExtendedPointSet = { + points: [number, number][]; + indexMap: number[]; +}; + +export type TIrregularDelaunayGrid = { + triangleIndices: Uint32Array; + positions: Float32Array; + latLonValues: Float32Array; + latitudes: Float32Array; + longitudes: Float32Array; +}; + +function createExtendedPointSet( + latitudes: Float32Array, + longitudes: Float32Array, + wrapThreshold: number +): TExtendedPointSet { + const points: [number, number][] = []; + const indexMap: number[] = []; + for (let index = 0; index < latitudes.length; index++) { + const longitude = longitudes[index]; + const latitude = latitudes[index]; + points.push([longitude, latitude]); + indexMap.push(index); + if (longitude < -wrapThreshold) { + points.push([longitude + 360, latitude]); + indexMap.push(index); + } + if (longitude > wrapThreshold) { + points.push([longitude - 360, latitude]); + indexMap.push(index); + } + } + return { points, indexMap }; +} + +function isValidTriangle( + points: [number, number][], + pointIndices: [number, number, number], + maxEdgeLength: number +) { + const [point0, point1, point2] = pointIndices.map( + (pointIndex) => points[pointIndex] + ); + if ( + (point0[0] < -180 || point0[0] > 180) && + (point1[0] < -180 || point1[0] > 180) && + (point2[0] < -180 || point2[0] > 180) + ) { + return false; + } + const edges = [ + [point0, point1], + [point1, point2], + [point2, point0], + ] as const; + for (const [[longitude0, latitude0], [longitude1, latitude1]] of edges) { + const longitudeDelta = longitude1 - longitude0; + const latitudeDelta = latitude1 - latitude0; + if ( + Math.sqrt( + longitudeDelta * longitudeDelta + latitudeDelta * latitudeDelta + ) > maxEdgeLength + ) { + return false; + } + } + return true; +} + +export function computeDelaunayTriangulation( + latitudes: Float32Array, + longitudes: Float32Array +): Uint32Array { + const { points, indexMap } = createExtendedPointSet( + latitudes, + longitudes, + 150 + ); + const triangles = Delaunay.from(points).triangles; + const validTriangles: number[] = []; + for ( + let triangleOffset = 0; + triangleOffset < triangles.length; + triangleOffset += 3 + ) { + const pointIndices: [number, number, number] = [ + triangles[triangleOffset], + triangles[triangleOffset + 1], + triangles[triangleOffset + 2], + ]; + if (isValidTriangle(points, pointIndices, 30)) { + validTriangles.push( + indexMap[pointIndices[0]], + indexMap[pointIndices[1]], + indexMap[pointIndices[2]] + ); + } + } + return new Uint32Array(validTriangles); +} + +export function computeTriangleAverage( + value0: number, + value1: number, + value2: number +) { + const valid0 = Number.isFinite(value0); + const valid1 = Number.isFinite(value1); + const valid2 = Number.isFinite(value2); + if (!valid0 && !valid1 && !valid2) { + return NaN; + } + let sum = 0; + let count = 0; + if (valid0) { + sum += value0; + count++; + } + if (valid1) { + sum += value1; + count++; + } + if (valid2) { + sum += value2; + count++; + } + return sum / count; +} + +function shouldFlipTriangle( + positions: Float32Array, + index0: number, + index1: number, + index2: number +) { + const offset0 = index0 * 3; + const offset1 = index1 * 3; + const offset2 = index2 * 3; + const x0 = positions[offset0]; + const y0 = positions[offset0 + 1]; + const z0 = positions[offset0 + 2]; + const edgeAx = positions[offset1] - x0; + const edgeAy = positions[offset1 + 1] - y0; + const edgeAz = positions[offset1 + 2] - z0; + const edgeBx = positions[offset2] - x0; + const edgeBy = positions[offset2 + 1] - y0; + const edgeBz = positions[offset2 + 2] - z0; + const normalX = edgeAy * edgeBz - edgeAz * edgeBy; + const normalY = edgeAz * edgeBx - edgeAx * edgeBz; + const normalZ = edgeAx * edgeBy - edgeAy * edgeBx; + const centerX = (x0 + positions[offset1] + positions[offset2]) / 3; + const centerY = (y0 + positions[offset1 + 1] + positions[offset2 + 1]) / 3; + const centerZ = (z0 + positions[offset1 + 2] + positions[offset2 + 2]) / 3; + return normalX * centerX + normalY * centerY + normalZ * centerZ < 0; +} + +export function buildIrregularDelaunayGrid( + latitudes: Float32Array, + longitudes: Float32Array, + latitudeShape: number[], + longitudeShape: number[], + dataLength: number +): TIrregularDelaunayGrid { + const coordinates = reconcileCoordinateArrays( + latitudes, + longitudes, + latitudeShape, + longitudeShape, + dataLength + ); + const positions = new Float32Array(dataLength * 3); + const latLonValues = new Float32Array(dataLength * 2); + const projection = new ProjectionHelper( + PROJECTION_TYPES.NEARSIDE_PERSPECTIVE, + { lat: 0, lon: 0 } + ); + for (let index = 0; index < dataLength; index++) { + projection.projectLatLonToArrays( + coordinates.latitudes[index], + coordinates.longitudes[index], + positions, + index * 3, + latLonValues, + index * 2 + ); + } + return { + ...coordinates, + positions, + latLonValues, + triangleIndices: computeDelaunayTriangulation( + coordinates.latitudes, + coordinates.longitudes + ), + }; +} + +export function getIrregularDelaunayBatchCount( + grid: TIrregularDelaunayGrid, + batchSize: number +) { + return Math.ceil(grid.triangleIndices.length / 3 / batchSize); +} + +function getTriangleAverage( + grid: TIrregularDelaunayGrid, + data: Float32Array, + triangleIndex: number +) { + const offset = triangleIndex * 3; + const index0 = grid.triangleIndices[offset]; + const index1 = grid.triangleIndices[offset + 1]; + const index2 = grid.triangleIndices[offset + 2]; + return computeTriangleAverage(data[index0], data[index1], data[index2]); +} + +export function buildIrregularDelaunayGeometryBatch( + grid: TIrregularDelaunayGrid, + data: Float32Array, + batchIndex: number, + batchSize: number +): TGridPointBatch { + const start = batchIndex * batchSize; + const end = Math.min(start + batchSize, grid.triangleIndices.length / 3); + const vertexCount = (end - start) * 3; + const positionValues = new Float32Array(vertexCount * 3); + const latLonValues = new Float32Array(vertexCount * 2); + const dataValues = new Float32Array(vertexCount); + for (let triangleIndex = start; triangleIndex < end; triangleIndex++) { + const triangleOffset = triangleIndex * 3; + const index0 = grid.triangleIndices[triangleOffset]; + let index1 = grid.triangleIndices[triangleOffset + 1]; + let index2 = grid.triangleIndices[triangleOffset + 2]; + const average = computeTriangleAverage( + data[index0], + data[index1], + data[index2] + ); + if (shouldFlipTriangle(grid.positions, index0, index1, index2)) { + [index1, index2] = [index2, index1]; + } + const batchTriangleIndex = triangleIndex - start; + for (let vertexIndex = 0; vertexIndex < 3; vertexIndex++) { + const sourceIndex = + vertexIndex === 0 ? index0 : vertexIndex === 1 ? index1 : index2; + const destinationIndex = batchTriangleIndex * 3 + vertexIndex; + const sourcePositionOffset = sourceIndex * 3; + const destinationPositionOffset = destinationIndex * 3; + positionValues[destinationPositionOffset] = + grid.positions[sourcePositionOffset]; + positionValues[destinationPositionOffset + 1] = + grid.positions[sourcePositionOffset + 1]; + positionValues[destinationPositionOffset + 2] = + grid.positions[sourcePositionOffset + 2]; + const sourceLatLonOffset = sourceIndex * 2; + const destinationLatLonOffset = destinationIndex * 2; + latLonValues[destinationLatLonOffset] = + grid.latLonValues[sourceLatLonOffset]; + latLonValues[destinationLatLonOffset + 1] = + grid.latLonValues[sourceLatLonOffset + 1]; + dataValues[destinationIndex] = average; + } + } + return { batchIndex, positionValues, latLonValues, dataValues }; +} + +export function buildIrregularDelaunayDataBatch( + grid: TIrregularDelaunayGrid, + data: Float32Array, + batchIndex: number, + batchSize: number +): TGridDataValueBatch { + const start = batchIndex * batchSize; + const end = Math.min(start + batchSize, grid.triangleIndices.length / 3); + const dataValues = new Float32Array((end - start) * 3); + for (let triangleIndex = start; triangleIndex < end; triangleIndex++) { + const average = getTriangleAverage(grid, data, triangleIndex); + const batchOffset = (triangleIndex - start) * 3; + dataValues[batchOffset] = average; + dataValues[batchOffset + 1] = average; + dataValues[batchOffset + 2] = average; + } + return { batchIndex, dataValues }; +} + +export function buildIrregularDelaunayHoverIndexData( + grid: TIrregularDelaunayGrid, + data: Float32Array +): TSerializedGeoSampleIndexData { + return buildSerializedGeoSampleIndexData( + Float64Array.from(grid.latitudes), + Float64Array.from(grid.longitudes), + data.slice() + ); +} diff --git a/src/lib/grids/irregularDelaunayWorkerClient.ts b/src/lib/grids/irregularDelaunayWorkerClient.ts new file mode 100644 index 00000000..6db6feaf --- /dev/null +++ b/src/lib/grids/irregularDelaunayWorkerClient.ts @@ -0,0 +1,61 @@ +import { + copyGridWorkerArray, + createGridGeometryWorkerClient, +} from "./gridGeometryWorkerClient.ts"; +import { GridGeometryWorkerMessageType } from "./gridGeometryWorkerProtocol.ts"; +import type { + TIrregularDelaunayBatch, + TIrregularDelaunayWorkerMetadata, + TIrregularDelaunayWorkerRequest, +} from "./irregularDelaunayWorkerProtocol.ts"; + +type TIrregularDelaunayBuildRequest = { + latitudes: Float32Array; + longitudes: Float32Array; + latitudeShape: number[]; + longitudeShape: number[]; + data: Float32Array; + batchSize: number; + forceGeometryRebuild: boolean; +}; + +const client = createGridGeometryWorkerClient< + TIrregularDelaunayWorkerRequest, + TIrregularDelaunayWorkerMetadata, + TIrregularDelaunayBatch +>( + () => + new Worker(new URL("./irregularDelaunay.worker.ts", import.meta.url), { + type: "module", + }), + true +); + +export function buildIrregularDelaunayGrid( + request: TIrregularDelaunayBuildRequest, + callbacks: { + onMetadata: (metadata: TIrregularDelaunayWorkerMetadata) => void; + onBatch: (batch: TIrregularDelaunayBatch) => void; + } +) { + return client.build((requestId) => { + const latitudes = copyGridWorkerArray(request.latitudes); + const longitudes = copyGridWorkerArray(request.longitudes); + const data = copyGridWorkerArray(request.data); + return { + message: { + ...request, + requestId, + type: GridGeometryWorkerMessageType.BUILD, + latitudes, + longitudes, + data, + latitudeShape: [...request.latitudeShape], + longitudeShape: [...request.longitudeShape], + }, + transfer: [latitudes.buffer, longitudes.buffer, data.buffer], + }; + }, callbacks); +} + +export const terminateIrregularDelaunayWorker = client.terminate; diff --git a/src/lib/grids/irregularDelaunayWorkerProtocol.ts b/src/lib/grids/irregularDelaunayWorkerProtocol.ts new file mode 100644 index 00000000..f60691a5 --- /dev/null +++ b/src/lib/grids/irregularDelaunayWorkerProtocol.ts @@ -0,0 +1,29 @@ +import type { TGridGeometryWorkerResponse } from "./gridGeometryWorkerProtocol.ts"; +import type { + TGridDataValueBatch, + TGridPointBatch, +} from "./gridWorkerTypes.ts"; + +export type TIrregularDelaunayBatch = TGridPointBatch | TGridDataValueBatch; + +export type TIrregularDelaunayWorkerRequest = { + requestId: number; + type: "build"; + latitudes: Float32Array; + longitudes: Float32Array; + latitudeShape: number[]; + longitudeShape: number[]; + data: Float32Array; + batchSize: number; + forceGeometryRebuild: boolean; +}; + +export type TIrregularDelaunayWorkerMetadata = { + totalBatches: number; + rebuildGeometry: boolean; +}; + +export type TIrregularDelaunayWorkerResponse = TGridGeometryWorkerResponse< + TIrregularDelaunayWorkerMetadata, + TIrregularDelaunayBatch +>; diff --git a/src/lib/grids/irregularWorkerClient.ts b/src/lib/grids/irregularWorkerClient.ts new file mode 100644 index 00000000..02794cf2 --- /dev/null +++ b/src/lib/grids/irregularWorkerClient.ts @@ -0,0 +1,70 @@ +import { + copyGridWorkerArray, + createGridGeometryWorkerClient, +} from "./gridGeometryWorkerClient.ts"; +import { GridGeometryWorkerMessageType } from "./gridGeometryWorkerProtocol.ts"; +import type { TGridPointBatch } from "./gridWorkerTypes.ts"; +import type { + TIrregularWorkerMetadata, + TIrregularWorkerRequest, +} from "./irregularWorkerProtocol.ts"; + +import type { + TProjectionCenter, + TProjectionType, +} from "@/lib/projection/projectionUtils.ts"; + +type TIrregularBuildRequest = { + latitudes: Float32Array; + longitudes: Float32Array; + latitudeShape: number[]; + longitudeShape: number[]; + data: Float32Array; + batchSize: number; + projectionType: TProjectionType; + projectionCenter: TProjectionCenter; +}; + +const client = createGridGeometryWorkerClient< + TIrregularWorkerRequest, + TIrregularWorkerMetadata, + TGridPointBatch +>( + () => + new Worker(new URL("./irregular.worker.ts", import.meta.url), { + type: "module", + }) +); + +export function buildIrregularGrid( + request: TIrregularBuildRequest, + callbacks: { + onMetadata: (metadata: TIrregularWorkerMetadata) => void; + onBatch: (batch: TGridPointBatch) => void; + } +) { + return client.build((requestId) => { + const latitudes = copyGridWorkerArray(request.latitudes); + const longitudes = copyGridWorkerArray(request.longitudes); + const data = copyGridWorkerArray(request.data); + return { + message: { + ...request, + requestId, + type: GridGeometryWorkerMessageType.BUILD, + latitudes, + longitudes, + data, + latitudeShape: [...request.latitudeShape], + longitudeShape: [...request.longitudeShape], + projectionCenter: { + lat: request.projectionCenter.lat, + lon: request.projectionCenter.lon, + }, + }, + transfer: [latitudes.buffer, longitudes.buffer, data.buffer], + }; + }, callbacks); +} + +export const terminateIrregularWorker = client.terminate; diff --git a/src/lib/grids/irregularWorkerProtocol.ts b/src/lib/grids/irregularWorkerProtocol.ts new file mode 100644 index 00000000..b20d67d2 --- /dev/null +++ b/src/lib/grids/irregularWorkerProtocol.ts @@ -0,0 +1,30 @@ +import type { TGridGeometryWorkerResponse } from "./gridGeometryWorkerProtocol.ts"; +import type { TGridPointBatch } from "./gridWorkerTypes.ts"; + +import type { + TProjectionCenter, + TProjectionType, +} from "@/lib/projection/projectionUtils.ts"; + +export type TIrregularWorkerRequest = { + requestId: number; + type: "build"; + latitudes: Float32Array; + longitudes: Float32Array; + latitudeShape: number[]; + longitudeShape: number[]; + data: Float32Array; + batchSize: number; + projectionType: TProjectionType; + projectionCenter: TProjectionCenter; +}; + +export type TIrregularWorkerMetadata = { + totalBatches: number; + estimatedSpacing: number; +}; + +export type TIrregularWorkerResponse = TGridGeometryWorkerResponse< + TIrregularWorkerMetadata, + TGridPointBatch +>; diff --git a/src/lib/grids/regularGridDataWorkerProtocol.ts b/src/lib/grids/regularGridDataWorkerProtocol.ts deleted file mode 100644 index acf6b7f2..00000000 --- a/src/lib/grids/regularGridDataWorkerProtocol.ts +++ /dev/null @@ -1,36 +0,0 @@ -import type * as zarr from "zarrita"; - -import type { TDataSource, TZarrFormat } from "@/lib/types/GlobeTypes.ts"; - -export const RegularGridDataWorkerMessageType = { - GET_DATA: "getData", - RESULT: "result", - ERROR: "error", -} as const; - -type TRegularGridDataWorkerMessageType = - (typeof RegularGridDataWorkerMessageType)[keyof typeof RegularGridDataWorkerMessageType]; - -export type TRegularGridDataWorkerRequest = { - requestId: number; - type: typeof RegularGridDataWorkerMessageType.GET_DATA; - source: Pick; - variable: string; - format: TZarrFormat; - selection: (number | null | zarr.Slice)[]; -}; - -type TRegularGridDataWorkerResponseBase = { - requestId: number; - type: TRegularGridDataWorkerMessageType; -}; - -export type TRegularGridDataWorkerResponse = - | (TRegularGridDataWorkerResponseBase & { - type: typeof RegularGridDataWorkerMessageType.RESULT; - data: zarr.TypedArray; - }) - | (TRegularGridDataWorkerResponseBase & { - type: typeof RegularGridDataWorkerMessageType.ERROR; - message: string; - }); diff --git a/src/lib/grids/serializedGeoSampleIndex.ts b/src/lib/grids/serializedGeoSampleIndex.ts new file mode 100644 index 00000000..e927a548 --- /dev/null +++ b/src/lib/grids/serializedGeoSampleIndex.ts @@ -0,0 +1,37 @@ +import { around } from "geokdbush"; +import KDBush from "kdbush"; + +import type { TSerializedGeoSampleIndexData } from "./gridWorkerTypes.ts"; + +import { ProjectionHelper } from "@/lib/projection/projectionUtils.ts"; + +type TSerializedGeoSampleIndex = { + findNearest( + lat: number, + lon: number + ): { lat: number; lon: number; value: number } | null; +}; + +export function createSerializedGeoSampleIndex( + data: TSerializedGeoSampleIndexData +): TSerializedGeoSampleIndex { + const index = KDBush.from(data.indexData); + return { + findNearest(lat, lon) { + const nearest = around( + index, + ProjectionHelper.normalizeLongitude(lon), + Math.max(-90, Math.min(90, lat)), + 1 + )[0]; + if (nearest === undefined) { + return null; + } + return { + lat: data.latitudes[nearest], + lon: data.longitudes[nearest], + value: data.values[nearest], + }; + }, + }; +} diff --git a/src/ui/grids/Curvilinear.vue b/src/ui/grids/Curvilinear.vue old mode 100755 new mode 100644 index 848fbed4..e4e44fd5 --- a/src/ui/grids/Curvilinear.vue +++ b/src/ui/grids/Curvilinear.vue @@ -1,13 +1,10 @@ diff --git a/src/ui/grids/GaussianReduced.vue b/src/ui/grids/GaussianReduced.vue index e41e7147..79563d70 100755 --- a/src/ui/grids/GaussianReduced.vue +++ b/src/ui/grids/GaussianReduced.vue @@ -1,13 +1,10 @@