home *** CD-ROM | disk | FTP | other *** search
- \ Original Date: November 4, 1985
- \ Last Modified: January 2, 1989
- \ Author: Jack W. Brown
- \ File name: L4P16.SEQ
- \ Function: Computes Area of a Polygon given the x,y
- \ coordinates of its verticies
-
- \ The following mathematical algorithm is often used to
- \ determine the area of cross-section provided it can be
- \ represented adequately by a finite number of straight line
- \ segments (this is almost always possible). The technique
- \ can also be applied to cross-sections with holes by moving
- \ around the hole in a counter clockwise direction and traversing
- \ to and from the hole along the same path.
-
- \ The general algorithm.
-
- \ p1 /---------\ p2 p1 = ( x1,y1 )
- \ / \ p2 = ( x2,y2 )
- \ / \ p3 p3 = ( x3,y3 )
- \ / / p4 = ( x4,y4 )
- \ p5 /--------------/ p4 p5 = ( x5,y5 )
- \
- \ AREA OF THE POLYGON =
- \ [(x1y5-x5y1)+(x2y1-x1y2)+(x3y2-x2y3)+(x4y3-x3y4)+(x5y4-x4y5)]/2
- \ In general:
- \ i=n
- \ AREA = 0.5*SUM [ x(i)y(i-1) - x(i-1)y(i) ]
- \ i=1
- \ where we define x0 to be x5 and y0 to be y5.
-
- \ Example without a hole.
- \ X Not drawn to scale!!
- \ | p1 = ( 8,4 )
- \ | p2 = ( 6,1 )
- \ | p4 ----------- p1 p3 = ( 2,1 )
- \ | / / p4 = ( 5,4 )
- \ | / /
- \ | / /
- \ | p3 ----------- p2
- \ |-----------------------Y
-
- \ A = [(8*4-5*4)+(6*4-8*1)+(2*1-6*1)+(5*1-2*4)]/2 = 10.5
-
-
- \ Example of a polygon with a hole removed
- \ Sorry but the diagram below is not to scale. units= centimeters
- \ Y
- \ p9 ---|-------------------------------- p1 p1 = (6,5)
- \ \ | p5 p4 / p2 = (2,0) = p8
- \ \ | +----------+ / p3 = (3,3) = p7
- \ \| |cut out | / p4 = (3,4)
- \ \ +----------+ / p5 = (1,4)
- \ |\ p6 p7,p2 / p6 = (1,3)
- \ | \ / p7 = (3,3)
- \ | \ / p8 = (2,0)
- \ | \ / p9 = (-1,5)
- \ | \ /
- \ | \ /
- \ | \/ p8,p2
- \ ---+-------------------------- X
- \
- \ Traverse outside clockwise and the cut out counter clockwise.
- \ A = [(6*5-(-1)*5)+(2*5-6*0)+(3*0-2*3)+(3*3-3*4)+(1*4-3*4)
- \ +(1*4-1*3)+(3*3-1*3)+(2*3-3*0)+(-1*0-2*5)]/2
- \ A = 15.5 sq cm
-
-
- \ This should be replaced by the bullet proof version in the
- \ file JBINPUT.SEQ
- : #IN ( -- n )
- QUERY INTERPRET ;
-
- CREATE X 102 ALLOT \ Array for x coordinates
- CREATE Y 102 ALLOT \ Array for y coordinates
-
- VARIABLE #POINTS \ Number of points in polygon
- VARIABLE AREA \ Sum of the x(i)y(i-1) - x(i-1)y(i)
-
- \ Fetch ith x component.
- : X@ ( i x{i} ) 2* X + @ ;
-
- \ Fetch ith y component.
- : Y@ ( i y{i} ) 2* Y + @ ;
-
- \ Store ith x component.
- : X! ( x i -- ) 2* X + ! ;
-
- \ Store ith y component.
- : Y! ( y i -- ) 2* Y + ! ;
-
- \ Move to the next tab stop.
- : TAB ( -- )
- BEGIN #OUT @ 8 MOD
- IF SPACE ELSE EXIT THEN
- AGAIN ;
-
- \ Get number from keyboard.
- : GET# ( -- n )
- ASCII > EMIT SPACE #IN ;
-
- \ Prompt and fetch number of data points.
- : GET_#POINTS ( -- )
- BEGIN
- CR ." Enter number of data points. "
- GET# DUP 3 <
- WHILE CR ." You need at least 3 data points!"
- REPEAT 50 MIN #POINTS ! ;
-
-
- \ Prompt and fetch all data points.
- : GET_DATA ( -- )
- CR CR ." Point " TAB ." X" TAB ." Y"
- #POINTS @ 1+ 1
- DO CR I 3 .R TAB GET# I X!
- TAB GET# I Y! LOOP
- #POINTS @ DUP X@ 0 X! Y@ 0 Y! ; \ Store last point in 0th slot
-
- \ Sum data points.
- : FIND_AREA ( -- )
- 0 AREA !
- #POINTS @ 1+ 1 ( n+1 so we loop n times )
- DO I X@ I 1- Y@ * ( X{i}*Y{i-1} )
- I 1- X@ I Y@ * ( X{i-1}*Y{i} )
- - AREA +!
- LOOP ;
-
-
- \ Display computed area.
- : PUT_AREA ( -- )
- AREA @ 2 /MOD
- CR ." AREA = " 6 .R ASCII . EMIT
- IF ASCII 5 EMIT ELSE ASCII 0 EMIT THEN SPACE ;
-
- \ Compute area of polygon.
- : POLY ( -- )
- GET_#POINTS
- GET_DATA
- FIND_AREA
- PUT_AREA ;
-
-
-
-