fix: resolve 32-bit truncation overflow in BigInt coordinate math
All checks were successful
npm-publish / publish (push) Successful in 4s
All checks were successful
npm-publish / publish (push) Successful in 4s
- 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
|
.npmrc
|
||||||
*.log
|
*.log
|
||||||
.DS_Store
|
.DS_Store
|
||||||
|
dist
|
||||||
4
package-lock.json
generated
4
package-lock.json
generated
@@ -1,12 +1,12 @@
|
|||||||
{
|
{
|
||||||
"name": "@gkucmierz/healpixjs-bigint",
|
"name": "@gkucmierz/healpixjs-bigint",
|
||||||
"version": "1.0.0",
|
"version": "1.0.2",
|
||||||
"lockfileVersion": 3,
|
"lockfileVersion": 3,
|
||||||
"requires": true,
|
"requires": true,
|
||||||
"packages": {
|
"packages": {
|
||||||
"": {
|
"": {
|
||||||
"name": "@gkucmierz/healpixjs-bigint",
|
"name": "@gkucmierz/healpixjs-bigint",
|
||||||
"version": "1.0.0",
|
"version": "1.0.2",
|
||||||
"license": "MIT"
|
"license": "MIT"
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|||||||
@@ -1,6 +1,6 @@
|
|||||||
{
|
{
|
||||||
"name": "@gkucmierz/healpixjs-bigint",
|
"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",
|
"description": "True BigInt port of healpixjs supporting orders beyond the 32-bit (<=14) limit",
|
||||||
"keywords": [
|
"keywords": [
|
||||||
"healpix",
|
"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
|
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. */
|
for the W, N, E and S neighbors), its entry is set to -1. */
|
||||||
neighbours(ipix) {
|
neighbours(ipix) {
|
||||||
let result = new Int32Array(8);
|
let result = new BigInt64Array(8);
|
||||||
let xyf = this.nest2xyf(ipix);
|
let xyf = this.nest2xyf(ipix);
|
||||||
let ix = xyf.ix;
|
let ix = xyf.ix;
|
||||||
let iy = xyf.iy;
|
let iy = xyf.iy;
|
||||||
let face_num = xyf.face;
|
let face_num = xyf.face;
|
||||||
var nsm1 = this.nside - 1;
|
var nsm1 = this.nside - 1;
|
||||||
if ((ix > 0) && (ix < nsm1) && (iy > 0) && (iy < nsm1)) {
|
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 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 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 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[0] = fpix + pxm + py0;
|
||||||
result[1] = fpix + pxm + pyp;
|
result[1] = fpix + pxm + pyp;
|
||||||
result[2] = fpix + px0 + pyp;
|
result[2] = fpix + px0 + pyp;
|
||||||
@@ -272,10 +272,10 @@ export class Healpix {
|
|||||||
x = y;
|
x = y;
|
||||||
y = tint;
|
y = tint;
|
||||||
}
|
}
|
||||||
result[i] = this.xyf2nest(x, y, f);
|
result[i] = this.xyf2nest(BigInt(x), BigInt(y), BigInt(f));
|
||||||
}
|
}
|
||||||
else {
|
else {
|
||||||
result[i] = -1;
|
result[i] = -1n;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@@ -286,6 +286,76 @@ export class Healpix {
|
|||||||
return ((nside & (nside - 1)) != 0) ? -1 : Math.log2(nside);
|
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) {
|
nest2xyf(ipix) {
|
||||||
ipix = BigInt(ipix);
|
ipix = BigInt(ipix);
|
||||||
let pix = ipix & (this.npface - 1n);
|
let pix = ipix & (this.npface - 1n);
|
||||||
@@ -364,7 +434,7 @@ export class Healpix {
|
|||||||
pix2loc(pix) {
|
pix2loc(pix) {
|
||||||
let loc = new Hploc(undefined);
|
let loc = new Hploc(undefined);
|
||||||
let xyf = this.nest2xyf(pix);
|
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;
|
let nr;
|
||||||
if (jr < this.nside) {
|
if (jr < this.nside) {
|
||||||
nr = jr;
|
nr = jr;
|
||||||
@@ -437,17 +507,24 @@ export class Healpix {
|
|||||||
;
|
;
|
||||||
compress_bits(v) {
|
compress_bits(v) {
|
||||||
v = BigInt(v);
|
v = BigInt(v);
|
||||||
let raw = (v & 0x5555n) | ((v & 0x55550000n) >> 15n);
|
let compressed = 0n;
|
||||||
let compressed = BigInt(this.ctab[Number(raw & 0xffn)]) | (BigInt(this.ctab[Number(raw >> 8n)]) << 4n);
|
for (let i = 0n; i < 32n; i++) {
|
||||||
|
if ((v & (1n << (2n * i))) !== 0n) {
|
||||||
|
compressed |= (1n << i);
|
||||||
|
}
|
||||||
|
}
|
||||||
return compressed;
|
return compressed;
|
||||||
}
|
}
|
||||||
;
|
;
|
||||||
spread_bits(v) {
|
spread_bits(v) {
|
||||||
v = BigInt(v);
|
v = BigInt(v);
|
||||||
return BigInt(this.utab[Number(v & 0xffn)])
|
let spread = 0n;
|
||||||
| (BigInt(this.utab[Number((v >> 8n) & 0xffn)]) << 16n)
|
for (let i = 0n; i < 32n; i++) {
|
||||||
| (BigInt(this.utab[Number((v >> 16n) & 0xffn)]) << 32n)
|
if ((v & (1n << i)) !== 0n) {
|
||||||
| (BigInt(this.utab[Number((v >> 24n) & 0xffn)]) << 48n);
|
spread |= (1n << (2n * i));
|
||||||
|
}
|
||||||
|
}
|
||||||
|
return spread;
|
||||||
}
|
}
|
||||||
;
|
;
|
||||||
/**
|
/**
|
||||||
@@ -777,5 +854,72 @@ export class Healpix {
|
|||||||
}
|
}
|
||||||
return pixset;
|
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
|
//# sourceMappingURL=Healpix.js.map
|
||||||
Reference in New Issue
Block a user