Skip to content

Instantly share code, notes, and snippets.

@sto3psl

sto3psl/after-esbuild.js

Last active March 22, 2021 10:53
Show Gist options
  • Save sto3psl/8bca4562096050559ec51544483edc33 to your computer and use it in GitHub Desktop.
Save sto3psl/8bca4562096050559ec51544483edc33 to your computer and use it in GitHub Desktop.
Download ZIP
esbuild renaming problem
Raw
after-esbuild.js
'use strict';
function AVSurface(coordList, radiusList, indexList) {
var nAtoms = radiusList.length;
var x = new Float32Array(nAtoms);
var y = new Float32Array(nAtoms);
var z = new Float32Array(nAtoms);
for (var i = 0; i < nAtoms; i++) {
var ci = 3 * i;
x[i] = coordList[ci];
y[i] = coordList[ci + 1];
z[i] = coordList[ci + 2];
}
var bbox = computeBoundingBox(coordList);
if (coordList.length === 0) {
bbox[0].set([0, 0, 0]);
bbox[1].set([0, 0, 0]);
}
var min = bbox[0];
var max = bbox[1];
var r, r2;
var maxRadius2;
var probeRadius, scaleFactor, setAtomID, probePositions;
var lastClip = -1;
var dim, matrix2, grid, atomIndex;
var gridx, gridy, gridz;
var sinTable, cosTable;
var hash;
var neighbours;
var atob2 = new Float32Array([0, 0, 0]);
var mid = new Float32Array([0, 0, 0]);
var n1 = new Float32Array([0, 0, 0]);
var n2 = new Float32Array([0, 0, 0]);
var ngTorus;
function init3(_probeRadius, _scaleFactor, _setAtomID, _probePositions) {
probeRadius = defaults(_probeRadius, 1.4);
scaleFactor = defaults(_scaleFactor, 2);
setAtomID = defaults(_setAtomID, true);
probePositions = defaults(_probePositions, 30);
r = new Float32Array(nAtoms);
r2 = new Float32Array(nAtoms);
for (var i2 = 0; i2 < r.length; ++i2) {
var rExt = radiusList[i2] + probeRadius;
r[i2] = rExt;
r2[i2] = rExt * rExt;
}
maxRadius2 = 0;
for (var j = 0; j < r.length; ++j) {
if (r[j] > maxRadius2) {
maxRadius2 = r[j];
}
}
initializeGrid();
initializeAngleTables();
initializeHash();
lastClip = -1;
}
__name(init3, "init");
function fillGridDim(a, start, step) {
for (var i2 = 0; i2 < a.length; i2++) {
a[i2] = start + step * i2;
}
}
__name(fillGridDim, "fillGridDim");
function initializeGrid() {
var surfGrid = getSurfaceGrid(min, max, maxRadius2, scaleFactor, 0);
scaleFactor = surfGrid.scaleFactor;
dim = surfGrid.dim;
matrix2 = surfGrid.matrix;
ngTorus = Math.max(5, 2 + Math.floor(probeRadius * scaleFactor));
grid = uniformArray(dim[0] * dim[1] * dim[2], -1001);
atomIndex = new Int32Array(grid.length);
gridx = new Float32Array(dim[0]);
gridy = new Float32Array(dim[1]);
gridz = new Float32Array(dim[2]);
fillGridDim(gridx, min[0], 1 / scaleFactor);
fillGridDim(gridy, min[1], 1 / scaleFactor);
fillGridDim(gridz, min[2], 1 / scaleFactor);
}
__name(initializeGrid, "initializeGrid");
function initializeAngleTables() {
var theta = 0;
var step = 2 * Math.PI / probePositions;
cosTable = new Float32Array(probePositions);
sinTable = new Float32Array(probePositions);
for (var i2 = 0; i2 < probePositions; i2++) {
cosTable[i2] = Math.cos(theta);
sinTable[i2] = Math.sin(theta);
theta += step;
}
}
__name(initializeAngleTables, "initializeAngleTables");
function initializeHash() {
hash = makeAVHash(x, y, z, r, min, max, 2.01 * maxRadius2);
neighbours = new Int32Array(hash.neighbourListLength);
}
__name(initializeHash, "initializeHash");
function obscured(x2, y2, z2, a, b) {
var ai;
if (lastClip !== -1) {
ai = lastClip;
if (ai !== a && ai !== b && singleAtomObscures(ai, x2, y2, z2)) {
return ai;
} else {
lastClip = -1;
}
}
var ni = 0;
ai = neighbours[ni];
while (ai >= 0) {
if (ai !== a && ai !== b && singleAtomObscures(ai, x2, y2, z2)) {
lastClip = ai;
return ai;
}
ai = neighbours[++ni];
}
lastClip = -1;
return -1;
}
__name(obscured, "obscured");
function singleAtomObscures(ai, x2, y2, z2) {
var ci2 = 3 * ai;
var ra2 = r2[ai];
var dx = coordList[ci2] - x2;
var dy = coordList[ci2 + 1] - y2;
var dz = coordList[ci2 + 2] - z2;
var d2 = dx * dx + dy * dy + dz * dz;
return d2 < ra2;
}
__name(singleAtomObscures, "singleAtomObscures");
function projectPoints() {
for (var i2 = 0; i2 < nAtoms; i2++) {
var ax = x[i2];
var ay = y[i2];
var az = z[i2];
var ar = r[i2];
var ar2 = r2[i2];
hash.withinRadii(ax, ay, az, ar, neighbours);
var ng = Math.ceil(ar * scaleFactor);
var iax = Math.floor(scaleFactor * (ax - min[0]));
var iay = Math.floor(scaleFactor * (ay - min[1]));
var iaz = Math.floor(scaleFactor * (az - min[2]));
var minx = Math.max(0, iax - ng);
var miny = Math.max(0, iay - ng);
var minz = Math.max(0, iaz - ng);
var maxx = Math.min(dim[0], iax + ng + 2);
var maxy = Math.min(dim[1], iay + ng + 2);
var maxz = Math.min(dim[2], iaz + ng + 2);
for (var ix = minx; ix < maxx; ix++) {
var dx = gridx[ix] - ax;
var xoffset = dim[1] * dim[2] * ix;
for (var iy = miny; iy < maxy; iy++) {
var dy = gridy[iy] - ay;
var dxy2 = dx * dx + dy * dy;
var xyoffset = xoffset + dim[2] * iy;
for (var iz = minz; iz < maxz; iz++) {
var dz = gridz[iz] - az;
var d2 = dxy2 + dz * dz;
if (d2 < ar2) {
var idx = iz + xyoffset;
if (grid[idx] < 0) {
grid[idx] = -grid[idx];
}
var d = Math.sqrt(d2);
var ap = ar / d;
var spx = dx * ap;
var spy = dy * ap;
var spz = dz * ap;
spx += ax;
spy += ay;
spz += az;
if (obscured(spx, spy, spz, i2, -1) === -1) {
var dd = ar - d;
if (dd < grid[idx]) {
grid[idx] = dd;
if (setAtomID) {
atomIndex[idx] = i2;
}
}
}
}
}
}
}
}
}
__name(projectPoints, "projectPoints");
function projectTorii() {
for (var i2 = 0; i2 < nAtoms; i2++) {
hash.withinRadii(x[i2], y[i2], z[i2], r[i2], neighbours);
var ia = 0;
var ni = neighbours[ia];
while (ni >= 0) {
if (i2 < ni) {
projectTorus(i2, ni);
}
ni = neighbours[++ia];
}
}
}
__name(projectTorii, "projectTorii");
function projectTorus(a, b) {
var r1 = r[a];
var r22 = r[b];
var dx = atob2[0] = x[b] - x[a];
var dy = atob2[1] = y[b] - y[a];
var dz = atob2[2] = z[b] - z[a];
var d2 = dx * dx + dy * dy + dz * dz;
var d = Math.sqrt(d2);
var cosA = (r1 * r1 + d * d - r22 * r22) / (2 * r1 * d);
var dmp = r1 * cosA;
v3normalize(atob2, atob2);
normalToLine(n1, atob2);
v3normalize(n1, n1);
v3cross(n2, atob2, n1);
v3normalize(n2, n2);
var rInt = Math.sqrt(r1 * r1 - dmp * dmp);
v3multiplyScalar(n1, n1, rInt);
v3multiplyScalar(n2, n2, rInt);
v3multiplyScalar(atob2, atob2, dmp);
mid[0] = atob2[0] + x[a];
mid[1] = atob2[1] + y[a];
mid[2] = atob2[2] + z[a];
lastClip = -1;
var ng = ngTorus;
for (var i2 = 0; i2 < probePositions; i2++) {
var cost = cosTable[i2];
var sint = sinTable[i2];
var px2 = mid[0] + cost * n1[0] + sint * n2[0];
var py2 = mid[1] + cost * n1[1] + sint * n2[1];
var pz2 = mid[2] + cost * n1[2] + sint * n2[2];
if (obscured(px2, py2, pz2, a, b) === -1) {
var iax = Math.floor(scaleFactor * (px2 - min[0]));
var iay = Math.floor(scaleFactor * (py2 - min[1]));
var iaz = Math.floor(scaleFactor * (pz2 - min[2]));
var minx = Math.max(0, iax - ng);
var miny = Math.max(0, iay - ng);
var minz = Math.max(0, iaz - ng);
var maxx = Math.min(dim[0], iax + ng + 2);
var maxy = Math.min(dim[1], iay + ng + 2);
var maxz = Math.min(dim[2], iaz + ng + 2);
for (var ix = minx; ix < maxx; ix++) {
dx = px2 - gridx[ix];
var xoffset = dim[1] * dim[2] * ix;
for (var iy = miny; iy < maxy; iy++) {
dy = py2 - gridy[iy];
var dxy2 = dx * dx + dy * dy;
var xyoffset = xoffset + dim[2] * iy;
for (var iz = minz; iz < maxz; iz++) {
dz = pz2 - gridz[iz];
d2 = dxy2 + dz * dz;
var idx = iz + xyoffset;
var current = grid[idx];
if (current > 0 && d2 < current * current) {
grid[idx] = Math.sqrt(d2);
if (setAtomID) {
var dp = dx * atob2[0] + dy * atob2[1] + dz * atob2[2];
atomIndex[idx] = dp < 0 ? b : a;
}
}
}
}
}
}
}
}
__name(projectTorus, "projectTorus");
function normalToLine(out, p) {
out[0] = out[1] = out[2] = 1;
if (p[0] !== 0) {
out[0] = (p[1] + p[2]) / -p[0];
} else if (p[1] !== 0) {
out[1] = (p[0] + p[2]) / -p[1];
} else if (p[2] !== 0) {
out[2] = (p[0] + p[1]) / -p[2];
}
return out;
}
__name(normalToLine, "normalToLine");
function fixNegatives() {
for (var i2 = 0; i2 < grid.length; i2++) {
if (grid[i2] < 0) {
grid[i2] = 0;
}
}
}
__name(fixNegatives, "fixNegatives");
function fixAtomIDs() {
for (var i2 = 0; i2 < atomIndex.length; i2++) {
atomIndex[i2] = indexList[atomIndex[i2]];
}
}
__name(fixAtomIDs, "fixAtomIDs");
function getVolume(probeRadius2, scaleFactor2, setAtomID2) {
console.time("AVSurface.getVolume");
console.time("AVSurface.init");
init3(probeRadius2, scaleFactor2, setAtomID2);
console.timeEnd("AVSurface.init");
console.time("AVSurface.projectPoints");
projectPoints();
console.timeEnd("AVSurface.projectPoints");
console.time("AVSurface.projectTorii");
projectTorii();
console.timeEnd("AVSurface.projectTorii");
fixNegatives();
fixAtomIDs();
console.timeEnd("AVSurface.getVolume");
}
__name(getVolume, "getVolume");
this.getSurface = function(type, probeRadius2, scaleFactor2, cutoff, setAtomID2, smooth, contour) {
getVolume(probeRadius2, scaleFactor2, setAtomID2);
var volsurf = new VolumeSurface(grid, dim[2], dim[1], dim[0], atomIndex);
return volsurf.getSurface(probeRadius2, false, void 0, matrix2, contour);
};
}
function EDTSurface(coordList, radiusList, indexList) {
var radiusDict = getRadiusDict(radiusList);
var bbox = computeBoundingBox(coordList);
if (coordList.length === 0) {
bbox[0].set([0, 0, 0]);
bbox[1].set([0, 0, 0]);
}
var min = bbox[0];
var max = bbox[1];
var probeRadius, scaleFactor, cutoff;
var pLength, pWidth, pHeight;
var matrix2, ptran;
var depty, widxz;
var cutRadius;
var setAtomID;
var vpBits, vpDistance, vpAtomID;
function init3(btype, _probeRadius, _scaleFactor, _cutoff, _setAtomID) {
probeRadius = _probeRadius || 1.4;
scaleFactor = _scaleFactor || 2;
setAtomID = _setAtomID || true;
var maxRadius2 = 0;
for (var radius in radiusDict) {
maxRadius2 = Math.max(maxRadius2, radius);
}
var grid = getSurfaceGrid(min, max, maxRadius2, scaleFactor, btype ? probeRadius : 0);
pLength = grid.dim[0];
pWidth = grid.dim[1];
pHeight = grid.dim[2];
matrix2 = grid.matrix;
ptran = grid.tran;
scaleFactor = grid.scaleFactor;
depty = {};
widxz = {};
boundingatom(btype);
cutRadius = probeRadius * scaleFactor;
if (_cutoff) {
cutoff = _cutoff;
} else {
cutoff = probeRadius / scaleFactor;
}
vpBits = new Uint8Array(pLength * pWidth * pHeight);
if (btype) {
vpDistance = new Float64Array(pLength * pWidth * pHeight);
}
if (setAtomID) {
vpAtomID = new Int32Array(pLength * pWidth * pHeight);
}
}
__name(init3, "init");
var INOUT = 1;
var ISDONE = 2;
var ISBOUND = 4;
var nb = [
new Int32Array([1, 0, 0]),
new Int32Array([-1, 0, 0]),
new Int32Array([0, 1, 0]),
new Int32Array([0, -1, 0]),
new Int32Array([0, 0, 1]),
new Int32Array([0, 0, -1]),
new Int32Array([1, 1, 0]),
new Int32Array([1, -1, 0]),
new Int32Array([-1, 1, 0]),
new Int32Array([-1, -1, 0]),
new Int32Array([1, 0, 1]),
new Int32Array([1, 0, -1]),
new Int32Array([-1, 0, 1]),
new Int32Array([-1, 0, -1]),
new Int32Array([0, 1, 1]),
new Int32Array([0, 1, -1]),
new Int32Array([0, -1, 1]),
new Int32Array([0, -1, -1]),
new Int32Array([1, 1, 1]),
new Int32Array([1, 1, -1]),
new Int32Array([1, -1, 1]),
new Int32Array([-1, 1, 1]),
new Int32Array([1, -1, -1]),
new Int32Array([-1, -1, 1]),
new Int32Array([-1, 1, -1]),
new Int32Array([-1, -1, -1])
];
this.getVolume = function(type, probeRadius2, scaleFactor2, cutoff2, setAtomID2) {
console.time("EDTSurface.getVolume");
var btype = type !== "vws";
init3(btype, probeRadius2, scaleFactor2, cutoff2, setAtomID2);
fillvoxels(btype);
buildboundary();
if (type === "ms" || type === "ses") {
fastdistancemap();
}
if (type === "ses") {
boundingatom(false);
fillvoxelswaals();
}
marchingcubeinit(type);
for (var i = 0, il = vpAtomID.length; i < il; ++i) {
vpAtomID[i] = indexList[vpAtomID[i]];
}
console.timeEnd("EDTSurface.getVolume");
return {
data: vpBits,
nx: pHeight,
ny: pWidth,
nz: pLength,
atomindex: vpAtomID
};
};
this.getSurface = function(type, probeRadius2, scaleFactor2, cutoff2, setAtomID2, smooth, contour) {
var vd = this.getVolume(type, probeRadius2, scaleFactor2, cutoff2, setAtomID2);
var volsurf = new VolumeSurface(vd.data, vd.nx, vd.ny, vd.nz, vd.atomindex);
return volsurf.getSurface(1, smooth, void 0, matrix2, contour);
};
function boundingatom(btype) {
var r;
var j;
var k;
var txz;
var tdept;
var sradius;
var tradius;
var widxzR;
var deptyName;
var indx;
for (var name in radiusDict) {
r = parseFloat(name);
if (depty[name]) {
continue;
}
if (!btype) {
tradius = r * scaleFactor + 0.5;
} else {
tradius = (r + probeRadius) * scaleFactor + 0.5;
}
sradius = tradius * tradius;
widxzR = Math.floor(tradius) + 1;
deptyName = new Int32Array(widxzR * widxzR);
indx = 0;
for (j = 0; j < widxzR; ++j) {
for (k = 0; k < widxzR; ++k) {
txz = j * j + k * k;
if (txz > sradius) {
deptyName[indx] = -1;
} else {
tdept = Math.sqrt(sradius - txz);
deptyName[indx] = Math.floor(tdept);
}
++indx;
}
}
widxz[name] = widxzR;
depty[name] = deptyName;
}
}
__name(boundingatom, "boundingatom");
function fillatom(idx) {
var ci = idx * 3;
var ri = idx;
var cx, cy, cz, ox, oy, oz, mi, mj, mk, i, j, k, si, sj, sk;
var ii, jj, kk;
cx = Math.floor(0.5 + scaleFactor * (coordList[ci] + ptran[0]));
cy = Math.floor(0.5 + scaleFactor * (coordList[ci + 1] + ptran[1]));
cz = Math.floor(0.5 + scaleFactor * (coordList[ci + 2] + ptran[2]));
var at = radiusList[ri];
var deptyAt = depty[at];
var nind = 0;
var pWH = pWidth * pHeight;
var n = widxz[at];
var deptyAtNind;
for (i = 0; i < n; ++i) {
for (j = 0; j < n; ++j) {
deptyAtNind = deptyAt[nind];
if (deptyAtNind !== -1) {
for (ii = -1; ii < 2; ++ii) {
for (jj = -1; jj < 2; ++jj) {
for (kk = -1; kk < 2; ++kk) {
if (ii !== 0 && jj !== 0 && kk !== 0) {
mi = ii * i;
mk = kk * j;
for (k = 0; k <= deptyAtNind; ++k) {
mj = k * jj;
si = cx + mi;
sj = cy + mj;
sk = cz + mk;
if (si < 0 || sj < 0 || sk < 0 || si >= pLength || sj >= pWidth || sk >= pHeight) {
continue;
}
var index = si * pWH + sj * pHeight + sk;
if (!setAtomID) {
vpBits[index] |= INOUT;
} else {
if (!(vpBits[index] & INOUT)) {
vpBits[index] |= INOUT;
vpAtomID[index] = idx;
} else if (vpBits[index] & INOUT) {
var ci2 = vpAtomID[index];
if (ci2 !== ci) {
ox = cx + mi - Math.floor(0.5 + scaleFactor * (coordList[ci2] + ptran[0]));
oy = cy + mj - Math.floor(0.5 + scaleFactor * (coordList[ci2 + 1] + ptran[1]));
oz = cz + mk - Math.floor(0.5 + scaleFactor * (coordList[ci2 + 2] + ptran[2]));
if (mi * mi + mj * mj + mk * mk < ox * ox + oy * oy + oz * oz) {
vpAtomID[index] = idx;
}
}
}
}
}
}
}
}
}
}
nind++;
}
}
}
__name(fillatom, "fillatom");
function fillvoxels(btype) {
console.time("EDTSurface fillvoxels");
var i, il;
for (i = 0, il = vpBits.length; i < il; ++i) {
vpBits[i] = 0;
if (btype) {
vpDistance[i] = -1;
}
if (setAtomID) {
vpAtomID[i] = -1;
}
}
for (i = 0, il = coordList.length / 3; i < il; ++i) {
fillatom(i);
}
for (i = 0, il = vpBits.length; i < il; ++i) {
if (vpBits[i] & INOUT) {
vpBits[i] |= ISDONE;
}
}
console.timeEnd("EDTSurface fillvoxels");
}
__name(fillvoxels, "fillvoxels");
function fillAtomWaals(idx) {
var ci = idx * 3;
var ri = idx;
var cx;
var cy;
var cz;
var ox;
var oy;
var oz;
var nind = 0;
var mi;
var mj;
var mk;
var si;
var sj;
var sk;
var i;
var j;
var k;
var ii;
var jj;
var kk;
var n;
cx = Math.floor(0.5 + scaleFactor * (coordList[ci] + ptran[0]));
cy = Math.floor(0.5 + scaleFactor * (coordList[ci + 1] + ptran[1]));
cz = Math.floor(0.5 + scaleFactor * (coordList[ci + 2] + ptran[2]));
var at = radiusList[ri];
var pWH = pWidth * pHeight;
for (i = 0, n = widxz[at]; i < n; ++i) {
for (j = 0; j < n; ++j) {
if (depty[at][nind] !== -1) {
for (ii = -1; ii < 2; ++ii) {
for (jj = -1; jj < 2; ++jj) {
for (kk = -1; kk < 2; ++kk) {
if (ii !== 0 && jj !== 0 && kk !== 0) {
mi = ii * i;
mk = kk * j;
for (k = 0; k <= depty[at][nind]; ++k) {
mj = k * jj;
si = cx + mi;
sj = cy + mj;
sk = cz + mk;
if (si < 0 || sj < 0 || sk < 0 || si >= pLength || sj >= pWidth || sk >= pHeight) {
continue;
}
var index = si * pWH + sj * pHeight + sk;
if (!(vpBits[index] & ISDONE)) {
vpBits[index] |= ISDONE;
if (setAtomID) {
vpAtomID[index] = idx;
}
} else if (setAtomID) {
var ci2 = vpAtomID[index];
ox = Math.floor(0.5 + scaleFactor * (coordList[ci2] + ptran[0]));
oy = Math.floor(0.5 + scaleFactor * (coordList[ci2 + 1] + ptran[1]));
oz = Math.floor(0.5 + scaleFactor * (coordList[ci2 + 2] + ptran[2]));
if (mi * mi + mj * mj + mk * mk < ox * ox + oy * oy + oz * oz) {
vpAtomID[index] = idx;
}
}
}
}
}
}
}
}
nind++;
}
}
}
__name(fillAtomWaals, "fillAtomWaals");
function fillvoxelswaals() {
var i, il;
for (i = 0, il = vpBits.length; i < il; ++i) {
vpBits[i] &= ~ISDONE;
}
for (i = 0, il = coordList.length / 3; i < il; ++i) {
fillAtomWaals(i);
}
}
__name(fillvoxelswaals, "fillvoxelswaals");
function buildboundary() {
var i, j, k;
var pWH = pWidth * pHeight;
for (i = 0; i < pLength; ++i) {
for (j = 0; j < pHeight; ++j) {
for (k = 0; k < pWidth; ++k) {
var index = i * pWH + k * pHeight + j;
if (vpBits[index] & INOUT) {
var ii = 0;
while (ii < 26) {
var ti = i + nb[ii][0];
var tj = j + nb[ii][2];
var tk = k + nb[ii][1];
if (ti > -1 && ti < pLength && tk > -1 && tk < pWidth && tj > -1 && tj < pHeight && !(vpBits[ti * pWH + tk * pHeight + tj] & INOUT)) {
vpBits[index] |= ISBOUND;
break;
} else {
ii++;
}
}
}
}
}
}
}
__name(buildboundary, "buildboundary");
function fastdistancemap() {
console.time("EDTSurface fastdistancemap");
var i, j, k, n;
var boundPoint = makeGrid(pLength, pWidth, pHeight, Uint16Array, 3);
var pWH = pWidth * pHeight;
var cutRSq = cutRadius * cutRadius;
var totalsurfacevox = 0;
var index;
for (i = 0; i < pLength; ++i) {
for (j = 0; j < pWidth; ++j) {
for (k = 0; k < pHeight; ++k) {
index = i * pWH + j * pHeight + k;
vpBits[index] &= ~ISDONE;
if (vpBits[index] & INOUT) {
if (vpBits[index] & ISBOUND) {
boundPoint.set(i, j, k, i, j, k);
vpDistance[index] = 0;
vpBits[index] |= ISDONE;
totalsurfacevox += 1;
}
}
}
}
}
var inarray = new Int32Array(3 * totalsurfacevox);
var positin = 0;
var outarray = new Int32Array(3 * totalsurfacevox);
var positout = 0;
for (i = 0; i < pLength; ++i) {
for (j = 0; j < pWidth; ++j) {
for (k = 0; k < pHeight; ++k) {
index = i * pWH + j * pHeight + k;
if (vpBits[index] & ISBOUND) {
inarray[positin] = i;
inarray[positin + 1] = j;
inarray[positin + 2] = k;
positin += 3;
vpBits[index] &= ~ISBOUND;
}
}
}
}
do {
positout = fastoneshell(inarray, boundPoint, positin, outarray);
positin = 0;
for (i = 0, n = positout; i < n; i += 3) {
index = pWH * outarray[i] + pHeight * outarray[i + 1] + outarray[i + 2];
vpBits[index] &= ~ISBOUND;
if (vpDistance[index] <= 1.0404 * cutRSq) {
inarray[positin] = outarray[i];
inarray[positin + 1] = outarray[i + 1];
inarray[positin + 2] = outarray[i + 2];
positin += 3;
}
}
} while (positin > 0);
var cutoffSq = cutoff * cutoff;
var index2;
var bp = new Uint16Array(3);
for (i = 0; i < pLength; ++i) {
for (j = 0; j < pWidth; ++j) {
for (k = 0; k < pHeight; ++k) {
index = i * pWH + j * pHeight + k;
vpBits[index] &= ~ISBOUND;
if (vpBits[index] & INOUT) {
if (!(vpBits[index] & ISDONE) || vpBits[index] & ISDONE && vpDistance[index] >= cutoffSq) {
vpBits[index] |= ISBOUND;
if (setAtomID && vpBits[index] & ISDONE) {
boundPoint.toArray(i, j, k, bp);
index2 = bp[0] * pWH + bp[1] * pHeight + bp[2];
vpAtomID[index] = vpAtomID[index2];
}
}
}
}
}
}
console.timeEnd("EDTSurface fastdistancemap");
}
__name(fastdistancemap, "fastdistancemap");
function fastoneshell(inarray, boundPoint, positin, outarray) {
var tx, ty, tz;
var dx, dy, dz;
var i, j, n;
var square;
var index;
var nbj;
var bp = new Uint16Array(3);
var positout = 0;
if (positin === 0) {
return positout;
}
var tnvix = -1;
var tnviy = -1;
var tnviz = -1;
var pWH = pWidth * pHeight;
for (i = 0, n = positin; i < n; i += 3) {
tx = inarray[i];
ty = inarray[i + 1];
tz = inarray[i + 2];
boundPoint.toArray(tx, ty, tz, bp);
for (j = 0; j < 6; ++j) {
nbj = nb[j];
tnvix = tx + nbj[0];
tnviy = ty + nbj[1];
tnviz = tz + nbj[2];
if (tnvix < pLength && tnvix > -1 && tnviy < pWidth && tnviy > -1 && tnviz < pHeight && tnviz > -1) {
index = tnvix * pWH + pHeight * tnviy + tnviz;
if (vpBits[index] & INOUT && !(vpBits[index] & ISDONE)) {
boundPoint.fromArray(tnvix, tnviy, tnviz, bp);
dx = tnvix - bp[0];
dy = tnviy - bp[1];
dz = tnviz - bp[2];
square = dx * dx + dy * dy + dz * dz;
vpDistance[index] = square;
vpBits[index] |= ISDONE;
vpBits[index] |= ISBOUND;
outarray[positout] = tnvix;
outarray[positout + 1] = tnviy;
outarray[positout + 2] = tnviz;
positout += 3;
} else if (vpBits[index] & INOUT && vpBits[index] & ISDONE) {
dx = tnvix - bp[0];
dy = tnviy - bp[1];
dz = tnviz - bp[2];
square = dx * dx + dy * dy + dz * dz;
if (square < vpDistance[index]) {
boundPoint.fromArray(tnvix, tnviy, tnviz, bp);
vpDistance[index] = square;
if (!(vpBits[index] & ISBOUND)) {
vpBits[index] |= ISBOUND;
outarray[positout] = tnvix;
outarray[positout + 1] = tnviy;
outarray[positout + 2] = tnviz;
positout += 3;
}
}
}
}
}
}
for (i = 0, n = positin; i < n; i += 3) {
tx = inarray[i];
ty = inarray[i + 1];
tz = inarray[i + 2];
boundPoint.toArray(tx, ty, tz, bp);
for (j = 6; j < 18; j++) {
nbj = nb[j];
tnvix = tx + nbj[0];
tnviy = ty + nbj[1];
tnviz = tz + nbj[2];
if (tnvix < pLength && tnvix > -1 && tnviy < pWidth && tnviy > -1 && tnviz < pHeight && tnviz > -1) {
index = tnvix * pWH + pHeight * tnviy + tnviz;
if (vpBits[index] & INOUT && !(vpBits[index] & ISDONE)) {
boundPoint.fromArray(tnvix, tnviy, tnviz, bp);
dx = tnvix - bp[0];
dy = tnviy - bp[1];
dz = tnviz - bp[2];
square = dx * dx + dy * dy + dz * dz;
vpDistance[index] = square;
vpBits[index] |= ISDONE;
vpBits[index] |= ISBOUND;
outarray[positout] = tnvix;
outarray[positout + 1] = tnviy;
outarray[positout + 2] = tnviz;
positout += 3;
} else if (vpBits[index] & INOUT && vpBits[index] & ISDONE) {
dx = tnvix - bp[0];
dy = tnviy - bp[1];
dz = tnviz - bp[2];
square = dx * dx + dy * dy + dz * dz;
if (square < vpDistance[index]) {
boundPoint.fromArray(tnvix, tnviy, tnviz, bp);
vpDistance[index] = square;
if (!(vpBits[index] & ISBOUND)) {
vpBits[index] |= ISBOUND;
outarray[positout] = tnvix;
outarray[positout + 1] = tnviy;
outarray[positout + 2] = tnviz;
positout += 3;
}
}
}
}
}
}
for (i = 0, n = positin; i < n; i += 3) {
tx = inarray[i];
ty = inarray[i + 1];
tz = inarray[i + 2];
boundPoint.toArray(tx, ty, tz, bp);
for (j = 18; j < 26; j++) {
nbj = nb[j];
tnvix = tx + nbj[0];
tnviy = ty + nbj[1];
tnviz = tz + nbj[2];
if (tnvix < pLength && tnvix > -1 && tnviy < pWidth && tnviy > -1 && tnviz < pHeight && tnviz > -1) {
index = tnvix * pWH + pHeight * tnviy + tnviz;
if (vpBits[index] & INOUT && !(vpBits[index] & ISDONE)) {
boundPoint.fromArray(tnvix, tnviy, tnviz, bp);
dx = tnvix - bp[0];
dy = tnviy - bp[1];
dz = tnviz - bp[2];
square = dx * dx + dy * dy + dz * dz;
vpDistance[index] = square;
vpBits[index] |= ISDONE;
vpBits[index] |= ISBOUND;
outarray[positout] = tnvix;
outarray[positout + 1] = tnviy;
outarray[positout + 2] = tnviz;
positout += 3;
} else if (vpBits[index] & INOUT && vpBits[index] & ISDONE) {
dx = tnvix - bp[0];
dy = tnviy - bp[1];
dz = tnviz - bp[2];
square = dx * dx + dy * dy + dz * dz;
if (square < vpDistance[index]) {
boundPoint.fromArray(tnvix, tnviy, tnviz, bp);
vpDistance[index] = square;
if (!(vpBits[index] & ISBOUND)) {
vpBits[index] |= ISBOUND;
outarray[positout] = tnvix;
outarray[positout + 1] = tnviy;
outarray[positout + 2] = tnviz;
positout += 3;
}
}
}
}
}
}
return positout;
}
__name(fastoneshell, "fastoneshell");
function marchingcubeinit(stype) {
var i;
var n = vpBits.length;
if (stype === "vws") {
for (i = 0; i < n; ++i) {
vpBits[i] &= ~ISBOUND;
vpBits[i] = vpBits[i] & ISDONE ? 1 : 0;
}
} else if (stype === "ms") {
for (i = 0; i < n; ++i) {
vpBits[i] &= ~ISDONE;
if (vpBits[i] & ISBOUND) {
vpBits[i] |= ISDONE;
}
vpBits[i] &= ~ISBOUND;
vpBits[i] = vpBits[i] & ISDONE ? 1 : 0;
}
} else if (stype === "ses") {
for (i = 0; i < n; ++i) {
if (vpBits[i] & ISBOUND && vpBits[i] & ISDONE) {
vpBits[i] &= ~ISBOUND;
} else if (vpBits[i] & ISBOUND && !(vpBits[i] & ISDONE)) {
vpBits[i] |= ISDONE;
}
vpBits[i] = vpBits[i] & ISDONE ? 1 : 0;
}
} else if (stype === "sas") {
for (i = 0; i < n; ++i) {
vpBits[i] &= ~ISBOUND;
vpBits[i] = vpBits[i] & ISDONE ? 1 : 0;
}
}
}
__name(marchingcubeinit, "marchingcubeinit");
}
function MarchingCubes(field, nx, ny, nz, atomindex) {
var isolevel = 0;
var noNormals = false;
var contour = false;
var wrap = false;
var isNegativeIso = false;
var normalFactor = -1;
var n = nx * ny * nz;
var yd = nx;
var zd = nx * ny;
var normalCache, vertexIndex;
var count, icount;
var ilist = new Int32Array(12);
var positionArray = [];
var normalArray = [];
var indexArray = [];
var atomindexArray = [];
var edgeTable = getEdgeTable();
var triTable = getTriTable();
var allowedContours = getAllowedContours();
var mx, my, mz;
this.triangulate = function(_isolevel, _noNormals, _box2, _contour, _wrap) {
isolevel = _isolevel;
isNegativeIso = isolevel < 0;
contour = _contour;
wrap = _wrap;
noNormals = _noNormals || contour;
if (!noNormals) {
normalFactor = isolevel > 0 ? -1 : 1;
if (!normalCache) {
normalCache = new Float32Array(n * 3);
}
}
var vIndexLength = n * 3;
if (!vertexIndex || vertexIndex.length !== vIndexLength) {
vertexIndex = new Int32Array(vIndexLength);
}
count = 0;
icount = 0;
if (_box2 !== void 0) {
var min = _box2[0].map(Math.round);
var max = _box2[1].map(Math.round);
mx = nx * Math.ceil(Math.abs(min[0]) / nx);
my = ny * Math.ceil(Math.abs(min[1]) / ny);
mz = nz * Math.ceil(Math.abs(min[2]) / nz);
triangulate(min[0], min[1], min[2], max[0], max[1], max[2]);
} else {
mx = my = mz = 0;
triangulate();
}
positionArray.length = count * 3;
if (!noNormals) {
normalArray.length = count * 3;
}
indexArray.length = icount;
if (atomindex) {
atomindexArray.length = count;
}
return {
position: new Float32Array(positionArray),
normal: noNormals ? void 0 : new Float32Array(normalArray),
index: getUintArray(indexArray, positionArray.length / 3),
atomindex: atomindex ? new Int32Array(atomindexArray) : void 0,
contour
};
};
function lerp2(a, b, t) {
return a + (b - a) * t;
}
__name(lerp2, "lerp");
function index(x, y, z) {
x = (x + mx) % nx;
y = (y + my) % ny;
z = (z + mz) % nz;
return zd * z + yd * y + x;
}
__name(index, "index");
function VIntX(q, offset, x, y, z, valp1, valp2) {
var _q = 3 * q;
if (vertexIndex[_q] < 0) {
var mu = (isolevel - valp1) / (valp2 - valp1);
var nc = normalCache;
var c = count * 3;
positionArray[c + 0] = x + mu;
positionArray[c + 1] = y;
positionArray[c + 2] = z;
if (!noNormals) {
var q3 = q * 3;
normalArray[c] = normalFactor * lerp2(nc[q3], nc[q3 + 3], mu);
normalArray[c + 1] = normalFactor * lerp2(nc[q3 + 1], nc[q3 + 4], mu);
normalArray[c + 2] = normalFactor * lerp2(nc[q3 + 2], nc[q3 + 5], mu);
}
if (atomindex) {
atomindexArray[count] = atomindex[q + Math.round(mu)];
}
vertexIndex[_q] = count;
ilist[offset] = count;
count += 1;
} else {
ilist[offset] = vertexIndex[_q];
}
}
__name(VIntX, "VIntX");
function VIntY(q, offset, x, y, z, valp1, valp2) {
var _q = 3 * q + 1;
if (vertexIndex[_q] < 0) {
var mu = (isolevel - valp1) / (valp2 - valp1);
var nc = normalCache;
var c = count * 3;
positionArray[c] = x;
positionArray[c + 1] = y + mu;
positionArray[c + 2] = z;
if (!noNormals) {
var q3 = q * 3;
var q6 = q3 + yd * 3;
normalArray[c] = normalFactor * lerp2(nc[q3], nc[q6], mu);
normalArray[c + 1] = normalFactor * lerp2(nc[q3 + 1], nc[q6 + 1], mu);
normalArray[c + 2] = normalFactor * lerp2(nc[q3 + 2], nc[q6 + 2], mu);
}
if (atomindex) {
atomindexArray[count] = atomindex[q + Math.round(mu) * yd];
}
vertexIndex[_q] = count;
ilist[offset] = count;
count += 1;
} else {
ilist[offset] = vertexIndex[_q];
}
}
__name(VIntY, "VIntY");
function VIntZ(q, offset, x, y, z, valp1, valp2) {
var _q = 3 * q + 2;
if (vertexIndex[_q] < 0) {
var mu = (isolevel - valp1) / (valp2 - valp1);
var nc = normalCache;
var c = count * 3;
positionArray[c] = x;
positionArray[c + 1] = y;
positionArray[c + 2] = z + mu;
if (!noNormals) {
var q3 = q * 3;
var q6 = q3 + zd * 3;
normalArray[c] = normalFactor * lerp2(nc[q3], nc[q6], mu);
normalArray[c + 1] = normalFactor * lerp2(nc[q3 + 1], nc[q6 + 1], mu);
normalArray[c + 2] = normalFactor * lerp2(nc[q3 + 2], nc[q6 + 2], mu);
}
if (atomindex) {
atomindexArray[count] = atomindex[q + Math.round(mu) * zd];
}
vertexIndex[_q] = count;
ilist[offset] = count;
count += 1;
} else {
ilist[offset] = vertexIndex[_q];
}
}
__name(VIntZ, "VIntZ");
function compNorm(q) {
var q3 = q * 3;
if (normalCache[q3] === 0) {
normalCache[q3] = field[(q - 1 + n) % n] - field[(q + 1) % n];
normalCache[q3 + 1] = field[(q - yd + n) % n] - field[(q + yd) % n];
normalCache[q3 + 2] = field[(q - zd + n) % n] - field[(q + zd) % n];
}
}
__name(compNorm, "compNorm");
function polygonize(fx, fy, fz, q, edgeFilter) {
var q1;
var qy;
var qz;
var q1y;
var q1z;
var qyz;
var q1yz;
if (wrap) {
q = index(fx, fy, fz);
q1 = index(fx + 1, fy, fz);
qy = index(fx, fy + 1, fz);
qz = index(fx, fy, fz + 1);
q1y = index(fx + 1, fy + 1, fz);
q1z = index(fx + 1, fy, fz + 1);
qyz = index(fx, fy + 1, fz + 1);
q1yz = index(fx + 1, fy + 1, fz + 1);
} else {
q1 = q + 1;
qy = q + yd;
qz = q + zd;
q1y = qy + 1;
q1z = qz + 1;
qyz = qy + zd;
q1yz = qyz + 1;
}
var cubeindex = 0;
var field0 = field[q];
var field1 = field[q1];
var field2 = field[qy];
var field3 = field[q1y];
var field4 = field[qz];
var field5 = field[q1z];
var field6 = field[qyz];
var field7 = field[q1yz];
if (field0 < isolevel) {
cubeindex |= 1;
}
if (field1 < isolevel) {
cubeindex |= 2;
}
if (field2 < isolevel) {
cubeindex |= 8;
}
if (field3 < isolevel) {
cubeindex |= 4;
}
if (field4 < isolevel) {
cubeindex |= 16;
}
if (field5 < isolevel) {
cubeindex |= 32;
}
if (field6 < isolevel) {
cubeindex |= 128;
}
if (field7 < isolevel) {
cubeindex |= 64;
}
var bits = edgeTable[cubeindex];
if (bits === 0) {
return 0;
}
var fx2 = fx + 1;
var fy2 = fy + 1;
var fz2 = fz + 1;
if (bits & 1) {
if (!noNormals) {
compNorm(q);
compNorm(q1);
}
VIntX(q, 0, fx, fy, fz, field0, field1);
}
if (bits & 2) {
if (!noNormals) {
compNorm(q1);
compNorm(q1y);
}
VIntY(q1, 1, fx2, fy, fz, field1, field3);
}
if (bits & 4) {
if (!noNormals) {
compNorm(qy);
compNorm(q1y);
}
VIntX(qy, 2, fx, fy2, fz, field2, field3);
}
if (bits & 8) {
if (!noNormals) {
compNorm(q);
compNorm(qy);
}
VIntY(q, 3, fx, fy, fz, field0, field2);
}
if (bits & 16) {
if (!noNormals) {
compNorm(qz);
compNorm(q1z);
}
VIntX(qz, 4, fx, fy, fz2, field4, field5);
}
if (bits & 32) {
if (!noNormals) {
compNorm(q1z);
compNorm(q1yz);
}
VIntY(q1z, 5, fx2, fy, fz2, field5, field7);
}
if (bits & 64) {
if (!noNormals) {
compNorm(qyz);
compNorm(q1yz);
}
VIntX(qyz, 6, fx, fy2, fz2, field6, field7);
}
if (bits & 128) {
if (!noNormals) {
compNorm(qz);
compNorm(qyz);
}
VIntY(qz, 7, fx, fy, fz2, field4, field6);
}
if (bits & 256) {
if (!noNormals) {
compNorm(q);
compNorm(qz);
}
VIntZ(q, 8, fx, fy, fz, field0, field4);
}
if (bits & 512) {
if (!noNormals) {
compNorm(q1);
compNorm(q1z);
}
VIntZ(q1, 9, fx2, fy, fz, field1, field5);
}
if (bits & 1024) {
if (!noNormals) {
compNorm(q1y);
compNorm(q1yz);
}
VIntZ(q1y, 10, fx2, fy2, fz, field3, field7);
}
if (bits & 2048) {
if (!noNormals) {
compNorm(qy);
compNorm(qyz);
}
VIntZ(qy, 11, fx, fy2, fz, field2, field6);
}
var triIndex = cubeindex << 4;
var e1;
var e2;
var e3;
var i = 0;
while (triTable[triIndex + i] !== -1) {
e1 = triTable[triIndex + i];
e2 = triTable[triIndex + i + 1];
e3 = triTable[triIndex + i + 2];
if (contour) {
if (allowedContours[e1][e2] & edgeFilter) {
indexArray[icount++] = ilist[e1];
indexArray[icount++] = ilist[e2];
}
if (allowedContours[e2][e3] & edgeFilter) {
indexArray[icount++] = ilist[e2];
indexArray[icount++] = ilist[e3];
}
if (allowedContours[e1][e3] & edgeFilter) {
indexArray[icount++] = ilist[e1];
indexArray[icount++] = ilist[e3];
}
} else {
indexArray[icount++] = ilist[isNegativeIso ? e1 : e2];
indexArray[icount++] = ilist[isNegativeIso ? e2 : e1];
indexArray[icount++] = ilist[e3];
}
i += 3;
}
}
__name(polygonize, "polygonize");
function triangulate(xBeg, yBeg, zBeg, xEnd, yEnd, zEnd) {
var q;
var q3;
var x;
var y;
var z;
var yOffset;
var zOffset;
xBeg = xBeg !== void 0 ? xBeg : 0;
yBeg = yBeg !== void 0 ? yBeg : 0;
zBeg = zBeg !== void 0 ? zBeg : 0;
xEnd = xEnd !== void 0 ? xEnd : nx - 1;
yEnd = yEnd !== void 0 ? yEnd : ny - 1;
zEnd = zEnd !== void 0 ? zEnd : nz - 1;
if (!wrap) {
if (noNormals) {
xBeg = Math.max(0, xBeg);
yBeg = Math.max(0, yBeg);
zBeg = Math.max(0, zBeg);
xEnd = Math.min(nx - 1, xEnd);
yEnd = Math.min(ny - 1, yEnd);
zEnd = Math.min(nz - 1, zEnd);
} else {
xBeg = Math.max(1, xBeg);
yBeg = Math.max(1, yBeg);
zBeg = Math.max(1, zBeg);
xEnd = Math.min(nx - 2, xEnd);
yEnd = Math.min(ny - 2, yEnd);
zEnd = Math.min(nz - 2, zEnd);
}
}
var xBeg2, yBeg2, zBeg2, xEnd2, yEnd2, zEnd2;
if (!wrap) {
xBeg2 = Math.max(0, xBeg - 2);
yBeg2 = Math.max(0, yBeg - 2);
zBeg2 = Math.max(0, zBeg - 2);
xEnd2 = Math.min(nx, xEnd + 2);
yEnd2 = Math.min(ny, yEnd + 2);
zEnd2 = Math.min(nz, zEnd + 2);
for (z = zBeg2; z < zEnd2; ++z) {
zOffset = zd * z;
for (y = yBeg2; y < yEnd2; ++y) {
yOffset = zOffset + yd * y;
for (x = xBeg2; x < xEnd2; ++x) {
q = 3 * (yOffset + x);
vertexIndex[q] = -1;
vertexIndex[q + 1] = -1;
vertexIndex[q + 2] = -1;
}
}
}
} else {
xBeg2 = xBeg - 2;
yBeg2 = yBeg - 2;
zBeg2 = zBeg - 2;
xEnd2 = xEnd + 2;
yEnd2 = yEnd + 2;
zEnd2 = zEnd + 2;
for (z = zBeg2; z < zEnd2; ++z) {
for (y = yBeg2; y < yEnd2; ++y) {
for (x = xBeg2; x < xEnd2; ++x) {
q3 = index(x, y, z) * 3;
vertexIndex[q3] = -1;
vertexIndex[q3 + 1] = -1;
vertexIndex[q3 + 2] = -1;
}
}
}
}
if (!wrap) {
var __break;
var __xBeg = xBeg;
var __yBeg = yBeg;
var __zBeg = zBeg;
var __xEnd = xEnd;
var __yEnd = yEnd;
var __zEnd = zEnd;
__break = false;
for (z = zBeg; z < zEnd; ++z) {
for (y = yBeg; y < yEnd; ++y) {
for (x = xBeg; x < xEnd; ++x) {
q = nx * ny * z + nx * y + x;
if (field[q] >= isolevel) {
__zBeg = z;
__break = true;
break;
}
}
if (__break) {
break;
}
}
if (__break) {
break;
}
}
__break = false;
for (y = yBeg; y < yEnd; ++y) {
for (z = __zBeg; z < zEnd; ++z) {
for (x = xBeg; x < xEnd; ++x) {
q = nx * ny * z + nx * y + x;
if (field[q] >= isolevel) {
__yBeg = y;
__break = true;
break;
}
}
if (__break) {
break;
}
}
if (__break) {
break;
}
}
__break = false;
for (x = xBeg; x < xEnd; ++x) {
for (y = __yBeg; y < yEnd; ++y) {
for (z = __zBeg; z < zEnd; ++z) {
q = nx * ny * z + nx * y + x;
if (field[q] >= isolevel) {
__xBeg = x;
__break = true;
break;
}
}
if (__break) {
break;
}
}
if (__break) {
break;
}
}
__break = false;
for (z = zEnd; z >= zBeg; --z) {
for (y = yEnd; y >= yBeg; --y) {
for (x = xEnd; x >= xBeg; --x) {
q = nx * ny * z + nx * y + x;
if (field[q] >= isolevel) {
__zEnd = z;
__break = true;
break;
}
}
if (__break) {
break;
}
}
if (__break) {
break;
}
}
__break = false;
for (y = yEnd; y >= yBeg; --y) {
for (z = __zEnd; z >= zBeg; --z) {
for (x = xEnd; x >= xBeg; --x) {
q = nx * ny * z + nx * y + x;
if (field[q] >= isolevel) {
__yEnd = y;
__break = true;
break;
}
}
if (__break) {
break;
}
}
if (__break) {
break;
}
}
__break = false;
for (x = xEnd; x >= xBeg; --x) {
for (y = __yEnd; y >= yBeg; --y) {
for (z = __zEnd; z >= zBeg; --z) {
q = nx * ny * z + nx * y + x;
if (field[q] >= isolevel) {
__xEnd = x;
__break = true;
break;
}
}
if (__break) {
break;
}
}
if (__break) {
break;
}
}
if (noNormals) {
xBeg = Math.max(0, __xBeg - 1);
yBeg = Math.max(0, __yBeg - 1);
zBeg = Math.max(0, __zBeg - 1);
xEnd = Math.min(nx - 1, __xEnd + 1);
yEnd = Math.min(ny - 1, __yEnd + 1);
zEnd = Math.min(nz - 1, __zEnd + 1);
} else {
xBeg = Math.max(1, __xBeg - 1);
yBeg = Math.max(1, __yBeg - 1);
zBeg = Math.max(1, __zBeg - 1);
xEnd = Math.min(nx - 2, __xEnd + 1);
yEnd = Math.min(ny - 2, __yEnd + 1);
zEnd = Math.min(nz - 2, __zEnd + 1);
}
}
var edgeFilter = 15;
for (z = zBeg; z < zEnd; ++z, edgeFilter &= ~4) {
zOffset = zd * z;
edgeFilter |= 2;
for (y = yBeg; y < yEnd; ++y, edgeFilter &= ~2) {
yOffset = zOffset + yd * y;
edgeFilter |= 1;
for (x = xBeg; x < xEnd; ++x, edgeFilter &= ~1) {
q = yOffset + x;
polygonize(x, y, z, q, edgeFilter);
}
}
}
}
__name(triangulate, "triangulate");
}
function VolumeSurface(data, nx, ny, nz, atomindex) {
var mc = new MarchingCubes(data, nx, ny, nz, atomindex);
function getSurface3(isolevel, smooth, box, matrix2, contour, wrap) {
if (wrap === void 0)
wrap = false;
var sd = mc.triangulate(isolevel, smooth, box, contour, wrap);
if (smooth && !contour) {
laplacianSmooth(sd.position, sd.index, smooth, true);
sd.normal = computeVertexNormals(sd.position, sd.index);
}
if (matrix2) {
applyMatrix4toVector3array(matrix2, sd.position);
if (sd.normal) {
var normalMatrix2 = m3new();
m3makeNormal(normalMatrix2, matrix2);
applyMatrix3toVector3array(normalMatrix2, sd.normal);
}
}
return sd;
}
__name(getSurface3, "getSurface");
this.getSurface = getSurface3;
}
function applyMatrix3toVector3array(m, a) {
for (var i = 0, il = a.length; i < il; i += 3) {
var x = a[i];
var y = a[i + 1];
var z = a[i + 2];
a[i] = m[0] * x + m[3] * y + m[6] * z;
a[i + 1] = m[1] * x + m[4] * y + m[7] * z;
a[i + 2] = m[2] * x + m[5] * y + m[8] * z;
}
}
function applyMatrix4toVector3array(m, a) {
for (var i = 0, il = a.length; i < il; i += 3) {
var x = a[i];
var y = a[i + 1];
var z = a[i + 2];
a[i] = m[0] * x + m[4] * y + m[8] * z + m[12];
a[i + 1] = m[1] * x + m[5] * y + m[9] * z + m[13];
a[i + 2] = m[2] * x + m[6] * y + m[10] * z + m[14];
}
}
function computeBoundingBox(array) {
var minX = Infinity;
var minY = Infinity;
var minZ = Infinity;
var maxX = -Infinity;
var maxY = -Infinity;
var maxZ = -Infinity;
for (var i = 0, l = array.length; i < l; i += 3) {
var x = array[i];
var y = array[i + 1];
var z = array[i + 2];
if (x < minX) {
minX = x;
}
if (y < minY) {
minY = y;
}
if (z < minZ) {
minZ = z;
}
if (x > maxX) {
maxX = x;
}
if (y > maxY) {
maxY = y;
}
if (z > maxZ) {
maxZ = z;
}
}
return [
v3new([minX, minY, minZ]),
v3new([maxX, maxY, maxZ])
];
}
function computeVertexNormals(position, index, normal) {
var i, il;
if (normal === void 0) {
normal = new Float32Array(position.length);
} else {
for (i = 0, il = normal.length; i < il; i++) {
normal[i] = 0;
}
}
var a = new Float32Array(3);
var b = new Float32Array(3);
var c = new Float32Array(3);
var cb = new Float32Array(3);
var ab = new Float32Array(3);
if (index) {
for (i = 0, il = index.length; i < il; i += 3) {
var ai = index[i] * 3;
var bi = index[i + 1] * 3;
var ci = index[i + 2] * 3;
v3fromArray(a, position, ai);
v3fromArray(b, position, bi);
v3fromArray(c, position, ci);
v3sub(cb, c, b);
v3sub(ab, a, b);
v3cross(cb, cb, ab);
normal[ai] += cb[0];
normal[ai + 1] += cb[1];
normal[ai + 2] += cb[2];
normal[bi] += cb[0];
normal[bi + 1] += cb[1];
normal[bi + 2] += cb[2];
normal[ci] += cb[0];
normal[ci + 1] += cb[1];
normal[ci + 2] += cb[2];
}
} else {
for (i = 0, il = position.length; i < il; i += 9) {
v3fromArray(a, position, i);
v3fromArray(b, position, i + 3);
v3fromArray(c, position, i + 6);
v3sub(cb, c, b);
v3sub(ab, a, b);
v3cross(cb, cb, ab);
normal[i] = cb[0];
normal[i + 1] = cb[1];
normal[i + 2] = cb[2];
normal[i + 3] = cb[0];
normal[i + 4] = cb[1];
normal[i + 5] = cb[2];
normal[i + 6] = cb[0];
normal[i + 7] = cb[1];
normal[i + 8] = cb[2];
}
}
normalizeVector3array(normal);
return normal;
}
function defaults(value2, defaultValue) {
return value2 !== void 0 ? value2 : defaultValue;
}
function degToRad(deg) {
return deg * 0.01745;
}
function getAllowedContours() {
return [
[0, 4, 4, 4, 2, 0, 0, 0, 2, 2, 0, 0],
[4, 0, 4, 4, 0, 8, 0, 0, 0, 8, 8, 0],
[4, 4, 0, 4, 0, 0, 8, 0, 0, 0, 8, 8],
[4, 4, 4, 0, 0, 0, 0, 1, 1, 0, 0, 1],
[2, 0, 0, 0, 0, 8, 8, 8, 2, 2, 0, 0],
[0, 8, 0, 0, 8, 0, 8, 8, 0, 8, 8, 0],
[0, 0, 8, 0, 8, 8, 0, 8, 0, 0, 8, 8],
[0, 0, 0, 1, 8, 8, 8, 0, 1, 0, 0, 1],
[2, 0, 0, 1, 2, 0, 0, 1, 0, 2, 0, 1],
[2, 8, 0, 0, 2, 8, 0, 0, 2, 0, 8, 0],
[0, 8, 8, 0, 0, 8, 8, 0, 0, 8, 0, 8],
[0, 0, 8, 1, 0, 0, 8, 1, 1, 0, 8, 0]
];
}
function getEdgeTable() {
return new Uint32Array([
0,
265,
515,
778,
1030,
1295,
1541,
1804,
2060,
2309,
2575,
2822,
3082,
3331,
3593,
3840,
400,
153,
915,
666,
1430,
1183,
1941,
1692,
2460,
2197,
2975,
2710,
3482,
3219,
3993,
3728,
560,
825,
51,
314,
1590,
1855,
1077,
1340,
2620,
2869,
2111,
2358,
3642,
3891,
3129,
3376,
928,
681,
419,
170,
1958,
1711,
1445,
1196,
2988,
2725,
2479,
2214,
4010,
3747,
3497,
3232,
1120,
1385,
1635,
1898,
102,
367,
613,
876,
3180,
3429,
3695,
3942,
2154,
2403,
2665,
2912,
1520,
1273,
2035,
1786,
502,
255,
1013,
764,
3580,
3317,
4095,
3830,
2554,
2291,
3065,
2800,
1616,
1881,
1107,
1370,
598,
863,
85,
348,
3676,
3925,
3167,
3414,
2650,
2899,
2137,
2384,
1984,
1737,
1475,
1226,
966,
719,
453,
204,
4044,
3781,
3535,
3270,
3018,
2755,
2505,
2240,
2240,
2505,
2755,
3018,
3270,
3535,
3781,
4044,
204,
453,
719,
966,
1226,
1475,
1737,
1984,
2384,
2137,
2899,
2650,
3414,
3167,
3925,
3676,
348,
85,
863,
598,
1370,
1107,
1881,
1616,
2800,
3065,
2291,
2554,
3830,
4095,
3317,
3580,
764,
1013,
255,
502,
1786,
2035,
1273,
1520,
2912,
2665,
2403,
2154,
3942,
3695,
3429,
3180,
876,
613,
367,
102,
1898,
1635,
1385,
1120,
3232,
3497,
3747,
4010,
2214,
2479,
2725,
2988,
1196,
1445,
1711,
1958,
170,
419,
681,
928,
3376,
3129,
3891,
3642,
2358,
2111,
2869,
2620,
1340,
1077,
1855,
1590,
314,
51,
825,
560,
3728,
3993,
3219,
3482,
2710,
2975,
2197,
2460,
1692,
1941,
1183,
1430,
666,
915,
153,
400,
3840,
3593,
3331,
3082,
2822,
2575,
2309,
2060,
1804,
1541,
1295,
1030,
778,
515,
265,
0
]);
}
function getRadiusDict(radiusList) {
var radiusDict = {};
for (var i = 0, il = radiusList.length; i < il; ++i) {
radiusDict[radiusList[i]] = true;
}
return radiusDict;
}
function getSurfaceGrid(min, max, maxRadius2, scaleFactor, extraMargin) {
var margin = 1 / scaleFactor * 3;
margin += maxRadius2;
v3subScalar(min, min, extraMargin + margin);
v3addScalar(max, max, extraMargin + margin);
v3multiplyScalar(min, min, scaleFactor);
v3floor(min, min);
v3divideScalar(min, min, scaleFactor);
v3multiplyScalar(max, max, scaleFactor);
v3ceil(max, max);
v3divideScalar(max, max, scaleFactor);
var dim = new Float32Array(3);
v3sub(dim, max, min);
v3multiplyScalar(dim, dim, scaleFactor);
v3ceil(dim, dim);
v3addScalar(dim, dim, 1);
var maxSize = Math.pow(10, 6) * 256;
var tmpSize = dim[0] * dim[1] * dim[2] * 3;
if (maxSize <= tmpSize) {
scaleFactor *= Math.pow(maxSize / tmpSize, 1 / 3);
v3multiplyScalar(min, min, scaleFactor);
v3floor(min, min);
v3divideScalar(min, min, scaleFactor);
v3multiplyScalar(max, max, scaleFactor);
v3ceil(max, max);
v3divideScalar(max, max, scaleFactor);
v3sub(dim, max, min);
v3multiplyScalar(dim, dim, scaleFactor);
v3ceil(dim, dim);
v3addScalar(dim, dim, 1);
}
var tran = new Float32Array(min);
v3negate(tran, tran);
var matrix2 = m4new();
var mroty = m4new();
m4makeRotationY(mroty, degToRad(90));
m4multiply(matrix2, matrix2, mroty);
var mscale = m4new();
m4makeScale(mscale, -1 / scaleFactor, 1 / scaleFactor, 1 / scaleFactor);
m4multiply(matrix2, matrix2, mscale);
var mtrans = m4new();
m4makeTranslation(mtrans, -scaleFactor * tran[2], -scaleFactor * tran[1], -scaleFactor * tran[0]);
m4multiply(matrix2, matrix2, mtrans);
return {
dim,
tran,
matrix: matrix2,
scaleFactor
};
}
function getTriTable() {
return new Int32Array([
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
1,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
8,
3,
9,
8,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
3,
1,
2,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
2,
10,
0,
2,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
8,
3,
2,
10,
8,
10,
9,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
11,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
11,
2,
8,
11,
0,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
9,
0,
2,
3,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
11,
2,
1,
9,
11,
9,
8,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
10,
1,
11,
10,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
10,
1,
0,
8,
10,
8,
11,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
9,
0,
3,
11,
9,
11,
10,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
8,
10,
10,
8,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
7,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
3,
0,
7,
3,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
1,
9,
8,
4,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
1,
9,
4,
7,
1,
7,
3,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
10,
8,
4,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
4,
7,
3,
0,
4,
1,
2,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
2,
10,
9,
0,
2,
8,
4,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
10,
9,
2,
9,
7,
2,
7,
3,
7,
9,
4,
-1,
-1,
-1,
-1,
8,
4,
7,
3,
11,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
11,
4,
7,
11,
2,
4,
2,
0,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
0,
1,
8,
4,
7,
2,
3,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
7,
11,
9,
4,
11,
9,
11,
2,
9,
2,
1,
-1,
-1,
-1,
-1,
3,
10,
1,
3,
11,
10,
7,
8,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
11,
10,
1,
4,
11,
1,
0,
4,
7,
11,
4,
-1,
-1,
-1,
-1,
4,
7,
8,
9,
0,
11,
9,
11,
10,
11,
0,
3,
-1,
-1,
-1,
-1,
4,
7,
11,
4,
11,
9,
9,
11,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
5,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
5,
4,
0,
8,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
5,
4,
1,
5,
0,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
8,
5,
4,
8,
3,
5,
3,
1,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
10,
9,
5,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
0,
8,
1,
2,
10,
4,
9,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
2,
10,
5,
4,
2,
4,
0,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
10,
5,
3,
2,
5,
3,
5,
4,
3,
4,
8,
-1,
-1,
-1,
-1,
9,
5,
4,
2,
3,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
11,
2,
0,
8,
11,
4,
9,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
5,
4,
0,
1,
5,
2,
3,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
1,
5,
2,
5,
8,
2,
8,
11,
4,
8,
5,
-1,
-1,
-1,
-1,
10,
3,
11,
10,
1,
3,
9,
5,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
9,
5,
0,
8,
1,
8,
10,
1,
8,
11,
10,
-1,
-1,
-1,
-1,
5,
4,
0,
5,
0,
11,
5,
11,
10,
11,
0,
3,
-1,
-1,
-1,
-1,
5,
4,
8,
5,
8,
10,
10,
8,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
7,
8,
5,
7,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
3,
0,
9,
5,
3,
5,
7,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
7,
8,
0,
1,
7,
1,
5,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
5,
3,
3,
5,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
7,
8,
9,
5,
7,
10,
1,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
1,
2,
9,
5,
0,
5,
3,
0,
5,
7,
3,
-1,
-1,
-1,
-1,
8,
0,
2,
8,
2,
5,
8,
5,
7,
10,
5,
2,
-1,
-1,
-1,
-1,
2,
10,
5,
2,
5,
3,
3,
5,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
7,
9,
5,
7,
8,
9,
3,
11,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
5,
7,
9,
7,
2,
9,
2,
0,
2,
7,
11,
-1,
-1,
-1,
-1,
2,
3,
11,
0,
1,
8,
1,
7,
8,
1,
5,
7,
-1,
-1,
-1,
-1,
11,
2,
1,
11,
1,
7,
7,
1,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
5,
8,
8,
5,
7,
10,
1,
3,
10,
3,
11,
-1,
-1,
-1,
-1,
5,
7,
0,
5,
0,
9,
7,
11,
0,
1,
0,
10,
11,
10,
0,
-1,
11,
10,
0,
11,
0,
3,
10,
5,
0,
8,
0,
7,
5,
7,
0,
-1,
11,
10,
5,
7,
11,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
6,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
3,
5,
10,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
0,
1,
5,
10,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
8,
3,
1,
9,
8,
5,
10,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
6,
5,
2,
6,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
6,
5,
1,
2,
6,
3,
0,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
6,
5,
9,
0,
6,
0,
2,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
9,
8,
5,
8,
2,
5,
2,
6,
3,
2,
8,
-1,
-1,
-1,
-1,
2,
3,
11,
10,
6,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
11,
0,
8,
11,
2,
0,
10,
6,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
1,
9,
2,
3,
11,
5,
10,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
10,
6,
1,
9,
2,
9,
11,
2,
9,
8,
11,
-1,
-1,
-1,
-1,
6,
3,
11,
6,
5,
3,
5,
1,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
11,
0,
11,
5,
0,
5,
1,
5,
11,
6,
-1,
-1,
-1,
-1,
3,
11,
6,
0,
3,
6,
0,
6,
5,
0,
5,
9,
-1,
-1,
-1,
-1,
6,
5,
9,
6,
9,
11,
11,
9,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
10,
6,
4,
7,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
3,
0,
4,
7,
3,
6,
5,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
9,
0,
5,
10,
6,
8,
4,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
6,
5,
1,
9,
7,
1,
7,
3,
7,
9,
4,
-1,
-1,
-1,
-1,
6,
1,
2,
6,
5,
1,
4,
7,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
5,
5,
2,
6,
3,
0,
4,
3,
4,
7,
-1,
-1,
-1,
-1,
8,
4,
7,
9,
0,
5,
0,
6,
5,
0,
2,
6,
-1,
-1,
-1,
-1,
7,
3,
9,
7,
9,
4,
3,
2,
9,
5,
9,
6,
2,
6,
9,
-1,
3,
11,
2,
7,
8,
4,
10,
6,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
10,
6,
4,
7,
2,
4,
2,
0,
2,
7,
11,
-1,
-1,
-1,
-1,
0,
1,
9,
4,
7,
8,
2,
3,
11,
5,
10,
6,
-1,
-1,
-1,
-1,
9,
2,
1,
9,
11,
2,
9,
4,
11,
7,
11,
4,
5,
10,
6,
-1,
8,
4,
7,
3,
11,
5,
3,
5,
1,
5,
11,
6,
-1,
-1,
-1,
-1,
5,
1,
11,
5,
11,
6,
1,
0,
11,
7,
11,
4,
0,
4,
11,
-1,
0,
5,
9,
0,
6,
5,
0,
3,
6,
11,
6,
3,
8,
4,
7,
-1,
6,
5,
9,
6,
9,
11,
4,
7,
9,
7,
11,
9,
-1,
-1,
-1,
-1,
10,
4,
9,
6,
4,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
10,
6,
4,
9,
10,
0,
8,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
0,
1,
10,
6,
0,
6,
4,
0,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
8,
3,
1,
8,
1,
6,
8,
6,
4,
6,
1,
10,
-1,
-1,
-1,
-1,
1,
4,
9,
1,
2,
4,
2,
6,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
0,
8,
1,
2,
9,
2,
4,
9,
2,
6,
4,
-1,
-1,
-1,
-1,
0,
2,
4,
4,
2,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
8,
3,
2,
8,
2,
4,
4,
2,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
4,
9,
10,
6,
4,
11,
2,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
2,
2,
8,
11,
4,
9,
10,
4,
10,
6,
-1,
-1,
-1,
-1,
3,
11,
2,
0,
1,
6,
0,
6,
4,
6,
1,
10,
-1,
-1,
-1,
-1,
6,
4,
1,
6,
1,
10,
4,
8,
1,
2,
1,
11,
8,
11,
1,
-1,
9,
6,
4,
9,
3,
6,
9,
1,
3,
11,
6,
3,
-1,
-1,
-1,
-1,
8,
11,
1,
8,
1,
0,
11,
6,
1,
9,
1,
4,
6,
4,
1,
-1,
3,
11,
6,
3,
6,
0,
0,
6,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
6,
4,
8,
11,
6,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
7,
10,
6,
7,
8,
10,
8,
9,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
7,
3,
0,
10,
7,
0,
9,
10,
6,
7,
10,
-1,
-1,
-1,
-1,
10,
6,
7,
1,
10,
7,
1,
7,
8,
1,
8,
0,
-1,
-1,
-1,
-1,
10,
6,
7,
10,
7,
1,
1,
7,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
6,
1,
6,
8,
1,
8,
9,
8,
6,
7,
-1,
-1,
-1,
-1,
2,
6,
9,
2,
9,
1,
6,
7,
9,
0,
9,
3,
7,
3,
9,
-1,
7,
8,
0,
7,
0,
6,
6,
0,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
7,
3,
2,
6,
7,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
3,
11,
10,
6,
8,
10,
8,
9,
8,
6,
7,
-1,
-1,
-1,
-1,
2,
0,
7,
2,
7,
11,
0,
9,
7,
6,
7,
10,
9,
10,
7,
-1,
1,
8,
0,
1,
7,
8,
1,
10,
7,
6,
7,
10,
2,
3,
11,
-1,
11,
2,
1,
11,
1,
7,
10,
6,
1,
6,
7,
1,
-1,
-1,
-1,
-1,
8,
9,
6,
8,
6,
7,
9,
1,
6,
11,
6,
3,
1,
3,
6,
-1,
0,
9,
1,
11,
6,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
7,
8,
0,
7,
0,
6,
3,
11,
0,
11,
6,
0,
-1,
-1,
-1,
-1,
7,
11,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
7,
6,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
0,
8,
11,
7,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
1,
9,
11,
7,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
8,
1,
9,
8,
3,
1,
11,
7,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
1,
2,
6,
11,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
10,
3,
0,
8,
6,
11,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
9,
0,
2,
10,
9,
6,
11,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
6,
11,
7,
2,
10,
3,
10,
8,
3,
10,
9,
8,
-1,
-1,
-1,
-1,
7,
2,
3,
6,
2,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
7,
0,
8,
7,
6,
0,
6,
2,
0,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
7,
6,
2,
3,
7,
0,
1,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
6,
2,
1,
8,
6,
1,
9,
8,
8,
7,
6,
-1,
-1,
-1,
-1,
10,
7,
6,
10,
1,
7,
1,
3,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
7,
6,
1,
7,
10,
1,
8,
7,
1,
0,
8,
-1,
-1,
-1,
-1,
0,
3,
7,
0,
7,
10,
0,
10,
9,
6,
10,
7,
-1,
-1,
-1,
-1,
7,
6,
10,
7,
10,
8,
8,
10,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
6,
8,
4,
11,
8,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
6,
11,
3,
0,
6,
0,
4,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
8,
6,
11,
8,
4,
6,
9,
0,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
4,
6,
9,
6,
3,
9,
3,
1,
11,
3,
6,
-1,
-1,
-1,
-1,
6,
8,
4,
6,
11,
8,
2,
10,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
10,
3,
0,
11,
0,
6,
11,
0,
4,
6,
-1,
-1,
-1,
-1,
4,
11,
8,
4,
6,
11,
0,
2,
9,
2,
10,
9,
-1,
-1,
-1,
-1,
10,
9,
3,
10,
3,
2,
9,
4,
3,
11,
3,
6,
4,
6,
3,
-1,
8,
2,
3,
8,
4,
2,
4,
6,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
4,
2,
4,
6,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
9,
0,
2,
3,
4,
2,
4,
6,
4,
3,
8,
-1,
-1,
-1,
-1,
1,
9,
4,
1,
4,
2,
2,
4,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
8,
1,
3,
8,
6,
1,
8,
4,
6,
6,
10,
1,
-1,
-1,
-1,
-1,
10,
1,
0,
10,
0,
6,
6,
0,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
6,
3,
4,
3,
8,
6,
10,
3,
0,
3,
9,
10,
9,
3,
-1,
10,
9,
4,
6,
10,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
9,
5,
7,
6,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
3,
4,
9,
5,
11,
7,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
0,
1,
5,
4,
0,
7,
6,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
11,
7,
6,
8,
3,
4,
3,
5,
4,
3,
1,
5,
-1,
-1,
-1,
-1,
9,
5,
4,
10,
1,
2,
7,
6,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
6,
11,
7,
1,
2,
10,
0,
8,
3,
4,
9,
5,
-1,
-1,
-1,
-1,
7,
6,
11,
5,
4,
10,
4,
2,
10,
4,
0,
2,
-1,
-1,
-1,
-1,
3,
4,
8,
3,
5,
4,
3,
2,
5,
10,
5,
2,
11,
7,
6,
-1,
7,
2,
3,
7,
6,
2,
5,
4,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
5,
4,
0,
8,
6,
0,
6,
2,
6,
8,
7,
-1,
-1,
-1,
-1,
3,
6,
2,
3,
7,
6,
1,
5,
0,
5,
4,
0,
-1,
-1,
-1,
-1,
6,
2,
8,
6,
8,
7,
2,
1,
8,
4,
8,
5,
1,
5,
8,
-1,
9,
5,
4,
10,
1,
6,
1,
7,
6,
1,
3,
7,
-1,
-1,
-1,
-1,
1,
6,
10,
1,
7,
6,
1,
0,
7,
8,
7,
0,
9,
5,
4,
-1,
4,
0,
10,
4,
10,
5,
0,
3,
10,
6,
10,
7,
3,
7,
10,
-1,
7,
6,
10,
7,
10,
8,
5,
4,
10,
4,
8,
10,
-1,
-1,
-1,
-1,
6,
9,
5,
6,
11,
9,
11,
8,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
6,
11,
0,
6,
3,
0,
5,
6,
0,
9,
5,
-1,
-1,
-1,
-1,
0,
11,
8,
0,
5,
11,
0,
1,
5,
5,
6,
11,
-1,
-1,
-1,
-1,
6,
11,
3,
6,
3,
5,
5,
3,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
10,
9,
5,
11,
9,
11,
8,
11,
5,
6,
-1,
-1,
-1,
-1,
0,
11,
3,
0,
6,
11,
0,
9,
6,
5,
6,
9,
1,
2,
10,
-1,
11,
8,
5,
11,
5,
6,
8,
0,
5,
10,
5,
2,
0,
2,
5,
-1,
6,
11,
3,
6,
3,
5,
2,
10,
3,
10,
5,
3,
-1,
-1,
-1,
-1,
5,
8,
9,
5,
2,
8,
5,
6,
2,
3,
8,
2,
-1,
-1,
-1,
-1,
9,
5,
6,
9,
6,
0,
0,
6,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
5,
8,
1,
8,
0,
5,
6,
8,
3,
8,
2,
6,
2,
8,
-1,
1,
5,
6,
2,
1,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
3,
6,
1,
6,
10,
3,
8,
6,
5,
6,
9,
8,
9,
6,
-1,
10,
1,
0,
10,
0,
6,
9,
5,
0,
5,
6,
0,
-1,
-1,
-1,
-1,
0,
3,
8,
5,
6,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
5,
6,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
11,
5,
10,
7,
5,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
11,
5,
10,
11,
7,
5,
8,
3,
0,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
11,
7,
5,
10,
11,
1,
9,
0,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
10,
7,
5,
10,
11,
7,
9,
8,
1,
8,
3,
1,
-1,
-1,
-1,
-1,
11,
1,
2,
11,
7,
1,
7,
5,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
3,
1,
2,
7,
1,
7,
5,
7,
2,
11,
-1,
-1,
-1,
-1,
9,
7,
5,
9,
2,
7,
9,
0,
2,
2,
11,
7,
-1,
-1,
-1,
-1,
7,
5,
2,
7,
2,
11,
5,
9,
2,
3,
2,
8,
9,
8,
2,
-1,
2,
5,
10,
2,
3,
5,
3,
7,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
8,
2,
0,
8,
5,
2,
8,
7,
5,
10,
2,
5,
-1,
-1,
-1,
-1,
9,
0,
1,
5,
10,
3,
5,
3,
7,
3,
10,
2,
-1,
-1,
-1,
-1,
9,
8,
2,
9,
2,
1,
8,
7,
2,
10,
2,
5,
7,
5,
2,
-1,
1,
3,
5,
3,
7,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
7,
0,
7,
1,
1,
7,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
0,
3,
9,
3,
5,
5,
3,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
8,
7,
5,
9,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
8,
4,
5,
10,
8,
10,
11,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
5,
0,
4,
5,
11,
0,
5,
10,
11,
11,
3,
0,
-1,
-1,
-1,
-1,
0,
1,
9,
8,
4,
10,
8,
10,
11,
10,
4,
5,
-1,
-1,
-1,
-1,
10,
11,
4,
10,
4,
5,
11,
3,
4,
9,
4,
1,
3,
1,
4,
-1,
2,
5,
1,
2,
8,
5,
2,
11,
8,
4,
5,
8,
-1,
-1,
-1,
-1,
0,
4,
11,
0,
11,
3,
4,
5,
11,
2,
11,
1,
5,
1,
11,
-1,
0,
2,
5,
0,
5,
9,
2,
11,
5,
4,
5,
8,
11,
8,
5,
-1,
9,
4,
5,
2,
11,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
5,
10,
3,
5,
2,
3,
4,
5,
3,
8,
4,
-1,
-1,
-1,
-1,
5,
10,
2,
5,
2,
4,
4,
2,
0,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
10,
2,
3,
5,
10,
3,
8,
5,
4,
5,
8,
0,
1,
9,
-1,
5,
10,
2,
5,
2,
4,
1,
9,
2,
9,
4,
2,
-1,
-1,
-1,
-1,
8,
4,
5,
8,
5,
3,
3,
5,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
4,
5,
1,
0,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
8,
4,
5,
8,
5,
3,
9,
0,
5,
0,
3,
5,
-1,
-1,
-1,
-1,
9,
4,
5,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
11,
7,
4,
9,
11,
9,
10,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
8,
3,
4,
9,
7,
9,
11,
7,
9,
10,
11,
-1,
-1,
-1,
-1,
1,
10,
11,
1,
11,
4,
1,
4,
0,
7,
4,
11,
-1,
-1,
-1,
-1,
3,
1,
4,
3,
4,
8,
1,
10,
4,
7,
4,
11,
10,
11,
4,
-1,
4,
11,
7,
9,
11,
4,
9,
2,
11,
9,
1,
2,
-1,
-1,
-1,
-1,
9,
7,
4,
9,
11,
7,
9,
1,
11,
2,
11,
1,
0,
8,
3,
-1,
11,
7,
4,
11,
4,
2,
2,
4,
0,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
11,
7,
4,
11,
4,
2,
8,
3,
4,
3,
2,
4,
-1,
-1,
-1,
-1,
2,
9,
10,
2,
7,
9,
2,
3,
7,
7,
4,
9,
-1,
-1,
-1,
-1,
9,
10,
7,
9,
7,
4,
10,
2,
7,
8,
7,
0,
2,
0,
7,
-1,
3,
7,
10,
3,
10,
2,
7,
4,
10,
1,
10,
0,
4,
0,
10,
-1,
1,
10,
2,
8,
7,
4,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
9,
1,
4,
1,
7,
7,
1,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
9,
1,
4,
1,
7,
0,
8,
1,
8,
7,
1,
-1,
-1,
-1,
-1,
4,
0,
3,
7,
4,
3,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
4,
8,
7,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
10,
8,
10,
11,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
0,
9,
3,
9,
11,
11,
9,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
1,
10,
0,
10,
8,
8,
10,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
1,
10,
11,
3,
10,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
2,
11,
1,
11,
9,
9,
11,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
0,
9,
3,
9,
11,
1,
2,
9,
2,
11,
9,
-1,
-1,
-1,
-1,
0,
2,
11,
8,
0,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
3,
2,
11,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
3,
8,
2,
8,
10,
10,
8,
9,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
9,
10,
2,
0,
9,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
2,
3,
8,
2,
8,
10,
0,
1,
8,
1,
10,
8,
-1,
-1,
-1,
-1,
1,
10,
2,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
1,
3,
8,
9,
1,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
9,
1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
0,
3,
8,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1,
-1
]);
}
function getUintArray(sizeOrArray, maxUint) {
var TypedArray = maxUint > 65535 ? Uint32Array : Uint16Array;
return new TypedArray(sizeOrArray);
}
function laplacianSmooth(verts, faces, numiter, inflate2) {
numiter = numiter || 1;
inflate2 = inflate2 || true;
var nv = verts.length / 3;
var nf = faces.length / 3;
var norms = void 0;
if (inflate2) {
norms = new Float32Array(nv * 3);
}
var tps = new Float32Array(nv * 3);
var i;
var ndeg = 20;
var vertdeg = new Array(ndeg);
for (i = 0; i < ndeg; ++i) {
vertdeg[i] = new Uint32Array(nv);
}
for (i = 0; i < nv; ++i) {
vertdeg[0][i] = 0;
}
var j, jl;
var flagvert;
for (i = 0; i < nf; ++i) {
var ao = i * 3;
var bo = i * 3 + 1;
var co = i * 3 + 2;
flagvert = true;
for (j = 0, jl = vertdeg[0][faces[ao]]; j < jl; ++j) {
if (faces[bo] === vertdeg[j + 1][faces[ao]]) {
flagvert = false;
break;
}
}
if (flagvert) {
vertdeg[0][faces[ao]]++;
vertdeg[vertdeg[0][faces[ao]]][faces[ao]] = faces[bo];
}
flagvert = true;
for (j = 0, jl = vertdeg[0][faces[ao]]; j < jl; ++j) {
if (faces[co] === vertdeg[j + 1][faces[ao]]) {
flagvert = false;
break;
}
}
if (flagvert) {
vertdeg[0][faces[ao]]++;
vertdeg[vertdeg[0][faces[ao]]][faces[ao]] = faces[co];
}
flagvert = true;
for (j = 0, jl = vertdeg[0][faces[bo]]; j < jl; ++j) {
if (faces[ao] === vertdeg[j + 1][faces[bo]]) {
flagvert = false;
break;
}
}
if (flagvert) {
vertdeg[0][faces[bo]]++;
vertdeg[vertdeg[0][faces[bo]]][faces[bo]] = faces[ao];
}
flagvert = true;
for (j = 0, jl = vertdeg[0][faces[bo]]; j < jl; ++j) {
if (faces[co] === vertdeg[j + 1][faces[bo]]) {
flagvert = false;
break;
}
}
if (flagvert) {
vertdeg[0][faces[bo]]++;
vertdeg[vertdeg[0][faces[bo]]][faces[bo]] = faces[co];
}
flagvert = true;
for (j = 0; j < vertdeg[0][faces[co]]; ++j) {
if (faces[ao] === vertdeg[j + 1][faces[co]]) {
flagvert = false;
break;
}
}
if (flagvert) {
vertdeg[0][faces[co]]++;
vertdeg[vertdeg[0][faces[co]]][faces[co]] = faces[ao];
}
flagvert = true;
for (j = 0, jl = vertdeg[0][faces[co]]; j < jl; ++j) {
if (faces[bo] === vertdeg[j + 1][faces[co]]) {
flagvert = false;
break;
}
}
if (flagvert) {
vertdeg[0][faces[co]]++;
vertdeg[vertdeg[0][faces[co]]][faces[co]] = faces[bo];
}
}
var wt = 1;
var wt2 = 0.5;
var i3, vi3, vdi, wtvi, wt2vi;
var ssign = -1;
var scaleFactor = 1;
var outwt = 0.75 / (scaleFactor + 3.5);
for (var k = 0; k < numiter; ++k) {
for (i = 0; i < nv; ++i) {
i3 = i * 3;
vdi = vertdeg[0][i];
if (vdi < 3) {
tps[i3] = verts[i3];
tps[i3 + 1] = verts[i3 + 1];
tps[i3 + 2] = verts[i3 + 2];
} else if (vdi === 3 || vdi === 4) {
tps[i3] = 0;
tps[i3 + 1] = 0;
tps[i3 + 2] = 0;
for (j = 0; j < vdi; ++j) {
vi3 = vertdeg[j + 1][i] * 3;
tps[i3] += verts[vi3];
tps[i3 + 1] += verts[vi3 + 1];
tps[i3 + 2] += verts[vi3 + 2];
}
tps[i3] += wt2 * verts[i3];
tps[i3 + 1] += wt2 * verts[i3 + 1];
tps[i3 + 2] += wt2 * verts[i3 + 2];
wt2vi = wt2 + vdi;
tps[i3] /= wt2vi;
tps[i3 + 1] /= wt2vi;
tps[i3 + 2] /= wt2vi;
} else {
tps[i3] = 0;
tps[i3 + 1] = 0;
tps[i3 + 2] = 0;
for (j = 0; j < vdi; ++j) {
vi3 = vertdeg[j + 1][i] * 3;
tps[i3] += verts[vi3];
tps[i3 + 1] += verts[vi3 + 1];
tps[i3 + 2] += verts[vi3 + 2];
}
tps[i3] += wt * verts[i3];
tps[i3 + 1] += wt * verts[i3 + 1];
tps[i3 + 2] += wt * verts[i3 + 2];
wtvi = wt + vdi;
tps[i3] /= wtvi;
tps[i3 + 1] /= wtvi;
tps[i3 + 2] /= wtvi;
}
}
verts.set(tps);
if (inflate2) {
computeVertexNormals(verts, faces, norms);
var nv3 = nv * 3;
for (i3 = 0; i3 < nv3; i3 += 3) {
verts[i3] += ssign * outwt * norms[i3];
verts[i3 + 1] += ssign * outwt * norms[i3 + 1];
verts[i3 + 2] += ssign * outwt * norms[i3 + 2];
}
}
}
}
function m3makeNormal(out, m4) {
var r0 = v3new([m4[0], m4[1], m4[2]]);
var r1 = v3new([m4[4], m4[5], m4[6]]);
var r2 = v3new([m4[8], m4[9], m4[10]]);
var cp = v3new();
v3cross(cp, r1, r2);
out[0] = cp[0];
out[1] = cp[1];
out[2] = cp[2];
v3cross(cp, r2, r0);
out[3] = cp[0];
out[4] = cp[1];
out[5] = cp[2];
v3cross(cp, r0, r1);
out[6] = cp[0];
out[7] = cp[1];
out[8] = cp[2];
}
function m3new() {
return new Float32Array([
1,
0,
0,
0,
1,
0,
0,
0,
1
]);
}
function m4makeRotationY(out, theta) {
var c = Math.cos(theta);
var s = Math.sin(theta);
m4set(out, c, 0, s, 0, 0, 1, 0, 0, -s, 0, c, 0, 0, 0, 0, 1);
}
function m4makeScale(out, x, y, z) {
m4set(out, x, 0, 0, 0, 0, y, 0, 0, 0, 0, z, 0, 0, 0, 0, 1);
}
function m4makeTranslation(out, x, y, z) {
m4set(out, 1, 0, 0, x, 0, 1, 0, y, 0, 0, 1, z, 0, 0, 0, 1);
}
function m4multiply(out, a, b) {
var a11 = a[0];
var a12 = a[4];
var a13 = a[8];
var a14 = a[12];
var a21 = a[1];
var a22 = a[5];
var a23 = a[9];
var a24 = a[13];
var a31 = a[2];
var a32 = a[6];
var a33 = a[10];
var a34 = a[14];
var a41 = a[3];
var a42 = a[7];
var a43 = a[11];
var a44 = a[15];
var b11 = b[0];
var b12 = b[4];
var b13 = b[8];
var b14 = b[12];
var b21 = b[1];
var b22 = b[5];
var b23 = b[9];
var b24 = b[13];
var b31 = b[2];
var b32 = b[6];
var b33 = b[10];
var b34 = b[14];
var b41 = b[3];
var b42 = b[7];
var b43 = b[11];
var b44 = b[15];
out[0] = a11 * b11 + a12 * b21 + a13 * b31 + a14 * b41;
out[4] = a11 * b12 + a12 * b22 + a13 * b32 + a14 * b42;
out[8] = a11 * b13 + a12 * b23 + a13 * b33 + a14 * b43;
out[12] = a11 * b14 + a12 * b24 + a13 * b34 + a14 * b44;
out[1] = a21 * b11 + a22 * b21 + a23 * b31 + a24 * b41;
out[5] = a21 * b12 + a22 * b22 + a23 * b32 + a24 * b42;
out[9] = a21 * b13 + a22 * b23 + a23 * b33 + a24 * b43;
out[13] = a21 * b14 + a22 * b24 + a23 * b34 + a24 * b44;
out[2] = a31 * b11 + a32 * b21 + a33 * b31 + a34 * b41;
out[6] = a31 * b12 + a32 * b22 + a33 * b32 + a34 * b42;
out[10] = a31 * b13 + a32 * b23 + a33 * b33 + a34 * b43;
out[14] = a31 * b14 + a32 * b24 + a33 * b34 + a34 * b44;
out[3] = a41 * b11 + a42 * b21 + a43 * b31 + a44 * b41;
out[7] = a41 * b12 + a42 * b22 + a43 * b32 + a44 * b42;
out[11] = a41 * b13 + a42 * b23 + a43 * b33 + a44 * b43;
out[15] = a41 * b14 + a42 * b24 + a43 * b34 + a44 * b44;
}
function m4new() {
return new Float32Array([
1,
0,
0,
0,
0,
1,
0,
0,
0,
0,
1,
0,
0,
0,
0,
1
]);
}
function m4set(out, n11, n12, n13, n14, n21, n22, n23, n24, n31, n32, n33, n34, n41, n42, n43, n44) {
out[0] = n11;
out[4] = n12;
out[8] = n13;
out[12] = n14;
out[1] = n21;
out[5] = n22;
out[9] = n23;
out[13] = n24;
out[2] = n31;
out[6] = n32;
out[10] = n33;
out[14] = n34;
out[3] = n41;
out[7] = n42;
out[11] = n43;
out[15] = n44;
}
function makeAVHash(atomsX, atomsY, atomsZ, atomsR, min, max, maxDistance) {
maxDistance = Math.max(0.1, maxDistance);
var nAtoms = atomsX.length;
var minX = min[0];
var minY = min[1];
var minZ = min[2];
var maxX = max[0];
var maxY = max[1];
var maxZ = max[2];
function hashFunc(w, minW) {
return Math.floor((w - minW) / maxDistance);
}
__name(hashFunc, "hashFunc");
var iDim = hashFunc(maxX, minX) + 1;
var jDim = hashFunc(maxY, minY) + 1;
var kDim = hashFunc(maxZ, minZ) + 1;
var nCells = iDim * jDim * kDim;
var jkDim = jDim * kDim;
var cellID = /* @__PURE__ */ __name(function(x, y, z) {
return (hashFunc(x, minX) * jDim + hashFunc(y, minY)) * kDim + hashFunc(z, minZ);
}, "cellID");
var preHash = [];
for (var i = 0; i < nAtoms; i++) {
var cid = cellID(atomsX[i], atomsY[i], atomsZ[i]);
if (preHash[cid] === void 0) {
preHash[cid] = [i];
} else {
preHash[cid].push(i);
}
}
var cellOffsets = new Uint32Array(nCells);
var cellLengths = new Uint16Array(nCells);
var data = new Uint32Array(nAtoms);
var offset = 0;
var maxCellLength = 0;
for (i = 0; i < nCells; i++) {
var start = cellOffsets[i] = offset;
var subArray = preHash[i];
if (subArray !== void 0) {
for (var j = 0; j < subArray.length; j++) {
data[offset] = subArray[j];
offset++;
}
}
var cellLength = offset - start;
cellLengths[i] = cellLength;
if (cellLength > maxCellLength) {
maxCellLength = cellLength;
}
}
var neighbourListLength = 27 * maxCellLength + 1;
var withinRadii = /* @__PURE__ */ __name(function(x, y, z, rExtra, out) {
var outIdx = 0;
var nearI = hashFunc(x, minX);
var nearJ = hashFunc(y, minY);
var nearK = hashFunc(z, minZ);
var loI = Math.max(0, nearI - 1);
var loJ = Math.max(0, nearJ - 1);
var loK = Math.max(0, nearK - 1);
var hiI = Math.min(iDim, nearI + 2);
var hiJ = Math.min(jDim, nearJ + 2);
var hiK = Math.min(kDim, nearK + 2);
for (var i2 = loI; i2 < hiI; ++i2) {
var iOffset = i2 * jkDim;
for (var j2 = loJ; j2 < hiJ; ++j2) {
var jOffset = j2 * kDim;
for (var k = loK; k < hiK; ++k) {
var cid2 = iOffset + jOffset + k;
var cellStart = cellOffsets[cid2];
var cellEnd = cellStart + cellLengths[cid2];
for (var dataIndex = cellStart; dataIndex < cellEnd; dataIndex++) {
var atomIndex = data[dataIndex];
var dx = atomsX[atomIndex] - x;
var dy = atomsY[atomIndex] - y;
var dz = atomsZ[atomIndex] - z;
var rSum = atomsR[atomIndex] + rExtra;
if (dx * dx + dy * dy + dz * dz <= rSum * rSum) {
out[outIdx++] = data[dataIndex];
}
}
}
}
}
out[outIdx] = -1;
}, "withinRadii");
return {
neighbourListLength,
withinRadii
};
}
function makeGrid(length2, width, height, DataCtor, elemSize) {
DataCtor = DataCtor || Int32Array;
elemSize = elemSize || 1;
var data = new DataCtor(length2 * width * height * elemSize);
function index(x, y, z) {
return ((x * width + y) * height + z) * elemSize;
}
__name(index, "index");
function set3(x, y, z) {
var args = [], len = arguments.length - 3;
while (len-- > 0)
args[len] = arguments[len + 3];
var i = index(x, y, z);
for (var j = 0; j < elemSize; ++j) {
data[i + j] = args[j];
}
}
__name(set3, "set");
function toArray3(x, y, z, array, offset) {
if (array === void 0)
array = [];
if (offset === void 0)
offset = 0;
var i = index(x, y, z);
for (var j = 0; j < elemSize; ++j) {
array[offset + j] = data[i + j];
}
}
__name(toArray3, "toArray");
function fromArray(x, y, z, array, offset) {
if (offset === void 0)
offset = 0;
var i = index(x, y, z);
for (var j = 0; j < elemSize; ++j) {
data[i + j] = array[offset + j];
}
}
__name(fromArray, "fromArray");
function copy(grid) {
data.set(grid.data);
}
__name(copy, "copy");
return {data, index, set: set3, toArray: toArray3, fromArray, copy};
}
function normalizeVector3array(a) {
for (var i = 0, il = a.length; i < il; i += 3) {
var x = a[i];
var y = a[i + 1];
var z = a[i + 2];
var len2 = x * x + y * y + z * z;
if (len2 > 0) {
var s = 1 / Math.sqrt(len2);
a[i] = x * s;
a[i + 1] = y * s;
a[i + 2] = z * s;
}
}
}
function uniformArray(n, a, optionalTarget) {
var array = optionalTarget || new Float32Array(n);
for (var i = 0; i < n; ++i) {
array[i] = a;
}
return array;
}
function v3addScalar(out, a, s) {
out[0] = a[0] + s;
out[1] = a[1] + s;
out[2] = a[2] + s;
}
function v3ceil(out, a) {
out[0] = Math.ceil(a[0]);
out[1] = Math.ceil(a[1]);
out[2] = Math.ceil(a[2]);
}
function v3cross(out, a, b) {
var ax = a[0];
var ay = a[1];
var az = a[2];
var bx = b[0];
var by = b[1];
var bz = b[2];
out[0] = ay * bz - az * by;
out[1] = az * bx - ax * bz;
out[2] = ax * by - ay * bx;
}
function v3divideScalar(out, a, s) {
v3multiplyScalar(out, a, 1 / s);
}
function v3floor(out, a) {
out[0] = Math.floor(a[0]);
out[1] = Math.floor(a[1]);
out[2] = Math.floor(a[2]);
}
function v3fromArray(out, array, offset) {
if (offset === void 0)
offset = 0;
out[0] = array[offset];
out[1] = array[offset + 1];
out[2] = array[offset + 2];
}
function v3length2(a) {
return a[0] * a[0] + a[1] * a[1] + a[2] * a[2];
}
function v3multiplyScalar(out, a, s) {
out[0] = a[0] * s;
out[1] = a[1] * s;
out[2] = a[2] * s;
}
function v3negate(out, a) {
out[0] = -a[0];
out[1] = -a[1];
out[2] = -a[2];
}
function v3new(array) {
return new Float32Array(array || 3);
}
function v3normalize(out, a) {
var length2 = v3length2(a);
if (length2 == 0) {
out[0] = a[0];
out[1] = a[1];
out[2] = a[2];
} else {
v3multiplyScalar(out, a, 1 / Math.sqrt(length2));
}
}
function v3sub(out, a, b) {
out[0] = a[0] - b[0];
out[1] = a[1] - b[1];
out[2] = a[2] - b[2];
}
function v3subScalar(out, a, s) {
out[0] = a[0] - s;
out[1] = a[1] - s;
out[2] = a[2] - s;
}
self.func = function func2(e, callback) {
var a = e.data.args;
var p = e.data.params;
if (a && p) {
var SurfClass = p.type === "av" ? AVSurface : EDTSurface;
var surf = new SurfClass(a.coordList, a.radiusList, a.indexList);
var sd = surf.getSurface(p.type, p.probeRadius, p.scaleFactor, p.cutoff, true, p.smooth, p.contour);
var transferList = [sd.position.buffer, sd.index.buffer];
if (sd.normal) {
transferList.push(sd.normal.buffer);
}
if (sd.atomindex) {
transferList.push(sd.atomindex.buffer);
}
var data = {
sd,
p
};
callback(data, transferList);
}
};
self.onmessage = function onmessage(e) {
var name = e.data.__name;
var postId = e.data.__postId;
if (name === void 0) {
console.error("message __name undefined");
} else if (self.func === void 0) {
console.error("worker func undefined", name);
} else {
var callback = /* @__PURE__ */ __name(function(aMessage, transferList) {
aMessage = aMessage || {};
if (postId !== void 0) {
aMessage.__postId = postId;
}
try {
self.postMessage(aMessage, transferList);
} catch (error) {
console.error("self.postMessage:", error);
self.postMessage(aMessage);
}
}, "callback");
self.func(e, callback);
}
};
Raw
av-surface.ts
/**
* @file AV Surface
* @author Fred Ludlow <fred.ludlow@gmail.com>
* @private
*/
import { getSurfaceGrid } from './surface-utils'
import { VolumeSurface } from './volume'
import { uniformArray } from '../math/array-utils'
import {
computeBoundingBox, v3multiplyScalar, v3cross, v3normalize
} from '../math/vector-utils'
import { defaults } from '../utils'
import { NumberArray } from '../types';
/**
* Modifed from SpatialHash
*
* Main differences are:
* - Optimized grid size to ensure we only ever need to look +/-1 cell
* - Aware of atomic radii and will only output atoms within rAtom + rExtra
* (see withinRadii method)
*
* (Uses rounding rather than bitshifting as consequence of arbitrary grid size)
* @class
* @param {Float32Array} atomsX - x coordinates
* @param {Float32Array} atomsY - y coordinates
* @param {Float32Array} atomsZ - z coordinates
* @param {Float32Array} atomsR - atom radii
* @param {Float32Array} min - xyz min coordinates
* @param {Float32Array} max - xyz max coordinates
* @param {Float} maxDistance - max distance
*/
export interface iAVHash {
neighbourListLength: number
withinRadii: (x: number, y: number, z: number, rExtra: number, out: Int32Array) => void
}
function makeAVHash (atomsX: Float32Array, atomsY: Float32Array, atomsZ: Float32Array, atomsR: Float32Array, min: Float32Array, max: Float32Array, maxDistance: number): iAVHash {
maxDistance = Math.max(0.1, maxDistance) // Avoid maxDistance of zero, see #802
var nAtoms = atomsX.length
var minX = min[ 0 ]
var minY = min[ 1 ]
var minZ = min[ 2 ]
var maxX = max[ 0 ]
var maxY = max[ 1 ]
var maxZ = max[ 2 ]
function hashFunc (w: number, minW: number) {
return Math.floor((w - minW) / maxDistance)
}
var iDim = hashFunc(maxX, minX) + 1
var jDim = hashFunc(maxY, minY) + 1
var kDim = hashFunc(maxZ, minZ) + 1
var nCells = iDim * jDim * kDim
var jkDim = jDim * kDim
/* Get cellID for cartesian x,y,z */
var cellID = function (x: number, y: number, z: number) {
return (((hashFunc(x, minX) * jDim) + hashFunc(y, minY)) * kDim) + hashFunc(z, minZ)
}
/* Initial building, could probably be optimized further */
var preHash = [] // preHash[ cellID ] = [ atomId1, atomId2 ];
for (var i = 0; i < nAtoms; i++) {
var cid = cellID(atomsX[ i ], atomsY[ i ], atomsZ[ i ])
if (preHash[ cid ] === undefined) {
preHash[ cid ] = [ i ]
} else {
preHash[ cid ].push(i)
}
}
var cellOffsets = new Uint32Array(nCells)
var cellLengths = new Uint16Array(nCells)
var data = new Uint32Array(nAtoms)
var offset = 0
var maxCellLength = 0
for (i = 0; i < nCells; i++) {
var start = cellOffsets[ i ] = offset
var subArray = preHash[ i ]
if (subArray !== undefined) {
for (var j = 0; j < subArray.length; j++) {
data[ offset ] = subArray[ j ]
offset++
}
}
var cellLength = offset - start
cellLengths[ i ] = cellLength
if (cellLength > maxCellLength) { maxCellLength = cellLength }
}
// Maximum number of neighbours we could ever produce (27 adjacent cells of equal population)
const neighbourListLength = (27 * maxCellLength) + 1
/**
* Populate the supplied out array with atom indices that are within rAtom + rExtra
* of x,y,z
*
* -1 in out array indicates the end of the list
*
* @param {Float} x - x coordinate
* @param {Float} y - y coordinate
* @param {Float} z - z coordinate
* @param {Float} rExtra - additional radius
* @param {Float32Array} out - pre-allocated output array
* @return {undefined}
*/
const withinRadii = function (x: number, y: number, z: number, rExtra: number, out: Int32Array) {
var outIdx = 0
var nearI = hashFunc(x, minX)
var nearJ = hashFunc(y, minY)
var nearK = hashFunc(z, minZ)
var loI = Math.max(0, nearI - 1)
var loJ = Math.max(0, nearJ - 1)
var loK = Math.max(0, nearK - 1)
var hiI = Math.min(iDim, nearI + 2)
var hiJ = Math.min(jDim, nearJ + 2)
var hiK = Math.min(kDim, nearK + 2)
for (var i = loI; i < hiI; ++i) {
var iOffset = i * jkDim
for (var j = loJ; j < hiJ; ++j) {
var jOffset = j * kDim
for (var k = loK; k < hiK; ++k) {
var cid = iOffset + jOffset + k
var cellStart = cellOffsets[ cid ]
var cellEnd = cellStart + cellLengths[ cid ]
for (var dataIndex = cellStart; dataIndex < cellEnd; dataIndex++) {
var atomIndex = data[ dataIndex ]
var dx = atomsX[ atomIndex ] - x
var dy = atomsY[ atomIndex ] - y
var dz = atomsZ[ atomIndex ] - z
var rSum = atomsR[ atomIndex ] + rExtra
if ((dx * dx + dy * dy + dz * dz) <= (rSum * rSum)) {
out[ outIdx++ ] = data[ dataIndex ]
}
}
}
}
}
// Add terminator
out[ outIdx ] = -1
}
return {
neighbourListLength: neighbourListLength,
withinRadii: withinRadii
}
}
interface AVSurface {
getSurface: (type: string, probeRadius: number, scaleFactor: number, cutoff: number, setAtomID: boolean, smooth: number, contour: boolean) => any
}
function AVSurface (this: AVSurface, coordList: Float32Array, radiusList: Float32Array, indexList: Uint16Array|Uint32Array) {
// Field generation method adapted from AstexViewer (Mike Hartshorn)
// by Fred Ludlow.
// Other parts based heavily on NGL (Alexander Rose) EDT Surface class
//
// Should work as a drop-in alternative to EDTSurface (though some of
// the EDT paramters are not relevant in this method).
const nAtoms = radiusList.length
const x = new Float32Array(nAtoms)
const y = new Float32Array(nAtoms)
const z = new Float32Array(nAtoms)
for (let i = 0; i < nAtoms; i++) {
const ci = 3 * i
x[ i ] = coordList[ ci ]
y[ i ] = coordList[ ci + 1 ]
z[ i ] = coordList[ ci + 2 ]
}
let bbox = computeBoundingBox(coordList)
if (coordList.length === 0) {
bbox[ 0 ].set([ 0, 0, 0 ])
bbox[ 1 ].set([ 0, 0, 0 ])
}
const min = bbox[0]
const max = bbox[1]
let r: Float32Array, r2: Float32Array // Atom positions, expanded radii (squared)
let maxRadius: number
// Parameters
let probeRadius: number, scaleFactor: number, setAtomID: boolean, probePositions: number
// Cache last value for obscured test
let lastClip = -1
// Grid params
let dim: Float32Array, matrix: Float32Array, grid: NumberArray, atomIndex: Int32Array
// grid indices -> xyz coords
let gridx: Float32Array, gridy: Float32Array, gridz: Float32Array
// Lookup tables:
let sinTable: Float32Array, cosTable: Float32Array
// Spatial Hash
let hash: iAVHash
// Neighbour array to be filled by hash
let neighbours: Int32Array
// Vectors for Torus Projection
const atob = new Float32Array([ 0.0, 0.0, 0.0 ])
const mid = new Float32Array([ 0.0, 0.0, 0.0 ])
const n1 = new Float32Array([ 0.0, 0.0, 0.0 ])
const n2 = new Float32Array([ 0.0, 0.0, 0.0 ])
let ngTorus: number
function init (_probeRadius?: number, _scaleFactor?: number, _setAtomID?: boolean, _probePositions?: number) {
probeRadius = defaults(_probeRadius, 1.4)
scaleFactor = defaults(_scaleFactor, 2.0)
setAtomID = defaults(_setAtomID, true)
probePositions = defaults(_probePositions, 30)
r = new Float32Array(nAtoms)
r2 = new Float32Array(nAtoms)
for (let i = 0; i < r.length; ++i) {
var rExt = radiusList[ i ] + probeRadius
r[ i ] = rExt
r2[ i ] = rExt * rExt
}
maxRadius = 0
for (let j = 0; j < r.length; ++j) {
if (r[ j ] > maxRadius) maxRadius = r[ j ]
}
initializeGrid()
initializeAngleTables()
initializeHash()
lastClip = -1
}
function fillGridDim (a: Float32Array, start: number, step: number) {
for (let i = 0; i < a.length; i++) {
a[i] = start + (step * i)
}
}
function initializeGrid () {
const surfGrid = getSurfaceGrid(
min, max, maxRadius, scaleFactor, 0.0
)
scaleFactor = surfGrid.scaleFactor
dim = surfGrid.dim
matrix = surfGrid.matrix
ngTorus = Math.max(5, 2 + Math.floor(probeRadius * scaleFactor))
grid = uniformArray(dim[0] * dim[1] * dim[2], -1001.0)
atomIndex = new Int32Array(grid.length)
gridx = new Float32Array(dim[0])
gridy = new Float32Array(dim[1])
gridz = new Float32Array(dim[2])
fillGridDim(gridx, min[0], 1 / scaleFactor)
fillGridDim(gridy, min[1], 1 / scaleFactor)
fillGridDim(gridz, min[2], 1 / scaleFactor)
}
function initializeAngleTables () {
var theta = 0.0
var step = 2 * Math.PI / probePositions
cosTable = new Float32Array(probePositions)
sinTable = new Float32Array(probePositions)
for (var i = 0; i < probePositions; i++) {
cosTable[ i ] = Math.cos(theta)
sinTable[ i ] = Math.sin(theta)
theta += step
}
}
function initializeHash () {
hash = makeAVHash(x, y, z, r, min, max, 2.01 * maxRadius)
neighbours = new Int32Array(hash.neighbourListLength)
}
function obscured (x: number, y: number, z: number, a: number, b: number) {
// Is the point at x,y,z obscured by any of the atoms
// specifeid by indices in neighbours. Ignore indices
// a and b (these are the relevant atoms in projectPoints/Torii)
// Cache the last clipped atom (as very often the same one in
// subsequent calls)
let ai: number
if (lastClip !== -1) {
ai = lastClip
if (ai !== a && ai !== b && singleAtomObscures(ai, x, y, z)) {
return ai
} else {
lastClip = -1
}
}
var ni = 0
ai = neighbours[ ni ]
while (ai >= 0) {
if (ai !== a && ai !== b && singleAtomObscures(ai, x, y, z)) {
lastClip = ai
return ai
}
ai = neighbours[ ++ni ]
}
lastClip = -1
return -1
}
function singleAtomObscures (ai: number, x: number, y: number, z: number) {
var ci = 3 * ai
var ra2 = r2[ ai ]
var dx = coordList[ ci ] - x
var dy = coordList[ ci + 1 ] - y
var dz = coordList[ ci + 2 ] - z
var d2 = dx * dx + dy * dy + dz * dz
return d2 < ra2
}
function projectPoints () {
// For each atom:
// Iterate over a subsection of the grid, for each point:
// If current value < 0.0, unvisited, set positive
//
// In any case: Project this point onto surface of the atomic sphere
// If this projected point is not obscured by any other atom
// Calcualte delta distance and set grid value to minimum of
// itself and delta
// Should we alias frequently accessed closure variables??
// Assume JS engine capable of optimizing this
// anyway...
for (var i = 0; i < nAtoms; i++) {
var ax = x[ i ]
var ay = y[ i ]
var az = z[ i ]
var ar = r[ i ]
var ar2 = r2[ i ]
hash.withinRadii(ax, ay, az, ar, neighbours)
// Number of grid points, round this up...
var ng = Math.ceil(ar * scaleFactor)
// Center of the atom, mapped to grid points (take floor)
var iax = Math.floor(scaleFactor * (ax - min[ 0 ]))
var iay = Math.floor(scaleFactor * (ay - min[ 1 ]))
var iaz = Math.floor(scaleFactor * (az - min[ 2 ]))
// Extents of grid to consider for this atom
var minx = Math.max(0, iax - ng)
var miny = Math.max(0, iay - ng)
var minz = Math.max(0, iaz - ng)
// Add two to these points:
// - iax are floor'd values so this ensures coverage
// - these are loop limits (exclusive)
var maxx = Math.min(dim[ 0 ], iax + ng + 2)
var maxy = Math.min(dim[ 1 ], iay + ng + 2)
var maxz = Math.min(dim[ 2 ], iaz + ng + 2)
for (var ix = minx; ix < maxx; ix++) {
var dx = gridx[ ix ] - ax
var xoffset = dim[ 1 ] * dim[ 2 ] * ix
for (var iy = miny; iy < maxy; iy++) {
var dy = gridy[ iy ] - ay
var dxy2 = dx * dx + dy * dy
var xyoffset = xoffset + dim[ 2 ] * iy
for (var iz = minz; iz < maxz; iz++) {
var dz = gridz[ iz ] - az
var d2 = dxy2 + dz * dz
if (d2 < ar2) {
var idx = iz + xyoffset
if (grid[idx] < 0.0) {
// Unvisited, make positive
grid[ idx ] = -grid[ idx ]
}
// Project on to the surface of the sphere
// sp is the projected point ( dx, dy, dz ) * ( ra / d )
var d = Math.sqrt(d2)
var ap = ar / d
var spx = dx * ap
var spy = dy * ap
var spz = dz * ap
spx += ax
spy += ay
spz += az
if (obscured(spx, spy, spz, i, -1) === -1) {
var dd = ar - d
if (dd < grid[ idx ]) {
grid[ idx ] = dd
if (setAtomID) atomIndex[ idx ] = i
}
}
}
}
}
}
}
}
function projectTorii () {
for (var i = 0; i < nAtoms; i++) {
hash.withinRadii(x[ i ], y[ i ], z[ i ], r[ i ], neighbours)
var ia = 0
var ni = neighbours[ ia ]
while (ni >= 0) {
if (i < ni) {
projectTorus(i, ni)
}
ni = neighbours[ ++ia ]
}
}
}
function projectTorus (a: number, b: number) {
var r1 = r[ a ]
var r2 = r[ b ]
var dx = atob[ 0 ] = x[ b ] - x[ a ]
var dy = atob[ 1 ] = y[ b ] - y[ a ]
var dz = atob[ 2 ] = z[ b ] - z[ a ]
var d2 = dx * dx + dy * dy + dz * dz
// This check now redundant as already done in AVHash.withinRadii
// if( d2 > (( r1 + r2 ) * ( r1 + r2 )) ){ return; }
var d = Math.sqrt(d2)
// Find angle between a->b vector and the circle
// of their intersection by cosine rule
var cosA = (r1 * r1 + d * d - r2 * r2) / (2.0 * r1 * d)
// distance along a->b at intersection
var dmp = r1 * cosA
v3normalize(atob, atob)
// Create normal to line
normalToLine(n1 as any, atob)
v3normalize(n1, n1)
// Cross together for second normal vector
v3cross(n2, atob, n1)
v3normalize(n2, n2)
// r is radius of circle of intersection
var rInt = Math.sqrt(r1 * r1 - dmp * dmp)
v3multiplyScalar(n1, n1, rInt)
v3multiplyScalar(n2, n2, rInt)
v3multiplyScalar(atob, atob, dmp)
mid[ 0 ] = atob[ 0 ] + x[ a ]
mid[ 1 ] = atob[ 1 ] + y[ a ]
mid[ 2 ] = atob[ 2 ] + z[ a ]
lastClip = -1
var ng = ngTorus
for (var i = 0; i < probePositions; i++) {
var cost = cosTable[ i ]
var sint = sinTable[ i ]
var px = mid[ 0 ] + cost * n1[ 0 ] + sint * n2[ 0 ]
var py = mid[ 1 ] + cost * n1[ 1 ] + sint * n2[ 1 ]
var pz = mid[ 2 ] + cost * n1[ 2 ] + sint * n2[ 2 ]
if (obscured(px, py, pz, a, b) === -1) {
// As above, iterate over our grid...
// px, py, pz in grid coords
var iax = Math.floor(scaleFactor * (px - min[ 0 ]))
var iay = Math.floor(scaleFactor * (py - min[ 1 ]))
var iaz = Math.floor(scaleFactor * (pz - min[ 2 ]))
var minx = Math.max(0, iax - ng)
var miny = Math.max(0, iay - ng)
var minz = Math.max(0, iaz - ng)
var maxx = Math.min(dim[ 0 ], iax + ng + 2)
var maxy = Math.min(dim[ 1 ], iay + ng + 2)
var maxz = Math.min(dim[ 2 ], iaz + ng + 2)
for (var ix = minx; ix < maxx; ix++) {
dx = px - gridx[ ix ]
var xoffset = dim[ 1 ] * dim[ 2 ] * ix
for (var iy = miny; iy < maxy; iy++) {
dy = py - gridy[ iy ]
var dxy2 = dx * dx + dy * dy
var xyoffset = xoffset + dim[ 2 ] * iy
for (var iz = minz; iz < maxz; iz++) {
dz = pz - gridz[ iz ]
d2 = dxy2 + dz * dz
var idx = iz + xyoffset
var current = grid[ idx ]
if (current > 0.0 && d2 < (current * current)) {
grid[ idx ] = Math.sqrt(d2)
if (setAtomID) {
// Is this grid point closer to a or b?
// Take dot product of atob and gridpoint->p (dx, dy, dz)
const dp = dx * atob[ 0 ] + dy * atob [ 1 ] + dz * atob[ 2 ]
atomIndex[ idx ] = dp < 0.0 ? b : a
}
}
}
}
}
}
}
}
function normalToLine (out: Int32Array, p: Float32Array) {
out[ 0 ] = out[ 1 ] = out[ 2 ] = 1.0
if (p[ 0 ] !== 0) {
out[ 0 ] = (p[ 1 ] + p[ 2 ]) / -p[ 0 ]
} else if (p[ 1 ] !== 0) {
out[ 1 ] = (p[ 0 ] + p[ 2 ]) / -p[ 1 ]
} else if (p[ 2 ] !== 0) {
out[ 2 ] = (p[ 0 ] + p[ 1 ]) / -p[ 2 ]
}
return out
}
function fixNegatives () {
for (var i = 0; i < grid.length; i++) {
if (grid[ i ] < 0) grid[ i ] = 0
}
}
function fixAtomIDs () {
for (var i = 0; i < atomIndex.length; i++) {
atomIndex[ i ] = indexList[ atomIndex[ i ] ]
}
}
function getVolume (probeRadius: number, scaleFactor: number, setAtomID: boolean) {
// Basic steps are:
// 1) Initialize
// 2) Project points
// 3) Project torii
console.time('AVSurface.getVolume')
console.time('AVSurface.init')
init(probeRadius, scaleFactor, setAtomID)
console.timeEnd('AVSurface.init')
console.time('AVSurface.projectPoints')
projectPoints()
console.timeEnd('AVSurface.projectPoints')
console.time('AVSurface.projectTorii')
projectTorii()
console.timeEnd('AVSurface.projectTorii')
fixNegatives()
fixAtomIDs()
console.timeEnd('AVSurface.getVolume')
}
this.getSurface = function (type: string, probeRadius: number, scaleFactor: number, cutoff: number, setAtomID: boolean, smooth: number, contour: boolean) {
// type and cutoff left in for compatibility with EDTSurface.getSurface
// function signature
getVolume(probeRadius, scaleFactor, setAtomID)
var volsurf = new (VolumeSurface as any)(
grid, dim[ 2 ], dim[ 1 ], dim[ 0 ], atomIndex
) as VolumeSurface
return volsurf.getSurface!(probeRadius, false, undefined, matrix, contour)
}
}
Object.assign(AVSurface, {__deps: [
getSurfaceGrid, VolumeSurface, uniformArray, computeBoundingBox,
v3multiplyScalar, v3cross, v3normalize,
makeAVHash,
defaults
]})
export { AVSurface, makeAVHash }
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment