poly2.c

关于网格剖分的

  facetT *facet;
  boolT waserror= False;
  realT dist;

  if (facetlist == qh facet_list)
    zzval_(Zflippedfacets)= 0;
  FORALLfacet_(facetlist) {
    if (facet->normal && !qh_checkflipped (facet, &dist, !qh_ALL)) {
      fprintf(qh ferr, "qhull precision error: facet f%d is flipped, distance= %6.12g\n",
	      facet->id, dist);
      if (!qh FORCEoutput) {
	qh_errprint("ERRONEOUS", facet, NULL, NULL, NULL);
	waserror= True;
      }
    }
  }
  if (waserror) {
    fprintf (qh ferr, "\n\
A flipped facet occurs when its distance to the interior point is\n\
greater than %2.2g, the maximum roundoff error.\n", -qh DISTround);
    qh_errexit(qh_ERRprec, NULL, NULL);
  }
} /* checkflipped_all */

/*-<a                             href="qh-c.htm#poly"
  >-------------------------------</a><a name="checkpolygon">-</a>
  
  qh_checkpolygon( facetlist )
    checks the correctness of the structure

  notes:
    call with either qh.facet_list or qh.newfacet_list
    checks num_facets and num_vertices if qh.facet_list

  design:
    for each facet
      checks facet and outside set
    initializes vertexlist
    for each facet
      checks vertex set
    if checking all facets (qh.facetlist)
      check facet count
      if qh.VERTEXneighbors
        check vertex neighbors and count
      check vertex count
*/
void qh_checkpolygon(facetT *facetlist) {
  facetT *facet;
  vertexT *vertex, **vertexp, *vertexlist;
  int numfacets= 0, numvertices= 0, numridges= 0;
  int totvneighbors= 0, totvertices= 0;
  boolT waserror= False, nextseen= False, visibleseen= False;
  
  trace1((qh ferr, "qh_checkpolygon: check all facets from f%d\n", facetlist->id));
  if (facetlist != qh facet_list || qh ONLYgood)
    nextseen= True;
  FORALLfacet_(facetlist) {
    if (facet == qh visible_list)
      visibleseen= True;
    if (!facet->visible) {
      if (!nextseen) {
	if (facet == qh facet_next)
	  nextseen= True;
	else if (qh_setsize (facet->outsideset)
        && !(qh NARROWhull && qh CHECKfrequently)) {
	  fprintf (qh ferr, "qhull internal error (qh_checkpolygon): f%d has outside set before qh facet_next\n",
		   facet->id);
	  qh_errexit (qh_ERRqhull, facet, NULL);
	}
      }
      numfacets++;
      qh_checkfacet(facet, False, &waserror);
    }
  }
  if (qh visible_list && !visibleseen && facetlist == qh facet_list) {
    fprintf (qh ferr, "qhull internal error (qh_checkpolygon): visible list f%d no longer on facet list\n", qh visible_list->id);
    qh_printlists();
    qh_errexit (qh_ERRqhull, qh visible_list, NULL);
  }
  if (facetlist == qh facet_list)
    vertexlist= qh vertex_list;
  else if (facetlist == qh newfacet_list)
    vertexlist= qh newvertex_list;
  else
    vertexlist= NULL;
  FORALLvertex_(vertexlist) {
    vertex->seen= False;
    vertex->visitid= 0;
  }  
  FORALLfacet_(facetlist) {
    if (facet->visible)
      continue;
    if (facet->simplicial)
      numridges += qh hull_dim;
    else
      numridges += qh_setsize (facet->ridges);
    FOREACHvertex_(facet->vertices) {
      vertex->visitid++;
      if (!vertex->seen) {
	vertex->seen= True;
	numvertices++;
	if (qh_pointid (vertex->point) == -1) {
	  fprintf (qh ferr, "qhull internal error (qh_checkpolygon): unknown point %p for vertex v%d first_point %p\n",
		   vertex->point, vertex->id, qh first_point);
	  waserror= True;
	}
      }
    }
  }
  qh vertex_visit += numfacets;
  if (facetlist == qh facet_list) {
    if (numfacets != qh num_facets - qh num_visible) {
      fprintf(qh ferr, "qhull internal error (qh_checkpolygon): actual number of facets is %d, cumulative facet count is %d\n",
	      numfacets, qh num_facets- qh num_visible);
      waserror= True;
    }
    qh vertex_visit++;
    if (qh VERTEXneighbors) {
      FORALLvertices {
	qh_setcheck (vertex->neighbors, "neighbors for v", vertex->id);
	if (vertex->deleted)
	  continue;
	totvneighbors += qh_setsize (vertex->neighbors);
      }
      FORALLfacet_(facetlist)
	totvertices += qh_setsize (facet->vertices);
      if (totvneighbors != totvertices) {
	fprintf(qh ferr, "qhull internal error (qh_checkpolygon): vertex neighbors inconsistent.  Totvneighbors %d, totvertices %d\n",
		totvneighbors, totvertices);
	waserror= True;
      }
    }
    if (numvertices != qh num_vertices - qh_setsize(qh del_vertices)) {
      fprintf(qh ferr, "qhull internal error (qh_checkpolygon): actual number of vertices is %d, cumulative vertex count is %d\n",
	      numvertices, qh num_vertices - qh_setsize(qh del_vertices));
      waserror= True;
    }
    if (qh hull_dim == 2 && numvertices != numfacets) {
      fprintf (qh ferr, "qhull internal error (qh_checkpolygon): #vertices %d != #facets %d\n",
        numvertices, numfacets);
      waserror= True;
    }
    if (qh hull_dim == 3 && numvertices + numfacets - numridges/2 != 2) {
      fprintf (qh ferr, "qhull internal error (qh_checkpolygon): #vertices %d + #facets %d - #edges %d != 2\n",
        numvertices, numfacets, numridges/2);
      waserror= True;
    }
  }
  if (waserror) 
    qh_errexit(qh_ERRqhull, NULL, NULL);
} /* checkpolygon */


/*-<a                             href="qh-c.htm#poly"
  >-------------------------------</a><a name="checkvertex">-</a>
  
  qh_checkvertex( vertex )
    check vertex for consistency
    checks vertex->neighbors

  notes:
    neighbors checked efficiently in checkpolygon
*/
void qh_checkvertex (vertexT *vertex) {
  boolT waserror= False;
  facetT *neighbor, **neighborp, *errfacet=NULL;

  if (qh_pointid (vertex->point) == -1) {
    fprintf (qh ferr, "qhull internal error (qh_checkvertex): unknown point id %p\n", vertex->point);
    waserror= True;
  }
  if (vertex->id >= qh vertex_id) {
    fprintf (qh ferr, "qhull internal error (qh_checkvertex): unknown vertex id %d\n", vertex->id);
    waserror= True;
  }
  if (!waserror && !vertex->deleted) {
    if (qh_setsize (vertex->neighbors)) {
      FOREACHneighbor_(vertex) {
        if (!qh_setin (neighbor->vertices, vertex)) {
          fprintf (qh ferr, "qhull internal error (qh_checkvertex): neighbor f%d does not contain v%d\n", neighbor->id, vertex->id);
	  errfacet= neighbor;
	  waserror= True;
	}
      }
    }
  }
  if (waserror) {
    qh_errprint ("ERRONEOUS", NULL, NULL, NULL, vertex);
    qh_errexit (qh_ERRqhull, errfacet, NULL);
  }
} /* checkvertex */
  
/*-<a                             href="qh-c.htm#poly"
  >-------------------------------</a><a name="clearcenters">-</a>
  
  qh_clearcenters( type )
    clear old data from facet->center

  notes:
    sets new centertype
    nop if CENTERtype is the same
*/
void qh_clearcenters (qh_CENTER type) {
  facetT *facet;
  
  if (qh CENTERtype != type) {
    FORALLfacets {
      if (qh CENTERtype == qh_ASvoronoi){
        if (facet->center) {
          qh_memfree (facet->center, qh center_size);
          facet->center= NULL;
        }
      }else /* qh CENTERtype == qh_AScentrum */ {
        if (facet->center) {
          qh_memfree (facet->center, qh normal_size);
	  facet->center= NULL;
        }
      }
    }
    qh CENTERtype= type;
  }
  trace2((qh ferr, "qh_clearcenters: switched to center type %d\n", type));
} /* clearcenters */

/*-<a                             href="qh-c.htm#poly"
  >-------------------------------</a><a name="createsimplex">-</a>
  
  qh_createsimplex( vertices )
    creates a simplex from a set of vertices

  returns:
    initializes qh.facet_list to the simplex
    initializes qh.newfacet_list, .facet_tail
    initializes qh.vertex_list, .newvertex_list, .vertex_tail

  design:
    initializes lists
    for each vertex
      create a new facet
    for each new facet
      create its neighbor set
*/
void qh_createsimplex(setT *vertices) {
  facetT *facet= NULL, *newfacet;
  boolT toporient= True;
  int vertex_i, vertex_n, nth;
  setT *newfacets= qh_settemp (qh hull_dim+1);
  vertexT *vertex;
  
  qh facet_list= qh newfacet_list= qh facet_tail= qh_newfacet();
  qh num_facets= qh num_vertices= 0;
  qh vertex_list= qh newvertex_list= qh vertex_tail= qh_newvertex(NULL);
  FOREACHvertex_i_(vertices) {
    newfacet= qh_newfacet();
    newfacet->vertices= qh_setnew_delnthsorted (vertices, vertex_n,
						vertex_i, 0);
    newfacet->toporient= toporient;
    qh_appendfacet(newfacet);
    newfacet->newfacet= True;
    qh_appendvertex (vertex);
    qh_setappend (&newfacets, newfacet);
    toporient ^= True;
  }
  FORALLnew_facets {
    nth= 0;
    FORALLfacet_(qh newfacet_list) {
      if (facet != newfacet) 
        SETelem_(newfacet->neighbors, nth++)= facet;
    }
    qh_settruncate (newfacet->neighbors, qh hull_dim);
  }
  qh_settempfree (&newfacets);
  trace1((qh ferr, "qh_createsimplex: created simplex\n"));
} /* createsimplex */

/*-<a                             href="qh-c.htm#poly"
  >-------------------------------</a><a name="delridge">-</a>
  
  qh_delridge( ridge )
    deletes ridge from data structures it belongs to
    frees up its memory

  notes:
    in merge.c, caller sets vertex->delridge for each vertex
    ridges also freed in qh_freeqhull
*/
void qh_delridge(ridgeT *ridge) {
  void **freelistp; /* used !qh_NOmem */
  
  qh_setdel(ridge->top->ridges, ridge);
  qh_setdel(ridge->bottom->ridges, ridge);
  qh_setfree(&(ridge->vertices));
  qh_memfree_(ridge, sizeof(ridgeT), freelistp);
} /* delridge */


/*-<a                             href="qh-c.htm#poly"
  >-------------------------------</a><a name="delvertex">-</a>
  
  qh_delvertex( vertex )
    deletes a vertex and frees its memory

  notes:
    assumes vertex->adjacencies have been updated if needed
    unlinks from vertex_list
*/
void qh_delvertex (vertexT *vertex) {

  if (vertex == qh tracevertex)
    qh tracevertex= NULL;
  qh_removevertex (vertex);
  qh_setfree (&vertex->neighbors);
  qh_memfree(vertex, sizeof(vertexT));
} /* delvertex */


/*-<a                             href="qh-c.htm#poly"
  >-------------------------------</a><a name="facet3vertex">-</a>
  
  qh_facet3vertex(  )
    return temporary set of 3-d vertices in qh_ORIENTclock order

  design:
    if simplicial facet
      build set from facet->vertices with facet->toporient
    else
      for each ridge in order
        build set from ridge's vertices
*/
setT *qh_facet3vertex (facetT *facet) {
  ridgeT *ridge, *firstridge;
  vertexT *vertex;
  int cntvertices, cntprojected=0;
  setT *vertices;

  cntvertices= qh_setsize(facet->vertices);
  vertices= qh_settemp (cntvertices);
  if (facet->simplicial) {
    if (cntvertices != 3) {
      fprintf (qh ferr, "qhull internal error (qh_facet3vertex): only %d vertices for simplicial facet f%d\n", 
                  cntvertices, facet->id);
      qh_errexit(qh_ERRqhull, facet, NULL);
    }
    qh_setappend (&vertices, SETfirst_(facet->vertices));
    if (facet->toporient ^ qh_ORIENTclock)
      qh_setappend (&vertices, SETsecond_(facet->vertices));
    else
      qh_setaddnth (&vertices, 0, SETsecond_(facet->vertices));
    qh_setappend (&vertices, SETelem_(facet->vertices, 2));
  }else {
    ridge= firstridge= SETfirstt_(facet->ridges, ridgeT);   /* no infinite */
    while ((ridge= qh_nextridge3d (ridge, facet, &vertex))) {
      qh_setappend (&vertices, vertex);
      if (++cntprojected > cntvertices || ridge == firstridge)
        break;
    }
    if (!ridge || cntprojected != cntvertices) {
      fprintf (qh ferr, "qhull internal error (qh_facet3vertex): ridges for facet %d don't match up.  got at least %d\n", 
                  facet->id, cntprojected);
      qh_errexit(qh_ERRqhull, facet, ridge);
    }
  }
  return vertices;
} /* facet3vertex */

/*-<a                             href="qh-c.htm#poly"
  >-------------------------------</a><a name="findbestfacet">-</a>
  
  qh_findbestfacet( point, bestoutside, bestdist, isoutside )
    find facet that is furthest below a point 

    for Delaunay triangulations, 
      Use qh_setdelaunay() to lift point to paraboloid and scale by 'Qbb' if needed
      Do not use options 'Qbk', 'QBk', or 'QbB' since they scale the coordinates. 

  returns: