poly2.c

这是一个新的用于图像处理的工具箱

  be sure to stop on any 2nd visit
*/
ridgeT *qh_nextridge3d (ridgeT *atridge, facetT *facet, vertexT **vertexp) {
  vertexT *atvertex, *vertex, *othervertex;
  ridgeT *ridge, **ridgep;

  if ((atridge->top == facet) ^ qh_ORIENTclock)
    atvertex= SETsecond_(atridge->vertices);
  else
    atvertex= SETfirst_(atridge->vertices);
  FOREACHridge_(facet->ridges) {
    if (ridge == atridge)
      continue;
    if ((ridge->top == facet) ^ qh_ORIENTclock) {
      othervertex= SETsecond_(ridge->vertices);
      vertex= SETfirst_(ridge->vertices);
    }else {
      vertex= SETsecond_(ridge->vertices);
      othervertex= SETfirst_(ridge->vertices);
    }
    if (vertex == atvertex) {
      if (vertexp)
        *vertexp= othervertex;
      return ridge;
    }
  }
  return NULL;
} /* nextridge3d */
#else /* qh_NOmerge */
void qh_matchduplicates (facetT *atfacet, int atskip, int hashsize, int *hashcount) {
}
ridgeT *qh_nextridge3d (ridgeT *atridge, facetT *facet, vertexT **vertexp) {

  return NULL;
}
#endif /* qh_NOmerge */
  
/*----------------------------------------
-nonupper- return first facet without upperdelaunay or flipped
*/
facetT *qh_nonupper (facetT *facetlist) {
  facetT *facet;

  FORALLfacet_(facetlist) {
    if (!facet->upperdelaunay && !facet->flipped)
      return facet;
  }
  fprintf(qh ferr, "qhull errror (qh_nonupper): all facets from f%d are upper-Delaunay or flipped\n",
    getid_(facetlist)); 
  qh_errexit( qh_ERRqhull, facetlist, NULL);
  return NULL; /* avoid warning */
} /* nonupper */

/*-------------------------------------------------
-point- return point for a point id, or NULL if unknown
alternative code:
    return ((pointT *)((unsigned long)qh first_point
           + (unsigned long)((id)*qh normal_size)));
*/
pointT *qh_point (int id) {

  if (id < 0)
    return NULL;
  if (id < qh num_points)
    return qh first_point + id * qh hull_dim;
  id -= qh num_points;
  if (id < qh_setsize (qh other_points))
    return SETelem_(qh other_points, id);
  return NULL;
} /* point */
  
/*-------------------------------------------------
-point_add- access function for pointfacet and pointvertex
*/
void qh_point_add (setT *set, pointT *point, void *elem) {
  int id, size;

  SETreturnsize_(set, size);
  if ((id= qh_pointid(point)) < 0)
    fprintf (qh ferr, "qhull internal warning (point_add): unknown point %p id %d\n", 
      point, id);
  else if (id >= size) {
    fprintf (qh ferr, "qhull internal errror (point_add): point p%d is out of bounds (%d)\n",
	     id, size);
    qh_errexit (qh_ERRqhull, NULL, NULL);
  }else
    SETelem_(set, id)= elem;
} /* point_add */


/*-------------------------------------------------
-pointfacet- return temporary set of facets indexed by point id
  for vertices, coplanarset, and outsideset
  access with FOREACHfacet_i_(facets) and SETelem_(facets, i)
  NULL if no facet for point (inside)
     this will include qh GOODpointp
*/
setT *qh_pointfacet (void /*qh facet_list*/) {
  int numpoints= qh num_points + qh_setsize (qh other_points);
  setT *facets;
  facetT *facet;
  vertexT *vertex, **vertexp;
  pointT *point, **pointp;
  
  facets= qh_settemp (numpoints);
  qh_setzero (facets, 0, numpoints);
  qh vertex_visit++;
  FORALLfacets {
    FOREACHvertex_(facet->vertices) {
      if (vertex->visitid != qh vertex_visit) {
        vertex->visitid= qh vertex_visit;
        qh_point_add (facets, vertex->point, facet);
      }
    }
    FOREACHpoint_(facet->coplanarset) 
      qh_point_add (facets, point, facet);
    FOREACHpoint_(facet->outsideset) 
      qh_point_add (facets, point, facet);
  }
  return facets;
} /* pointfacet */

/*-------------------------------------------------
-pointvertex- return temporary set of vertices indexed by point id
  access with FOREACHvertex_i_(vertices) and SETelem_(vertices, i)
  NULL if no vertex for point
     this will include qh GOODpointp
*/
setT *qh_pointvertex (void /*qh facet_list*/) {
  int numpoints= qh num_points + qh_setsize (qh other_points);
  setT *vertices;
  vertexT *vertex;
  
  vertices= qh_settemp (numpoints);
  qh_setzero (vertices, 0, numpoints);
  FORALLvertices 
    qh_point_add (vertices, vertex->point, vertex);
  return vertices;
} /* pointvertex */


/*-------------------------------------------------
-prependfacet- prepend facet to the start of a facetlist
  increments qh numfacets
  updates facetlist, qh facet_list, facet_next
notes:
  be careful of prepending since it can lose a pointer.
  e.g., can lose _next by deleting and then prepending before _next
*/
void qh_prependfacet(facetT *facet, facetT **facetlist) {
  facetT *prevfacet, *list= *facetlist;
  

  trace4((qh ferr, "qh_prependfacet: prepend f%d before f%d\n",
	  facet->id, list->id));
  prevfacet= list->previous;
  facet->previous= prevfacet;
  if (prevfacet)
    prevfacet->next= facet;
  list->previous= facet;
  facet->next= *facetlist;
  if (qh facet_list == list)  /* this may change *facetlist */
    qh facet_list= facet;
  if (qh facet_next == list)
    qh facet_next= facet;
  *facetlist= facet;
  qh num_facets++;
} /* prependfacet */


/*-----------------------------------------
-printhashtable- print hash table
  not in I/O to avoid bringing io.c in
*/
void qh_printhashtable(FILE *fp) {
  facetT *facet, *neighbor;
  int id, facet_i, facet_n, neighbor_i= 0, neighbor_n= 0;
  vertexT *vertex, **vertexp;

  FOREACHfacet_i_(qh hash_table) {
    if (facet) {
      FOREACHneighbor_i_(facet) {
        if (!neighbor || neighbor == qh_MERGEridge || neighbor == qh_DUPLICATEridge) 
          break;
      }
      if (neighbor_i == neighbor_n)
        continue;
      fprintf (fp, "hash %d f%d ", facet_i, facet->id);
      FOREACHvertex_(facet->vertices)
        fprintf (fp, "v%d ", vertex->id);
      fprintf (fp, "\n neighbors:");
      FOREACHneighbor_i_(facet) {
	if (neighbor == qh_MERGEridge)
	  id= -3;
	else if (neighbor == qh_DUPLICATEridge)
	  id= -2;
	else
	  id= getid_(neighbor);
        fprintf (fp, " %d", id);
      }
      fprintf (fp, "\n");
    }
  }
} /* printhashtable */
     

/*-------------------------------------------------
-printlists- print out facet and vertex list for debugging (without 'f/v' tags)
*/
void qh_printlists (void) {
  facetT *facet;
  vertexT *vertex;
  
  fprintf (qh ferr, "qh_printlists: facets:");
  FORALLfacets 
    fprintf (qh ferr, " %d", facet->id);
  fprintf (qh ferr, "\n  new facets %d visible facets %d next facet for addpoint %d\n  vertices (new %d):",
     getid_(qh newfacet_list), getid_(qh visible_list), getid_(qh facet_next),
     getid_(qh newvertex_list));
  FORALLvertices
    fprintf (qh ferr, " %d", vertex->id);
  fprintf (qh ferr, "\n");
} /* printlists */
  
/*-------------------------------------------------
-resetlists- reset newvertex_list, newfacet_list, visible_list
  clears num_visible
  if stats, maintains statistics
  visible_list is empty if qh_deletevisible was called
*/
void qh_resetlists (boolT stats /*qh newvertex_list newfacet_list visible_list*/) {
  vertexT *vertex;
  facetT *newfacet, *visible;
  int totnew=0, totver=0;
  
  if (stats) {
    FORALLvertex_(qh newvertex_list)
      totver++;
    FORALLnew_facets 
      totnew++;
    zadd_(Zvisvertextot, totver);
    zmax_(Zvisvertexmax, totver);
    zadd_(Znewfacettot, totnew);
    zmax_(Znewfacetmax, totnew);
  }
  FORALLvertex_(qh newvertex_list)
    vertex->newlist= False;
  qh newvertex_list= NULL;
  FORALLnew_facets
    newfacet->newfacet= False;
  qh newfacet_list= NULL;
  FORALLvisible_facets {
    visible->f.replace= NULL;
    visible->visible= False;
  }
  qh visible_list= NULL;
  qh num_visible= 0;
  qh NEWfacets= False;
} /* resetlists */

/*-------------------------------------------------
-setvoronoi_all- compute Voronoi centers for all facets
returns:
  facet->center is the Voronoi center
  includes upperDelaunay facets if !qh ATinfinity (i.e., 'Qu')
notes:
  unused/untested code: please email [e-mail] if this works ok for you
  Use FORALLvertices to locate the vertex for a point.  
  Use FOREACHneighbor_(vertex) to visit the Voronoi centers for a Voronoi cell.
*/
void qh_setvoronoi_all (void) {
  facetT *facet;

  qh_clearcenters (qh_ASvoronoi);
  qh_vertexneighbors();
  
  FORALLfacets {
    if (!facet->normal || !facet->upperdelaunay || !qh ATinfinity) {
      if (!facet->center)
        facet->center= qh_facetcenter (facet->vertices);
    }
  }
} /* setvoronoi_all */

/*-------------------------------------------------
-vertexintersect- intersects two vertex sets (inverse id ordered)
  temporary set vertexsetA is replaced by the intersection
     must be at top of stack
  could overwrite vertexsetA if currently too slow
*/
void qh_vertexintersect(setT **vertexsetA,setT *vertexsetB) {
  setT *intersection;

  intersection= qh_vertexintersect_new (*vertexsetA, vertexsetB);
  qh_settempfree (vertexsetA);
  *vertexsetA= intersection;
  qh_settemppush (intersection);
} /* vertexintersect */

/*-------------------------------------------------
-vertexintersect_new- intersects two vertex sets (inverse id ordered)
returns:
  a new set
*/
setT *qh_vertexintersect_new (setT *vertexsetA,setT *vertexsetB) {
  setT *intersection= qh_setnew (qh hull_dim - 1);
  vertexT **vertexA= SETaddr_(vertexsetA, vertexT); 
  vertexT **vertexB= SETaddr_(vertexsetB, vertexT); 

  while (*vertexA && *vertexB) {
    if (*vertexA  == *vertexB) {
      qh_setappend(&intersection, *vertexA);
      vertexA++; vertexB++;
    }else {
      if ((*vertexA)->id > (*vertexB)->id)
        vertexA++;
      else
        vertexB++;
    }
  }
  return intersection;
} /* vertexintersect_new */

/*-------------------------------------------
-vertexneighhbors- for each vertex in hull, determine facet neighbors
  nop if VERTEXneighbors
  assumes all vertex->neighbors are NULL
returns:
  sets qh VERTEXneighbors, qh_addpoint() will maintain them
*/
void qh_vertexneighbors (void /*qh facet_list*/) {
  facetT *facet;
  vertexT *vertex, **vertexp;

  if (qh VERTEXneighbors)
    return;
  trace1((qh ferr, "qh_vertexneighbors: determing neighboring facets for each vertex\n"));
  qh vertex_visit++;
  FORALLfacets {
    if (facet->visible)
      continue;
    FOREACHvertex_(facet->vertices) {
      if (vertex->visitid != qh vertex_visit) {
        vertex->visitid= qh vertex_visit;
        vertex->neighbors= qh_setnew (qh hull_dim);
      }
      qh_setappend (&vertex->neighbors, facet);
    }
  }
  qh VERTEXneighbors= True;
} /* vertexneighbors */

/*-------------------------------------------------
-vertexsubset- returns True if vertexsetA is a subset of vertexsetB, False
  otherwise; relies on vertexsets being sorted;
an empty set is a subset of any other set
*/
boolT qh_vertexsubset(setT *vertexsetA, setT *vertexsetB) {
  vertexT **vertexA= (vertexT **) SETaddr_(vertexsetA, vertexT);
  vertexT **vertexB= (vertexT **) SETaddr_(vertexsetB, vertexT);

  while (True) {
    if (!*vertexA)
      return True;
    if (!*vertexB)
      return False;
    if ((*vertexA)->id > (*vertexB)->id)
      return False;
    if (*vertexA  == *vertexB)
      vertexA++;
    vertexB++; 
  }
} /* vertexsubset */