utils/TerrainWorker.worker.jsjavascript
"use strict";
//
//Terrain worker
//

function binarySearchFast(arr, x) {
    let start = 0,
        end = arr.length - 1;
    while (start <= end) {
        let k = Math.floor((start + end) * 0.5);
        if (Math.abs(arr[k] - x) < 1e-3) return k;
        else if (arr[k] < x) start = k + 1;
        else end = k - 1;
    }
    return -1;
}

function checkNoDataValue(noDataValues, value) {
    return binarySearchFast(noDataValues, value) !== -1;
}

var Vec3 = function (x, y, z) {
    this.x = x;
    this.y = y;
    this.z = z;
};

var doubleToTwoFloats = function (v, high, low) {
    let x = v.x,
        y = v.y,
        z = v.z;

    var doubleHigh;

    if (x >= 0.0) {
        doubleHigh = Math.floor(x / 65536.0) * 65536.0;
        high.x = Math.fround(doubleHigh);
        low.x = Math.fround(x - doubleHigh);
    } else {
        doubleHigh = Math.floor(-x / 65536.0) * 65536.0;
        high.x = Math.fround(-doubleHigh);
        low.x = Math.fround(x + doubleHigh);
    }

    if (y >= 0.0) {
        doubleHigh = Math.floor(y / 65536.0) * 65536.0;
        high.y = Math.fround(doubleHigh);
        low.y = Math.fround(y - doubleHigh);
    } else {
        doubleHigh = Math.floor(-y / 65536.0) * 65536.0;
        high.y = Math.fround(-doubleHigh);
        low.y = Math.fround(y + doubleHigh);
    }

    if (z >= 0.0) {
        doubleHigh = Math.floor(z / 65536.0) * 65536.0;
        high.z = Math.fround(doubleHigh);
        low.z = Math.fround(z - doubleHigh);
    } else {
        doubleHigh = Math.floor(-z / 65536.0) * 65536.0;
        high.z = Math.fround(-doubleHigh);
        low.z = Math.fround(z + doubleHigh);
    }
};

Vec3.prototype.sub = function (v) {
    return new Vec3(this.x - v.x, this.y - v.y, this.z - v.z);
};

Vec3.prototype.add = function (v) {
    return new Vec3(this.x + v.x, this.y + v.y, this.z + v.z);
};

Vec3.prototype.cross = function (v) {
    return new Vec3(this.y * v.z - this.z * v.y, this.z * v.x - this.x * v.z, this.x * v.y - this.y * v.x);
};

Vec3.prototype.normalize = function () {
    var x = this.x,
        y = this.y,
        z = this.z;
    var length = 1.0 / Math.sqrt(x * x + y * y + z * z);
    this.x = x * length;
    this.y = y * length;
    this.z = z * length;
    return this;
};

Vec3.prototype.distance = function (v) {
    return this.sub(v).length();
};

Vec3.prototype.length = function () {
    return Math.sqrt(this.x * this.x + this.y * this.y + this.z * this.z);
};

var blerp = function (x, y, fQ11, fQ21, fQ12, fQ22) {
    return fQ11 * (1.0 - x) * (1.0 - y) + fQ21 * x * (1.0 - y) + fQ12 * (1.0 - x) * y + fQ22 * x * y;
};

var _tempVec = new Vec3(0.0, 0.0, 0.0);
var _rtcTempVec = new Vec3(0.0, 0.0, 0.0);

var _tempHigh = new Vec3(0.0, 0.0, 0.0),
    _tempLow = new Vec3(0.0, 0.0, 0.0);

self.onmessage = function (e) {
    var elevations = e.data.elevations,
        this_plainVertices = e.data.this_plainVertices,
        this_plainNormals = e.data.this_plainNormals,
        this_normalMapVertices = e.data.this_normalMapVertices,
        this_normalMapNormals = e.data.this_normalMapNormals,
        heightFactor = e.data.heightFactor,
        gridSize = e.data.gridSize,
        noDataValues = e.data.noDataValues,
        id = e.data.id;

    let rtc_x = e.data.relativeCenter[0],
        rtc_y = e.data.relativeCenter[1],
        rtc_z = e.data.relativeCenter[2];

    var xmin = 549755748352.0,
        xmax = -549755748352.0,
        ymin = 549755748352.0,
        ymax = -549755748352.0,
        zmin = 549755748352.0,
        zmax = -549755748352.0;

    // Squared geocentric distance of the same vertices. The axis aligned bounds above cannot tell
    // relief from segment size, because the box of a tilted patch is far larger than the patch, so
    // the radius range is accumulated here instead. Kept squared to avoid a root per vertex.
    var r2min = Infinity,
        r2max = -Infinity;

    const fileGridSize = Math.sqrt(elevations.length) - 1;

    const fileGridSize_one = fileGridSize + 1;
    const fileGridSize_one_x2 = fileGridSize_one * fileGridSize_one;
    const tgs = gridSize;
    const dg = fileGridSize / tgs;
    const gs = tgs + 1;
    const hf = heightFactor;

    var vInd = 0,
        noDataInd = 0;

    var gsgs3 = gs * gs * 3;

    var terrainVertices = new Float64Array(gsgs3),
        terrainVerticesHigh = new Float32Array(gsgs3),
        terrainVerticesLow = new Float32Array(gsgs3),
        noDataVertices = new Uint8Array(gs * gs);

    var normalMapNormals, normalMapVertices;

    var nv = this_normalMapVertices,
        nn = this_normalMapNormals;

    if (fileGridSize >= tgs) {
        normalMapNormals = new Float32Array(fileGridSize_one_x2 * 3);
        normalMapVertices = new Float64Array(fileGridSize_one_x2 * 3);

        for (let k = 0; k < fileGridSize_one_x2; k++) {
            var j = k % fileGridSize_one,
                i = ~~(k / fileGridSize_one);

            //
            // V0
            //
            var hInd0 = k;
            var vInd0 = hInd0 * 3;
            var currElv = elevations[hInd0];
            if (checkNoDataValue(noDataValues, currElv)) {
                currElv = 0.0;
            }
            var h0 = hf * currElv;
            var v0 = new Vec3(
                rtc_x + nv[vInd0] + h0 * nn[vInd0],
                rtc_y + nv[vInd0 + 1] + h0 * nn[vInd0 + 1],
                rtc_z + nv[vInd0 + 2] + h0 * nn[vInd0 + 2]
            );

            normalMapVertices[vInd0] = v0.x - rtc_x;
            normalMapVertices[vInd0 + 1] = v0.y - rtc_y;
            normalMapVertices[vInd0 + 2] = v0.z - rtc_z;

            //
            // The vertex goes into screen buffer
            if (i % dg === 0 && j % dg === 0) {
                let currVert = new Vec3(rtc_x + nv[vInd0], rtc_y + nv[vInd0 + 1], rtc_z + nv[vInd0 + 2]);
                let nextVert = new Vec3(rtc_x + nv[vInd0 + 3], rtc_y + nv[vInd0 + 4], rtc_z + nv[vInd0 + 5]);

                let nextElv = elevations[hInd0 + 1];
                if (checkNoDataValue(noDataValues, nextElv)) {
                    nextElv = 0.0;
                }

                let eps = false;
                if (noDataValues.length === 0) {
                    let step = currVert.distance(nextVert);
                    let deltaElv = Math.abs(currElv - nextElv);
                    eps = deltaElv / step > 10.0 || currElv < -5000;
                }

                if (eps) {
                    noDataVertices[noDataInd] = 1;
                } else {
                    noDataVertices[noDataInd] = 0;
                    if (v0.x < xmin) xmin = v0.x;
                    if (v0.x > xmax) xmax = v0.x;
                    if (v0.y < ymin) ymin = v0.y;
                    if (v0.y > ymax) ymax = v0.y;
                    if (v0.z < zmin) zmin = v0.z;
                    if (v0.z > zmax) zmax = v0.z;

                    let r2_v0 = v0.x * v0.x + v0.y * v0.y + v0.z * v0.z;
                    if (r2_v0 < r2min) r2min = r2_v0;
                    if (r2_v0 > r2max) r2max = r2_v0;
                }

                let rtc_v0 = new Vec3(v0.x - rtc_x, v0.y - rtc_y, v0.z - rtc_z);
                doubleToTwoFloats(rtc_v0, _tempHigh, _tempLow);

                terrainVerticesHigh[vInd] = _tempHigh.x;
                terrainVerticesLow[vInd] = _tempLow.x;
                terrainVertices[vInd++] = rtc_v0.x;

                terrainVerticesHigh[vInd] = _tempHigh.y;
                terrainVerticesLow[vInd] = _tempLow.y;
                terrainVertices[vInd++] = rtc_v0.y;

                terrainVerticesHigh[vInd] = _tempHigh.z;
                terrainVerticesLow[vInd] = _tempLow.z;
                terrainVertices[vInd++] = rtc_v0.z;

                noDataInd++;
            }

            if (i !== fileGridSize && j !== fileGridSize) {
                //
                //  V1
                //
                var hInd1 = k + 1;
                var vInd1 = hInd1 * 3;
                var elv = elevations[hInd1];
                if (checkNoDataValue(noDataValues, elv)) {
                    elv = 0.0;
                }
                var h1 = hf * elv;
                var v1 = new Vec3(
                    rtc_x + nv[vInd1] + h1 * nn[vInd1],
                    rtc_y + nv[vInd1 + 1] + h1 * nn[vInd1 + 1],
                    rtc_z + nv[vInd1 + 2] + h1 * nn[vInd1 + 2]
                );

                normalMapVertices[vInd1] = v1.x - rtc_x;
                normalMapVertices[vInd1 + 1] = v1.y - rtc_y;
                normalMapVertices[vInd1 + 2] = v1.z - rtc_z;

                //
                //  V2
                //
                var hInd2 = k + fileGridSize_one;
                var vInd2 = hInd2 * 3;
                elv = elevations[hInd2];
                if (checkNoDataValue(noDataValues, elv)) {
                    elv = 0.0;
                }
                var h2 = hf * elv;
                var v2 = new Vec3(
                    rtc_x + nv[vInd2] + h2 * nn[vInd2],
                    rtc_y + nv[vInd2 + 1] + h2 * nn[vInd2 + 1],
                    rtc_z + nv[vInd2 + 2] + h2 * nn[vInd2 + 2]
                );

                normalMapVertices[vInd2] = v2.x - rtc_x;
                normalMapVertices[vInd2 + 1] = v2.y - rtc_y;
                normalMapVertices[vInd2 + 2] = v2.z - rtc_z;

                //
                //  V3
                //
                var hInd3 = k + fileGridSize_one + 1;
                var vInd3 = hInd3 * 3;
                elv = elevations[hInd3];
                if (checkNoDataValue(noDataValues, elv)) {
                    elv = 0.0;
                }
                var h3 = hf * elv;
                var v3 = new Vec3(
                    rtc_x + nv[vInd3] + h3 * nn[vInd3],
                    rtc_y + nv[vInd3 + 1] + h3 * nn[vInd3 + 1],
                    rtc_z + nv[vInd3 + 2] + h3 * nn[vInd3 + 2]
                );

                normalMapVertices[vInd3] = v3.x - rtc_x;
                normalMapVertices[vInd3 + 1] = v3.y - rtc_y;
                normalMapVertices[vInd3 + 2] = v3.z - rtc_z;

                //
                // Normal
                //
                var e10 = v1.sub(v0),
                    e20 = v2.sub(v0),
                    e30 = v3.sub(v0);

                var sw = e20.cross(e30).normalize();
                var ne = e30.cross(e10).normalize();
                var n0 = ne.add(sw).normalize();

                normalMapNormals[vInd0] += n0.x;
                normalMapNormals[vInd0 + 1] += n0.y;
                normalMapNormals[vInd0 + 2] += n0.z;

                normalMapNormals[vInd1] += ne.x;
                normalMapNormals[vInd1 + 1] += ne.y;
                normalMapNormals[vInd1 + 2] += ne.z;

                normalMapNormals[vInd2] += sw.x;
                normalMapNormals[vInd2 + 1] += sw.y;
                normalMapNormals[vInd2 + 2] += sw.z;

                normalMapNormals[vInd3] += n0.x;
                normalMapNormals[vInd3 + 1] += n0.y;
                normalMapNormals[vInd3 + 2] += n0.z;
            }
        }
    } else {
        normalMapNormals = new Float32Array(gsgs3);
        normalMapVertices = new Float64Array(gsgs3);

        var oneSize = tgs / fileGridSize;
        var gsgs = gsgs3 / 3;
        var fgsOne = fileGridSize + 1;

        for (let i = 0; i < gsgs; i++) {
            let ii = Math.floor(i / gs),
                ij = i % gs;

            let qii = ii % oneSize,
                qij = ij % oneSize;

            let hlt_ind = Math.floor(ii / oneSize) * fgsOne + Math.floor(ij / oneSize);

            if (ij === tgs) {
                hlt_ind -= 1;
                qij = oneSize;
            }

            if (ii === tgs) {
                hlt_ind -= fgsOne;
                qii = oneSize;
            }

            let hrt_ind = hlt_ind + 1,
                hlb_ind = hlt_ind + fgsOne,
                hrb_ind = hlb_ind + 1;

            let h_lt = elevations[hlt_ind],
                h_rt = elevations[hrt_ind],
                h_lb = elevations[hlb_ind],
                h_rb = elevations[hrb_ind];

            if (checkNoDataValue(noDataValues, h_lt)) {
                h_lt = 0.0;
            }

            if (checkNoDataValue(noDataValues, h_rt)) {
                h_rt = 0.0;
            }

            if (checkNoDataValue(noDataValues, h_lb)) {
                h_lb = 0.0;
            }

            if (checkNoDataValue(noDataValues, h_rb)) {
                h_rb = 0.0;
            }

            let hi = blerp(qij / oneSize, qii / oneSize, h_lt, h_rt, h_lb, h_rb);

            let i3 = i * 3;

            _rtcTempVec.x = this_plainVertices[i3] + hi * this_plainNormals[i3];
            _rtcTempVec.y = this_plainVertices[i3 + 1] + hi * this_plainNormals[i3 + 1];
            _rtcTempVec.z = this_plainVertices[i3 + 2] + hi * this_plainNormals[i3 + 2];

            doubleToTwoFloats(_rtcTempVec, _tempHigh, _tempLow);

            terrainVertices[i3] = _rtcTempVec.x;
            terrainVertices[i3 + 1] = _rtcTempVec.y;
            terrainVertices[i3 + 2] = _rtcTempVec.z;

            terrainVerticesHigh[i3] = _tempHigh.x;
            terrainVerticesHigh[i3 + 1] = _tempHigh.y;
            terrainVerticesHigh[i3 + 2] = _tempHigh.z;

            terrainVerticesLow[i3] = _tempLow.x;
            terrainVerticesLow[i3 + 1] = _tempLow.y;
            terrainVerticesLow[i3 + 2] = _tempLow.z;

            _tempVec.x = _rtcTempVec.x + rtc_x;
            _tempVec.y = _rtcTempVec.y + rtc_y;
            _tempVec.z = _rtcTempVec.z + rtc_z;

            if (_tempVec.x < xmin) xmin = _tempVec.x;
            if (_tempVec.x > xmax) xmax = _tempVec.x;
            if (_tempVec.y < ymin) ymin = _tempVec.y;
            if (_tempVec.y > ymax) ymax = _tempVec.y;
            if (_tempVec.z < zmin) zmin = _tempVec.z;
            if (_tempVec.z > zmax) zmax = _tempVec.z;

            let r2_tmp = _tempVec.x * _tempVec.x + _tempVec.y * _tempVec.y + _tempVec.z * _tempVec.z;
            if (r2_tmp < r2min) r2min = r2_tmp;
            if (r2_tmp > r2max) r2max = r2_tmp;
        }

        normalMapVertices.set(terrainVertices);

        for (let k = 0; k < gsgs; k++) {
            let j = k % gs,
                i = ~~(k / gs);

            if (i !== tgs && j !== tgs) {
                let v0ind = k * 3,
                    v1ind = v0ind + 3,
                    v2ind = v0ind + gs * 3,
                    v3ind = v2ind + 3;

                let v0 = new Vec3(
                        terrainVertices[v0ind] + rtc_x,
                        terrainVertices[v0ind + 1] + rtc_y,
                        terrainVertices[v0ind + 2] + rtc_z
                    ),
                    v1 = new Vec3(
                        terrainVertices[v1ind] + rtc_x,
                        terrainVertices[v1ind + 1] + rtc_y,
                        terrainVertices[v1ind + 2] + rtc_z
                    ),
                    v2 = new Vec3(
                        terrainVertices[v2ind] + rtc_x,
                        terrainVertices[v2ind + 1] + rtc_y,
                        terrainVertices[v2ind + 2] + rtc_z
                    ),
                    v3 = new Vec3(
                        terrainVertices[v3ind] + rtc_x,
                        terrainVertices[v3ind + 1] + rtc_y,
                        terrainVertices[v3ind + 2] + rtc_z
                    );

                let e10 = v1.sub(v0).normalize(),
                    e20 = v2.sub(v0).normalize(),
                    e30 = v3.sub(v0).normalize();

                let sw = e20.cross(e30).normalize();
                let ne = e30.cross(e10).normalize();
                let n0 = ne.add(sw).normalize();

                normalMapNormals[v0ind] += n0.x;
                normalMapNormals[v0ind + 1] += n0.y;
                normalMapNormals[v0ind + 2] += n0.z;

                normalMapNormals[v1ind] += ne.x;
                normalMapNormals[v1ind + 1] += ne.y;
                normalMapNormals[v1ind + 2] += ne.z;

                normalMapNormals[v2ind] += sw.x;
                normalMapNormals[v2ind + 1] += sw.y;
                normalMapNormals[v2ind + 2] += sw.z;

                normalMapNormals[v3ind] += n0.x;
                normalMapNormals[v3ind + 1] += n0.y;
                normalMapNormals[v3ind + 2] += n0.z;
            }
        }
    }

    self.postMessage(
        {
            id: id,
            normalMapNormals: normalMapNormals,
            normalMapVertices: normalMapVertices,
            terrainVertices: terrainVertices,
            terrainVerticesHigh: terrainVerticesHigh,
            terrainVerticesLow: terrainVerticesLow,
            noDataVertices: noDataVertices,
            bounds: [xmin, ymin, zmin, xmax, ymax, zmax],
            radiusRange: r2max >= r2min ? [Math.sqrt(r2min), Math.sqrt(r2max)] : null,
            relativeCenter: e.data.relativeCenter
        },
        [
            normalMapNormals.buffer,
            normalMapVertices.buffer,
            terrainVertices.buffer,
            terrainVerticesHigh.buffer,
            terrainVerticesLow.buffer,
            noDataVertices.buffer
        ]
    );
};