3D空间共面四边形给定角度下对边正交距离计算及FORTRAN实现

核心思路 #

所有点、待求线段都在同一个平面上，完全没必要在3D空间硬写直线求交、距离计算的逻辑——这类3D实现不仅代码量大，还很容易出数值精度bug。最省事的方案是做正交投影把整个问题转成你已经跑通的2D场景计算，正交投影不会改变平面内的长度、角度、垂直/平行关系，算出来的结果和3D直接求解完全一致，还能100%复用你已有的2D代码，开发量极小。

具体计算步骤 #

• 步骤1：计算平面单位法向量
  取平面内两个不共线的向量（例如v1 = P1 - P0、v2 = P3 - P0），通过叉乘得到法向量，归一化得到单位法向量n：

  n = cross(v1, v2)
  n = n / sqrt(dot(n, n))
  

  如果叉乘结果模长接近0，说明选取的两个向量共线，换一组平面内向量重新计算即可。
• 步骤2：构建平面内局部正交2D坐标系
  取平面内第一个单位基向量u沿P0指向P1的方向：

  u = v1 / sqrt(dot(v1, v1))
  

  第二个单位基向量w通过法向量和u叉乘得到，自然满足在平面内、与u正交的要求：

  w = cross(n, u)
  

  局部坐标系原点取P0即可。
• 步骤3：将所有3D点转换为局部2D坐标
  对任意点Pi，其相对P0的位移向量为vi = Pi - P0，对应的2D坐标为：

  xi = dot(vi, u)
  yi = dot(vi, w)
  

• 步骤4：复用2D逻辑计算线段长度
  得到所有点的2D坐标后，直接传入你已经验证过的2D计算函数，结合已知角度α，即可得到红色、橙色线段的长度，结果无需额外换算，就是3D场景下的真实长度。

Fortran实现代码 #

module vec_op
    implicit none
    integer, parameter :: dp = kind(1.0d0)
contains
    ! 向量点乘
    function dot(a, b) result(res)
        real(dp), intent(in) :: a(3), b(3)
        real(dp) :: res
        res = a(1)*b(1) + a(2)*b(2) + a(3)*b(3)
    end function dot

    ! 向量叉乘
    function cross(a, b) result(res)
        real(dp), intent(in) :: a(3), b(3)
        real(dp) :: res(3)
        res(1) = a(2)*b(3) - a(3)*b(2)
        res(2) = a(3)*b(1) - a(1)*b(3)
        res(3) = a(1)*b(2) - a(2)*b(1)
    end function cross

    ! 向量模长
    function norm(a) result(res)
        real(dp), intent(in) :: a(3)
        real(dp) :: res
        res = sqrt(dot(a,a))
    end function norm
end module vec_op

program quad_dist
    use vec_op
    implicit none
    ! 输入：P0-P4的3D坐标，alpha为弧度制角度
    real(dp) :: P0(3), P1(3), P2(3), P3(3), P4(3), alpha
    ! 输出：红线段长度len_red，橙线段长度len_orange
    real(dp) :: len_red, len_orange
    ! 内部变量
    real(dp) :: v1(3), v2(3), n(3), u(3), w(3)
    real(dp) :: p0_2d(2), p1_2d(2), p2_2d(2), p3_2d(2), p4_2d(2)

    ! --------------- 替换为实际输入赋值逻辑 ---------------
    ! P0 = [x0,y0,z0]
    ! ... 其余点坐标依次赋值
    ! alpha = 角度值 * acos(-1.0_dp)/180.0_dp  ! 角度转弧度用
    ! -----------------------------------------------------------

    ! 1. 计算平面单位法向量
    v1 = P1 - P0
    v2 = P3 - P0
    n = cross(v1, v2)
    if (norm(n) < 1e-12_dp) then
        ! 换一组向量避免共线
        v1 = P2 - P1
        v2 = P3 - P1
        n = cross(v1, v2)
    end if
    n = n / norm(n)

    ! 2. 构建局部正交基
    u = v1 / norm(v1)
    w = cross(n, u)

    ! 3. 所有点转2D局部坐标
    p0_2d = [dot(P0-P0, u), dot(P0-P0, w)]
    p1_2d = [dot(P1-P0, u), dot(P1-P0, w)]
    p2_2d = [dot(P2-P0, u), dot(P2-P0, w)]
    p3_2d = [dot(P3-P0, u), dot(P3-P0, w)]
    p4_2d = [dot(P4-P0, u), dot(P4-P0, w)]

    ! 4. 直接调用已有的2D计算函数，传入p0_2d~p4_2d和alpha即可
    ! 示例：call calc_2d_dist(p0_2d,p1_2d,p2_2d,p3_2d,p4_2d,alpha,len_red,len_orange)

    print *, "红色线段长度：", len_red
    print *, "橙色线段长度：", len_orange

end program quad_dist

注意事项 #

• Fortran内置三角函数默认接收弧度值输入，如果你的α是角度单位，一定要先转成弧度再传入2D计算逻辑，避免结果错误。
• 基向量的选取不影响最终长度计算结果，只要保证两个基向量正交、都在平面内、为单位长度即可。
• 如果点坐标存在测量/计算误差导致不完全共面，该方法会自动将点正交投影到拟合平面上，鲁棒性优于直接3D求交。

内容的提问来源于stack exchange，提问作者learner123