有向面积&有向体积

平面三角形有向面积

对于2D平面上的一个三角形，其有向面积的计算公式如下

$A = \frac{1}{2} \left|\begin{matrix} x_{0} & y_{0} & 1 \\ x_{1} & y_{1} & 1 \\ x_{2} & y_{2} & 1 \end{matrix}\right|$

C++代码如下：

double signedAreaOfTriangle(const Eigen::Vector2d &p0, const Eigen::Vector2d &p1, const Eigen::Vector2d &p2)
{
    double v12 = (p1(0)-p0(0)) * (p2(1)-p0(1));
    double v21 = (p2(0)-p0(0)) * (p1(1)-p0(1));
    return 0.5 * (v12 - v21);
}

如果三角形的点按照逆时针排布，其面积为正，如果按照顺时针排布，其面积为负。

空间四面体有向体积

对于3D空间中的一个四面体，其有向体积的计算公式如下

$V = \frac{1}{6} \left|\begin{matrix} x_{0} & y_{0} & z_{0} & 1 \\ x_{1} & y_{1} & z_{1} & 1 \\ x_{2} & y_{2} & z_{2} & 1 \\ x_{3} & y_{3} & z_{3} & 1 \end{matrix}\right|$

C++代码如下：

double signedVolumeOfTetrahedron(const Eigen::Vector3d &p0, const Eigen::Vector3d &p1, const Eigen::Vector3d &p2, const Eigen::Vector3d &p3) 
{
    Eigen::Vector3d p01 = p1 - p0;
    Eigen::Vector3d p02 = p2 - p0;
    Eigen::Vector3d p03 = p3 - p0;
    double v123 = p01(0)*p02(1)*p03(2);
    double v231 = p02(0)*p03(1)*p01(2);
    double v312 = p03(0)*p01(1)*p02(2);
    double v321 = p03(0)*p02(1)*p01(2);
    double v213 = p02(0)*p01(1)*p03(2);
    double v132 = p01(0)*p03(1)*p02(2);
    return (1.0/6.0)*(v123 - v231 - v312 - v321 + v213 + v132);
}

应用

平面多边形面积计算

对于二维平面上的一个封闭多边形，其面积可以通过计算多边形的每条边与原点构成的三角形的有向面积之和得到，注意多边形的顶点必须有序，即为顺时针或者逆时针。

double areaOfPolygon(const vector<Eigen::Vector2d> &points)
{
    // points must be oriented, clock-wise or anti-clock-wise
    double areas = 0.0;
    int num_pts = points.size();
    if (num_pts < 3)
        return areas;
    Eigen::Vector2d o(0.0, 0.0);  // origin in 2d plane
    for (int i = 0; i < num_pts-1; ++i)
    {
      areas += signedAreaOfTriangle(o, points[i], points[i+1]);
    }
    areas += signedAreaOfTriangle(o, points[num_pts-1], points[0]);
    return abs(areas);
}

平面多边形顶点顺序

判断平面多边形的有序顶点是顺时针还是逆时针，可以通过计算其有向面积得到，如果面积为正即为逆时针，如面积为负即为顺时针。

空间网格的体积计算

对于三维空间中的一个封闭的三角网格，通过计算网格的每个三角面片与原点构成的四面体的有向体积之和，便可以得到该三角网格的有向体积。

double volumeOfMesh(const vector<Eigen::Vector3d> &points, const vector<Eigen::Vector3i> &triangles) 
{
    // the normals of the mesh must be coninsistent
    double vols = 0.0;
    if (points.size() < 4)
        return vols;
    Eigen::Vector3d o(0.0, 0.0, 0.0);  // origin in 3d space
    Eigen::Vector3d p1, p2, p3;
    for (const auto tri& : triangles)
    {
        p1 = points[tri(0)];
        p2 = points[tri(1)];
        p3 = points[tri(2)];
        vols += signedVolumeOfTetrahedron(o, p1, p2, p3);
    }
    return abs(vols);
}