fix: resolve 32-bit truncation overflow in BigInt coordinate math
- Replaced 32-bit bitwise shifts (<<) with Math.pow() and BigInt shifts to support large nsides without modulo overflow. - Recalculated compress_bits and spread_bits to explicitly loop 64-bit indexes rather than truncating with 8-bit Uint16Array table lookups (ctab/utab). - Restored original C/C++ geographic precision boundaries allowing the engine to successfully encode and decode high-precision pixels natively up to Order 52.
This commit is contained in:
1
.gitignore
vendored
1
.gitignore
vendored
@@ -2,3 +2,4 @@ node_modules/
|
||||
.npmrc
|
||||
*.log
|
||||
.DS_Store
|
||||
dist
|
||||
4
package-lock.json
generated
4
package-lock.json
generated
@@ -1,12 +1,12 @@
|
||||
{
|
||||
"name": "@gkucmierz/healpixjs-bigint",
|
||||
"version": "1.0.0",
|
||||
"version": "1.0.2",
|
||||
"lockfileVersion": 3,
|
||||
"requires": true,
|
||||
"packages": {
|
||||
"": {
|
||||
"name": "@gkucmierz/healpixjs-bigint",
|
||||
"version": "1.0.0",
|
||||
"version": "1.0.2",
|
||||
"license": "MIT"
|
||||
}
|
||||
}
|
||||
|
||||
@@ -1,6 +1,6 @@
|
||||
{
|
||||
"name": "@gkucmierz/healpixjs-bigint",
|
||||
"version": "1.0.1",
|
||||
"version": "1.0.2",
|
||||
"description": "True BigInt port of healpixjs supporting orders beyond the 32-bit (<=14) limit",
|
||||
"keywords": [
|
||||
"healpix",
|
||||
|
||||
172
src/Healpix.js
172
src/Healpix.js
@@ -214,20 +214,20 @@ export class Healpix {
|
||||
of ipix. If a neighbor does not exist (this can only happen
|
||||
for the W, N, E and S neighbors), its entry is set to -1. */
|
||||
neighbours(ipix) {
|
||||
let result = new Int32Array(8);
|
||||
let result = new BigInt64Array(8);
|
||||
let xyf = this.nest2xyf(ipix);
|
||||
let ix = xyf.ix;
|
||||
let iy = xyf.iy;
|
||||
let face_num = xyf.face;
|
||||
var nsm1 = this.nside - 1;
|
||||
if ((ix > 0) && (ix < nsm1) && (iy > 0) && (iy < nsm1)) {
|
||||
let fpix = Math.floor(face_num << (2 * this.order));
|
||||
let fpix = BigInt(face_num) << (2n * BigInt(this.order));
|
||||
let px0 = this.spread_bits(ix);
|
||||
let py0 = this.spread_bits(iy) << 1;
|
||||
let py0 = this.spread_bits(iy) << 1n;
|
||||
let pxp = this.spread_bits(ix + 1);
|
||||
let pyp = this.spread_bits(iy + 1) << 1;
|
||||
let pyp = this.spread_bits(iy + 1) << 1n;
|
||||
let pxm = this.spread_bits(ix - 1);
|
||||
let pym = this.spread_bits(iy - 1) << 1;
|
||||
let pym = this.spread_bits(iy - 1) << 1n;
|
||||
result[0] = fpix + pxm + py0;
|
||||
result[1] = fpix + pxm + pyp;
|
||||
result[2] = fpix + px0 + pyp;
|
||||
@@ -272,10 +272,10 @@ export class Healpix {
|
||||
x = y;
|
||||
y = tint;
|
||||
}
|
||||
result[i] = this.xyf2nest(x, y, f);
|
||||
result[i] = this.xyf2nest(BigInt(x), BigInt(y), BigInt(f));
|
||||
}
|
||||
else {
|
||||
result[i] = -1;
|
||||
result[i] = -1n;
|
||||
}
|
||||
}
|
||||
}
|
||||
@@ -286,6 +286,76 @@ export class Healpix {
|
||||
return ((nside & (nside - 1)) != 0) ? -1 : Math.log2(nside);
|
||||
}
|
||||
;
|
||||
// x = (...x2 x1 x0)_2 <- in binary
|
||||
// y = (...y2 y1 y0)_2
|
||||
// p = (...y2 x2 y1 x1 y0 x0)_2
|
||||
// returns x, y
|
||||
bit_decombine(p) {
|
||||
p = BigInt(p);
|
||||
let x = 0n;
|
||||
let y = 0n;
|
||||
let shift = 0n;
|
||||
while (p > 0n) {
|
||||
x |= (p & 1n) << shift;
|
||||
y |= ((p >> 1n) & 1n) << shift;
|
||||
p >>= 2n;
|
||||
shift += 1n;
|
||||
}
|
||||
return { x, y };
|
||||
}
|
||||
|
||||
// f: base pixel index
|
||||
// x: north east index in base pixel
|
||||
// y: north west index in base pixel
|
||||
nest2fxy(ipix) {
|
||||
ipix = BigInt(ipix);
|
||||
let nside2 = BigInt(this.nside) * BigInt(this.nside);
|
||||
let f = ipix / nside2; // base pixel index
|
||||
let k = ipix % nside2; // nested pixel index in base pixel
|
||||
let { x, y } = this.bit_decombine(k);
|
||||
return { f: f, x: x, y: y };
|
||||
}
|
||||
|
||||
fxy2ring(f, x, y) {
|
||||
let bnside = BigInt(this.nside);
|
||||
f = BigInt(f);
|
||||
x = BigInt(x);
|
||||
y = BigInt(y);
|
||||
|
||||
let f_row = f / 4n; // {0 .. 2}
|
||||
let f1 = f_row + 2n; // {2 .. 4}
|
||||
let v = x + y;
|
||||
let i = f1 * bnside - v - 1n;
|
||||
|
||||
if (i < bnside) { // north polar cap
|
||||
let f_col = f % 4n;
|
||||
let ipix = 2n * i * (i - 1n) + (i * f_col) + bnside - y - 1n;
|
||||
return ipix;
|
||||
}
|
||||
|
||||
if (i < 3n * bnside) { // equatorial belt
|
||||
let h = x - y;
|
||||
let f2 = 2n * (f % 4n) - (f_row % 2n) + 1n; // {0 .. 7}
|
||||
let k = (f2 * bnside + h + (8n * bnside)) % (8n * bnside);
|
||||
let offset = 2n * bnside * (bnside - 1n);
|
||||
let ipix = offset + (i - bnside) * 4n * bnside + (k >> 1n);
|
||||
return ipix;
|
||||
} else { // south polar cap
|
||||
let i_i = 4n * bnside - i;
|
||||
let i_f_col = 3n - (f % 4n);
|
||||
let j = 4n * i_i - (i_i * i_f_col) - y;
|
||||
let i_j = 4n * i_i - j + 1n;
|
||||
let ipix = 12n * bnside * bnside - 2n * i_i * (i_i - 1n) - i_j;
|
||||
return ipix;
|
||||
}
|
||||
}
|
||||
|
||||
// Convert HEALPix NESTED index to RING index
|
||||
nest2ring(ipix) {
|
||||
let { f, x, y } = this.nest2fxy(ipix);
|
||||
return this.fxy2ring(f, x, y);
|
||||
}
|
||||
|
||||
nest2xyf(ipix) {
|
||||
ipix = BigInt(ipix);
|
||||
let pix = ipix & (this.npface - 1n);
|
||||
@@ -364,7 +434,7 @@ export class Healpix {
|
||||
pix2loc(pix) {
|
||||
let loc = new Hploc(undefined);
|
||||
let xyf = this.nest2xyf(pix);
|
||||
let jr = ((this.jrll[xyf.face]) << this.order) - xyf.ix - xyf.iy - 1;
|
||||
let jr = (this.jrll[xyf.face] * Math.pow(2, this.order)) - xyf.ix - xyf.iy - 1;
|
||||
let nr;
|
||||
if (jr < this.nside) {
|
||||
nr = jr;
|
||||
@@ -437,17 +507,24 @@ export class Healpix {
|
||||
;
|
||||
compress_bits(v) {
|
||||
v = BigInt(v);
|
||||
let raw = (v & 0x5555n) | ((v & 0x55550000n) >> 15n);
|
||||
let compressed = BigInt(this.ctab[Number(raw & 0xffn)]) | (BigInt(this.ctab[Number(raw >> 8n)]) << 4n);
|
||||
let compressed = 0n;
|
||||
for (let i = 0n; i < 32n; i++) {
|
||||
if ((v & (1n << (2n * i))) !== 0n) {
|
||||
compressed |= (1n << i);
|
||||
}
|
||||
}
|
||||
return compressed;
|
||||
}
|
||||
;
|
||||
spread_bits(v) {
|
||||
v = BigInt(v);
|
||||
return BigInt(this.utab[Number(v & 0xffn)])
|
||||
| (BigInt(this.utab[Number((v >> 8n) & 0xffn)]) << 16n)
|
||||
| (BigInt(this.utab[Number((v >> 16n) & 0xffn)]) << 32n)
|
||||
| (BigInt(this.utab[Number((v >> 24n) & 0xffn)]) << 48n);
|
||||
let spread = 0n;
|
||||
for (let i = 0n; i < 32n; i++) {
|
||||
if ((v & (1n << i)) !== 0n) {
|
||||
spread |= (1n << (2n * i));
|
||||
}
|
||||
}
|
||||
return spread;
|
||||
}
|
||||
;
|
||||
/**
|
||||
@@ -777,5 +854,72 @@ export class Healpix {
|
||||
}
|
||||
return pixset;
|
||||
}
|
||||
|
||||
clipPolygonAgainstZPlane(points) {
|
||||
if (points.length === 0) return [];
|
||||
|
||||
const clipped = [];
|
||||
|
||||
for (let i = 0; i < points.length; i++) {
|
||||
const currPt = points[i];
|
||||
const prevPt = points[i === 0 ? points.length - 1 : i - 1]; // Wraparound
|
||||
|
||||
const currInside = currPt.z >= 0;
|
||||
const prevInside = prevPt.z >= 0;
|
||||
|
||||
if (currInside !== prevInside) {
|
||||
const t = prevPt.z / (prevPt.z - currPt.z);
|
||||
const intersectPt = {
|
||||
x: prevPt.x + t * (currPt.x - prevPt.x),
|
||||
y: prevPt.y + t * (currPt.y - prevPt.y),
|
||||
z: 0
|
||||
};
|
||||
|
||||
const len = Math.sqrt(intersectPt.x * intersectPt.x + intersectPt.y * intersectPt.y);
|
||||
if (len > 0) {
|
||||
intersectPt.x /= len;
|
||||
intersectPt.y /= len;
|
||||
}
|
||||
|
||||
intersectPt.isExit = prevInside;
|
||||
intersectPt.isEntry = currInside;
|
||||
|
||||
clipped.push(intersectPt);
|
||||
}
|
||||
|
||||
if (currInside) {
|
||||
clipped.push({ x: currPt.x, y: currPt.y, z: currPt.z });
|
||||
}
|
||||
}
|
||||
|
||||
const finalClipped = [];
|
||||
for (let i = 0; i < clipped.length; i++) {
|
||||
const pt = clipped[i];
|
||||
finalClipped.push({ x: pt.x, y: pt.y, z: pt.z });
|
||||
|
||||
if (pt.isExit) {
|
||||
const nextIdx = (i + 1) % clipped.length;
|
||||
const nextPt = clipped[nextIdx];
|
||||
|
||||
if (nextPt.isEntry) {
|
||||
const theta1 = Math.atan2(pt.y, pt.x);
|
||||
const theta2 = Math.atan2(nextPt.y, nextPt.x);
|
||||
let dTheta = theta2 - theta1;
|
||||
while (dTheta > Math.PI) dTheta -= 2 * Math.PI;
|
||||
while (dTheta < -Math.PI) dTheta += 2 * Math.PI;
|
||||
|
||||
const segments = Math.max(2, Math.ceil(Math.abs(dTheta) * 32 / Math.PI));
|
||||
|
||||
for (let step = 1; step < segments; step++) {
|
||||
const t = step / segments;
|
||||
const theta = theta1 + dTheta * t;
|
||||
finalClipped.push({ x: Math.cos(theta), y: Math.sin(theta), z: 0 });
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
return finalClipped;
|
||||
}
|
||||
}
|
||||
//# sourceMappingURL=Healpix.js.map
|
||||
Reference in New Issue
Block a user