176 lines
6.7 KiB
TypeScript
176 lines
6.7 KiB
TypeScript
// Quick test script to verify DEM tile decoding accuracy
|
|
// Tests the exact coordinates from the user's screenshot
|
|
import * as fs from 'fs';
|
|
import * as path from 'path';
|
|
import * as zlib from 'zlib';
|
|
import * as https from 'https';
|
|
|
|
const DEM_DIR = process.env.DEM_DATA_DIR || path.join(process.cwd(), 'infrastructure/osm-data/dem_tiles');
|
|
|
|
function decodePngRgba(buffer: Buffer): { width: number; height: number; pixels: Uint8Array } | null {
|
|
if (buffer.length < 8 || buffer[0] !== 0x89 || buffer[1] !== 0x50) return null;
|
|
|
|
let offset = 8;
|
|
let width = 256, height = 256, bitDepth = 8, colorType = 2;
|
|
const idatChunks: Buffer[] = [];
|
|
|
|
while (offset < buffer.length) {
|
|
const length = buffer.readUInt32BE(offset);
|
|
const type = buffer.toString('ascii', offset + 4, offset + 8);
|
|
if (type === 'IHDR') {
|
|
width = buffer.readUInt32BE(offset + 8);
|
|
height = buffer.readUInt32BE(offset + 12);
|
|
bitDepth = buffer[offset + 16];
|
|
colorType = buffer[offset + 17];
|
|
} else if (type === 'IDAT') {
|
|
idatChunks.push(buffer.subarray(offset + 8, offset + 8 + length));
|
|
} else if (type === 'IEND') break;
|
|
offset += 12 + length;
|
|
}
|
|
|
|
if (idatChunks.length === 0) return null;
|
|
const decompressed = zlib.inflateSync(Buffer.concat(idatChunks));
|
|
const channels = colorType === 6 ? 4 : colorType === 2 ? 3 : 1;
|
|
const bytesPerPixel = channels;
|
|
const scanlineLength = 1 + width * bytesPerPixel;
|
|
const rawPixels = new Uint8Array(width * height * 4);
|
|
let prevRow = new Uint8Array(width * bytesPerPixel);
|
|
|
|
for (let y = 0; y < height; y++) {
|
|
const rowStart = y * scanlineLength;
|
|
const filterType = decompressed[rowStart];
|
|
const currentRow = new Uint8Array(width * bytesPerPixel);
|
|
|
|
for (let x = 0; x < width * bytesPerPixel; x++) {
|
|
const rawByte = decompressed[rowStart + 1 + x];
|
|
const a = x >= bytesPerPixel ? currentRow[x - bytesPerPixel] : 0;
|
|
const b = prevRow[x];
|
|
const c = x >= bytesPerPixel ? prevRow[x - bytesPerPixel] : 0;
|
|
|
|
let val = rawByte;
|
|
if (filterType === 1) val = (rawByte + a) & 0xff;
|
|
else if (filterType === 2) val = (rawByte + b) & 0xff;
|
|
else if (filterType === 3) val = (rawByte + Math.floor((a + b) / 2)) & 0xff;
|
|
else if (filterType === 4) {
|
|
// CORRECTED Paeth per PNG spec (using <=)
|
|
const p = a + b - c;
|
|
const pa = Math.abs(p - a);
|
|
const pb = Math.abs(p - b);
|
|
const pc = Math.abs(p - c);
|
|
let pr: number;
|
|
if (pa <= pb && pa <= pc) pr = a;
|
|
else if (pb <= pc) pr = b;
|
|
else pr = c;
|
|
val = (rawByte + pr) & 0xff;
|
|
}
|
|
currentRow[x] = val;
|
|
}
|
|
|
|
for (let px = 0; px < width; px++) {
|
|
const srcIdx = px * bytesPerPixel;
|
|
const dstIdx = (y * width + px) * 4;
|
|
rawPixels[dstIdx] = currentRow[srcIdx];
|
|
rawPixels[dstIdx + 1] = currentRow[srcIdx + 1];
|
|
rawPixels[dstIdx + 2] = currentRow[srcIdx + 2];
|
|
rawPixels[dstIdx + 3] = channels === 4 ? currentRow[srcIdx + 3] : 255;
|
|
}
|
|
prevRow = currentRow;
|
|
}
|
|
|
|
return { width, height, pixels: rawPixels };
|
|
}
|
|
|
|
function coordToTile(lat: number, lng: number, zoom: number) {
|
|
const n = Math.pow(2, zoom);
|
|
const xVal = ((lng + 180.0) / 360.0) * n;
|
|
const latRad = lat * (Math.PI / 180.0);
|
|
const yVal = (1.0 - Math.log(Math.tan(latRad) + 1.0 / Math.cos(latRad)) / Math.PI) / 2.0 * n;
|
|
return {
|
|
tileX: Math.floor(xVal),
|
|
tileY: Math.floor(yVal),
|
|
subX: (xVal - Math.floor(xVal)) * 256.0,
|
|
subY: (yVal - Math.floor(yVal)) * 256.0,
|
|
};
|
|
}
|
|
|
|
function sampleElevation(pixels: Uint8Array, width: number, subX: number, subY: number): number {
|
|
const x0 = Math.min(Math.floor(subX), width - 1);
|
|
const y0 = Math.min(Math.floor(subY), width - 1);
|
|
const idx = (y0 * width + x0) * 4;
|
|
const r = pixels[idx], g = pixels[idx + 1], b = pixels[idx + 2];
|
|
return (r * 256.0 + g + b / 256.0) - 32768.0;
|
|
}
|
|
|
|
// Test points from the screenshot:
|
|
// Observer A: 35.9106, 31.9539
|
|
// Target B: 36.1032, 32.0608
|
|
// Path should be ~500-700m elevation throughout (contour lines show 540m-600m+ area)
|
|
|
|
const testPoints = [
|
|
{ name: 'Observer A', lat: 31.9539, lng: 35.9106 },
|
|
{ name: 'Target B', lat: 32.0608, lng: 36.1032 },
|
|
// Midpoints along the path
|
|
{ name: 'Mid 25%', lat: 31.9539 + (32.0608-31.9539)*0.25, lng: 35.9106 + (36.1032-35.9106)*0.25 },
|
|
{ name: 'Mid 50%', lat: 31.9539 + (32.0608-31.9539)*0.50, lng: 35.9106 + (36.1032-35.9106)*0.50 },
|
|
{ name: 'Mid 75%', lat: 31.9539 + (32.0608-31.9539)*0.75, lng: 35.9106 + (36.1032-35.9106)*0.75 },
|
|
];
|
|
|
|
async function fetchTile(zoom: number, x: number, y: number): Promise<Buffer | null> {
|
|
const localPath = path.join(DEM_DIR, `${zoom}_${x}_${y}.png`);
|
|
if (fs.existsSync(localPath)) {
|
|
console.log(` -> Tile ${zoom}/${x}/${y} loaded from disk cache`);
|
|
return fs.readFileSync(localPath);
|
|
}
|
|
|
|
return new Promise((resolve) => {
|
|
const url = `https://s3.amazonaws.com/elevation-tiles-prod/terrarium/${zoom}/${x}/${y}.png`;
|
|
console.log(` -> Fetching tile from ${url}`);
|
|
https.get(url, (res) => {
|
|
const chunks: Buffer[] = [];
|
|
res.on('data', (c: Buffer) => chunks.push(c));
|
|
res.on('end', () => {
|
|
const buf = Buffer.concat(chunks);
|
|
// Save to disk
|
|
if (!fs.existsSync(DEM_DIR)) fs.mkdirSync(DEM_DIR, { recursive: true });
|
|
fs.writeFileSync(localPath, buf);
|
|
resolve(buf);
|
|
});
|
|
}).on('error', () => resolve(null));
|
|
});
|
|
}
|
|
|
|
async function main() {
|
|
console.log('=== DEM Tile Elevation Verification ===\n');
|
|
console.log(`DEM Directory: ${DEM_DIR}\n`);
|
|
|
|
for (const pt of testPoints) {
|
|
const { tileX, tileY, subX, subY } = coordToTile(pt.lat, pt.lng, 13);
|
|
console.log(`${pt.name}: (${pt.lat.toFixed(4)}, ${pt.lng.toFixed(4)}) -> Tile 13/${tileX}/${tileY}, subpixel (${subX.toFixed(1)}, ${subY.toFixed(1)})`);
|
|
|
|
const buf = await fetchTile(13, tileX, tileY);
|
|
if (!buf) {
|
|
console.log(` FAILED to load tile\n`);
|
|
continue;
|
|
}
|
|
|
|
const decoded = decodePngRgba(buf);
|
|
if (!decoded) {
|
|
console.log(` FAILED to decode PNG\n`);
|
|
continue;
|
|
}
|
|
|
|
const elev = sampleElevation(decoded.pixels, decoded.width, subX, subY);
|
|
console.log(` PNG: ${decoded.width}x${decoded.height}, colorType channels, size=${buf.length} bytes`);
|
|
console.log(` ELEVATION = ${elev.toFixed(1)}m`);
|
|
|
|
// Also check raw pixel values
|
|
const x0 = Math.min(Math.floor(subX), decoded.width - 1);
|
|
const y0 = Math.min(Math.floor(subY), decoded.height - 1);
|
|
const idx = (y0 * decoded.width + x0) * 4;
|
|
console.log(` Raw RGB: R=${decoded.pixels[idx]}, G=${decoded.pixels[idx+1]}, B=${decoded.pixels[idx+2]}`);
|
|
console.log(` Formula: (${decoded.pixels[idx]}*256 + ${decoded.pixels[idx+1]} + ${decoded.pixels[idx+2]}/256) - 32768 = ${elev.toFixed(1)}\n`);
|
|
}
|
|
}
|
|
|
|
main().catch(console.error);
|