78 lines
1.7 KiB
Fortran
78 lines
1.7 KiB
Fortran
module Ray_Casting_Algo
|
|
use Polygons
|
|
implicit none
|
|
|
|
real, parameter, private :: eps = 0.00001
|
|
private :: ray_intersects_seg
|
|
|
|
contains
|
|
|
|
function ray_intersects_seg(p0, a0, b0) result(intersect)
|
|
type(point), intent(in) :: p0, a0, b0
|
|
logical :: intersect
|
|
|
|
type(point) :: a, b, p
|
|
real :: m_red, m_blue
|
|
|
|
p = p0
|
|
! let variable "a" be the point with smallest y coordinate
|
|
if ( a0%y > b0%y ) then
|
|
b = a0
|
|
a = b0
|
|
else
|
|
a = a0
|
|
b = b0
|
|
end if
|
|
|
|
if ( (p%y == a%y) .or. (p%y == b%y) ) p%y = p%y + eps
|
|
|
|
intersect = .false.
|
|
|
|
if ( (p%y > b%y) .or. (p%y < a%y) ) return
|
|
if ( p%x > max(a%x, b%x) ) return
|
|
|
|
if ( p%x < min(a%x, b%x) ) then
|
|
intersect = .true.
|
|
else
|
|
if ( abs(a%x - b%x) > tiny(a%x) ) then
|
|
m_red = (b%y - a%y) / (b%x - a%x)
|
|
else
|
|
m_red = huge(m_red)
|
|
end if
|
|
if ( abs(a%x - p%x) > tiny(a%x) ) then
|
|
m_blue = (p%y - a%y) / (p%x - a%x)
|
|
else
|
|
m_blue = huge(m_blue)
|
|
end if
|
|
if ( m_blue >= m_red ) then
|
|
intersect = .true.
|
|
else
|
|
intersect = .false.
|
|
end if
|
|
end if
|
|
|
|
end function ray_intersects_seg
|
|
|
|
function point_is_inside(p, pol) result(inside)
|
|
logical :: inside
|
|
type(point), intent(in) :: p
|
|
type(polygon), intent(in) :: pol
|
|
|
|
integer :: i, cnt, pa, pb
|
|
|
|
cnt = 0
|
|
do i = lbound(pol%vertices,1), ubound(pol%vertices,1), 2
|
|
pa = pol%vertices(i)
|
|
pb = pol%vertices(i+1)
|
|
if ( ray_intersects_seg(p, pol%points(pa), pol%points(pb)) ) cnt = cnt + 1
|
|
end do
|
|
|
|
inside = .true.
|
|
if ( mod(cnt, 2) == 0 ) then
|
|
inside = .false.
|
|
end if
|
|
|
|
end function point_is_inside
|
|
|
|
end module Ray_Casting_Algo
|