[UP]


Manual Reference Pages  - locpt (3)

NAME

locpt(3f) - [M_math:geometry] find if a point is inside a polygonal path

CONTENTS

Synopsis
Description
Example

SYNOPSIS

Usage:

   subroutine locpt (x0,y0,x,y,n,l,m)
   real, intent(in)     :: x0, y0, x(:), y(:)
   integer, intent(in)  :: n
   integer, intent(out) :: l, m

DESCRIPTION

Given a polygonal line connecting the vertices (X(I),Y(I)) (I = 1,...,N) taken in this order. it is assumed that the polygonal path is a loop, where (X(N),Y(N)) = (X(1),Y(1)) or there is an arc from (X(N),Y(N)) to (X(1),Y(1)). N.B. The polygon may cross itself any number of times.

(X0,Y0) is an arbitrary point and l and m are variables. On output, L and M are assigned the following values ...

     L = -1   If (X0,Y0) is outside the polygonal path
     L =  0   If (X0,Y0) lies on the polygonal path
     L =  1   If (X0,Y0) is inside the polygonal path

M = 0 if (X0,Y0) is on or outside the path. If (X0,Y0) is inside the path then M is the winding number of the path around the point (X0,Y0).
o Fortran 66 version by A.H. Morris
o Converted to ELF90 compatibility by Alan Miller, 15 February 1997 saved from url=(0050)http://users.bigpond.net.au/amiller/NSWC/locpt.f90

EXAMPLE

Draw a polygon and the envelope of the polygon, and find the area of each polygon. Also place a number of small circles in the plot area colored according to whether they fall within the border of the original polygon.

   program demo_envelope
   use M_draw
   use M_math,     only : envelope        ! Find vertices (in clockwise order) of a polygon enclosing the points
   use M_math,     only : locpt           ! find if a point is inside a polygonal path
   use M_math,     only : polyarea        ! compute the area bounded by a closed polygonal curve
   implicit none
   integer,parameter :: n=6
   !  3--------------4
   !   \           /
   !     \       /
   !       \   /
   !         X 2,5
   !       /  \
   !     /      \
   !   /          \
   !  1--------------6
   real,parameter    :: x(n)=[-5.0, 0.0,-5.0, 5.0, 0.0, 5.0]
   real,parameter    :: y(n)=[-5.0, 0.0, 5.0, 5.0, 0.0,-5.0]
   real              :: xy(2,n)
   integer           :: vertex(n)
   integer           :: nvert
   integer           :: i
   integer           :: idum
      xy(1,:)=x
      xy(2,:)=y
      call vinit(’ ’)
      call page(-10.0,10.0,-10.0,10.0)
      call color(D_BLACK) ! set current color to black
      call clear()        ! clear to current color
      call polyfill(.true.)
      call color(D_BLUE)  ! we want to draw polygon in this color
      call poly2(n,xy)    ! draw filled polygon using points given
      idum=getkey()       ! pause for some input
      call color(D_CYAN)
      call polyhatch(.true.)
      call envelope(x, y, n, vertex, nvert)   ! calculate envelope
      call poly2(nvert,xy(:,vertex(1:nvert))) ! draw hatched envelope
      idum=getkey()       ! pause for some input
      call polyhatch(.false.)
      call linewidth(50)
      call color(D_WHITE)
      call poly2(n,xy)    ! draw line along original points
      idum=getkey()       ! pause for some input
      call random_seed()
      do i=1,70
         call pickrandom()
      enddo
      idum=getkey()       ! pause for some input
      call vexit()        ! wrap up and exit graphics mode
      write(*,*)’polyarea=’,polyarea(x,y)
      write(*,*)’polyarea=’,polyarea( xy(1,vertex(1:nvert)), xy(2,vertex(1:nvert)))
   contains
   subroutine pickrandom()
   ! randomly pick a point in the plot area and color it according to whether it is inside
   ! the original polygon
   real :: pointx, pointy
   integer :: l, m
      call random_number(pointx)
      call random_number(pointy)
      pointx=pointx*20.0-10.0
      pointy=pointy*20.0-10.0
      call locpt(pointx,pointy,x,y,n,l,m)
      select case(l)
       case(-1)
         call color(D_RED)
       case(0)
         call color(D_YELLOW)
       case(1)
         call color(D_GREEN)
       case default
         write(*,*)’*pickrandom* internal error: L value unknown’
         call color(D_WHITE)
      end select
      call circle(pointx,pointy,0.2)
   end subroutine pickrandom
   end program demo_envelope


locpt (3) March 11, 2021
Generated by manServer 1.08 from f356b7cb-34d0-4434-9900-9e9388b3354c using man macros.