import {asin, atan2, cos, degrees, epsilon, epsilon2, radians, sin, sqrt} from "./math.js"; import noop from "./noop.js"; import stream from "./stream.js"; var W0, W1, X0, Y0, Z0, X1, Y1, Z1, X2, Y2, Z2, lambda00, phi00, // first point x0, y0, z0; // previous point var centroidStream = { sphere: noop, point: centroidPoint, lineStart: centroidLineStart, lineEnd: centroidLineEnd, polygonStart: function() { centroidStream.lineStart = centroidRingStart; centroidStream.lineEnd = centroidRingEnd; }, polygonEnd: function() { centroidStream.lineStart = centroidLineStart; centroidStream.lineEnd = centroidLineEnd; } }; // Arithmetic mean of Cartesian vectors. function centroidPoint(lambda, phi) { lambda *= radians, phi *= radians; var cosPhi = cos(phi); centroidPointCartesian(cosPhi * cos(lambda), cosPhi * sin(lambda), sin(phi)); } function centroidPointCartesian(x, y, z) { ++W0; X0 += (x - X0) / W0; Y0 += (y - Y0) / W0; Z0 += (z - Z0) / W0; } function centroidLineStart() { centroidStream.point = centroidLinePointFirst; } function centroidLinePointFirst(lambda, phi) { lambda *= radians, phi *= radians; var cosPhi = cos(phi); x0 = cosPhi * cos(lambda); y0 = cosPhi * sin(lambda); z0 = sin(phi); centroidStream.point = centroidLinePoint; centroidPointCartesian(x0, y0, z0); } function centroidLinePoint(lambda, phi) { lambda *= radians, phi *= radians; var cosPhi = cos(phi), x = cosPhi * cos(lambda), y = cosPhi * sin(lambda), z = sin(phi), w = atan2(sqrt((w = y0 * z - z0 * y) * w + (w = z0 * x - x0 * z) * w + (w = x0 * y - y0 * x) * w), x0 * x + y0 * y + z0 * z); W1 += w; X1 += w * (x0 + (x0 = x)); Y1 += w * (y0 + (y0 = y)); Z1 += w * (z0 + (z0 = z)); centroidPointCartesian(x0, y0, z0); } function centroidLineEnd() { centroidStream.point = centroidPoint; } // See J. E. Brock, The Inertia Tensor for a Spherical Triangle, // J. Applied Mechanics 42, 239 (1975). function centroidRingStart() { centroidStream.point = centroidRingPointFirst; } function centroidRingEnd() { centroidRingPoint(lambda00, phi00); centroidStream.point = centroidPoint; } function centroidRingPointFirst(lambda, phi) { lambda00 = lambda, phi00 = phi; lambda *= radians, phi *= radians; centroidStream.point = centroidRingPoint; var cosPhi = cos(phi); x0 = cosPhi * cos(lambda); y0 = cosPhi * sin(lambda); z0 = sin(phi); centroidPointCartesian(x0, y0, z0); } function centroidRingPoint(lambda, phi) { lambda *= radians, phi *= radians; var cosPhi = cos(phi), x = cosPhi * cos(lambda), y = cosPhi * sin(lambda), z = sin(phi), cx = y0 * z - z0 * y, cy = z0 * x - x0 * z, cz = x0 * y - y0 * x, m = sqrt(cx * cx + cy * cy + cz * cz), w = asin(m), // line weight = angle v = m && -w / m; // area weight multiplier X2 += v * cx; Y2 += v * cy; Z2 += v * cz; W1 += w; X1 += w * (x0 + (x0 = x)); Y1 += w * (y0 + (y0 = y)); Z1 += w * (z0 + (z0 = z)); centroidPointCartesian(x0, y0, z0); } export default function(object) { W0 = W1 = X0 = Y0 = Z0 = X1 = Y1 = Z1 = X2 = Y2 = Z2 = 0; stream(object, centroidStream); var x = X2, y = Y2, z = Z2, m = x * x + y * y + z * z; // If the area-weighted ccentroid is undefined, fall back to length-weighted ccentroid. if (m < epsilon2) { x = X1, y = Y1, z = Z1; // If the feature has zero length, fall back to arithmetic mean of point vectors. if (W1 < epsilon) x = X0, y = Y0, z = Z0; m = x * x + y * y + z * z; // If the feature still has an undefined ccentroid, then return. if (m < epsilon2) return [NaN, NaN]; } return [atan2(y, x) * degrees, asin(z / sqrt(m)) * degrees]; }