giovedì 17 agosto 2017

Convex Hull in 3D

Problema: Dato un insieme di $n$ punti nello spazio, trovare il più piccolo involucro convesso che contiene gli $n$ punti. Trovare cioè i punti, tra gli $n$ dell'insieme dato, che si trovano sui bordi del più piccolo poliedro convesso che contiene tutti gli altri.

Soluzione: Incrementale
Complessità: $O(n^2)$




Realizzare una implementazione robusta dell'involucro (o inviluppo) convesso in 3D richiede parecchio codice per gestire casi degeneri come punti tutti complanari, rimozione dei punti coincidenti o che giacciono sui bordi, ecc. Per questo motivo la cosa più importante, prima di iniziare, è quella di comprendere il problema con un esempio semplice ma concreto, che eviti tutti i casi degeneri e si concentri sugli aspetti più importanti. A quel punto gestire gli altri dettagli risulterà più agevole.
L'esempio in questo breve articolo è preso integralmente da [1], dove vengono spiegati bene anche la teoria e i dettagli implementativi. Si faccia riferimento ad esso per maggiori informazioni.

Si prenderà in esame la costruzione dell'involucro convesso di un insieme di vertici che rappresentano un cubo ed il risultato sarà il cubo stesso. Gli stessi vertici dati in input verranno restituiti in output. Se il risultato in sé non è molto interessante lo è, al contrario, il modo attraverso il quale l'algoritmo incrementale costruisce l'involucro convesso e restituisce il risultato. Per mantenere una certa coerenza nell'output e nelle strutture di dati utilizzate si richiederà che i vertici delle facce del poliedro convesso vengano conservati in senso antiorario (volendo anche in senso opposto va bene, l'importante è stabilire un ordine e mantenerlo).
L'input è costituito dai vertici di un cubo unitario:

$v_0 = (0, 0, 0)$, $v_1 = (0, 1, 0)$, $v_2 = (1, 1, 0)$, $v_3 = (1, 0, 0)$,
$v_4 = (0, 0, 1)$, $v_5 = (0, 1, 1)$, $v_6 = (1, 1, 1)$, $v_7 = (1, 0, 1)$

Prima di tutto sarebbe utile vedere, anche molto sommariamente, come sono organizzati i dati e le strutture dati che li mantengono. Si vedano i commenti alle strutture Vertex, Edge e Face nel codice fornito in fondo a questo articolo o, in alternativa, [1]. La cosa importante da tenere a mente è che l'involucro convesso verrà creato in modo incrementale tramite la costruzione di triangoli. Ogni triangolo ha due facce ma a parte il caso iniziale si è interessati solo alla faccia esterna dell'involucro convesso. Ogni lato di un triangolo ha due punti estremi ed è condiviso da due facce adiacenti (esterne ed appartenenti a triangoli differenti). Ogni faccia, invece, ha tre lati e tre vertici.

Il primo passaggio (metodo double_triangle) prevede di prendere i primi tre punti che non giacciono sulla stessa retta (metodo collinear) e costruire un triangolo con due facce (metodo make_face). Questo è l'unico caso in cui si prenderanno in considerazione entrambe le facce di un triangolo. Queste devono avere i vertici in ordine antiorario quindi, prendendo la figura A come esempio, la faccia $f_0$ (sotto il piano $xy$) avrà vertici $v_0, v_1, v_2$, mentre la faccia $f_1$ (sopra il piano $xy$) avrà vertici $v_2, v_1, v_0$. Costruite queste due facce l' "involucro convesso" (tra virgolette) per questi primi tre vertici è stato creato e non resta che individuare il successivo punto non complanare ai tre precedenti (metodo volume_sign) e procedere con la costruzione incrementale dell'involucro convesso per questi 4 punti. Nella figura B si può vedere che il punto $v_4$ è stato selezionato perché $v_3$ è complanare a $v_0, v_1, v_2$. Prima di vedere come procedere nella costruzione iterativa è utile illustrare un paio di dettagli importanti.
I metodi collinear e volume_sign sfruttano quanto visto in [2] (cioè la relazione tra moduli dei prodotti vettoriale e misto) per il calcolo, rispettivamente, dell'area del triangolo e del volume del poliedro individuato dai vertici inviati in input.
Il metodo make_face si occupa della creazione di una faccia del primo triangolo a due facce (infatti viene chiamato due volte), dei relativi lati e della registrazione, per ognuno di questi lati, di una delle due facce adiacenti al lato (in questo caso si tratta semplicemente della faccia stessa appena creata). Questo primo passaggio, gestito dal metodo double_triangle, si conclude con la registrazione dell'altra faccia adiacente a ogni lato. Essendoci per il momento solo due facce si tratta praticamente di un collegamento incrociato tra i sei lati e le due facce di questo primo triangolo.

Trovati 4 punti non complanari si procede alla costruzione incrementale del poliedro convesso con il metodo construct_hull. Si iterano tutti i punti rimasti e si controlla, per ognuno di essi, se è il caso di aggiungerlo (se esterno all'involucro convesso corrente) o rimuoverlo (se è interno) dalla lista dei vertici che fanno parte del nuovo involucro convesso.
Il metodo add_one (attraverso volume_sign e quanto visto in [2]) marca tutte le facce visibili dal punto nell'iterazione corrente. Se nessuna delle facce è visibile il punto corrente è evidentemente interno all'involucro convesso e può essere marcato per essere eliminato. Altrimenti ad essere marcati per la rimozione sono i lati che hanno entrambe le facce adiacenti visibili. Qualsiasi altro bordo con una sola faccia visibile verrà marcato come bordo di base per la costruzione di una nuova faccia del nuovo involucro convesso. Infatti (prendendo come esempio la figura B) il vertice $v_4$ vede i vertici della faccia $f_0$ in senso orario e quelli della faccia $f_1$ in senso antiorario, marcando quindi solo quest'ultima come visibile da $v_4$. Tutti e tre i lati di $f_1$ hanno infatti una faccia visibile ($f_1$ stessa) e l'altra no ($f_0$) e vengono usati come base per costruire le nuove facce dei triangoli del nuovo involucro convesso. Essendoci un singolo vertice all'apice delle nuove facce create viene a formarsi una sorta di cono.
E' infatti il metodo make_cone_face che si occupa di creare queste nuove facce che formano un cono. Questo metodo ha un compito importante e, allo stesso tempo, delicato all'interno dell'algoritmo. In primo luogo, come detto, ogni faccia deve avere tre lati ed ogni lato deve avere due facce adiacenti. La creazione di una nuova faccia deve evitare quindi che venga creato due volte lo stesso lato per le sue due facce adiacenti. Per fare un esempio concreto nella figura B quando make_cone_face crea una nuova faccia (mettiamo il caso che crei prima $v_0, v_2, v_4$) e successivamente va a creare la faccia adiacente (in questo caso $v_1, v_4, v_2$) non deve ricreare il lato (in questo caso $v_2 - v_4$) due volte. Il campo duplicate della struttura Vertex è usato a questo scopo e registra se il lato della nuova faccia che si sta creando è già stato creato. Guardando ancora la figura B si nota che una volta creata la faccia $v_0, v_2, v_4$ make_cone_face può controllare il campo duplicate di entrambi i vertici del lato $v_1 - v_2$ della nuova faccia $v_1, v_4, v_2$ e verificare che il bordo $v_2 - v_4$ è già stato creato. Avendo quindi già il lato base make_cone_face costruisce gli altri due lati della nuova faccia solo nel caso in cui entrambi non siano già stati creati in precedenza. In secondo luogo make_cone_face si occupa di mantenere la stessa orientazione (metodo make_ccw) della faccia marcata come visibile in add_one. Tornando per un momento alla figura A, $f_1$ è l'unica faccia visibile da $v_4$ i cui lati (tutti e tre) hanno solo una faccia adiacente visibile ($f_1$ stessa, $f_0$ non è visibile da $v_4$) e quindi andranno a costituire i lati base di tre nuove facce. Passando alla figura B, ognuna di queste tre nuove facce costruite da make_cone_ace manterrà la stessa orientazione dei vertici di $f_1$ (antioraria). A questo scopo, un controllo incrociato tra gli estremi del lato base ed i vertici della faccia visibile indica l'ordine corretto. Infine make_cone_face imposta la faccia adiacente di uno solo dei lati creati, anche quando i lati creati sono due. Questo perché l'ultima faccia del cono di triangoli creata avrà solo un lato con adiacenza da impostare. Ad esempio, guardando la figura B, dopo la costruzione delle facce $v_0, v_2, v_4$ e $v_1, v_4, v_2$, l'unico lato che necessita di impostare la faccia adiacente è uno solo tra $v_0 - v_4$ e $v_1 - v_4$.
Non resta che ripulire gli array che mantengono facce, lati e vertici da tutti gli elementi residui e ridondanti. Sembrerebbe una operazione banale (dato che sono stati marcati vertici interni, facce visibili e i bordi da eliminare) tuttavia è necessaria maggiore accuratezza perché, ad esempio, finora si sono marcati i vertici interni all'involucro convesso ma non quelli che si trovano sui bordi da cancellare. Inoltre si è impostata l'adiacenza dei due lati creati da make_cone_face ma non del lato base che, a causa della creazione della nuova faccia, ha ora una sola faccia adiacente valida. L'altra è stata marcata come visibile e verrà presto rimossa. Certo, sarebbe anche stato possibile impostare questa adiacenza insieme alle altre in make_cone_face ma a livello di responsabilità e competenza del codice sarebbe stata un piccola violazione. make_cone_face si occupa di costruire ed aggiungere dati alle strutture esistenti e non di ripulire quelle esistenti dai dati superflui. Per questi motivi è preferibile affidare questo compito al metodo clean_edge (tramite il campo newface della struttura Edge) che si occupa di cancellare i bordi. Tutto ciò però ha una conseguenza: si rivela necessario cancellare i bordi prima delle facce. Infatti per impostare la faccia adiacente dei lati base è necessario aver accesso a quest'ultime dato che è attraverso le facce adiacenti che si può verificare qual'è quella visibile da sostituire. Dopo aver impostato le facce adiacenti dei lati base è possibile rimuovere i bordi marcati (in add_one) come eliminabili (metodo clean_edges). Si passa quindi, alla rimozione delle facce marcate come visibili (metodo clean_faces). Infine c'è la rimozione dei vertici (metodo clean_vertices) che sfrutta i bordi non cancellati per marcare i vertici sull'involucro convesso e quindi cancellare solo i vertici finora processati ma che non sono sull'involucro convesso. Infine, sempre il metodo clean_vertices, reimposta i campi di Vertex pronti per una nuova iterazione in construct_hull.

Per chiarire brevemente quanto detto finora si prenda come esempio il passaggio dalla figura C alla figura D. L'involucro convesso costruito all'iterazione precedente è quella della figura C e comprende i vertici $v_0, v_1, v_2, v_4, v_5$. L'iterazione corrente seleziona $v_6$ e, dato che questo vertice vede solo la faccia individuata dai vertici $v_4, v_2, v_5$, vengono costruite tre facce con basi i lati di questa faccia e con stessa orientazione antioraria. Alla fine dell'iterazione che porta alla figura D la faccia visibile da $v_6$ è l'unico elemento che viene cancellato. Nessun bordo o vertice viene eliminato.

Per quanto riguarda l'implementazione in codice C/C++ è utile notare che nei metodi vertex_list_delete, vertex_list_delete e vertex_list_delete il puntatore al primo elemento punta all'elemento successivo se l'elemento da eliminare è il primo. Per questo motivo la gestione dell'eliminazione degli elementi da parte dei metodi che li utilizzano (in particolare clean_vertices, clean_edges e clean_faces) richiede particolare attenzione visto che la condizione di fine ciclo prevede che l'iterazione finisca quando si incontra nuovamente il primo elemento della lista. Ma se il primo elemento da eliminare è proprio il primo della lista il puntatore al primo elemento della lista punterà al secondo elemento e la condizione di fine ciclo verrà soddisfatta immediatamente perché all'iterazione successiva il secondo elemento sarà uguale al puntatore al primo elemento della lista (che ora, dopo l'eliminazione del primo elemento, è proprio il secondo). Per ovviare a questo inconveniente basta scorrere preliminarmente la lista per rimuovere tutti gli elementi marcati come eliminabili che si trovano all'inizio della stessa fino a quando non si incontra un elemento non marcato (il puntatore al primo elemento ora punterà proprio a questo elemento). Si salva un puntatore all'elemento successivo a quello non marcato e così si può procede come al solito sui restanti elementi, sicuri di aver iterato su tutti quelli della lista.

/*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 <cfloat>
#include <cmath>
 
typedef float Vertex_Coord[3];
struct Edge;
struct Face;
 
// Vertice
struct Vertex
{
 int vnum;          // indice del vertice
 Edge *duplicate;   // non null se lato che ha estremo questo vertice è già stato creato 
 Vertex_Coord pos;  // coordinate vertice
 Vertex *prev;
 Vertex *next;
 bool onhull;       // vero se vertice su involusco convesso
 bool mark;         // vero se vertice è stato processato dall'algoritmo
};
 
// Lato
struct Edge
{
 Face *adjface[2];    // puntatori alle due facce adiacenti
 Vertex *endpts[2];   // puntatori ai due estremi
 Face *newface;       // non null se lato è base di nuovo triangolo
 Edge *prev;
 Edge *next;
 bool remove;         // vero se lato si può cancellare; solo se entrambe le facce adiacenti sono visibili
};
 
// Faccia
struct Face
{
 Edge *edge[3];       // puntatori ai 3 lati
 Vertex *vertex[3];   // puntatori ai 3 vertici
 Face *prev;
 Face *next;
 bool visible;        // vero se è visibile da vertice processato all'iterazione corrente
};
 
Edge* make_null_edge(Edge*& elist);
Face* make_null_face(Face*& flist);
void double_triangle(Vertex*& vlistEdge*& elistFace*& flist);
Face* make_face(Vertexv0Vertexv1Vertexv2FacefoldEdge*& elistFace*& flist);
int volume_sign(Face*, Vertex*);
void construct_hull(Vertex*&, Edge*&, Face*&);
void add_one(VertexpEdge*& elistFace*& flist);
Face* make_cone_face(EdgeeVertexpEdge*& elistFace*& flist);
void make_ccw(Face*, Edge*, Vertex*);
void clean_up(Vertex*&, Edge*&, Face*&);
void clean_edges(Edge*&);
void clean_faces(Face*&);
void clean_vertices(Vertex*&, Edge*);
void print_faces(Face*);
void print_edges(Edge*);
bool collinear(const Vertex*, const Vertex*, const Vertex*);
bool AlmostEqual(floatfloatfloat);
 
void print_list(Vertex*);
 
inline float max_3(float afloat bfloat c)
{
 return ((a > b) ? ((a > c) ? a : c) : (b > c) ? b : c);
}
 
template<typename Tvoid list_insert(T*& headTp)
{
 if (head && head->next != nullptr)
 {
  p->next = head;
  p->prev = head->prev;
  head->prev = p;
  p->prev->next = p;
 }
 else
 {
  head = p;
  head->next = head->prev = p;
 }
}
 
template<typename Tvoid list_delete(T*& headTp)
{
 if (head)
 {
  if (head == head->next)
   head = nullptr;
  else if (p == head)
   head = head->next;
 
  p->next->prev = p->prev;
  p->prev->next = p->next;
 }
 
 delete p;
}
 
int main(int argccharargv[], charenv[])
{
 Vertex_Coord coord[] = {
 {0.0, 0.0, 0.0}, {0.0, 10.0, 0.0}, {10.0, 10.0, 0.0}, {10.0, 0.0, 0.0},
 {0.0, 0.0, 10.0}, {0.0, 10.0, 10.0}, {10.0, 10.0, 10.0}, {10.0, 0.0, 10.0} };
 int numv = sizeof(coord) / (sizeof(float) * 3);
 
 Vertex *vlist = nullptr;
 Edge *elist = nullptr;
 Face *flist = nullptr;
 
 for (int i = 0; i < numv; i++)
 {
  Vertex *tmpv = new Vertex;
  tmpv->vnum = i;
  tmpv->pos[0] = coord[i][0];
  tmpv->pos[1] = coord[i][1];
  tmpv->pos[2] = coord[i][2];
  tmpv->duplicate = nullptr;
  tmpv->onhull = false;
  tmpv->mark = false;
  tmpv->next = nullptr;
  tmpv->prev = nullptr;
 
  list_insert<Vertex>(vlist, tmpv);
 }
 
 double_triangle(vlist, elist, flist);
 construct_hull(vlist, elist, flist);
 print_faces(flist);
 
 return 0;
}
 
Edge* make_null_edge(Edge *&elist)
{
 Edge *tmpe = new Edge;
 tmpe->adjface[0] = tmpe->adjface[1] = tmpe->newface = nullptr;
 tmpe->endpts[0] = tmpe->endpts[1] = nullptr;
 tmpe->remove = false;
 list_insert<Edge>(elist, tmpe);
 return tmpe;
}
 
Face* make_null_face(Face *&flist)
{
 Face *tmpf = new Face;
 tmpf->edge[0] = tmpf->edge[1] = tmpf->edge[2] = nullptr;
 tmpf->vertex[0] = tmpf->vertex[1] = tmpf->vertex[2] = nullptr;
 tmpf->visible = false;
 list_insert<Face>(flist, tmpf);
 return tmpf;
}
 
void double_triangle(Vertex *&vlistEdge *&elistFace *&flist)
{
 Vertex *v0, *v1, *v2, *v3;
 Face *f0, *f1 = nullptr;
 int vol;
 
 v0 = vlist;
 
 // itera finchè non trova tre vertici non collineari
 while (collinear(v0, v0->next, v0->next->next))
  if ((v0 = v0->next) == vlist)
  {
   std::cout << "Tutti i punti sono collineari!" << std::endl;
   exit(1);
  }
 
 v1 = v0->next;
 v2 = v1->next;
 v0->mark = v1->mark = v2->mark = true;
 
 f0 = make_face(v0, v1, v2, f1, elistflist);
 f1 = make_face(v2, v1, v0, f0, elistflist);
 
 f0->edge[0]->adjface[1] = f1;
 f0->edge[1]->adjface[1] = f1;
 f0->edge[2]->adjface[1] = f1;
 f1->edge[0]->adjface[1] = f0;
 f1->edge[1]->adjface[1] = f0;
 f1->edge[2]->adjface[1] = f0;
 
 v3 = v2->next;
 vol = volume_sign(f0, v3);
 
 // itera finchè non trova quarto vertice non complanare ai precedenti 3
 while (!vol)
 {
  if ((v3 = v3->next) == v0)
  {
   std::cout << "Tutti i punti sono collineari!" << std::endl;
   exit(1);
  }
  vol = volume_sign(f0, v3);
 }
 
 // trovato vertice non complanare lo imposta come primo della lista
 vlist = v3;
}
 
Face* make_face(Vertex *v0Vertex *v1Vertex *v2Face *foldEdge *&elistFace *&flist)
{
 Face *f;
 Edge *e0, *e1, *e2;
 
 // se fold null crea nuovi lati 
 if (!fold)
 {
  e0 = make_null_edge(elist);
  e1 = make_null_edge(elist);
  e2 = make_null_edge(elist);
 }
 else // alrimenti usa quelli di fold in ordine contrario
 {
  e0 = fold->edge[2];
  e1 = fold->edge[1];
  e2 = fold->edge[0];
 }
 
 e0->endpts[0] = v0;
 e0->endpts[1] = v1;
 e1->endpts[0] = v1;
 e1->endpts[1] = v2;
 e2->endpts[0] = v2;
 e2->endpts[1] = v0;
 
 f = make_null_face(flist);
 f->edge[0] = e0; f->edge[1] = e1; f->edge[2] = e2;
 f->vertex[0] = v0; f->vertex[1] = v1; f->vertex[2] = v2;
 
 e0->adjface[0] = e1->adjface[0] = e2->adjface[0] = f;
 
 return f;
}
 
int volume_sign(Face *fVertex *p)
{
 float vol;
 float ax, ay, az, bx, by, bz, cx, cy, cz;
 
 ax = f->vertex[0]->pos[0] - p->pos[0];
 ay = f->vertex[0]->pos[1] - p->pos[1];
 az = f->vertex[0]->pos[2] - p->pos[2];
 bx = f->vertex[1]->pos[0] - p->pos[0];
 by = f->vertex[1]->pos[1] - p->pos[1];
 bz = f->vertex[1]->pos[2] - p->pos[2];
 cx = f->vertex[2]->pos[0] - p->pos[0];
 cy = f->vertex[2]->pos[1] - p->pos[1];
 cz = f->vertex[2]->pos[2] - p->pos[2];
 
 vol = ax * (by*cz - bz * cy) + ay * (bz*cx - bx * cz) + az * (bx*cy - by * cx);
 
 if (vol > 0.5)
  return 1;
 else if (vol < -0.5)
  return -1;
 else
  return 0;
}
 
void construct_hull(Vertex *&vlistEdge *&elistFace *&flist)
{
 Vertex *v, *vn;
 v = vlist;
 
 // itera su tutti i vertici
 do {
  vn = v->next;
 
  if (!(v->mark))
  {
   v->mark = true;
   add_one(v, elistflist);
   clean_up(vlistelistflist);
  }
  v = vn;
 } while (v != vlist);
}
 
void add_one(Vertex *pEdge *&elistFace *&flist)
{
 Face *f;
 Edge *e, *en;
 bool visible = false;
 
 f = flist;
 
 // marca facce visibili da vertice p
 do {
  if (volume_sign(f, p) < 0)
  {
   f->visible = true;
   visible = true;
  }
  f = f->next;
 } while (f != flist);
 
 // se nessuna faccia è visibile p è interno a involucro corrente
 if (!visible)
 {
  p->onhull = false;
  return;
 }
 
 e = elist;
 
 // marca lati per la rimozione solo se entrambe facce adiacenti visibili, altrimenti è lato base
 do {
  en = e->next;
  if (e->adjface[0]->visible && e->adjface[1]->visible)
   e->remove = true;
  else if (e->adjface[0]->visible || e->adjface[1]->visible)
  {
   e->newface = make_cone_face(e, pelistflist);
  }
 
  e = en;
 } while (e != elist);
}
 
Face* make_cone_face(Edge *eVertex *pEdge *&elistFace *&flist)
{
 Edge *new_edge[2];
 Face *new_face;
 
 for (int i = 0; i < 2; ++i)
 {
  // crea solo i lati che non sono già stati creati controllando /
  // campo duplicate di entrambi gli estremi di lato base
  if (!(new_edge[i] = e->endpts[i]->duplicate))
  {
   new_edge[i] = make_null_edge(elist);
   new_edge[i]->endpts[0] = e->endpts[i];
   new_edge[i]->endpts[1] = p;
   e->endpts[i]->duplicate = new_edge[i];
  }
 }
 
 // crea faccia del nuovo triangolo e ne imposta i lati
 new_face = make_null_face(flist);
 new_face->edge[0] = e;
 new_face->edge[1] = new_edge[0];
 new_face->edge[2] = new_edge[1];
 
 make_ccw(new_face, ep);
 
 // imposta al max una delle quattro facce adiacenti per i dei due lati /
 // che non sono il lato base
 for (int i = 0; i < 2; ++i)
  for (int j = 0; j < 2; ++j)
   if (!(new_edge[i]->adjface[j]))
   {
    new_edge[i]->adjface[j] = new_face;
    break;
   }
 
 return new_face;
}
 
void make_ccw(Face *fEdge *eVertex *p)
{
 Face *fv;
 Edge *s;
 int i;
 
 // imposta fv a faccia visibile da p
 if (e->adjface[0]->visible)
  fv = e->adjface[0];
 else
  fv = e->adjface[1];
 
 // scorre vertici di faccia visibile fino a trovare quello che coincide /
 // con primo estremo di lato base
 for (i = 0; fv->vertex[i] != e->endpts[0]; ++i)
  ;
 
 // se vertice successivo di faccia visibile != da estremo successivo lato base /
 // è neccessario invertire i vertici per mantenere orientazione di faccia visibile fv /
 // in nuova faccia f
 if (fv->vertex[(i + 1) % 3] != e->endpts[1])
 {
  f->vertex[0] = e->endpts[1];
  f->vertex[1] = e->endpts[0];
 }
 else
 {
  f->vertex[0] = e->endpts[0];
  f->vertex[1] = e->endpts[1];
 }
 
 f->vertex[2] = p;
}
 
void clean_up(Vertex *&vlistEdge *&elistFace *&flist)
{
 clean_edges(elist);
 clean_faces(flist);
 clean_vertices(vlistelist);
}
 
void clean_edges(Edge *&elist)
{
 Edge *e;
 Edge *t;
 
 e = elist;
 
 // imposta prima altro lato adiacente di lati base
 do {
  if (e->newface)
  {
   if (e->adjface[0]->visible)
    e->adjface[0] = e->newface;
   else
    e->adjface[1] = e->newface;
 
   e->newface = nullptr;
  }
  e = e->next;
 } while (e != elist);
 
 // rimuove lati marcati iniziando da quelli in cima alla lista /
 // fino a trovare il puntatore al primo non marcato che fungerà /
 // da primo della lista nella condizione del prossimo do-while
 while (elist && elist->remove)
 {
  e = elist;
  list_delete<Edge>(elist, e);
 }
 
 e = elist->next;
 
 do {
  if (e->remove)
  {
   t = e;
   e = e->next;
   list_delete<Edge>(elist, t);
  }
  else
   e = e->next;
 } while (e != elist);
}
 
void clean_faces(Face *&flist)
{
 Face *f;
 Face *t;
 
 while (flist && flist->visible)
 {
  f = flist;
  list_delete<Face>(flist, f);
 }
 
 f = flist->next;
 
 do {
  if (f->visible)
  {
   t = f;
   f = f->next;
   list_delete<Face>(flist, t);
  }
  else
   f = f->next;
 } while (f != flist);
}
 
void clean_vertices(Vertex *&vlistEdge *elist)
{
 Edge *e;
 Vertex *v, *t;
 
 e = elist;
 
 // marca vertici di lati non cancellati come vertici sul'involucro
 do {
  e->endpts[0]->onhull = e->endpts[1]->onhull = true;
  e = e->next;
 } while (e != elist);
 
 // cancella quelli marcati come processati e che non sono sull'involucro
 while (vlist && vlist->mark && !(vlist->onhull))
 {
  v = vlist;
  list_delete<Vertex>(vlist, v);
 }
 
 v = vlist->next;
 
 do {
  if (v->mark && !(v->onhull))
  {
   t = v;
   v = v->next;
   list_delete<Vertex>(vlist, t);
  }
  else
   v = v->next;
 } while (v != vlist);
 
 v = vlist;
 // reimposta campi di Vertex pronti per nuova iterazioni
 do {
  v->duplicate = nullptr;
  v->onhull = false;
  v = v->next;
 } while (v != vlist);
}
 
void print_faces(Face *flist)
{
 Face *f = flist;
 std::cout << "List of all faces:\t(vertex indices)\n";
 std::cout << "v0\tv1\tv2\n";
 do {
  std::cout << f->vertex[0]->vnum << "\t" << f->vertex[1]->vnum << "\t" << f->vertex[2]->vnum << std::endl;
  f = f->next;
 } while (f != flist);
}
 
void print_edges(Edge *elist)
{
 Edge *e = elist;
 int i = 0;
 do {
  std::cout << "Edge " << i << "(" << e->endpts[0] << ", " << e->endpts[1] << ")" << std::endl;
  e = e->next;
  ++i;
 } while (e != elist);
 --i;
 std::cout << "Total edges: " << i << std::endl;
}
 
void print_list(Vertex *head)
{
 Vertex *vp = head;
 do {
  std::cout << "(" << vp->pos[0] << ", " << vp->pos[1] << ", " << vp->pos[2] << ")" << std::endl;
  vp = vp->next;
 } while (vp != head);
}
 
bool AlmostEqual(float afloat bfloat epsilon = FLT_EPSILON)
{
 if (fabs(a - b) <= epsilon * max_3(1.0f, fabs(a), fabs(b)))
  return true;
 
 return false;
}
 
bool collinear(const Vertex *aconst Vertex *bconst Vertex *c)
{
 return 
  AlmostEqual(
   (c->pos[2] - a->pos[2]) * (b->pos[1] - a->pos[1]) - (b->pos[2] - a->pos[2]) * (c->pos[1] - a->pos[1]), 
   0.0f) 
  &&
  AlmostEqual(
   (b->pos[2] - a->pos[2]) * (c->pos[0] - a->pos[0]) - (b->pos[0] - a->pos[0]) * (c->pos[2] - a->pos[2]), 
   0.0f) 
  &&
  AlmostEqual(
   (b->pos[0] - a->pos[0]) * (c->pos[1] - a->pos[1]) - (b->pos[1] - a->pos[1]) * (c->pos[0] - a->pos[0]), 
   0.0f);
}


Riferimeti:
[1] Computational Geometry in C - O'Rourke
[2] Collinearità e Complanarità in 3D

Nessun commento:

Posta un commento