Skip to content

Commit

Permalink
feat: bilinear resampling method
Browse files Browse the repository at this point in the history
  • Loading branch information
hongfaqiu committed Jul 17, 2024
1 parent 082d6ae commit 2aea3b1
Show file tree
Hide file tree
Showing 10 changed files with 141 additions and 128 deletions.
3 changes: 2 additions & 1 deletion .gitignore
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
node_modules
/lerna-debug.log
/vite-example/dist
/.turbo
/.turbo
/vite-example/public/*.tif*
3 changes: 3 additions & 0 deletions .vscode/settings.json
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@
{
"commentTranslate.targetLanguage": "zh-CN"
}
1 change: 1 addition & 0 deletions packages/TIFFImageryProvider/README.md
Original file line number Diff line number Diff line change
Expand Up @@ -172,6 +172,7 @@ interface TIFFImageryProviderOptions {
cache?: number;
/** resample web worker pool size, defaults to the number of CPUs available. When this parameter is `null` or 0, then the resampling will be done in the main thread. */
workerPoolSize?: number;
/** resample method, defaults to nearest */
resampleMethod?: 'bilinear' | 'nearest';
}

Expand Down
1 change: 1 addition & 0 deletions packages/TIFFImageryProvider/README_CN.md
Original file line number Diff line number Diff line change
Expand Up @@ -168,6 +168,7 @@ interface TIFFImageryProviderOptions {
cache?: number;
/** 重新采样 Web Worker 工作池大小,默认为可用 CPU 数量。当该参数为null或 0,则重采样将在主线程中完成。 */
workerPoolSize?: number;
/** 重采样方法,默认为 nearest */
resampleMethod?: 'bilinear' | 'nearest'
}

Expand Down
60 changes: 34 additions & 26 deletions packages/TIFFImageryProvider/src/TIFFImageryProvider.ts
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ import { Event, GeographicTilingScheme, Credit, Rectangle, ImageryLayerFeatureIn
import GeoTIFF, { Pool, fromUrl, fromBlob, GeoTIFFImage, TypedArrayArrayWithDimensions } from 'geotiff';

import { addColorScale, plot } from './plotty'
import { getMinMax, generateColorScale, findAndSortBandNumbers, stringColorToRgba, ResampleDataOptions } from "./helpers/utils";
import { getMinMax, generateColorScale, findAndSortBandNumbers, stringColorToRgba } from "./helpers/utils";
import { ColorScaleNames, TypedArray } from "./plotty/typing";
import TIFFImageryProviderTilingScheme from "./TIFFImageryProviderTilingScheme";
import { BBox, reprojection } from "./helpers/reprojection";
Expand All @@ -11,6 +11,7 @@ import { GenerateImageOptions, generateImage } from "./helpers/generateImage";
import { reverseArray } from "./helpers/utils";
import { createCanavas } from "./helpers/createCanavas";
import WorkerPool from "./worker/pool";
import { ResampleDataOptions } from "./helpers/resample";

export interface SingleBandRenderOptions {
/** band index start from 1, defaults to 1 */
Expand Down Expand Up @@ -138,8 +139,12 @@ export interface TIFFImageryProviderOptions {
} | undefined;
/** cache survival time, defaults to 60 * 1000 ms */
cache?: number;
/** resample web worker pool size, defaults to the number of CPUs available. When this parameter is `null` or 0, then the resampling will be done in the main thread. */
/** resample web worker pool size, defaults to the number of CPUs available.
* When this parameter is `null` or 0,
* then the resampling will be done in the main thread.
* */
workerPoolSize?: number;
/** resample method, defaults to nearest */
resampleMethod?: ResampleDataOptions['method']
}

Expand Down Expand Up @@ -240,7 +245,7 @@ export class TIFFImageryProvider {

this._source = source;

// 获取空间范围
// get bounding box
this.origin = this._getOrigin(image);
this.bbox = image.getBoundingBox();
this.reverseY = this._checkIfReversed(image);
Expand Down Expand Up @@ -273,27 +278,27 @@ export class TIFFImageryProvider {
}

this.rectangle = this.tilingScheme.rectangle
// 处理跨180度经线的情况
// dealing with situations crossing the 180 degree meridian
// https://github.com/CesiumGS/cesium/blob/da00d26473f663db180cacd8e662ca4309e09560/packages/engine/Source/Core/TileAvailability.js#L195
if (this.rectangle.east < this.rectangle.west) {
this.rectangle.east += CesiumMath.TWO_PI;
}
this._imageCount = await source.getImageCount();
this.tileSize = this.tileWidth = tileSize || (this._isTiled ? image.getTileWidth() : image.getWidth()) || 256;
this.tileHeight = tileSize || (this._isTiled ? image.getTileHeight() : image.getHeight()) || 256;
// 获取合适的COG层级
// get the appropriate COG level
this.requestLevels = this._isTiled ? await this._getCogLevels() : [0];
this._images = new Array(this._imageCount).fill(null);

// 获取波段数
// Get the number of bands
const samples = image.getSamplesPerPixel();
this.samples = samples;
this.renderOptions = renderOptions ?? {}
// 获取nodata值
// Get nodata value
const noData = image.getGDALNoData();
this.noData = this.renderOptions.nodata ?? noData;

// 赋初值
// Assign initial value
if (samples < 3 && this.renderOptions.convertToRGB) {
const error = new DeveloperError('Can not render the image as RGB, please check the convertToRGB parameter')
throw error;
Expand Down Expand Up @@ -323,7 +328,7 @@ export class TIFFImageryProvider {
this.readSamples = findAndSortBandNumbers(single.expression);
}

// 获取波段最大最小值信息
// Get the maximum and minimum value information of the band
const bands: Record<number, {
min: number;
max: number;
Expand Down Expand Up @@ -363,7 +368,7 @@ export class TIFFImageryProvider {
}

if (!single?.expression && !bands[bandNum]) {
// 尝试获取波段最大最小值
// Try to get the maximum and minimum values ​​of the band
console.warn(`Can not get band${bandNum} min/max, try to calculate min/max values, or setting ${single ? 'domain' : 'min / max'}`)

const previewImage = await source.getImage(this.requestLevels[0])
Expand All @@ -377,7 +382,7 @@ export class TIFFImageryProvider {
}))
this.bands = bands;

// 如果是单通道渲染, 则构建plot对象
// If it is single-pass rendering, build the plot object
try {
if (this.renderOptions.single) {
const band = this.bands[single.band];
Expand Down Expand Up @@ -464,7 +469,7 @@ export class TIFFImageryProvider {
const height = image.getHeight();
const size = Math.max(width, height);

// 如果第一张瓦片的image tileSize大于256,则顺位后延,以减少请求量
// If the image tileSize of the first tile is greater than 256, the sequence will be delayed to reduce the amount of requests.
if (i === this._imageCount - 1) {
const firstImageLevel = Math.ceil((size - this.tileSize) / this.tileSize)
levels.push(...new Array(firstImageLevel).fill(i))
Expand All @@ -484,7 +489,7 @@ export class TIFFImageryProvider {
}

/**
* 获取瓦片数据
* Get tile data
* @param x
* @param y
* @param z
Expand Down Expand Up @@ -539,7 +544,9 @@ export class TIFFImageryProvider {
if (this.reverseY) {
window = [window[0], height - window[3], window[2], height - window[1]];
}
const windowWidth = window[2] - window[0], windowHeight = window[3] - window[1];
const buffer = 1;
window = [window[0] - buffer, window[1] - buffer, window[2] + buffer, window[3] + buffer]
const sourceWidth = window[2] - window[0], sourceHeight = window[3] - window[1];

const options = {
window,
Expand Down Expand Up @@ -570,14 +577,15 @@ export class TIFFImageryProvider {

const result: TypedArray[] = [];
for (let i = 0; i < res.length; i++) {
// TODO Buffer effects are not considered
const prjData = reprojection({
data: res[i] as any,
sourceWidth: windowWidth,
sourceHeight: windowHeight,
sourceWidth,
sourceHeight,
nodata: this.noData,
project: this._proj.project,
sourceBBox,
targetBBox
targetBBox,
})
result.push(prjData)
}
Expand All @@ -592,14 +600,14 @@ export class TIFFImageryProvider {
const y1 = y0 + step;

res = await Promise.all(res.map(async (data) => this.workerPool.resample(data, {
sourceWidth: windowWidth,
sourceHeight: windowHeight,
targetWidth: this.tileWidth,
targetHeight: this.tileHeight,
window: [x0, y0, x1, y1],
method: this.options.resampleMethod ?? 'nearest'
})
));
sourceWidth,
sourceHeight,
targetWidth: this.tileWidth,
targetHeight: this.tileHeight,
window: [x0, y0, x1, y1],
method: this.options.resampleMethod,
buffer,
})));

return {
data: res,
Expand Down Expand Up @@ -708,7 +716,7 @@ export class TIFFImageryProvider {
let posX: number, posY: number, window: number[];
const { west, south, north, width: lonWidth } = this.rectangle;
let lonGap = longitude - west;
// 处理跨180°经线的情况
// Handling cases across 180° longitude
if (longitude < west) {
lonGap += CesiumMath.TWO_PI;
}
Expand Down
2 changes: 1 addition & 1 deletion packages/TIFFImageryProvider/src/helpers/reprojection.ts
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
import { TypedArray } from "geotiff";
import { copyNewSize } from "./utils";
import { copyNewSize } from "./resample";

export type BBox = [minX: number, minY: number, maxX: number, maxY: number];

Expand Down
97 changes: 97 additions & 0 deletions packages/TIFFImageryProvider/src/helpers/resample.ts
Original file line number Diff line number Diff line change
@@ -0,0 +1,97 @@
import { TypedArray } from "geotiff";

export type ResampleDataOptions = {
sourceWidth: number;
sourceHeight: number;
targetWidth: number;
targetHeight: number;
/** start from 0 to 1, examples: [0, 0, 0.5, 0.5] */
window: [number, number, number, number];
method?: 'bilinear' | 'nearest';
buffer?: number;
}

export function resampleNearest(data: TypedArray, options: ResampleDataOptions) {
const { sourceWidth, sourceHeight, targetWidth, targetHeight, window, buffer = 0 } = options;
const [x0, y0, x1, y1] = window;

const effectiveSourceWidth = sourceWidth - 2 * buffer;
const effectiveSourceHeight = sourceHeight - 2 * buffer;

const resampledData = copyNewSize(data, targetWidth, targetHeight);

for (let y = 0; y < targetHeight; y++) {
for (let x = 0; x < targetWidth; x++) {
const col = buffer + ((effectiveSourceWidth * (x0 + x / targetWidth * (x1 - x0))) >>> 0);
const row = buffer + ((effectiveSourceHeight * (y0 + y / targetHeight * (y1 - y0))) >>> 0);
resampledData[y * targetWidth + x] = data[row * sourceWidth + col];
}
}

return resampledData;
}

export function resampleData(data: TypedArray, options: ResampleDataOptions) {
const { method = 'nearest' } = options;

switch (method) {
case "nearest":
return resampleNearest(data, options);
case "bilinear":
return resampleBilinear(data, options);
}
}

export function copyNewSize(array: TypedArray, width: number, height: number, samplesPerPixel = 1) {
return new (Object.getPrototypeOf(array).constructor)(width * height * samplesPerPixel) as typeof array;
}

function lerp(v0: number, v1: number, t: number) {
return ((1 - t) * v0) + (t * v1);
}

/**
* Resample the input array using bilinear interpolation.
* @param data The input arrays to resample
* @param options The options for resampling
* @returns The resampled rasters
*/
export function resampleBilinear(data: TypedArray, options: ResampleDataOptions) {
const { sourceWidth, sourceHeight, targetWidth, targetHeight, window, buffer = 0 } = options;
const [x0, y0, x1, y1] = window;

const windowWidth = x1 - x0;
const windowHeight = y1 - y0;

const newArray = copyNewSize(data, targetWidth, targetHeight);

const effectiveSourceWidth = sourceWidth - 2 * buffer;
const effectiveSourceHeight = sourceHeight - 2 * buffer;

for (let y = 0; y < targetHeight; y++) {
const rawY = effectiveSourceHeight * (y0 + y / targetHeight * windowHeight) + buffer;
const yl = Math.floor(rawY);
const yh = Math.min(Math.ceil(rawY), sourceHeight - 1);

for (let x = 0; x < targetWidth; x++) {
const rawX = effectiveSourceWidth * (x0 + x / targetWidth * windowWidth) + buffer;
const tx = rawX % 1;

const xl = Math.floor(rawX);
const xh = Math.min(Math.ceil(rawX), sourceWidth - 1);

const ll = data[(yl * sourceWidth) + xl];
const hl = data[(yl * sourceWidth) + xh];
const lh = data[(yh * sourceWidth) + xl];
const hh = data[(yh * sourceWidth) + xh];

const value = lerp(
lerp(ll, hl, tx),
lerp(lh, hh, tx),
rawY % 1,
);
newArray[(y * targetWidth) + x] = value;
}
}
return newArray;
}
Loading

0 comments on commit 2aea3b1

Please sign in to comment.