Fortran Coder

标题: 如何实现判别一个点是否在四边形内 [打印本页]

作者: hipeilei    时间: 2015-10-25 17:49
标题: 如何实现判别一个点是否在四边形内
如何实现使用Fortran判别一个点点是否在四边形内?能不能列举个简单的程序说明一下,比如叉乘判别法、面积判别法

作者: vvt    时间: 2015-10-25 19:58
本帖最后由 vvt 于 2015-10-25 20:00 编辑

John Burkardt 写过一个几何函数库,里面有很多类似的函数。2D的,3D的都有。
参考 http://fcode.cn/code_prof-11-1.html
我列几个函数给你,主要是调用 polygon_contains_point_2_2d ,其他函数是内部实用的。

这个函数可以判断一个二维的点,是否在一个(凸)多边形里。当多边形为四边形时,就是你想要的。

[Fortran] 纯文本查看 复制代码
Program www_fcode_cn
  Implicit None
  Logical :: ins
  Call polygon_contains_point_2_2d(4, [0.,1.,1.,0.], 1.1, [1.,1.,0.,0.], 0.5, ins)
  Write (*, *) ins
End Program www_fcode_cn
Subroutine polygon_contains_point_2_2d(n, xn, xval, yn, yval, inside)
  ! *******************************************************************************
  ! ! POLYGON_CONTAINS_POINT_2_2D finds if a point is inside a convex polygon in 2D.
  ! Modified:
  ! 06 February 1999
  ! Author:
  ! John Burkardt
  ! Parameters:
  ! Input, integer N, the number of nodes or vertices in the polygon.
  ! N must be at least 3.
  ! Input, real XN(N), the X coordinates of the vertices.
  ! Input, real XVAL, the X coordinate of the point to be tested.
  ! Input, real YN(N), the Y coordinates of the vertices.
  ! Input, real YVAL, the Y coordinate of the point to be tested.
  ! Output, logical INSIDE, is .TRUE. if ( X,Y) is inside
  ! the polygon or on its boundary, and .FALSE. otherwise.
  Implicit None
  Integer n,i
  Logical inside
  Real x1,x2,x3
  Real xn(n),xval,y1,y2,y3,yn(n),yval
  inside = .False.
  ! A point is inside a convex polygon if and only if it is inside
  ! one of the triangles formed by X(1),Y(1) and any two consecutive
  ! points on the polygon's circumference.
  x1 = xn(1)
  y1 = yn(1)
  Do i = 2, n - 1
    x2 = xn(i)
    y2 = yn(i)
    x3 = xn(i+1)
    y3 = yn(i+1)
    Call triangle_contains_point_1_2d(x1, y1, x2, y2, x3, y3, xval, yval, inside)
    If (inside) Then
      Return
    End If
  End Do
  Return
End Subroutine polygon_contains_point_2_2d
Subroutine triangle_contains_point_1_2d(x1, y1, x2, y2, x3, y3, x, y, inside)
  ! *******************************************************************************
  ! ! TRIANGLE_CONTAINS_POINT_1_2D finds if a point is inside a triangle in 2D.
  ! Modified:
  ! 16 June 2001
  ! Author:
  ! John Burkardt
  ! Parameters:
  ! Input, real X1, Y1, X2, Y2, X3, Y3, the triangle vertices.
  ! The vertices should be given in counter clockwise order.
  ! Input, real X, Y, the point to be checked.
  ! Output, logical INSIDE, is .TRUE. if (X,Y) is inside
  ! the triangle or on its boundary, and .FALSE. otherwise.
  Implicit None
  Real c(3)
  Logical inside
  Real x,x1,x2,x3
  Real y,y1,y2,y3
  Call triangle_barycentric_2d(x1, y1, x2, y2, x3, y3, x, y, c)
  inside = .not.(any(c(1:3)<0.0E+00))
End Subroutine triangle_contains_point_1_2d
Subroutine triangle_barycentric_2d(x1, y1, x2, y2, x3, y3, x, y, c)
  ! *******************************************************************************
  ! ! TRIANGLE_BARYCENTRIC_2D finds the barycentric coordinates of a point in 2D.
  ! Discussion:
  ! The barycentric coordinate of point X related to vertex A can be
  ! interpreted as the ratio of the area of the triangle with
  ! vertex A replaced by vertex X to the area of the original
  ! triangle.
  ! Modified:
  ! 20 October 2001
  ! Author:
  ! John Burkardt
  ! Parameters:
  ! Input, real X1, Y1, X2, Y2, X3, Y3, the triangle vertices.
  ! The vertices should be given in counter clockwise order.
  ! Input, real X, Y, the point to be checked.
  ! Output, real C(3), the barycentric coordinates of (X,Y) with respect
  ! to the triangle.
  Implicit None
  Integer, Parameter :: n = 2
  Integer, Parameter :: nrhs = 1
  Real a(n, n+nrhs)
  Real c(3)
  Integer info
  Real x,x1,x2,x3
  Real y,y1,y2,y3
  ! Set up the linear system
  ! ( X2-X1  X3-X1 ) C1  = X-X1
  ! ( Y2-Y1  Y3-Y1 ) C2    Y-Y1
  ! which is satisfied by the barycentric coordinates of (X,Y).
  a(1, 1) = x2 - x1
  a(1, 2) = x3 - x1
  a(1, 3) = x - x1
  a(2, 1) = y2 - y1
  a(2, 2) = y3 - y1
  a(2, 3) = y - y1
  ! Solve the linear system.
  Call rmat_solve(a, n, nrhs, info)
  If (info/=0) Then
    Write (*, '(a)') ' '
    Write (*, '(a)') 'TRIANGLE_BARYCENTRIC_2D - Fatal error!'
    Write (*, '(a)') '  The linear system is singular.'
    Write (*, '(a)') '  The input data does not form a proper triangle.'
    Stop
  End If
  c(1) = a(1, 3)
  c(2) = a(2, 3)
  c(3) = 1.0E+00 - c(1) - c(2)
End Subroutine triangle_barycentric_2d
Subroutine rmat_solve(a, n, nrhs, info)
  ! *******************************************************************************
  ! ! RMAT_SOLVE uses Gauss-Jordan elimination to solve an N by N linear system.
  ! Modified:
  ! 08 November 2000
  ! Author:
  ! John Burkardt
  ! Parameters:
  ! Input/output, real A(N,N+NRHS), contains in rows and columns 1
  ! to N the coefficient matrix, and in columns N+1 through
  ! N+NRHS, the right hand sides.  On output, the coefficient matrix
  ! area has been destroyed, while the right hand sides have
  ! been overwritten with the corresponding solutions.
  ! Input, integer NRHS, the number of right hand sides.  NRHS
  ! must be at least 0.
  ! Output, integer INFO, singularity flag.
  ! 0, the matrix was not singular, the solutions were computed;
  ! J, factorization failed on step J, and the solutions could not
  ! be computed.
  Implicit None
  Integer n
  Integer nrhs
  Real a(n, n+nrhs)
  Real apivot,factor,temp
  Integer i,j,k
  Integer info,ipivot
  info = 0
  Do j = 1, n
    ipivot = j
    apivot = a(j, j)
    Do i = j + 1, n
      If (abs(a(i,j))>abs(apivot)) Then
        apivot = a(i, j)
        ipivot = i
      End If
    End Do
    If (apivot==0.0E+00) Then
      info = j
      Return
    End If
    Do i = 1, n + nrhs
      Call r_swap(a(ipivot,i), a(j,i))
    End Do
    a(j, j) = 1.0E+00
    a(j, j+1:n+nrhs) = a(j, j+1:n+nrhs)/apivot
    Do i = 1, n
      If (i/=j) Then
        factor = a(i, j)
        a(i, j) = 0.0E+00
        a(i, j+1:n+nrhs) = a(i, j+1:n+nrhs) - factor*a(j, j+1:n+nrhs)
      End If
    End Do
  End Do
End Subroutine rmat_solve
Subroutine r_swap(x, y)
  Implicit None
  Real x,y,z
  z = x
  x = y
  y = z
End Subroutine r_swap


作者: hipeilei    时间: 2015-10-26 18:52
本帖最后由 hipeilei 于 2015-10-31 15:28 编辑
vvt 发表于 2015-10-25 19:58
John Burkardt 写过一个几何函数库,里面有很多类似的函数。2D的,3D的都有。
参考 http://fcode.cn/code_p ...

非常感谢!后来我用面积法做出来了。
作者: kerb    时间: 2015-11-9 23:43
四边形从对角线分成两个三角形,对于其中一个三角形,其内点与三个角点连线构成的三个小三角形面积之和等于这个三角的面积,如果其中一个小三角形面积等于零说明这个点在某条边上,如果三个小三角形的面积之和大于大三角形的面积,点在其外,已知三点坐标其面积有一个行列式公式计算其面积,编程应该容易
作者: aliouying    时间: 2015-11-10 12:43
另外还有一种思路:以这个点链接四个端点,构成四个三角形,计算这四给三角形的面积和,与这个四边形面积做对比:若相等,则在其内或者边上,若大于四边形面积,则在四边形之外。
作者: hipeilei    时间: 2015-11-12 10:11
kerb 发表于 2015-11-9 23:43
四边形从对角线分成两个三角形,对于其中一个三角形,其内点与三个角点连线构成的三个小三角形面积之和等于 ...

我就是这么做的,多谢
作者: hipeilei    时间: 2015-11-12 10:12
aliouying 发表于 2015-11-10 12:43
另外还有一种思路:以这个点链接四个端点,构成四个三角形,计算这四给三角形的面积和,与这个四边形面积做 ...

我正是这么做的,我觉得这样的方法就已经满足我的要求了,多谢!
作者: 珊瑚虫    时间: 2015-11-17 09:33
再取一个四边形内的点(如形心),现在有两个点,将这个两个点分别带入4条线的方程,如果两个点带入每条线的方程所得数的符号都一致,则点在多变形内。同样的方法可以判断,点在任意多变形内,或者多面体内。

作者: aliouying    时间: 2015-11-18 12:19
珊瑚虫 发表于 2015-11-17 09:33
再取一个四边形内的点(如形心),现在有两个点,将这个两个点分别带入4条线的方程,如果两个点带入每条线 ...

这个思想不错,不知道计算量怎么样,要不大家都写给小程序,一起测试测试?
作者: 珊瑚虫    时间: 2015-11-18 15:20
aliouying 发表于 2015-11-18 12:19
这个思想不错,不知道计算量怎么样,要不大家都写给小程序,一起测试测试? ...

速度快,实现简单,精度高,判断点与多边形、多面体关系时的不二之选,值得拥有
作者: 糖盒love玲珑    时间: 2015-12-28 19:13
珊瑚虫 的方法 思路简单, 就是需要找个合适的内点啊. 给一堆坐标, 找个内点也是不容易实现的吧??
作者: 糖盒love玲珑    时间: 2016-3-25 21:54
其实 这个问题,用射线法是最简单有效的算法了,你可以去知乎上看看。
作者: 深流水静水流深    时间: 2016-3-27 08:25
学习了!很好很给力!





欢迎光临 Fortran Coder (http://bbs.fcode.cn/) Powered by Discuz! X3.2