|
  
- UID
- 133
- 帖子
- 51
- 精华
- 1
- 积分
- 186
- 金币
- 55
- 威望
- 2
- 贡献
- 0

|
C语言中显示 点在多边形内 算法
本文是采用射线法判断点是否在多边形内的C语言程序。多年前,我自己实现了这样一个算法。但是随着时间的推移,我决定重写这个代码。参考周培德的《计算几何》一书,结合我的实践和经验,我相信,在这个算法的实现上,这是你迄今为止遇到的最优的代码。* M, J N% n L2 Q1 k& R
$ ?3 I- R! z0 B1 s" }6 ?$ |" ` 这是个C语言的小算法的实现程序,本来不想放到这里。可是,当我自己要实现这样一个算法的时候,想在网上找个现成的,考察下来竟然一个符合需要的也没有。我对自己大学读书时写的代码没有信心,所以,决定重新写一个,并把它放到这里,以飨读者。也增加一下BLOG的点击量。
+ X+ K# R1 J S/ w$ T: x# k3 L" |! t. \ X
首先定义点结构如下:
/ [" z" ]! T5 _2 S* h; v
9 J T+ \: _, C2 { i) f以下是引用片段:
& O+ e* g* u/ n, s4 m. ?1 y$ q* y% A /* Vertex structure */ ' v. ]) ~$ @) F: p' ]
typedef struct * |8 f" f( \7 f. u/ F
{
$ v; L& G8 {- G, u6 Q double x, y; " a/ _$ v$ B1 {% N) B: d
} vertex_t; + Y4 D/ f5 `6 f7 V( d" s& `
( ~( r# {5 k- w1 P4 t. d \* @8 |* H& ~. O! E9 l- d
本算法里所指的多边形,是指由一系列点序列组成的封闭简单多边形。它的首尾点可以是或不是同一个点(不强制要求首尾点是同一个点)。这样的多边形可以是任意形状的,包括多条边在一条绝对直线上。因此,定义多边形结构如下:; v1 w5 ]6 @. a* Y) e
! L# \: V# L7 N0 Q5 S; i
以下是引用片段:8 f! a! h: |" h& P8 d
/* Vertex list structure – polygon */ ! M* _; E- v( Y2 K" `
typedef struct ( S6 C7 i4 L( y
{
& z- L2 `' P7 I# h# `2 A1 a int num_vertices; /* Number of vertices in list */ 1 I- y L3 H+ }/ z9 B; r9 y
vertex_t *vertex; /* Vertex array pointer */
! R0 d. X* t8 @3 r# y; V } vertexlist_t;
( j j0 P' J" k, \& j9 x' [
* V0 `! C" c% y' q& R9 F t2 p; |+ d% i8 M/ R4 [$ m
为加快判别速度,首先计算多边形的外包矩形(rect_t),判断点是否落在外包矩形内,只有满足落在外包矩形内的条件的点,才进入下一步的计算。为此,引入外包矩形结构rect_t和求点集合的外包矩形内的方法vertices_get_extent,代码如下:
$ r7 U7 j1 K& h5 }7 l& `; E/ r
3 ?) z& v, Y. y3 m* ?9 J. t8 \以下是引用片段:. r/ z% f: b; n/ R& S+ {! ?
/* bounding rectangle type */
5 g- q" O' N& ^" S. t typedef struct
( W6 [ j9 C5 C9 ?5 j0 o { 9 u+ f( M2 Y3 j2 W1 F1 f0 h6 c- `
double min_x, min_y, max_x, max_y;
% F3 w1 m$ \2 e* Z; | } rect_t;
* K! }+ j* C' x! G5 _- N8 W /* gets extent of vertices */ . F$ r: h {5 Z! i1 c/ V) \4 C8 P1 S
void vertices_get_extent (const vertex_t* vl, int np, /* in vertices */ " @6 t6 X" Q3 j' B& g5 O
rect_t* rc /* out extent*/ )
' {: N* s( q2 z7 x7 f {
: K5 M- T) _- x" R) n; b9 U. b int i;
! w& p% c' x2 z+ G if (np > 0){
5 z( o5 w% O# e* r- ^# }2 ^' Z rc->min_x = rc->max_x = vl[0].x; rc->min_y = rc->max_y = vl[0].y;
( H9 }3 R) k, a- [/ T6 g }else{ 5 z: x' ?/ H8 |& t; y, k8 ^
rc->min_x = rc->min_y = rc->max_x = rc->max_y = 0; /* =0 ? no vertices at all */
9 V$ [' ]4 D! k" F7 Q }
' A) M: u4 i9 w* r8 ~ for(i=1; i
. b; s0 V( T' J( s8 w {
5 ^; }) ~( i# D# r" L if(vl.x < rc->min_x) rc->min_x = vl.x; 1 [6 S1 H, P7 N' t* u. q
if(vl.y < rc->min_y) rc->min_y = vl.y;
" _$ u, ]4 R3 P* s/ r4 r if(vl.x > rc->max_x) rc->max_x = vl.x; 3 D( e# E' S3 u4 i% A) K
if(vl.y > rc->max_y) rc->max_y = vl.y; 2 H: `( B7 U6 B/ ?6 W* a
}
`% `# k! V" e% w, ]9 ?' B4 U% S5 j } ) W& z2 C# _6 V! Z2 w
9 R. A" U C! I" ?! k' i0 |- h
% T- H5 [: d s( b9 b, }: Q5 g: | 当点满足落在多边形外包矩形内的条件,要进一步判断点(v)是否在多边形(vl:np)内。本程序采用射线法,由待测试点(v)水平引出一条射线B(v,w),计算B与vl边线的交点数目,记为c,根据奇内偶外原则(c为奇数说明v在vl内,否则v不在vl内)判断点是否在多边形内。! ?5 W0 o1 R) x) t% Q
9 R0 F( Z! h9 E
具体原理就不多说。为计算线段间是否存在交点,引入下面的函数:3 S- q1 r" K( m9 G( E; G3 I& c
4 u1 ?$ F1 a0 b. C$ n
(1)is_same判断2(p、q)个点是(1)否(0)在直线l(l_start,l_end)的同侧;) [7 t" e3 C" B3 F# N& G# _
9 [, I. `( S4 Y7 B (2)is_intersect用来判断2条线段(不是直线)s1、s2是(1)否(0)相交;
1 J9 V: a' d4 l7 g$ F5 \2 d
* F( S' `- t9 C! q1 E以下是引用片段: O+ a: ^% z8 S' ?& P( U3 K
/* p, q is on the same of line l */ " d2 Z# |" v9 c. N0 P
static int is_same(const vertex_t* l_start, const vertex_t* l_end, /* line l */
; l! A+ S, O8 f6 T' n7 X% E5 M+ | const vertex_t* p,
+ o& Y) P6 K6 D5 d const vertex_t* q)
- {; p. X9 `& ]: H( v3 M {
" T! J3 w1 | Y* o double dx = l_end->x - l_start->x; # D6 v( j3 r( {' V( n. g
double dy = l_end->y - l_start->y;
6 ~7 ?8 W& Z4 k1 `. p8 C1 } double dx1= p->x - l_start->x; $ s9 J* m; Z5 a2 p l% B
double dy1= p->y - l_start->y;
/ w4 K G2 n: k9 O# F5 d, A double dx2= q->x - l_end->x;
8 `4 g Y& F4 o1 ]- y7 z9 f double dy2= q->y - l_end->y; 7 @# i2 \1 B! R7 b( C
return ((dx*dy1-dy*dx1)*(dx*dy2-dy*dx2) > 0? 1 : 0);
: ]7 Y# Y2 M7 s9 ?( n7 s, H2 |" p# l } ' ?+ T' d6 O# x" B
/* 2 line segments (s1, s2) are intersect? */ ; { u* d- t5 j9 s* S7 W" D3 @
static int is_intersect(const vertex_t* s1_start, const vertex_t* s1_end, 2 |! W5 \5 Q5 |, ]
const vertex_t* s2_start, const vertex_t* s2_end)
. c3 U( ^' Z" w {
/ H) ~1 c, v: P) h return (is_same(s1_start, s1_end, s2_start, s2_end)==0 && - h( g- ~3 c. t: x9 C1 q% y4 B
is_same(s2_start, s2_end, s1_start, s1_end)==0)? 1: 0; , y( S+ h: ]/ c. h, V% }
} ( y# a j+ ~& x1 f+ z6 `
$ g0 P* X- X/ M1 v+ S7 R3 E9 v. I0 z* D# J" c% d4 B! \, C
下面的函数pt_in_poly就是判断点(v)是(1)否(0)在多边形(vl:np)内的程序:7 T* {- B" Q; j. ^
0 k: T) j% o. X& |- l% h
以下是引用片段:
" T, Q; H" S2 ?5 p! S/ Z- g/ Y int pt_in_poly ( const vertex_t* vl, int np, /* polygon vl with np vertices */ , e- m7 c5 ?5 `
const vertex_t* v)
0 Q; g0 v1 _, C. g7 @& \ {
1 i/ U+ \% C1 n* U int i, j, k1, k2, c; 4 n' e! s: F2 N$ i/ N
rect_t rc;
' n6 p9 f- W* [& e: ^$ Z vertex_t w; $ i7 t! @& n3 n# `+ A% `/ ~9 _3 \; `
if (np < 3) 2 Y) l% @( s5 S3 ^! t# c
return 0;
$ e, W9 X% ^. E8 _2 O0 @ vertices_get_extent(vl, np, &rc);
( R( B. D K% O* B' W& w% a8 K if (v->x < rc.min_x || v->x > rc.max_x || v->y < rc.min_y || v->y > rc.max_y) ( ~4 O5 u2 ~- C& o5 ^4 V
return 0;
7 v! y5 i, p/ e* ~8 d& h /* Set a horizontal beam l(*v, w) from v to the ultra right */ / A# e0 {) @5 \0 w
w.x = rc.max_x + DBL_EPSILON;
( \1 F5 E+ C( @( a g6 g w.y = v->y; 3 D! D; X0 I8 Z
c = 0; /* Intersection points counter */ - o. p1 ^: F2 a: O' }% k5 ]; {, G
for(i=0; i
+ M/ n# Y' \6 P4 O" p0 M% X {
9 E% O, X7 Y7 M. D1 K, o3 [ j = (i+1) % np; y1 g- G( x& ~( { Y
if(is_intersect(vl+i, vl+j, v, &w)) , q$ P' m# c6 J! }2 R
{
* V0 b# O$ S9 q C++;
i5 M" E, K3 S# Z }
4 u3 k9 d" V* q4 P else if(vl.y==w.y) 0 F5 d9 g8 c w) Q; z. S
{ ! n& P* W; R- _8 k" I
k1 = (np+i-1)%np; 2 E! n; K8 k# D3 W
while(k1!=i && vl[k1].y==w.y)
! d" z v3 C; k T, ], | k1 = (np+k1-1)%np;
( k7 M; F1 K' y& _1 V k2 = (i+1)%np;
+ q7 }- f" v: F6 H+ m while(k2!=i && vl[k2].y==w.y) 7 V& L( ^4 N% l) ^: I
k2 = (k2+1)%np;
! w9 ]; P9 S$ i, a6 `$ W4 G if(k1 != k2 && is_same(v, &w, vl+k1, vl+k2)==0)
% O; ^) H" o z C++;
5 C3 z2 v2 ~+ m" Y. v) T- R' X0 S if(k2 <= i)
2 J) q0 {7 c# w0 `8 P% J break;
7 q6 T' u* o* p: z* ? J, j( Z; |; | i = k2; % ^5 Y" W" Z$ C+ k
}
' N5 `; Q1 [9 ?5 z } 8 t3 T( y8 {! o8 K% H
return c%2; 5 ~; W- q1 w& n" E9 @
}
7 z$ Y7 Q# y3 s: w- g) K
# {3 v( G; a0 W) {, r' F+ I. u4 h k% A8 O/ M( z6 X( Q9 E6 C
本想配些插图说明问题,但是,CSDN的文章里放图片我还没用过。以后再试吧!实践证明,本程序算法的适应性极强。但是,对于点正好落在多边形边上的极端情形,有可能得出2种不同的结果。 |
|