IfcOpenShell/src/ifcviewer-full/Measurement.cpp

/********************************************************************************
 *                                                                              *
 * This file is part of IfcOpenShell.                                           *
 *                                                                              *
 * IfcOpenShell is free software: you can redistribute it and/or modify         *
 * it under the terms of the Lesser GNU General Public License as published by  *
 * the Free Software Foundation, either version 3.0 of the License, or          *
 * (at your option) any later version.                                          *
 *                                                                              *
 * IfcOpenShell is distributed in the hope that it will be useful,              *
 * but WITHOUT ANY WARRANTY; without even the implied warranty of               *
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the                 *
 * Lesser GNU General Public License for more details.                          *
 *                                                                              *
 * You should have received a copy of the Lesser GNU General Public License     *
 * along with this program. If not, see <http://www.gnu.org/licenses/>.         *
 *                                                                              *
 ********************************************************************************/

#include "Measurement.h"

#include "ViewportWindow.h"

#include <QtGlobal>
#include <Eigen/Dense>

#include <algorithm>
#include <cmath>
#include <cstring>
#include <queue>
#include <unordered_map>
#include <unordered_set>
#include <vector>

namespace {

double meshLocalVolume(const ViewportWindow::MeshTriangles& tris) {
    // Signed tetrahedra from the origin: V = sum( a · (b × c) ) / 6.
    // Absolute value at the end so winding convention doesn't matter.
    double sum = 0.0;
    const size_t n = tris.indices.size();
    for (size_t i = 0; i + 2 < n; i += 3) {
        const uint32_t ia = tris.indices[i + 0];
        const uint32_t ib = tris.indices[i + 1];
        const uint32_t ic = tris.indices[i + 2];
        const float* a = &tris.positions[3 * ia];
        const float* b = &tris.positions[3 * ib];
        const float* c = &tris.positions[3 * ic];
        const double cx = double(b[1]) * c[2] - double(b[2]) * c[1];
        const double cy = double(b[2]) * c[0] - double(b[0]) * c[2];
        const double cz = double(b[0]) * c[1] - double(b[1]) * c[0];
        sum += double(a[0]) * cx + double(a[1]) * cy + double(a[2]) * cz;
    }
    return std::abs(sum) / 6.0;
}

double det3(const float M[16]) {
    // Upper-left 3x3 of a column-major 4x4: M[col * 4 + row].
    const double m00 = M[0],  m10 = M[1],  m20 = M[2];
    const double m01 = M[4],  m11 = M[5],  m21 = M[6];
    const double m02 = M[8],  m12 = M[9],  m22 = M[10];
    return m00 * (m11 * m22 - m12 * m21)
         - m01 * (m10 * m22 - m12 * m20)
         + m02 * (m10 * m21 - m11 * m20);
}

} // namespace

double volumeOfObjects(ViewportWindow& vp,
                       const std::vector<uint32_t>& object_ids) {
    if (object_ids.empty()) return 0.0;

    // Group selected instances by (model_id, mesh_id) so each unique mesh
    // is read back at most once per call.  Each entry stores the |det| of
    // every instance of that mesh in the request.
    std::unordered_map<uint64_t, std::vector<double>> by_mesh;
    by_mesh.reserve(object_ids.size());
    for (uint32_t oid : object_ids) {
        ViewportWindow::InstanceLookup lk;
        if (!vp.findInstance(oid, lk)) continue;
        const uint64_t key = (uint64_t(lk.model_id) << 32) | lk.mesh_id;
        by_mesh[key].push_back(std::abs(det3(lk.placement_transformation)));
    }

    double total = 0.0;
    ViewportWindow::MeshTriangles tris;
    for (const auto& [key, dets] : by_mesh) {
        const uint32_t model_id = uint32_t(key >> 32);
        const uint32_t mesh_id  = uint32_t(key & 0xffffffffu);
        if (!vp.readbackMeshTriangles(model_id, mesh_id, tris)) continue;
        const double v = meshLocalVolume(tris);
        for (double d : dets) total += v * d;
    }
    return total;
}

namespace {

// edge_key: undirected edge between two mesh-local vertex indices.
uint64_t edgeKey(uint32_t a, uint32_t b) {
    if (a > b) std::swap(a, b);
    return (uint64_t(a) << 32) | uint64_t(b);
}

// Triangle area = 0.5 * |(b - a) × (c - a)|.  Also returns the unit normal
// (zeroed for degenerate tris).
double triAreaAndNormal(const float* a, const float* b, const float* c,
                        float n_out[3]) {
    const double bax = double(b[0]) - a[0];
    const double bay = double(b[1]) - a[1];
    const double baz = double(b[2]) - a[2];
    const double cax = double(c[0]) - a[0];
    const double cay = double(c[1]) - a[1];
    const double caz = double(c[2]) - a[2];
    const double nx = bay * caz - baz * cay;
    const double ny = baz * cax - bax * caz;
    const double nz = bax * cay - bay * cax;
    const double len = std::sqrt(nx * nx + ny * ny + nz * nz);
    if (len > 0.0) {
        n_out[0] = float(nx / len);
        n_out[1] = float(ny / len);
        n_out[2] = float(nz / len);
    } else {
        n_out[0] = n_out[1] = n_out[2] = 0.0f;
    }
    return 0.5 * len;
}

// Squared distance from `p` to triangle (a, b, c) — clipped to the
// triangle's interior or boundary, whichever is closest.  Standard
// implementation (Ericson, "Real-Time Collision Detection").
double pointTriangleDistSq(const float p[3],
                           const float a[3], const float b[3], const float c[3]) {
    auto sub = [](const float u[3], const float v[3], double r[3]) {
        r[0] = double(u[0]) - v[0];
        r[1] = double(u[1]) - v[1];
        r[2] = double(u[2]) - v[2];
    };
    auto dot = [](const double u[3], const double v[3]) {
        return u[0] * v[0] + u[1] * v[1] + u[2] * v[2];
    };
    double ab[3], ac[3], ap[3];
    sub(b, a, ab);
    sub(c, a, ac);
    sub(p, a, ap);
    const double d1 = dot(ab, ap);
    const double d2 = dot(ac, ap);
    if (d1 <= 0.0 && d2 <= 0.0) {
        return ap[0]*ap[0] + ap[1]*ap[1] + ap[2]*ap[2];
    }
    double bp[3];
    sub(p, b, bp);
    const double d3 = dot(ab, bp);
    const double d4 = dot(ac, bp);
    if (d3 >= 0.0 && d4 <= d3) {
        return bp[0]*bp[0] + bp[1]*bp[1] + bp[2]*bp[2];
    }
    const double vc = d1 * d4 - d3 * d2;
    if (vc <= 0.0 && d1 >= 0.0 && d3 <= 0.0) {
        const double v = d1 / (d1 - d3);
        const double qx = ap[0] - v * ab[0];
        const double qy = ap[1] - v * ab[1];
        const double qz = ap[2] - v * ab[2];
        return qx*qx + qy*qy + qz*qz;
    }
    double cp[3];
    sub(p, c, cp);
    const double d5 = dot(ab, cp);
    const double d6 = dot(ac, cp);
    if (d6 >= 0.0 && d5 <= d6) {
        return cp[0]*cp[0] + cp[1]*cp[1] + cp[2]*cp[2];
    }
    const double vb = d5 * d2 - d1 * d6;
    if (vb <= 0.0 && d2 >= 0.0 && d6 <= 0.0) {
        const double w = d2 / (d2 - d6);
        const double qx = ap[0] - w * ac[0];
        const double qy = ap[1] - w * ac[1];
        const double qz = ap[2] - w * ac[2];
        return qx*qx + qy*qy + qz*qz;
    }
    const double va = d3 * d6 - d5 * d4;
    if (va <= 0.0 && (d4 - d3) >= 0.0 && (d5 - d6) >= 0.0) {
        const double w = (d4 - d3) / ((d4 - d3) + (d5 - d6));
        const double qx = double(b[0]) + w * (double(c[0]) - b[0]) - p[0];
        const double qy = double(b[1]) + w * (double(c[1]) - b[1]) - p[1];
        const double qz = double(b[2]) + w * (double(c[2]) - b[2]) - p[2];
        return qx*qx + qy*qy + qz*qz;
    }
    // Inside the triangle — return perpendicular distance to its plane.
    const double denom = 1.0 / (va + vb + vc);
    const double v = vb * denom;
    const double w = vc * denom;
    const double qx = double(a[0]) + v * ab[0] + w * ac[0] - p[0];
    const double qy = double(a[1]) + v * ab[1] + w * ac[1] - p[1];
    const double qz = double(a[2]) + v * ab[2] + w * ac[2] - p[2];
    return qx*qx + qy*qy + qz*qz;
}

constexpr double kCoplanarDot = 0.9999;  // ~0.81° tolerance

} // namespace

AreaMeasurement::AreaMeasurement() = default;

void AreaMeasurement::clear(ViewportWindow& vp) {
    mesh_cache_.clear();
    selected_.clear();
    total_area_m2_ = 0.0;
    vp.setHighlightTriangles({}, 0, 0, 0, 0);
}

void AreaMeasurement::rebuildHighlight(ViewportWindow& vp) {
    // Push every selected triangle's three world-space vertices to the
    // overlay.  Mesh-local positions × per-instance composed transform.
    std::vector<float> world_xyz;
    world_xyz.reserve(selected_.size() * 9);
    for (const auto& [key, sel] : selected_) {
        const uint64_t cache_key = (uint64_t(sel.model_id) << 32)
                                 | uint64_t(sel.mesh_id);
        auto cit = mesh_cache_.find(cache_key);
        if (cit == mesh_cache_.end()) continue;
        const MeshCache& c = cit->second;
        if (size_t(sel.tri) * 3 + 2 >= c.indices.size()) continue;
        const float* M = sel.composed_transform;  // column-major
        for (int e = 0; e < 3; ++e) {
            const uint32_t vi = c.indices[3 * sel.tri + e];
            const float* p = &c.positions[3 * vi];
            // World = M * (p, 1).  Column-major: M[col*4 + row].
            const float wx = M[0]*p[0] + M[4]*p[1] + M[8]*p[2]  + M[12];
            const float wy = M[1]*p[0] + M[5]*p[1] + M[9]*p[2]  + M[13];
            const float wz = M[2]*p[0] + M[6]*p[1] + M[10]*p[2] + M[14];
            world_xyz.push_back(wx);
            world_xyz.push_back(wy);
            world_xyz.push_back(wz);
        }
    }
    // Translucent cyan-ish tint — readable on both light and dark surfaces.
    vp.setHighlightTriangles(world_xyz, 0.20f, 0.85f, 1.00f, 0.45f);
}

AreaMeasurement::MeshCache* AreaMeasurement::meshCache(ViewportWindow& vp,
                                                       uint32_t model_id,
                                                       uint32_t mesh_id) {
    const uint64_t key = (uint64_t(model_id) << 32) | uint64_t(mesh_id);
    auto it = mesh_cache_.find(key);
    if (it != mesh_cache_.end()) return &it->second;

    ViewportWindow::MeshTriangles tris;
    if (!vp.readbackMeshTriangles(model_id, mesh_id, tris)) return nullptr;

    MeshCache c;
    c.positions = std::move(tris.positions);
    c.indices   = std::move(tris.indices);
    const size_t n_tris = c.indices.size() / 3;
    c.tri_normals.resize(n_tris * 3);
    c.tri_areas.resize(n_tris);
    c.edges.reserve(n_tris * 3);
    for (size_t t = 0; t < n_tris; ++t) {
        const uint32_t ia = c.indices[3 * t + 0];
        const uint32_t ib = c.indices[3 * t + 1];
        const uint32_t ic = c.indices[3 * t + 2];
        const float* a = &c.positions[3 * ia];
        const float* b = &c.positions[3 * ib];
        const float* cc = &c.positions[3 * ic];
        float n[3];
        c.tri_areas[t] = triAreaAndNormal(a, b, cc, n);
        c.tri_normals[3 * t + 0] = n[0];
        c.tri_normals[3 * t + 1] = n[1];
        c.tri_normals[3 * t + 2] = n[2];
        c.edges[edgeKey(ia, ib)].push_back(uint32_t(t));
        c.edges[edgeKey(ib, ic)].push_back(uint32_t(t));
        c.edges[edgeKey(ic, ia)].push_back(uint32_t(t));
    }
    return &mesh_cache_.emplace(key, std::move(c)).first->second;
}

void AreaMeasurement::onPick(ViewportWindow& vp, int x, int y, bool alt) {
    ViewportWindow::MeshLocalPick pick;
    if (!vp.pickMeshLocalAt(x, y, pick)) return;

    MeshCache* cache = meshCache(vp, pick.model_id, pick.mesh_id);
    if (!cache) return;
    const size_t n_tris = cache->indices.size() / 3;
    if (n_tris == 0) return;

    // Find the seed triangle: the one whose interior (or boundary) is
    // closest to the pick's mesh-local point.
    uint32_t seed = 0;
    double best = std::numeric_limits<double>::infinity();
    for (size_t t = 0; t < n_tris; ++t) {
        const uint32_t ia = cache->indices[3 * t + 0];
        const uint32_t ib = cache->indices[3 * t + 1];
        const uint32_t ic = cache->indices[3 * t + 2];
        const double d = pointTriangleDistSq(pick.mesh_local,
                                             &cache->positions[3 * ia],
                                             &cache->positions[3 * ib],
                                             &cache->positions[3 * ic]);
        if (d < best) {
            best = d;
            seed = uint32_t(t);
        }
    }

    // Expand to coplanar patch (BFS over shared edges).  Alt skips it.
    std::vector<uint32_t> patch;
    if (alt) {
        patch.push_back(seed);
    } else {
        const float* sn = &cache->tri_normals[3 * seed];
        std::unordered_set<uint32_t> visited;
        visited.insert(seed);
        std::queue<uint32_t> frontier;
        frontier.push(seed);
        while (!frontier.empty()) {
            const uint32_t t = frontier.front(); frontier.pop();
            patch.push_back(t);
            for (int e = 0; e < 3; ++e) {
                const uint32_t ia = cache->indices[3 * t + e];
                const uint32_t ib = cache->indices[3 * t + (e + 1) % 3];
                auto it = cache->edges.find(edgeKey(ia, ib));
                if (it == cache->edges.end()) continue;
                for (uint32_t nt : it->second) {
                    if (nt == t || visited.count(nt)) continue;
                    const float* nn = &cache->tri_normals[3 * nt];
                    const double dot = double(sn[0]) * nn[0]
                                     + double(sn[1]) * nn[1]
                                     + double(sn[2]) * nn[2];
                    if (dot < kCoplanarDot) continue;
                    visited.insert(nt);
                    frontier.push(nt);
                }
            }
        }
    }

    // Toggle: if the seed was already in the set, remove the patch;
    // otherwise add it.
    const uint64_t seed_key = triKey(pick.object_id, seed);
    const bool removing = selected_.count(seed_key) > 0;
    double delta = 0.0;
    for (uint32_t t : patch) {
        const uint64_t k = triKey(pick.object_id, t);
        if (removing) {
            auto it = selected_.find(k);
            if (it != selected_.end()) {
                delta -= cache->tri_areas[t];
                selected_.erase(it);
            }
        } else {
            SelectedTri sel;
            sel.model_id = pick.model_id;
            sel.mesh_id  = pick.mesh_id;
            sel.tri      = t;
            std::memcpy(sel.composed_transform, pick.composed_transform,
                        sizeof(sel.composed_transform));
            if (selected_.emplace(k, sel).second) {
                delta += cache->tri_areas[t];
            }
        }
    }
    total_area_m2_ += delta;

    rebuildHighlight(vp);

    qInfo("Area %s%.6f m^2  (total: %.6f m^2, %zu tris)",
          delta >= 0.0 ? "+" : "", delta,
          total_area_m2_, selected_.size());
}

// ----- LengthMeasurement -----------------------------------------------------

namespace {

double dist3(const std::array<float, 3>& a, const std::array<float, 3>& b) {
    const double dx = double(b[0]) - a[0];
    const double dy = double(b[1]) - a[1];
    const double dz = double(b[2]) - a[2];
    return std::sqrt(dx*dx + dy*dy + dz*dz);
}

double triArea3(const std::array<float, 3>& a,
                const std::array<float, 3>& b,
                const std::array<float, 3>& c) {
    const double bax = double(b[0]) - a[0];
    const double bay = double(b[1]) - a[1];
    const double baz = double(b[2]) - a[2];
    const double cax = double(c[0]) - a[0];
    const double cay = double(c[1]) - a[1];
    const double caz = double(c[2]) - a[2];
    const double nx = bay * caz - baz * cay;
    const double ny = baz * cax - bax * caz;
    const double nz = bax * cay - bay * cax;
    return 0.5 * std::sqrt(nx*nx + ny*ny + nz*nz);
}

// Polygon area via best-fit plane + shoelace, falling back to fan
// triangulation when the points stray off the plane.  Returns the
// resulting area and a label naming which path was taken.
struct PolygonAreaResult {
    double      area_m2;
    const char* method;
};

PolygonAreaResult polygonArea(const std::vector<std::array<float, 3>>& pts) {
    using Vec3d = Eigen::Vector3d;
    using Mat3d = Eigen::Matrix3d;
    const size_t n = pts.size();

    // Centroid + bounding box (for the planarity threshold).
    Vec3d centroid = Vec3d::Zero();
    Vec3d bbox_min = Vec3d::Constant(std::numeric_limits<double>::infinity());
    Vec3d bbox_max = Vec3d::Constant(-std::numeric_limits<double>::infinity());
    for (const auto& p : pts) {
        const Vec3d v(p[0], p[1], p[2]);
        centroid += v;
        bbox_min = bbox_min.cwiseMin(v);
        bbox_max = bbox_max.cwiseMax(v);
    }
    centroid /= double(n);
    const double bbox_diag = (bbox_max - bbox_min).norm();

    // 3x3 covariance.  Smallest eigenvector of this is the plane normal.
    Mat3d cov = Mat3d::Zero();
    for (const auto& p : pts) {
        const Vec3d d = Vec3d(p[0], p[1], p[2]) - centroid;
        cov += d * d.transpose();
    }

    Eigen::SelfAdjointEigenSolver<Mat3d> es(cov);
    const Vec3d normal = es.eigenvectors().col(0);  // smallest eigenvalue

    // RMS plane distance, normalised against the bounding-box diagonal.
    double sq_sum = 0.0;
    for (const auto& p : pts) {
        const double d = (Vec3d(p[0], p[1], p[2]) - centroid).dot(normal);
        sq_sum += d * d;
    }
    const double rms = std::sqrt(sq_sum / double(n));
    const bool planar = bbox_diag > 0.0 && (rms / bbox_diag) < 1e-3;

    if (planar) {
        // Build an in-plane orthonormal basis.
        Vec3d u = normal.cross(Vec3d::UnitX());
        if (u.squaredNorm() < 1e-6) u = normal.cross(Vec3d::UnitY());
        u.normalize();
        const Vec3d v = normal.cross(u);

        // Project + shoelace.
        std::vector<std::array<double, 2>> uv(n);
        for (size_t i = 0; i < n; ++i) {
            const Vec3d d = Vec3d(pts[i][0], pts[i][1], pts[i][2]) - centroid;
            uv[i][0] = d.dot(u);
            uv[i][1] = d.dot(v);
        }
        double s = 0.0;
        for (size_t i = 0; i < n; ++i) {
            const auto& a = uv[i];
            const auto& b = uv[(i + 1) % n];
            s += a[0] * b[1] - b[0] * a[1];
        }
        return { 0.5 * std::abs(s), "planar" };
    }

    // Fan from p0.  Works for star-shaped polygons; for genuinely twisted
    // 3D point sets it's a heuristic — flagged in the method label.
    double area = 0.0;
    for (size_t i = 1; i + 1 < n; ++i) {
        area += triArea3(pts[0], pts[i], pts[i + 1]);
    }
    return { area, "fan-triangulated (non-planar)" };
}

} // namespace

LengthMeasurement::LengthMeasurement() = default;

void LengthMeasurement::clear(ViewportWindow& vp) {
    points_.clear();
    vp.setOverlayPoints({}, 0,0,0,0, 0,
                        0,0,0,0, 0);
    vp.setOverlayLines({},  0,0,0,0, 0,
                        0,0,0,0, 0);
    vp.setOverlayLabels({});
    vp.setHudText(QString());
}

void LengthMeasurement::onPick(ViewportWindow& vp, int x, int y, bool /*alt*/) {
    ViewportWindow::MeshLocalPick pick;
    if (!vp.pickMeshLocalAt(x, y, pick)) return;
    points_.push_back({pick.world_pos[0], pick.world_pos[1], pick.world_pos[2]});
    rebuildOverlay(vp);
}

void LengthMeasurement::removeLastPoint(ViewportWindow& vp) {
    if (points_.empty()) return;
    points_.pop_back();
    rebuildOverlay(vp);
}

void LengthMeasurement::rebuildOverlay(ViewportWindow& vp) {
    // Points: orange inner with thin black halo — readable on every
    // background.  Inner 8px disc + 2px halo each side.
    std::vector<float> pts_xyz;
    pts_xyz.reserve(points_.size() * 3);
    for (const auto& p : points_) {
        pts_xyz.push_back(p[0]);
        pts_xyz.push_back(p[1]);
        pts_xyz.push_back(p[2]);
    }
    vp.setOverlayPoints(pts_xyz,
                        /*inner*/  1.00f, 1.00f, 1.00f, 1.00f,
                        /*size*/   8.0f,
                        /*stroke*/ 0.00f, 0.00f, 0.00f, 0.85f,
                        /*extra*/  2.0f);

    // Connecting polyline.  For 4+ points also close the polygon since
    // that's the area-readout shape.  Same orange + halo treatment.
    std::vector<float> seg_xyz;
    std::vector<OverlayRenderer::Label> labels;

    if (points_.size() >= 2) {
        const size_t n = points_.size();
        seg_xyz.reserve(n * 6);
        labels.reserve(n);
        auto pushSegment = [&](const std::array<float, 3>& a,
                               const std::array<float, 3>& b) {
            seg_xyz.insert(seg_xyz.end(), a.begin(), a.end());
            seg_xyz.insert(seg_xyz.end(), b.begin(), b.end());
            OverlayRenderer::Label lbl;
            lbl.world_pos[0] = 0.5f * (a[0] + b[0]);
            lbl.world_pos[1] = 0.5f * (a[1] + b[1]);
            lbl.world_pos[2] = 0.5f * (a[2] + b[2]);
            lbl.text = QString::number(dist3(a, b), 'f', 3) + " m";
            labels.push_back(std::move(lbl));
        };
        for (size_t i = 0; i + 1 < n; ++i) pushSegment(points_[i], points_[i + 1]);
        if (n >= 4) pushSegment(points_[n - 1], points_[0]);
    }
    vp.setOverlayLines(seg_xyz,
                       /*inner*/  1.00f, 1.00f, 1.00f, 1.00f,
                       /*width*/  2.0f,
                       /*stroke*/ 0.00f, 0.00f, 0.00f, 0.85f,
                       /*extra*/  1.5f);
    vp.setOverlayLabels(labels);
    vp.setHudText(formatReadout());
}

QString LengthMeasurement::formatReadout() const {
    const size_t n = points_.size();
    if (n == 0) return QStringLiteral("Length tool: click first point");
    if (n == 1) return QStringLiteral("1 point  (click another)");

    if (n == 2) {
        const auto& a = points_[0];
        const auto& b = points_[1];
        const double d  = dist3(a, b);
        const double dx = std::abs(double(b[0]) - a[0]);
        const double dy = std::abs(double(b[1]) - a[1]);
        const double dz = std::abs(double(b[2]) - a[2]);
        return QString("Length: %1 m\nΔX: %2  ΔY: %3  ΔZ: %4 m")
            .arg(d, 0, 'f', 4)
            .arg(dx, 0, 'f', 4)
            .arg(dy, 0, 'f', 4)
            .arg(dz, 0, 'f', 4);
    }

    if (n == 3) {
        const auto& a = points_[0];
        const auto& b = points_[1];
        const auto& c = points_[2];
        // Angle at b (the middle-clicked vertex).
        const double bax = double(a[0]) - b[0];
        const double bay = double(a[1]) - b[1];
        const double baz = double(a[2]) - b[2];
        const double bcx = double(c[0]) - b[0];
        const double bcy = double(c[1]) - b[1];
        const double bcz = double(c[2]) - b[2];
        const double la = std::sqrt(bax*bax + bay*bay + baz*baz);
        const double lc = std::sqrt(bcx*bcx + bcy*bcy + bcz*bcz);
        double angle_deg = 0.0;
        if (la > 0.0 && lc > 0.0) {
            const double cosang = std::clamp(
                (bax*bcx + bay*bcy + baz*bcz) / (la * lc), -1.0, 1.0);
            angle_deg = std::acos(cosang) * 180.0 / M_PI;
        }
        return QString("Angle at pt 2: %1°\nTriangle area: %2 m²\nPerimeter: %3 m")
            .arg(angle_deg, 0, 'f', 2)
            .arg(triArea3(a, b, c), 0, 'f', 4)
            .arg(dist3(a, b) + dist3(b, c) + dist3(c, a), 0, 'f', 4);
    }

    // 4+ points: polygon area.
    const PolygonAreaResult r = polygonArea(points_);
    double perimeter = 0.0;
    for (size_t i = 0; i < n; ++i) {
        perimeter += dist3(points_[i], points_[(i + 1) % n]);
    }
    return QString("Polygon (%1 pts, %2)\nArea: %3 m²\nPerimeter: %4 m")
        .arg(n)
        .arg(r.method)
        .arg(r.area_m2, 0, 'f', 4)
        .arg(perimeter, 0, 'f', 4);
}