标题:
C语言中显示 点在多边形内 算法
[打印本页]
作者:
zw2004
时间:
2008-1-21 17:20
标题:
C语言中显示 点在多边形内 算法
本文是采用射线法判断点是否在多边形内的C语言程序。多年前，我自己实现了这样一个算法。但是随着时间的推移，我决定重写这个代码。参考周培德的《计算几何》一书，结合我的实践和经验，我相信，在这个算法的实现上，这是你迄今为止遇到的最优的代码。
- S: x. w$ L0 c' u
7 B" I. I# p3 p* N2 h7 p+ G: O
这是个C语言的小算法的实现程序，本来不想放到这里。可是，当我自己要实现这样一个算法的时候，想在网上找个现成的，考察下来竟然一个符合需要的也没有。我对自己大学读书时写的代码没有信心，所以，决定重新写一个，并把它放到这里，以飨读者。也增加一下BLOG的点击量。
0 ^ A1 u+ q( k) p9 s
; ]: ?% l$ Q$ d8 P% z
首先定义点结构如下：
S( s. S+ F Z c
. E' `! q8 V- V3 f7 m3 ]
以下是引用片段：
: E4 n3 a) E* ]" @# T1 X' l* }0 F1 m% b2 D
/* Vertex structure */
8 @" ~ o( E8 e6 \! p$ j' q
typedef struct
' T) m' x- Y( C1 O" |
{
0 d, P9 p/ x# X( ~5 Y" L) X& V4 @+ K$ \
double x, y;
) G4 T0 G2 U" i! a
} vertex_t;
. Q, R g9 w, T/ d+ C- C
7 V( t; j8 \' u8 c' a1 }$ k/ F
: z5 v0 [' B8 m$ l/ a. Q+ _
本算法里所指的多边形，是指由一系列点序列组成的封闭简单多边形。它的首尾点可以是或不是同一个点(不强制要求首尾点是同一个点)。这样的多边形可以是任意形状的，包括多条边在一条绝对直线上。因此，定义多边形结构如下：
) \8 x/ X! s1 L3 R& A& Q+ a% G
* z, ~% v; ~7 \( O2 R: n& I5 O
以下是引用片段：
' C, ~0 ^0 ]; z! D% j& E6 f6 f8 r- l
/* Vertex list structure – polygon */
& \& t5 E3 W+ q
typedef struct
+ s; H0 z; x% U* F: K1 `& W# D
{
( g" i/ \! P8 L$ F
int num_vertices; /* Number of vertices in list */
& x$ E6 a& I- Z; ?2 y7 A
vertex_t *vertex; /* Vertex array pointer */
* ~0 I w8 h; ~" h# P5 x6 u' N
} vertexlist_t;
9 M5 P: w& v% N d) D! g( L
* X3 w0 l1 |$ l0 Q2 c: M
- W. U% |* F1 G9 c( M
为加快判别速度，首先计算多边形的外包矩形(rect_t)，判断点是否落在外包矩形内，只有满足落在外包矩形内的条件的点，才进入下一步的计算。为此，引入外包矩形结构rect_t和求点集合的外包矩形内的方法vertices_get_extent，代码如下：
7 }) v; R9 _8 s/ U
% ^' s1 s5 G' ^. L Q& y
以下是引用片段：
0 Y6 M/ l6 }4 t
/* bounding rectangle type */
4 b# q+ W# q& q h; }4 J
typedef struct
, m ?" r. v3 H, C
{
; D \; r/ Z: y S8 Y$ _# ?: L
double min_x, min_y, max_x, max_y;
$ X0 P: l/ L& X. ?' V6 x
} rect_t;
( q" P. [% R5 X$ O9 M
/* gets extent of vertices */
$ w6 ]0 ]- h; T# H
void vertices_get_extent (const vertex_t* vl, int np, /* in vertices */
8 D. m; v- [! _9 i0 q# y
rect_t* rc /* out extent*/ )
5 K5 n& V: T6 g' W( ~ r* V8 k
{
! q% i# j: l' Q4 ] c. _
int i;
1 I7 S$ L$ f6 D) n0 h" L% C
if (np > 0){
% [ E B, [/ }& O- {. q
rc->min_x = rc->max_x = vl[0].x; rc->min_y = rc->max_y = vl[0].y;
) [7 |* N% b5 q% m4 _% F8 n8 x
}else{
" _& B; H8 z: D
rc->min_x = rc->min_y = rc->max_x = rc->max_y = 0; /* =0 ? no vertices at all */
% U: B" Z) S5 P e) |% [, m5 o
}
/ C! g- _& o E" A: D
for(i=1; i
* m. t4 m" q7 Q% x6 s
{
. w, f" {3 H8 N. J i
if(vl
.x < rc->min_x) rc->min_x = vl
.x;
, N3 a4 c9 v8 _) i. w
if(vl
.y < rc->min_y) rc->min_y = vl
.y;
: d: B6 C$ [2 _! o& h# V6 y
if(vl
.x > rc->max_x) rc->max_x = vl
.x;
! W2 J Q9 R1 V' O0 K7 R
if(vl
.y > rc->max_y) rc->max_y = vl
.y;
: `% {* Q7 v4 x: e* w1 i
}
8 w4 d1 m- r" C
}
6 Q6 {2 [0 U @
7 E$ `( c) x* W8 _5 V" s
3 G+ }2 A7 }( O# F& X
当点满足落在多边形外包矩形内的条件，要进一步判断点(v)是否在多边形(vl：np)内。本程序采用射线法，由待测试点(v)水平引出一条射线B(v，w)，计算B与vl边线的交点数目，记为c，根据奇内偶外原则(c为奇数说明v在vl内，否则v不在vl内)判断点是否在多边形内。
7 w4 x0 R0 x( n. L/ `& J
. {+ }9 O8 Q2 b4 C& j4 V. m" j
具体原理就不多说。为计算线段间是否存在交点，引入下面的函数：
$ }8 R* L) l& t L) P
" ?9 ~4 k8 G$ ]8 F8 }
(1)is_same判断2(p、q)个点是(1)否(0)在直线l(l_start，l_end)的同侧;
% C) r0 W' P- F; j& d: z+ ~* W
6 _, H# l8 M. y7 ~ J" k
(2)is_intersect用来判断2条线段(不是直线)s1、s2是(1)否(0)相交;
1 _: D# U" Z" p4 Z8 S1 n8 C! r( u4 F
5 M2 P4 p6 K5 s0 W5 [( h2 e
以下是引用片段：
9 ^1 i0 g4 p5 E# H
/* p, q is on the same of line l */
: B* [9 h: F+ p; q- F
static int is_same(const vertex_t* l_start, const vertex_t* l_end, /* line l */
4 }3 d- N- G" G( T) s$ p; o. N
const vertex_t* p,
/ @/ P/ y( L) h. }8 f
const vertex_t* q)
! _8 R" I" H7 ^% s. C
{
/ f1 k4 _& D* n2 q5 _; b, `
double dx = l_end->x - l_start->x;
: M2 R- `* ?7 _0 |
double dy = l_end->y - l_start->y;
) N4 Z2 w# |# U" \" {/ J' d* g3 ]
double dx1= p->x - l_start->x;
% F+ i% J+ h& I6 k0 A% A
double dy1= p->y - l_start->y;
/ |+ H8 _7 @6 S
double dx2= q->x - l_end->x;
" z$ X0 |' Q/ ~# ~
double dy2= q->y - l_end->y;
& {. J+ e; k( ?9 A2 M7 n- ]. c
return ((dx*dy1-dy*dx1)*(dx*dy2-dy*dx2) > 0? 1 : 0);
0 [ ]0 ?8 j" Z
}
6 T5 |) h, i+ H
/* 2 line segments (s1, s2) are intersect? */
- H1 k# ]" T0 x/ }1 [/ N7 _& J! {
static int is_intersect(const vertex_t* s1_start, const vertex_t* s1_end,
6 P6 T% z" z, Z, B& Q
const vertex_t* s2_start, const vertex_t* s2_end)
7 I" o ^- K9 E6 W+ n
{
3 J5 b3 C3 J+ M8 g' L2 x
return (is_same(s1_start, s1_end, s2_start, s2_end)==0 &&
0 Y: G; H* V7 X1 Z9 ~
is_same(s2_start, s2_end, s1_start, s1_end)==0)? 1: 0;
( c5 l; ]8 W& x# j! z3 Q
}
! p" i& p8 R" \5 p, I ^( F
+ s* Z3 a6 W8 m9 E! U
7 d: V' j: G% R: b8 o
下面的函数pt_in_poly就是判断点(v)是(1)否(0)在多边形(vl：np)内的程序：
# P' ~0 I0 K, x) }: x/ u, {$ J
8 o; E$ ?3 N# {* J+ r8 F* R8 Z' i/ j
以下是引用片段：
. I0 d: F' g }7 r
int pt_in_poly ( const vertex_t* vl, int np, /* polygon vl with np vertices */
) R6 W" {& M# f+ Q
const vertex_t* v)
# n; e1 @$ c( N2 {3 t
{
6 u f, X; l& e* |/ L
int i, j, k1, k2, c;
4 [/ {# g) D5 w- i6 b6 p" V
rect_t rc;
) h: ^! Y0 V6 l# M. y+ q' i: z3 G
vertex_t w;
+ `' l: ~ W& ~
if (np < 3)
% M. r0 I& F0 U% i) i+ L6 A0 v: u! g
return 0;
8 ~3 E$ v. ~0 ^ W& H6 \5 c2 v
vertices_get_extent(vl, np, &rc);
/ @8 t8 t K- A3 W: G) ?
if (v->x < rc.min_x || v->x > rc.max_x || v->y < rc.min_y || v->y > rc.max_y)
4 S8 A5 W. w% y, F2 p; x) L$ H
return 0;
" j c# W9 R7 L0 G
/* Set a horizontal beam l(*v, w) from v to the ultra right */
! S, @' \) e! I+ I
w.x = rc.max_x + DBL_EPSILON;
, @; ^: v; Z+ k8 a! R
w.y = v->y;
$ Z7 U5 h. |! p7 M
c = 0; /* Intersection points counter */
$ ?( d: T0 t6 h u7 M. d$ W# f! M
for(i=0; i
, z; `+ Y( c7 \; Q- q5 `; T. r
{
+ i; h" V& ], M: J
j = (i+1) % np;
% h, D, V* M$ v$ j1 W% X- ` f
if(is_intersect(vl+i, vl+j, v, &w))
$ Q" F( v t' ~8 s0 u
{
' n( J5 G( U% J5 S0 q
C++;
$ R( {; I% h8 u( r) z5 l
}
( F4 N3 L3 b- u6 q \& ?1 F- u
else if(vl
.y==w.y)
0 h5 ?( e+ }; @' i( x2 v
{
# o' a( D( t; _0 l) Y
k1 = (np+i-1)%np;
) C4 Y/ G( _# `0 o. |: K. v) o
while(k1!=i && vl[k1].y==w.y)
, ~4 w2 K5 t) _0 P
k1 = (np+k1-1)%np;
( v+ E+ r$ ?- C
k2 = (i+1)%np;
7 F& B/ f: P; @0 k9 t
while(k2!=i && vl[k2].y==w.y)
$ P* V8 G/ ^4 h: B
k2 = (k2+1)%np;
) e6 q/ \( ~% T$ n8 Z+ F1 w/ y
if(k1 != k2 && is_same(v, &w, vl+k1, vl+k2)==0)
) ~0 q* J- ^+ D2 m
C++;
2 v: a/ Q6 f8 j8 `/ s
if(k2 <= i)
9 C& X7 C3 Q3 B% M
break;
: x' c3 ?( W$ d$ S5 F% ~2 z6 z; z
i = k2;
7 u+ q* M( K" ~( R7 J5 m: X
}
8 R! S# B1 l" h: h# _2 X7 |
}
1 x, W: v# j' S$ U' W* Z
return c%2;
) F, [; t' `/ s1 U9 b
}
9 K6 Y$ U& e" N; x: `' M
4 z' y* A7 {4 ^/ l) r2 o
7 E( B n7 F3 F! M8 E; t
本想配些插图说明问题，但是，CSDN的文章里放图片我还没用过。以后再试吧!实践证明，本程序算法的适应性极强。但是，对于点正好落在多边形边上的极端情形，有可能得出2种不同的结果。
欢迎光临 捌玖网络工作室 (http://www.89w.org/)
Powered by Discuz! 7.2
