C语言中显示点在多边形内算法

本文是采用射线法判断点是否在多边形内的C语言程序。多年前，我自己实现了这样一个算法。但是随着时间的推移，我决定重写这个代码。参考周培德的《计算几何》一书，结合我的实践和经验，我相信，在这个算法的实现上，这是你迄今为止遇到的最优的代码。

　　这是个C语言的小算法的实现程序，本来不想放到这里。可是，当我自己要实现这样一个算法的时候，想在网上找个现成的，考察下来竟然一个符合需要的也没有。我对自己大学读书时写的代码没有信心，所以，决定重新写一个，并把它放到这里，以飨读者。也增加一下BLOG的点击量。6 N1 z, @1 v4 M. \
$ j! d7 J1 [! o7 m
　　首先定义点结构如下：
. B4 S/ ]4 ]9 p; f4 [  ~
以下是引用片段：
　　/* Vertex structure */ 5 |1 w2 r6 U/ c, ]/ o- f( v
　　typedef struct $ P1 }; ^  J( ?
　　{ + |3 S3 b# |) H1 k
　　double x, y; 7 z5 a% \7 h2 j3 P1 ^# g
　　} vertex_t;

" a  K& k7 U6 b4 Z' c8 a
　　本算法里所指的多边形，是指由一系列点序列组成的封闭简单多边形。它的首尾点可以是或不是同一个点(不强制要求首尾点是同一个点)。这样的多边形可以是任意形状的，包括多条边在一条绝对直线上。因此，定义多边形结构如下：  C. d6 f9 S# \! x: e; w' F
! L4 L& ^9 u7 r9 x0 M
以下是引用片段：, K, t; j, V# j4 [
　　/* Vertex list structure – polygon */
　　typedef struct . N; R, m) D" B& M9 N
　　{ ' ^1 I9 A$ O( S1 \1 A
　　int num_vertices; /* Number of vertices in list */
　　vertex_t *vertex; /* Vertex array pointer */
　　} vertexlist_t;
0 |$ g  s$ f: n7 W7 n# s9 p: m

　　为加快判别速度，首先计算多边形的外包矩形(rect_t)，判断点是否落在外包矩形内，只有满足落在外包矩形内的条件的点，才进入下一步的计算。为此，引入外包矩形结构rect_t和求点集合的外包矩形内的方法vertices_get_extent，代码如下：6 H( v* R+ V9 n
' G; V! a2 g- h4 F: _# D' j& U0 M
以下是引用片段：2 l" N  H6 W$ C6 c- Y3 z
　　/* bounding rectangle type */
　　typedef struct
　　{ 3 x0 S+ M! D9 q% o! v
　　double min_x, min_y, max_x, max_y;
　　} rect_t;   r! D. V. w. ~8 \/ U! Z
　　/* gets extent of vertices */ # R8 ]8 x, e7 E6 Z
　　void vertices_get_extent (const vertex_t* vl, int np, /* in vertices */
　　rect_t* rc /* out extent*/ ) - o4 o% ?$ H* L  W1 `$ v
　　{
　　int i;
　　if (np > 0){
　　rc->min_x = rc->max_x = vl[0].x; rc->min_y = rc->max_y = vl[0].y;
　　}else{ 3 Z0 E& P6 X- X8 C" K* ~+ o
　　rc->min_x = rc->min_y = rc->max_x = rc->max_y = 0; /* =0 ? no vertices at all */ 9 @# F7 i5 n% P8 }1 \
　　} / h. E) x) P  p5 Q, m0 z
　　for(i=1; i  0 I# Y2 [, u" v9 H% ?* [; o
　　{ 1 O) W" Y/ g& F
　　if(vl.x < rc->min_x) rc->min_x = vl.x; $ b7 \% ?# I) P: t6 X
　　if(vl.y < rc->min_y) rc->min_y = vl.y; 9 X! j* y, N! C
　　if(vl.x > rc->max_x) rc->max_x = vl.x;
2 y2 E( y/ n- w　　if(vl.y > rc->max_y) rc->max_y = vl.y;   Q: A7 }9 e. {
　　} 2 C: f% `, R# H4 R$ W: R2 ^
　　} % \% H# z2 D( h3 a# X7 J- [6 ^

0 f' V7 t, j. ~, j" M/ I, n8 z, y% _
　　当点满足落在多边形外包矩形内的条件，要进一步判断点(v)是否在多边形(vl：np)内。本程序采用射线法，由待测试点(v)水平引出一条射线B(v，w)，计算B与vl边线的交点数目，记为c，根据奇内偶外原则(c为奇数说明v在vl内，否则v不在vl内)判断点是否在多边形内。
3 x# H) ^+ s$ k9 i' o$ N- U( d. @
1 U+ U9 L5 {6 Z9 [+ ]1 B8 f- z4 {　　具体原理就不多说。为计算线段间是否存在交点，引入下面的函数：
6 b8 W+ s$ u& k% U$ t' y* _; ?% `! [2 z, X! ^! I# J7 o! Y+ i
　　(1)is_same判断2(p、q)个点是(1)否(0)在直线l(l_start，l_end)的同侧;4 D; @7 ]$ V8 y. [$ |

9 i( g& G8 K- A! Y　　(2)is_intersect用来判断2条线段(不是直线)s1、s2是(1)否(0)相交;9 O9 _1 B8 b1 j. J: U
* m6 P. K: X- A4 k7 R7 V
以下是引用片段：/ t( h8 J0 T6 z$ u9 w! ?! y0 I" K0 e
　　/* p, q is on the same of line l */ ! W3 Q! Y; A/ C( W
　　static int is_same(const vertex_t* l_start, const vertex_t* l_end, /* line l */
$ B2 H( J% B; U　　const vertex_t* p, + L. G# l/ n5 X( S; N% w5 V
　　const vertex_t* q) . b" w; _: E2 d' f- G% Y, m/ ~
　　{ 4 x" P" l: f6 p/ q1 q- f6 R
　　double dx = l_end->x - l_start->x;
8 G8 G; w2 s/ ?9 G　　double dy = l_end->y - l_start->y; 0 R  M" P' {3 h) C, u( j
　　double dx1= p->x - l_start->x; 7 y" j9 `5 H0 Z$ {- J5 t1 Y
　　double dy1= p->y - l_start->y;
5 F' U/ e2 X+ S% Z* t, H　　double dx2= q->x - l_end->x; ' M- c% r( y, j) k6 [" P
　　double dy2= q->y - l_end->y; $ c4 R- ~* G( C7 T' W
　　return ((dx*dy1-dy*dx1)*(dx*dy2-dy*dx2) > 0? 1 : 0);
6 S5 n+ r! o, i/ ~" `& _' [3 m　　} ) V4 R% N& }- W; {7 Z# @3 k* u. }
　　/* 2 line segments (s1, s2) are intersect? */
1 g, K# q7 ^/ \! O　　static int is_intersect(const vertex_t* s1_start, const vertex_t* s1_end, ; _2 f' L5 G! O, X* f; w
　　const vertex_t* s2_start, const vertex_t* s2_end) 4 r7 ~) |/ l" d" I- u  m
　　{
8 [1 M6 r8 @1 q7 G　　return (is_same(s1_start, s1_end, s2_start, s2_end)==0 && $ J% k% l. Z: r( |' Z1 m
　　is_same(s2_start, s2_end, s1_start, s1_end)==0)? 1: 0; " y; P+ \! b1 v
　　}
9 \: p, Q& |3 b" l) [
# v* v  m. \  D: g: x+ i7 k# ?
/ b0 a+ E9 D9 w　　下面的函数pt_in_poly就是判断点(v)是(1)否(0)在多边形(vl：np)内的程序：, c5 n& q+ D+ O# ?+ l/ p, J3 k
. y; [6 x% H% d3 i8 Q
以下是引用片段：
# @+ E3 a7 W5 P0 v3 L　　int pt_in_poly ( const vertex_t* vl, int np, /* polygon vl with np vertices */
7 n% [. J  P9 w/ C. S; @; |　　const vertex_t* v)
. w6 H4 z1 {: h' N9 i: ?+ X　　{
) e6 T; `( ^' a& X# f　　int i, j, k1, k2, c; 7 k! V2 B3 W' e3 [' O0 @, v
　　rect_t rc;
* ^  ?  f# N/ w6 s: z/ l+ @　　vertex_t w;
. b- r+ c+ {6 G& m. |7 ?( g& F* C　　if (np < 3)
( y! V  C$ w: I5 Z( K. J3 s% \　　return 0;
* f+ i( x  p! G　　vertices_get_extent(vl, np, &rc);
  w2 k! d! w7 e　　if (v->x < rc.min_x || v->x > rc.max_x || v->y < rc.min_y || v->y > rc.max_y)
# K1 u" s4 A/ W1 [0 J2 d) ~　　return 0; ) {& K3 R- d8 a# M" }
　　/* Set a horizontal beam l(*v, w) from v to the ultra right */
* R9 a  E- N0 V- _　　w.x = rc.max_x + DBL_EPSILON; , v. A& f) a6 v) g* h9 c
　　w.y = v->y;
9 O$ Y+ Y" b7 B0 m/ H( m3 \( ?　　c = 0; /* Intersection points counter */ ! S- J! |& p" h2 O4 ]" t
　　for(i=0; i  1 H3 W' r- L* p4 O3 E7 H' ?, G
　　{ . z6 Z, F" a6 k- Z
　　j = (i+1) % np; 6 c$ V$ A& q1 p( R, P$ e
　　if(is_intersect(vl+i, vl+j, v, &w)) 2 ?# W! m/ {! @- k% A
　　{
) ?. B) g, U7 z; C+ {  P　　C++; : _: r+ M( k$ a
　　}
/ p2 B& t6 e5 }2 {; z0 n　　else if(vl.y==w.y)
) c* t0 l& }% ]) n( U6 Z　　{ 3 g% e! Z( K; Q+ M
　　k1 = (np+i-1)%np; ) U0 n5 \. U) d) \
　　while(k1!=i && vl[k1].y==w.y) . G) N% B) P6 D1 s& K" t9 t9 B
　　k1 = (np+k1-1)%np;
* ~) R4 r* k0 {　　k2 = (i+1)%np;
, H2 F) |" Q# ^, `　　while(k2!=i && vl[k2].y==w.y)
  ~3 R, b) ], W+ B+ g# P　　k2 = (k2+1)%np; , |: n' Z# ^+ R8 j
　　if(k1 != k2 && is_same(v, &w, vl+k1, vl+k2)==0)
/ ]& v4 O: f4 ^- _6 ~% w　　C++; / w! M8 \4 d1 }
　　if(k2 <= i) 7 o! [1 J$ ?+ m, S. d; T( G* }- p. Q
　　break; 0 t% j* u5 `! B2 d; M, V" ^
　　i = k2;
& J) I) A& _6 N& |0 B# _　　} ) u* R7 B4 H. o
　　} 0 c; }; B% h& i
　　return c%2; " o4 ?. n( n, W6 c7 w3 l+ M
　　}
6 M# \) i0 B- p* }# o
! l; R( E; z+ P7 w: q3 {4 |" E3 p$ z  K( ^1 ~
　　本想配些插图说明问题，但是，CSDN的文章里放图片我还没用过。以后再试吧!实践证明，本程序算法的适应性极强。但是，对于点正好落在多边形边上的极端情形，有可能得出2种不同的结果。

返回列表

C语言中显示点在多边形内算法

[收藏此主题] [关注此主题的新回复]

[通过 QQ、MSN 分享给朋友]

C语言中显示 点在多边形内 算法

[收藏此主题] [关注此主题的新回复]

[通过 QQ、MSN 分享给朋友]

C语言中显示点在多边形内算法