giovedì 27 luglio 2017

Area di un poligono

Per calcolare l'area di un poligono la soluzione più logica e semplice sarebbe quella di triangolare il poligono e sommare le aree dei triangoli che compongono il poligono.


Questa soluzione risulta abbastanza intuitiva per poligoni convessi, meno per quelli concavi. Utilizzando però il giusto metodo si può sfruttare questa stessa soluzione a prescindere dal tipo di poligono.
Il classico calcolo dell'area di un triangolo $A=(b \times h)/2$, però, si rivela insufficiente in questo caso. Invece, dall'algebra lineare sappiamo che dati tre vertici $a, b, c$, sia il modulo del prodotto vettoriale ($(b - a)\times (c - a)$) che il determinante della matrice composta dalle coordinate (omogenee) dei tre vertici

\begin{equation}\label{eq:1}2A=\begin{vmatrix} a_{0} & a_{1} & 1\\ b_{0} & b_{1} & 1\\ c_{0} & c_{1} & 1\\ \end{vmatrix}=(b_0 - a_0)(c_1 - a_1) - (c_0 - a_0)(b_1 - a_1)\end{equation}

danno l'area (con segno nel caso del determinante) del parallelogramma composto dai vettori $(b - a)$ e $(c - a)$. Quindi, di fatto, il doppio dell'area del triangolo con vertici $a, b, c$. Per ulteriori dettagli si veda [3].
Usando quindi il determinante per il calcolo dell'area di un triangolo si ha che se i vertici sono in ordine antiorario l'area avrà segno positivo altrimenti sarà negativo (in caso di ordine orario).
Guardando la figura sotto si intuisce che il metodo di sommare le aree dei triangoli che compongono il poligono funziona anche per quelli concavi. L'area di $\triangle(v_0, v_1, v_2)$, ad esempio, avrà segno negativo e sarà sottratta all'area di $\triangle(v_0, v_2, v_3)$.

L'area di un poligono $P$ di $n$ vertici, convesso o concavo che sia, sarà quindi

\begin{equation}\label{eq:2}A_p=A_{\triangle}(v_0, v_1, v_2) + A_{\triangle}(v_0, v_2, v_3) + \dots + A_{\triangle}(v_0, v_{n-2}, v_{n-1})\end{equation}

Per l'implementazione in codice C/C++ risulta utile considerare i vertici del poligono come una lista doppiamente concatenata. A questo scopo viene utilizzata una struttura come nodo della lista per mantenere informazioni utili per ogni vertice (come l'indice del vertice corrente e le sue coordinate). Per il resto si segue semplicemente quanto visto per le formule \eqref{eq:1} e \eqref{eq:2}.

/*Per facilitare la stesura, la lettura e la comprensione del codice si è volutamente deciso di 
utilizzare un mix di C/C++ cercando di restare il più possibile all'interno del sottoinsieme 
C del C++ e di utilizzare solo quei costrutti e funzionalità di base del C++ che permettono 
una maggiore semplificazione.*/
 
#include <iostream>
#include <list>
 
typedef float Vertex_Coord[2];
 
struct vertex_node
{
 int vnum;
 Vertex_Coord pos;
};
 
float area2(const Vertex_Coordconst Vertex_Coordconst Vertex_Coord);
float area_poly2(std::list<vertex_node>&);
 
int main(int argccharargv[], charenv[])
{
 // poligono convesso
 Vertex_Coord poly_convex_coord[] = { {-1.0, -1.0}, {1.0, -1.0}, {1.0, 1.0}, {0.0, 2.0}, {-1.0, 1.0} };
 int numv = sizeof(poly_convex_coord) / (sizeof(float) * 2);
 
 std::list<vertex_node> vertex_list;
 
 for (int i = 0; i < numv; ++i)
 {
  vertex_list.push_back({ i, {poly_convex_coord[i][0], poly_convex_coord[i][1]} });
 }
 
 std::cout << "Area del poligono convesso: " << area_poly2(vertex_list) / 2 << std::endl;
 
 // poligono concavo
 Vertex_Coord poly_concave_coord[] = { {-1.0, -1.0}, {0.0, 1.0}, {1.0, 1.0}, {0.0, 2.0}, {-1.0, 1.0} };
 numv = sizeof(poly_concave_coord) / (sizeof(float) * 2);
 
 vertex_list.clear();
 
 for (int i = 0; i < numv; ++i)
 {
  vertex_list.push_back({ i, {poly_concave_coord[i][0], poly_concave_coord[i][1]} });
 }
 
 std::cout << "Area del poligono concavo: " << area_poly2(vertex_list) / 2 << std::endl;
 
 return 0;
}
 
float area2(const Vertex_Coord aconst Vertex_Coord bconst Vertex_Coord c)
{
 return (b[0] - a[0]) * (c[1] - a[1]) - (c[0] - a[0]) * (b[1] - a[1]);
}
 
float area_poly2(std::list<vertex_node> &vn)
{
 int sum = 0;
 std::list<vertex_node>::const_iterator c0 = vn.cbegin();
 
 for (std::list<vertex_node>::const_iterator ci = vn.cbegin(); ci != vn.cend(); ci++)
 {
  if (ci == c0)
  {
   ci++;
   Vertex_Coord vcp1 = { ci->pos[0], ci->pos[1] };
   ci++;
   Vertex_Coord vcp2 = { ci->pos[0], ci->pos[1] };
   sum += area2(c0->pos, vcp1, vcp2);
  }
  else
  {
   ci--;
   Vertex_Coord vcp1 = { ci->pos[0], ci->pos[1] };
   ci++;
   Vertex_Coord vcp2 = { ci->pos[0], ci->pos[1] };
   sum += area2(c0->pos, vcp1, vcp2);
  }
 }
 
 return sum;
}


Riferimenti:
[1] Computational Geometry in C - O'Rourke
[2] Elementi di Geometria Computazionale
[3] Collinearità di tre punti in 2D

Nessun commento:

Posta un commento