poly.c

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

  switch (size-firstindex) {
  case 1:
    hash= (ptr_intT)(*elemp) - (ptr_intT) skipelem;
    break;
  case 2:
    hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] - (ptr_intT) skipelem;
    break;
  case 3:
    hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] + (ptr_intT)elemp[2]
      - (ptr_intT) skipelem;
    break;
  case 4:
    hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] + (ptr_intT)elemp[2]
      + (ptr_intT)elemp[3] - (ptr_intT) skipelem;
    break;
  case 5:
    hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] + (ptr_intT)elemp[2]
      + (ptr_intT)elemp[3] + (ptr_intT)elemp[4] - (ptr_intT) skipelem;
    break;
  case 6:
    hash= (ptr_intT)(*elemp) + (ptr_intT)elemp[1] + (ptr_intT)elemp[2]
      + (ptr_intT)elemp[3] + (ptr_intT)elemp[4]+ (ptr_intT)elemp[5]
      - (ptr_intT) skipelem;
    break;
  default:
    hash= 0;
    i= 3;
    do {     /* this is about 10% in 10-d */
      if ((elem= (ptr_intT)*elemp++) != (ptr_intT)skipelem) {
        hash ^= (elem << i) + (elem >> (32-i));
	i += 3;
	if (i >= 32)
	  i -= 32;
      }
    }while(*elemp);
    break;
  }
  hash %= (ptr_intT) hashsize;
  /* hash= 0;   for debugging purposes */
  return hash;
} /* gethash */

/*-------------------------------------------------
-makenewfacet- creates a toporient? facet from vertices and apex
   modifies vertices 
returns:
    adds newfacet to qh facet_list 
       facet->neighbor= horizon, but not vice versa
    facet->vertices= vertices= apex+vertices
    newvertex_list updated
*/
facetT *qh_makenewfacet(setT *vertices, boolT toporient,facetT *horizon) {
  facetT *newfacet;
  vertexT *vertex, **vertexp;

  FOREACHvertex_(vertices) {
    if (!vertex->newlist) {
      qh_removevertex (vertex);
      qh_appendvertex (vertex);
    }
  }
  newfacet= qh_newfacet();
  newfacet->vertices= vertices;
  newfacet->toporient= toporient;
  qh_setappend(&(newfacet->neighbors), horizon);
  qh_appendfacet(newfacet);
  return(newfacet);
} /* makenewfacet */


/*--------------------------------------------
-makenewplanes- make new hyperplanes for facets
  ->f.samecycle defined for ->mergehorizon facets
returns:
  all facets have hyperplanes or are marked for merging
    doesn't create plane if horizon is coplanar (will merge)
*/
void qh_makenewplanes (void /* newfacet_list */) {
  facetT *newfacet;

  FORALLnew_facets {
    if (!newfacet->mergehorizon)
      qh_setfacetplane (newfacet);  
  }
} /* makenewplanes */

/*---------------------------------------------
-makenew_nonsimplicial- make new facets for ridges of visible facets
  qh visit_id if visible has already been seen
  attaches new facets if !qh ONLY good
  assumes all 'seen' flags false
returns:
  newfacet or NULL, bumps numnew as needed
  marks ridge neighbors for simplicial visible
  if (qh ONLYgood)
    ridges on newfacet, horizon, and visible
  else
    ridge and neighbors between newfacet and horizon
    visible facet's ridges are deleted    
*/
#ifndef qh_NOmerge
facetT *qh_makenew_nonsimplicial (facetT *visible, vertexT *apex, int *numnew) {
  void **freelistp;
  ridgeT *ridge, **ridgep;
  facetT *neighbor, *newfacet= NULL, *samecycle;
  setT *vertices;
  boolT toporient;

  FOREACHridge_(visible->ridges) {
    neighbor= otherfacet_(ridge, visible);
    if (neighbor->visible) {
      if (!qh ONLYgood) {
        if (neighbor->visitid == qh visit_id) {
          qh_setfree (&(ridge->vertices));  /* delete on 2nd visit */
	  qh_memfree_(ridge, sizeof(ridgeT), freelistp);
	}
      }
    }else {  /* neighbor is an horizon facet */
      toporient= (ridge->top == visible);
      vertices= qh_setnew (qh hull_dim); /* makes sure this is quick */
      qh_setappend (&vertices, apex);
      qh_setappend_set (&vertices, ridge->vertices);
      newfacet= qh_makenewfacet(vertices, toporient, neighbor);
      (*numnew)++;
      if (neighbor->coplanar) {
	newfacet->mergehorizon= True;
        if (!neighbor->seen) {
          newfacet->f.samecycle= newfacet;
          neighbor->f.newcycle= newfacet;
        }else {
          samecycle= neighbor->f.newcycle;
          newfacet->f.samecycle= samecycle->f.samecycle;
          samecycle->f.samecycle= newfacet;
	}
      }
      if (qh ONLYgood) {
        if (!neighbor->simplicial)
 	  qh_setappend(&(newfacet->ridges), ridge);
      }else {  /* qh_attachnewfacets */
        if (neighbor->seen) {
	  if (neighbor->simplicial) {
	    fprintf (qh ferr, "qhull internal error (qh_makenew_nonsimplicial): simplicial f%d sharing two ridges with f%d\n", 
	           neighbor->id, visible->id);
	    qh_errexit2 (qh_ERRqhull, neighbor, visible);
	  }
	  qh_setappend (&(neighbor->neighbors), newfacet);
	}else
          qh_setreplace (neighbor->neighbors, visible, newfacet);
        if (neighbor->simplicial) {
          qh_setdel (neighbor->ridges, ridge);
          qh_setfree (&(ridge->vertices)); 
	  qh_memfree (ridge, sizeof(ridgeT));
	}else {
 	  qh_setappend(&(newfacet->ridges), ridge);
 	  if (toporient)
 	    ridge->top= newfacet;
 	  else
 	    ridge->bottom= newfacet;
 	}
      trace4((qh ferr, "qh_makenew_nonsimplicial: created facet f%d from v%d and r%d of horizon f%d\n",
	    newfacet->id, apex->id, ridge->id, neighbor->id));
      }
    }
    neighbor->seen= True;        
  } /* for each ridge */
  if (!qh ONLYgood)
    SETfirst_(visible->ridges)= NULL;
  return newfacet;
} /* makenew_nonsimplicial */
#else /* qh_NOmerge */
facetT *qh_makenew_nonsimplicial (facetT *visible, vertexT *apex, int *numnew) {
  return NULL;
}
#endif /* qh_NOmerge */

/*---------------------------------------------
-makenew_simplicial- make new facets for simplicial facet
  nop if neighbor->seen or neighbor->visible (see makenew_nonsimplicial)
  attaches new facets if !qh ONLY good
returns:
  newfacet or NULL, bumps numnew as needed
  if (!qh ONLYgood)
    neighbors between newfacet and horizon
*/
facetT *qh_makenew_simplicial (facetT *visible, vertexT *apex, int *numnew) {
  facetT *neighbor, **neighborp, *newfacet= NULL;
  setT *vertices;
  boolT flip, toporient;
  int horizonskip, visibleskip;

  FOREACHneighbor_(visible) {
    if (!neighbor->seen && !neighbor->visible) {
      vertices= qh_facetintersect(neighbor,visible, &horizonskip, &visibleskip, 1);
      SETfirst_(vertices)= apex;
      flip= ((horizonskip & 0x1) ^ (visibleskip & 0x1));
      if (neighbor->toporient)         
	toporient= horizonskip & 0x1;
      else
	toporient= (horizonskip & 0x1) ^ 0x1;
      newfacet= qh_makenewfacet(vertices, toporient, neighbor);
      (*numnew)++;
      if (neighbor->coplanar && (qh PREmerge || qh MERGEexact)) {
#ifndef qh_NOmerge
	newfacet->f.samecycle= newfacet;
	newfacet->mergehorizon= True;
#endif
      }
      if (!qh ONLYgood)
        SETelem_(neighbor->neighbors, horizonskip)= newfacet;
      trace4((qh ferr, "qh_makenew_simplicial: create facet f%d top %d from v%d and horizon f%d skip %d top %d and visible f%d skip %d, flip? %d\n",
	    newfacet->id, toporient, apex->id, neighbor->id, horizonskip,
	      neighbor->toporient, visible->id, visibleskip, flip));
    }
  }
  return newfacet;
} /* makenew_simplicial */

/*-------------------------------------------------
-matchneighbor- match subridge of newfacet with neighbor or add to hash_table
  ridge is newfacet->vertices w/o newskip vertex
returns:
  duplicate ridges are unmatched and marked by qh_DUPLICATEridge
notes:
  do not allocate memory (need to free hash_table cleanly)
  similar to matchduplicates
  uses linear hash chains
*/
void qh_matchneighbor (facetT *newfacet, int newskip, int hashsize, int *hashcount) {
  boolT newfound= False;   /* True, if new facet is already in hash chain */
  boolT same, ismatch;
  unsigned hash, scan;
  facetT *facet, *matchfacet;
  int skip, matchskip;

  hash= qh_gethash (hashsize, newfacet->vertices, qh hull_dim, 1, 
                     SETelem_(newfacet->vertices, newskip));
  trace4((qh ferr, "qh_matchneighbor: newfacet f%d skip %d hash %d hashcount %d\n",
	  newfacet->id, newskip, hash, *hashcount));
  zinc_(Zhashlookup);
  for (scan= hash; (facet= SETelem_(qh hash_table, scan)); 
       scan= (++scan >= hashsize ? 0 : scan)) {
    if (facet == newfacet) {
      newfound= True;
      continue;
    }
    zinc_(Zhashtests);
    if (qh_matchvertices (1, newfacet->vertices, newskip, facet->vertices, &skip, &same)) {
      if (SETelem_(newfacet->vertices, newskip) == 
          SETelem_(facet->vertices, skip)) {
        fprintf (qh ferr, "qhull precision error: Vertex sets are the same for f%d and f%d.  Can not force output.\n",
          facet->id, newfacet->id);
        qh_errexit2 (qh_ERRprec, facet, newfacet);
      }
      ismatch= (same == (newfacet->toporient ^ facet->toporient));
      matchfacet= SETelem_(facet->neighbors, skip);
      if (ismatch && !matchfacet) {
        SETelem_(facet->neighbors, skip)= newfacet;
        SETelem_(newfacet->neighbors, newskip)= facet;
        (*hashcount)--;
        trace4((qh ferr, "qh_matchneighbor: f%d skip %d matched with new f%d skip %d\n",
           facet->id, skip, newfacet->id, newskip));
        return;
      }
      if (!qh PREmerge && !qh MERGEexact) {
	fprintf (qh ferr, "qhull precision error: facets f%d, f%d and f%d meet at a ridge with more than 2 neighbors.  Can not continue.\n",
		 facet->id, newfacet->id, getid_(matchfacet));
	qh_errexit2 (qh_ERRprec, facet, newfacet);
      }
      SETelem_(newfacet->neighbors, newskip)= qh_DUPLICATEridge;
      newfacet->dupridge= True;
      if (!newfacet->normal)
	qh_setfacetplane (newfacet);
      qh_addhash (newfacet, qh hash_table, hashsize, hash);
      (*hashcount)++;
      if (!facet->normal)
	qh_setfacetplane (facet);
      if (matchfacet != qh_DUPLICATEridge) {
	SETelem_(facet->neighbors, skip)= qh_DUPLICATEridge;
	facet->dupridge= True;
	if (!facet->normal)
	  qh_setfacetplane (facet);
	if (matchfacet) {
	  matchskip= qh_setindex (matchfacet->neighbors, facet);
	  SETelem_(matchfacet->neighbors, matchskip)= qh_DUPLICATEridge;
	  matchfacet->dupridge= True;
	  if (!matchfacet->normal)
	    qh_setfacetplane (matchfacet);
	  qh_addhash (matchfacet, qh hash_table, hashsize, hash);
	  *hashcount += 2;
	}
      }
      trace4((qh ferr, "qh_matchneighbor: new f%d skip %d duplicates ridge for f%d skip %d matching f%d ismatch %d at hash %d\n",
	   newfacet->id, newskip, facet->id, skip, 
	   (matchfacet == qh_DUPLICATEridge ? -2 : getid_(matchfacet)), 
	   ismatch, hash));
      return; /* end of duplicate ridge */
    }
  }
  if (!newfound) 
    SETelem_(qh hash_table, scan)= newfacet;  /* same as qh_addhash */
  (*hashcount)++;
  trace4((qh ferr, "qh_matchneighbor: no match for f%d skip %d at hash %d\n",
           newfacet->id, newskip, hash));
} /* matchneighbor */


/*-------------------------------------------------
-matchnewfacets- match newfacets in newfacet_list to their newfacet neighbors
  newfacets already have neighbor[0] (horizon facet)
  assumes qh hash_table is NULL
  vertex->neighbors has not been updated yet
returns:
  qh newfacet_list with full neighbor sets