lib/geometry/polygon.nit

   1 # This file is part of NIT (http://www.nitlanguage.org).
   2 #
   3 # Copyright 2015 Romain Chanoir <romain.chanoir@viacesi.fr>
   4 #
   5 # Licensed under the Apache License, Version 2.0 (the "License");
   6 # you may not use this file except in compliance with the License.
   7 # You may obtain a copy of the License at
   8 #
   9 #     http://www.apache.org/licenses/LICENSE-2.0
  10 #
  11 # Unless required by applicable law or agreed to in writing, software
  12 # distributed under the License is distributed on an "AS IS" BASIS,
  13 # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14 # See the License for the specific language governing permissions and
  15 # limitations under the License.
  16
  17 # Convex Polygons manipulations
  18 module polygon
  19
  20 import points_and_lines
  21
  22 # Abstraction of a Polygon
  23 abstract class APolygon
  24         # Vertices of this polygon
  25         var points: Array[Point[Float]]
  26
  27         init do assert points.length >= 3
  28
  29         # Get an array of the x coordinates of the vertices
  30         private fun x_coordinates: Array[Float] do
  31                 return [for p in points do p.x]
  32         end
  33
  34         # Get an array of the y coordinates of the vertices
  35         private fun y_coordinates: Array[Float] do
  36                 return [for p in points do p.y]
  37         end
  38
  39         # Get a matrice containing the coordinates of the vertices
  40         private fun vertices: Array[Array[Float]] do
  41                 var vertices = new Array[Array[Float]]
  42                 for i in [0..points.length[ do
  43                         var temp = new Array[Float]
  44                         temp.add(points[i].x)
  45                         temp.add(points[i].y)
  46                         vertices.add(temp)
  47                 end
  48                 return vertices
  49         end
  50
  51         # Returns the axes corresponding to the edges of the polygon, used for collision detection
  52         private fun axes: Array[Point[Float]] do
  53                 var axes = new Array[Point[Float]]
  54                 for i in [0..points.length[ do
  55                         var p1 = new Point[Float](points[i].x, points[i].y)
  56                         var p2 = new Point[Float](points[(i+1) % points.length].x, points[(i+1) % points.length].y)
  57                         var edge = new Point[Float](p1.x - p2.x, p1.y - p2.y)
  58                         var normal = new Point[Float](-edge.y, edge.x)
  59                         axes[i] = normal
  60                 end
  61                 return axes
  62         end
  63
  64         # Sort the vertices in counter clockwise order
  65         #
  66         # ~~~
  67         # var p1 = new Point[Float](0.0, 0.0)
  68         # var p2 = new Point[Float](5.0, 0.0)
  69         # var p3 = new Point[Float](0.0, 5.0)
  70         # var p4 = new Point[Float](5.0, 5.0)
  71         # var arr = new Array[Point[Float]].with_items(p1, p2, p3, p4)
  72         # var poly = new ConvexPolygon(arr)
  73         # poly.sort_ccw
  74         # assert poly.points == [p4, p2, p1, p3]
  75         # ~~~
  76         fun sort_ccw do
  77                 var sorter = new CounterClockWiseSort.with_center(vertices)
  78                 sorter.sort(points)
  79         end
  80
  81         # Sort the vertices in clockwise order
  82         #
  83         # ~~~
  84         # var p1 = new Point[Float](0.0, 0.0)
  85         # var p2 = new Point[Float](5.0, 0.0)
  86         # var p3 = new Point[Float](0.0, 5.0)
  87         # var p4 = new Point[Float](5.0, 5.0)
  88         # var arr = new Array[Point[Float]].with_items(p1, p2, p3, p4)
  89         # var poly = new ConvexPolygon(arr)
  90         # poly.sort_cw
  91         # assert poly.points == [p3, p1, p2, p4]
  92         # ~~~
  93         fun sort_cw do
  94                 var sorter = new ClockWiseSort.with_center(vertices)
  95                 sorter.sort(points)
  96         end
  97
  98         # Is `self` convex ?
  99         #
 100         # ~~~
 101         # var p1 = new Point[Float](0.0, 0.0)
 102         # var p2 = new Point[Float](5.0, 0.0)
 103         # var p3 = new Point[Float](0.0, 5.0)
 104         # var p4 = new Point[Float](5.0, 5.0)
 105         # var arr = new Array[Point[Float]].with_items(p1, p2, p3, p4)
 106         # var poly = new ConvexPolygon(arr)
 107         # poly.sort_ccw
 108         # assert poly.is_convex
 109         # ~~~
 110         fun is_convex: Bool do
 111                 var prev = points[points.length - 2]
 112                 var curr = points[points.length - 1]
 113                 var next = points[0]
 114                 var is_ccw = turn_left(prev, curr, next)
 115                 for i in [1..points.length[ do
 116                         prev = curr
 117                         curr= next
 118                         next = points[i]
 119                         if turn_left(prev ,curr, next) != is_ccw then return false
 120                 end
 121                 return true
 122         end
 123
 124         # Generate a projection of an edge of the polygon on a given axis
 125         private fun project(axis: Point[Float]): Projection do
 126                 var min = axis.x * points[0].x + axis.y * points[0].y
 127                 var max = min
 128                 for i in [0..points.length[ do
 129                         var p = axis.x * points[i].x + axis.y * points[i].y
 130                         if p < min then min = p
 131                         if p > max then max = p
 132                 end
 133                 var projection = new Projection(min, max)
 134                 return projection
 135         end
 136
 137         # Remove  `p` from the vertices, keeping at least 3 vertices
 138         fun delete_vertex(p: Point[Float]) do
 139                 assert points.length > 3
 140                 points.remove(p)
 141         end
 142
 143         # Does `self` intersects with `other`
 144         fun intersects(other: APolygon): Bool is abstract
 145
 146         # Is `p` contained in `self` ?
 147         fun contain(p: Point[Float]): Bool is abstract
 148
 149         # Add a vertex to `self`
 150         fun add_vertex(p: Point[Float]): Bool do
 151                 points.add(p)
 152                 return true
 153         end
 154 end
 155
 156 # A simple polygon
 157 class Polygon
 158         super APolygon
 159         # TODO: implement algorithms for concave polygons ?
 160 end
 161
 162 # Convex Polygon class
 163 class ConvexPolygon
 164         super APolygon
 165
 166         # Does this polygon intersects `other` ?
 167         #
 168         # ~~~
 169         # var p1 = new Point[Float](0.0, 0.0)
 170         # var p2 = new Point[Float](5.0, 0.0)
 171         # var p3 = new Point[Float](0.0, 5.0)
 172         # var p4 = new Point[Float](5.0, 5.0)
 173         # var arr = new Array[Point[Float]].with_items(p1, p2, p3, p4)
 174         # var poly = new ConvexPolygon(arr)
 175         # poly.sort_ccw
 176         # p1 = new Point[Float](2.5, 2.5)
 177         # p2 = new Point[Float](7.5, 2.5)
 178         # p3 = new Point[Float](2.5, 7.5)
 179         # p4 = new Point[Float](7.5, 7.5)
 180         # arr = new Array[Point[Float]].with_items(p1, p2, p3, p4)
 181         # var poly2 = new ConvexPolygon(arr)
 182         # poly2.sort_ccw
 183         # assert poly.intersects(poly2)
 184         # ~~~
 185         redef fun intersects(other) do
 186                 assert is_convex
 187
 188                 var axes1 = axes
 189                 var axes2 = other.axes
 190                 for axis in axes1 do
 191                         var project1 = project(axis)
 192                         var project2 = other.project(axis)
 193                         if not project1.overlap(project2) then return false
 194                 end
 195                 for axis in axes2 do
 196                         var project1 = project(axis)
 197                         var project2 = other.project(axis)
 198                         if not project1.overlap(project2) then return false
 199                 end
 200                 return true
 201         end
 202
 203         # ~~~
 204         # var p1 = new Point[Float](0.0, 0.0)
 205         # var p2 = new Point[Float](5.0, 0.0)
 206         # var p3 = new Point[Float](0.0, 5.0)
 207         # var p4 = new Point[Float](5.0, 5.0)
 208         # var p5 = new Point[Float](2.5, 2.5)
 209         # var arr = new Array[Point[Float]].with_items(p1, p2, p3, p4)
 210         # var poly = new ConvexPolygon(arr)
 211         # poly.sort_ccw
 212         # assert poly.contain(p5)
 213         # ~~~
 214         redef fun contain(p) do
 215                 var prev = points[points.length - 1]
 216                 var curr = p
 217                 var next = points[0]
 218                 var is_ccw = turn_left(prev, curr, next)
 219                 for i in [1..points.length[ do
 220                         prev = next
 221                         next = points[i]
 222                         if turn_left(prev, curr, next) != is_ccw then return false
 223                 end
 224                 return true
 225         end
 226
 227         # Check if the order of the points in the polygon is counter-clockwise
 228         # The vertices in the polygon need to be sorted
 229         #
 230         # ~~~
 231         # var p1 = new Point[Float](0.0, 0.0)
 232         # var p2 = new Point[Float](5.0, 0.0)
 233         # var p3 = new Point[Float](0.0, 5.0)
 234         # var p4 = new Point[Float](5.0, 5.0)
 235         # var arr = new Array[Point[Float]].with_items(p1, p2, p3, p4)
 236         # var poly = new ConvexPolygon(arr)
 237         # poly.sort_ccw
 238         # assert poly.is_ccw
 239         # ~~~
 240         fun is_ccw: Bool do
 241                 var min = points[0].y
 242                 var min_index = 0
 243                 for i in [1..points.length - 1[ do
 244                         if points[i].y < min then
 245                                 min = points[i].y
 246                                 min_index = i
 247                         end
 248                 end
 249                 var prev = points[(min_index - 1 + points.length) % points.length]
 250                 var next = points[(min_index + 1) % points.length]
 251                 return not turn_left(prev, points[min_index], next)
 252         end
 253
 254
 255
 256         # Add a vertex to the polygon
 257         #
 258         # ~~~
 259         # var p1 = new Point[Float](0.0, 0.0)
 260         # var p2 = new Point[Float](5.0, 0.0)
 261         # var p3 = new Point[Float](0.0, 5.0)
 262         # var p4 = new Point[Float](5.0, 5.0)
 263         # var arr = new Array[Point[Float]].with_items(p1, p2, p3, p4)
 264         # var poly = new ConvexPolygon(arr)
 265         # var p5 = new Point[Float](2.5, 7.5)
 266         # assert poly.add_vertex(p5)
 267         # ~~~
 268         redef fun add_vertex(p) do
 269                 assert points.length >= 3
 270                 var temp_list = points.clone
 271                 temp_list.add(p)
 272                 var temp_polygon = new ConvexPolygon(temp_list)
 273                 temp_polygon.sort_ccw
 274                 if temp_polygon.is_convex then
 275                         points = temp_polygon.points
 276                         return true
 277                 else
 278                         return false
 279                 end
 280         end
 281 end
 282
 283 # Projection of an edge of a `ConvexPolygon` used for intersection test
 284 private class Projection
 285         var min: Float is writable
 286         var max: Float is writable
 287
 288         fun overlap(other: Projection): Bool do return not (min > other.max or other.min > max)
 289 end
 290
 291 private class PointXCompare
 292         super Comparator
 293
 294         redef type COMPARED: Point[Float]
 295
 296         redef fun compare(a,b) do
 297                 if a.x == b.x then
 298                         if a.y == b.y then return 0
 299                         if a.y > b.y then return - 1
 300                         return 1
 301                 else
 302                         if a.x > b.x then return -1
 303                         return 1
 304                 end
 305         end
 306 end
 307
 308 # Sorter for polygon vertices
 309 private abstract class PolygonSorter
 310         super Comparator
 311
 312         redef type COMPARED: Point[Float]
 313
 314         # Center of the polygon's points
 315         var center: COMPARED
 316
 317         # init calculating the center
 318         init with_center(pts : Array[Array[Float]]) do
 319                 init(calc_center(pts))
 320         end
 321
 322         # Calculate the center
 323         fun calc_center(pts : Array[Array[Float]]): COMPARED do
 324                 var sumx = 0.0
 325                 var sumy = 0.0
 326                 for ap in pts do
 327                         sumx += ap[0]
 328                         sumy += ap[1]
 329                 end
 330                 return new Point[Float](sumx / pts.length.to_f, sumy / pts.length.to_f)
 331         end
 332 end
 333
 334 #TODO: CounterClockWise and ClockWise sorts are broken in some cases of concave polygons
 335
 336 # Sort the vertices of a polygon in counter clockwise order
 337 private class CounterClockWiseSort
 338         super PolygonSorter
 339
 340         redef fun compare(a,b) do
 341                 if a.x == b.x and a.y == b.y then return 0
 342                 if a.x - center.x >= 0.0 and b.x - center.x < 0.0 then return -1
 343                 if a.x - center.x < 0.0 and b.x - center.x >= 0.0 then return 1
 344                 if a.x - center.x == 0.0 and b.x - center.x == 0.0 then
 345                         if a.y - center.y >= 0.0 or b.y - center.y >= 0.0 then
 346                                 if a.y > b.y then return -1
 347                                 return 1
 348                         end
 349                         if b.y > a.y then return -1
 350                         return 1
 351                 end
 352
 353                 var det = (a.x - center.x) * (b.y - center.y) - (b.x - center.x) * (a.y - center.y)
 354                 if det > 0.0 then return 1
 355                 if det < 0.0 then return -1
 356
 357                 var d1 = (a.x - center.x) * (a.x - center.x) + (a.y - center.y) * (a.y - center.y)
 358                 var d2 = (b.x - center.x) * (b.x - center.x) + (b.y - center.y) * (b.y - center.y)
 359                 if d1 > d2 then return -1
 360                 return 1
 361         end
 362 end
 363
 364 # Sort the vertices of a polygon in clockwise order
 365 private class ClockWiseSort
 366         super PolygonSorter
 367
 368         redef fun compare(a,b) do
 369                 if a.x == b.x and a.y == b.y then return 0
 370                 if a.x - center.x >= 0.0 and b.x - center.x < 0.0 then return 1
 371                 if a.x - center.x < 0.0 and b.x - center.x >= 0.0 then return -1
 372                 if a.x - center.x == 0.0 and b.x - center.x == 0 then
 373                         if a.y - center.y >= 0.0 or b.y - center.y >= 0.0 then
 374                                 if a.y > b.y then return 1
 375                                 return -1
 376                         end
 377                         if b.y > a.y then return 1
 378                         return -1
 379                 end
 380
 381                 var det = (a.x - center.x) * (b.y - center.y) - (b.x - center.x) * (a.y - center.y)
 382                 if det > 0.0 then return -1
 383                 if det < 0.0 then return 1
 384
 385                 var d1 = (a.x - center.x) * (a.x - center.x) + (a.y - center.y) * (a.y - center.y)
 386                 var d2 = (b.x - center.x) * (b.x - center.x) + (b.y - center.y) * (b.y - center.y)
 387                 if d1 > d2 then return 1
 388                 return -1
 389         end
 390 end
 391
 392 # Get the convex hull of the set of `points`
 393 #
 394 # ~~~
 395 # var p1 = new Point[Float](0.0, 0.0)
 396 # var p2 = new Point[Float](5.0, 0.0)
 397 # var p3 = new Point[Float](0.0, 5.0)
 398 # var p4 = new Point[Float](5.0, 5.0)
 399 # var p5 = new Point[Float](2.5, 2.5)
 400 # var arr = new Array[Point[Float]].with_items(p1, p2, p3, p4, p5)
 401 # var poly = convex_hull(arr)
 402 # assert poly.points == [p4, p3, p1, p2]
 403 # ~~~
 404 fun convex_hull(points: Array[Point[Float]]): ConvexPolygon do
 405         var sorter = new PointXCompare
 406         sorter.sort(points)
 407         var l = points.length
 408
 409         var pl1 = new Array[Point[Float]]
 410         var pl2 = new Array[Point[Float]]
 411         for i in [0..l[ do
 412                 while pl2.length >= 2 and not turn_left(pl2[pl2.length - 2], pl2[pl2.length - 1], points[i]) do
 413                         pl2.remove(pl2.last)
 414                 end
 415                 pl2.add(points[i])
 416         end
 417         var i = l - 1
 418         while i >= 0 do
 419                 while pl1.length >= 2 and not turn_left(pl1[pl1.length - 2], pl1[pl1.length - 1], points[i]) do
 420                         pl1.remove(pl1.last)
 421                 end
 422                 pl1.add(points[i])
 423                 i-= 1
 424         end
 425         pl1.remove_at(pl1.length - 1)
 426         pl2.remove_at(pl2.length -1)
 427         pl2.add_all(pl1)
 428         return new ConvexPolygon(pl2)
 429 end
 430
 431 # Is the angle between [p1,p2] and [p2,p3] going left (counter clockwise) or right (clockwise) ?
 432 #
 433 # ~~~
 434 # var p1 = new Point[Float](0.0, 0.0)
 435 # var p2 = new Point[Float](5.0, 0.0)
 436 # var p3 = new Point[Float](0.0, 5.0)
 437 # assert turn_left(p1, p2, p3)
 438 # ~~~
 439 fun turn_left(p1, p2, p3: Point[Float]): Bool do
 440         return ((p2.x - p1.x) * (p3.y - p2.y) - (p3.x - p2.x) * (p2.y - p1.y)) > 0.0
 441 end
 442
 443 # Split a polygon into triangles
 444 # Useful for converting a concave polygon into multiple convex ones
 445 fun triangulate(pts: Array[Point[Float]], results: Array[ConvexPolygon]) do
 446         var poly = new Polygon(pts)
 447         pts = poly.points
 448         recursive_triangulate(pts, results)
 449 end
 450
 451 private fun recursive_triangulate(pts: Array[Point[Float]], results: Array[ConvexPolygon]) do
 452         if pts.length == 3 then
 453                 results.add(new ConvexPolygon(pts))
 454                 return
 455         end
 456         var prev = pts[pts.length - 1]
 457         var curr = pts[0]
 458         var next = pts[1]
 459         for i in [1..pts.length[ do
 460                 if turn_left(prev, curr, next) then
 461                         prev = pts[i-1]
 462                         curr = next
 463                         if i+1 == pts.length then next = pts[pts.length - 1] else next = pts[i+1]
 464                         continue
 465                 end
 466                 var contains = false
 467                 var triangle = new ConvexPolygon(new Array[Point[Float]].with_items(prev, curr, next))
 468                 for p in pts do
 469                         if p != prev and p != curr and p != next then
 470                                 if triangle.contain(p) then
 471                                         contains = true
 472                                         break
 473                                 end
 474                         end
 475                 end
 476                 if not contains then
 477                         results.add(triangle)
 478                         pts.remove(curr)
 479                         recursive_triangulate(pts, results)
 480                         break
 481                 end
 482                 prev = pts[i-1]
 483                 curr = next
 484                 if i+1 == pts.length then next = pts[pts.length - 1] else next = pts[i+1]
 485         end
 486 end