merge.c

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

    trace0((qh ferr, "qh_flippedmerges: merge flipped f%d into f%d dist %2.2g during p%d\n",
      facet1->id, neighbor->id, dist, qh furthest_id));
    qh_mergefacet (facet1, neighbor, &mindist, &maxdist, !qh_MERGEapex);
    nummerge++;
    if (qh PRINTstatistics) {
      zinc_(Zflipped);
      wadd_(Wflippedtot, dist);
      wmax_(Wflippedmax, dist);
    }
    qh_merge_degenredundant();
  }
  FOREACHmerge_(othermerges) {
    if (merge->facet1->visible || merge->facet2->visible)
      qh_memfree (merge, sizeof(mergeT));
    else
      qh_setappend (&qh facet_mergeset, merge);
  }
  qh_settempfree (&othermerges);
  if (nummerge)
    *wasmerge= True;
  trace1((qh ferr, "qh_flippedmerges: merged %d flipped facets into a good neighbor\n", nummerge));
} /* flippedmerges */


/*-------------------------------------------------
-forcedmerges- merge across duplicated ridges
  duplicate ridges on facet_mergeset
  assumes qh_merge_degenredundant has not be called
returns:
  no duplicate ridges
  wasmerge set if merge
  facet_mergeset includes non-forced merges (none for now)
  degen_mergeset includes degen/redun merges
notes: 
  duplicate ridges occur when the horizon is pinched,
        i.e. a subridge occurs in more than two horizon ridges.
     could rename vertices that pinch the horizon
  othermerges isn't needed since facet_mergeset is empty afterwards
    keep it in case of change
*/
void qh_forcedmerges(boolT *wasmerge) {
  facetT *facet1, *facet2;
  mergeT *merge, **mergep;
  realT dist1, dist2, mindist1, mindist2, maxdist1, maxdist2;
  setT *othermerges;
  int nummerge=0, numflip=0;

  if (qh TRACEmerge-1 == zzval_(Ztotmerge))
    qhmem.IStracing= qh IStracing= qh TRACElevel;
  trace4((qh ferr, "qh_forcedmerges: begin\n"));  
  othermerges= qh_settemppop(); /* was facet_mergeset */
  qh facet_mergeset= qh_settemp (qh TEMPsize);
  qh_settemppush (othermerges);
  FOREACHmerge_(othermerges) {
    if (merge->type != MRGridge) 
    	continue;
    facet1= merge->facet1;
    facet2= merge->facet2;
    while (facet1->visible)    	 /* must exist, no qh_merge_degenredunant */
      facet1= facet1->f.replace; /* previously merged facet */
    while (facet2->visible)
      facet2= facet2->f.replace; /* previously merged facet */
    if (facet1 == facet2)
      continue;
    if (!qh_setin (facet2->neighbors, facet1)) {
      fprintf (qh ferr, "qhull internal error (qh_forcedmerges): f%d and f%d had a duplicate ridge but as f%d and f%d they are no longer neighbors\n",
	       merge->facet1->id, merge->facet2->id, facet1->id, facet2->id);
      qh_errexit2 (qh_ERRqhull, facet1, facet2);
    }
    if (qh TRACEmerge-1 == zzval_(Ztotmerge))
      qhmem.IStracing= qh IStracing= qh TRACElevel;
    dist1= qh_getdistance (facet1, facet2, &mindist1, &maxdist1);
    dist2= qh_getdistance (facet2, facet1, &mindist2, &maxdist2);
    trace0((qh ferr, "qh_forcedmerges: duplicate ridge between f%d and f%d, dist %2.2g and reverse dist %2.2g during p%d\n",
	    facet1->id, facet2->id, dist1, dist2, qh furthest_id));
    if (dist1 < dist2) 
      qh_mergefacet (facet1, facet2, &mindist1, &maxdist1, !qh_MERGEapex);
    else {
      qh_mergefacet (facet2, facet1, &mindist2, &maxdist2, !qh_MERGEapex);
      dist1= dist2;
      facet1= facet2;
    }
    if (facet1->flipped) {
      zinc_(Zmergeflipdup);
      numflip++;
    }else
      nummerge++;
    if (qh PRINTstatistics) {
      zinc_(Zduplicate);
      wadd_(Wduplicatetot, dist1);
      wmax_(Wduplicatemax, dist1);
    }
  }
  FOREACHmerge_(othermerges) {
    if (merge->type == MRGridge)
      qh_memfree (merge, sizeof(mergeT));
    else
      qh_setappend (&qh facet_mergeset, merge);
  }
  qh_settempfree (&othermerges);
  if (nummerge)
    *wasmerge= True;
  trace1((qh ferr, "qh_forcedmerges: merged %d facets and %d flipped facets across duplicated ridges\n", 
                nummerge, numflip));
} /* forcedmerges */


/*-------------------------------------------------
-getmergeset- returns facet_mergeset of facet-neighbor pairs to be merged
    only tests !tested or nonconvex ridges of !tested facets
    assumes no nonconvex ridges with both facets tested
returns:
    sorted mergeset
    all ridges tested
notes:
  uses facet->tested/ridge->tested to prevent duplicate tests
  can not limit tests to modified ridges since the centrum changed
  uses qh visit_id
  see qh_getmergeset_initial
*/
void qh_getmergeset(facetT *facetlist) {
  facetT *facet, *neighbor, **neighborp;
  ridgeT *ridge, **ridgep;
  int nummerges;
  
  nummerges= qh_setsize (qh facet_mergeset);
  trace4((qh ferr, "qh_getmergeset: started.\n"));
  qh visit_id++;
  FORALLfacet_(facetlist) {
    if (facet->tested)
      continue;
    facet->visitid= qh visit_id;
    facet->tested= True;  /* must be non-simplicial due to merge */
    FOREACHneighbor_(facet)
      neighbor->seen= False;
    FOREACHridge_(facet->ridges) {
      if (ridge->tested && !ridge->nonconvex)
	continue;
      /* if tested & nonconvex, need to append merge */
      neighbor= otherfacet_(ridge, facet);
      if (neighbor->seen) {
	ridge->tested= True;
	ridge->nonconvex= False;
      }else if (neighbor->visitid != qh visit_id) {
        ridge->tested= True;
        ridge->nonconvex= False;
	neighbor->seen= True;      /* only one ridge is marked nonconvex */
	if (qh_test_appendmerge (facet, neighbor))
	  ridge->nonconvex= True;
      }
    }
  }
  nummerges= qh_setsize (qh facet_mergeset);
  if (qh ANGLEmerge)
    qsort(SETaddr_(qh facet_mergeset, mergeT), nummerges,sizeof(mergeT *),qh_compareangle);
  else
    qsort(SETaddr_(qh facet_mergeset, mergeT), nummerges,sizeof(mergeT *),qh_comparemerge);
  if (qh POSTmerging) {
    zadd_(Zmergesettot2, nummerges);
  }else {
    zadd_(Zmergesettot, nummerges);
    zmax_(Zmergesetmax, nummerges);
  }
  trace2((qh ferr, "qh_getmergeset: %d merges found\n", nummerges));
} /* getmergeset */


/*------------------------------------------------
-getmergeset_initial- initial mergeset for facets
  tests all facet/neighbor pairs on facetlist
  uses visit_id, assumes ridge->nonconvex is False
  facet_mergeset may have degen/redundant from flipped and forced merges
returns:
  sorted mergeset
  sets facet->tested, ridge->tested, and ridge->nonconvex
notes:
  see qh_getmergeset
*/
void qh_getmergeset_initial (facetT *facetlist) {
  facetT *facet, *neighbor, **neighborp;
  ridgeT *ridge, **ridgep;
  int nummerges;

  qh visit_id++;
  FORALLfacet_(facetlist) {
    facet->visitid= qh visit_id;
    facet->tested= True;
    FOREACHneighbor_(facet) {
      if (neighbor->visitid != qh visit_id) {
        if (qh_test_appendmerge (facet, neighbor)) {
          FOREACHridge_(neighbor->ridges) {
            if (facet == otherfacet_(ridge, neighbor)) {
              ridge->nonconvex= True;
              break;	/* only one ridge is marked nonconvex */
            }
          }
        }
      }
    }
    FOREACHridge_(facet->ridges)
      ridge->tested= True;
  }
  nummerges= qh_setsize (qh facet_mergeset);
  if (qh ANGLEmerge)
    qsort(SETaddr_(qh facet_mergeset, mergeT), nummerges,sizeof(mergeT *),qh_compareangle);
  else
    qsort(SETaddr_(qh facet_mergeset, mergeT), nummerges,sizeof(mergeT *),qh_comparemerge);
  if (qh POSTmerging) {
    zadd_(Zmergeinittot2, nummerges);
  }else {
    zadd_(Zmergeinittot, nummerges);
    zmax_(Zmergeinitmax, nummerges);
  }
  trace2((qh ferr, "qh_getmergeset_initial: %d merges found\n", nummerges));
} /* getmergeset_initial */


/*-----------------------------------------
-hashridge- add ridge to hashtable without oldvertex
  assumes hashtable is large enough
*/
void qh_hashridge (setT *hashtable, int hashsize, ridgeT *ridge, vertexT *oldvertex) {
  unsigned hash;
  ridgeT *ridgeA;

  hash= qh_gethash (hashsize, ridge->vertices, qh hull_dim-1, 0, oldvertex);
  while (True) {
    if (!(ridgeA= SETelem_(hashtable, hash))) {
      SETelem_(hashtable, hash)= ridge;
      break;
    }else if (ridgeA == ridge)
      break;
    if (++hash == hashsize)
      hash= 0;
  }
} /* hashridge */


/*-----------------------------------------
-hashridge_find- returns matching ridge in hashtable without oldvertex
  assumes hashtable is large enough
  can't match ridge to itself
  if oldvertex is NULL matches with one skip
returns:
  returns matching ridge;
  if no match,
    hashslot= -1 if ridge already in table
      else next NULL index
*/
ridgeT *qh_hashridge_find (setT *hashtable, int hashsize, ridgeT *ridge, 
              vertexT *vertex, vertexT *oldvertex, int *hashslot) {
  unsigned hash;
  ridgeT *ridgeA;

  *hashslot= 0;
  zinc_(Zhashridge);
  hash= qh_gethash (hashsize, ridge->vertices, qh hull_dim-1, 0, vertex);
  while ((ridgeA= SETelem_(hashtable, hash))) {
    if (ridgeA == ridge)
      *hashslot= -1;      
    else {
      zinc_(Zhashridgetest);
      if (qh_setequal_except (ridge->vertices, vertex, ridgeA->vertices, oldvertex))
        return ridgeA;
    }
    if (++hash == hashsize)
      hash= 0;
  }
  if (!*hashslot)
    *hashslot= hash;
  return NULL;
} /* hashridge_find */


/*-------------------------------------------------
-makeridges- creates explicit ridges between simplicial facets
  allows qh_MERGEridge flag
  uses existing ridges
returns:
  facet with ridges and without qh_MERGEridge
  ->simplicial is False
notes:
  similar to mergecycle_ridges
  duplicate neighbors ok if ridges already exist (qh_mergecycle_ridges)
*/
void qh_makeridges(facetT *facet) {
  facetT *neighbor, **neighborp;
  ridgeT *ridge, **ridgep;
  int neighbor_i, neighbor_n;
  boolT toporient, mergeridge= False;
  
  if (!facet->simplicial)
    return;
  trace4((qh ferr, "qh_makeridges: make ridges for f%d\n", facet->id));
  facet->simplicial= False;
  FOREACHneighbor_(facet) {
    if (neighbor == qh_MERGEridge)
      mergeridge= True;
    else
      neighbor->seen= False;
  }
  FOREACHridge_(facet->ridges)
    otherfacet_(ridge, facet)->seen= True;
  FOREACHneighbor_i_(facet) {
    if (neighbor == qh_MERGEridge)
      continue;  /* fixed by qh_mark_dupridges */
    else if (!neighbor->seen) {  /* no current ridges */
      ridge= qh_newridge();
      ridge->vertices= qh_setnew_delnthsorted (facet->vertices, qh hull_dim,
					                  neighbor_i, 0);
      toporient= facet->toporient ^ (neighbor_i & 0x1);
      if (toporient) {
        ridge->top= facet;
        ridge->bottom= neighbor;
      }else {
        ridge->top= neighbor;
        ridge->bottom= facet;
      }
#if 0 /* this also works */
      flip= (facet->toporient ^ neighbor->toporient)^(skip1 & 0x1) ^ (skip2 & 0x1);
      if (facet->toporient ^ (skip1 & 0x1) ^ flip) {
        ridge->top= neighbor;
        ridge->bottom= facet;
      }else {
        ridge->top= facet;
        ridge->bottom= neighbor;
      }
#endif
      qh_setappend(&(facet->ridges), ridge);
      qh_setappend(&(neighbor->ridges), ridge);
    }
  }
  if (mergeridge) {
    while (qh_setdel (facet->neighbors, qh_MERGEridge))
      ; /* delete each one */
  }
} /* makeridges */


/*-------------------------------------------------
-mark_dupridges- add duplicated ridges to facet_mergeset
  duplicate ridges marked by MERGEridge and both sides facet->dupridge
  uses qh visit_id
returns:
  duplicate ridges on facet_mergeset
  no MERGEridges in neighbor sets
  ->mergeridge/->mergeridge2 set and duplicated ridge on facet_mergeset
  all ->mergeridgex facets have ->normal 
  if ->mergeridgex facet in a samecycle, then another facet has a normal
    with flipped retested
notes: 
  duplicate ridges occur when the horizon is pinched,
        i.e. a subridge occurs in more than two horizon ridges.
     could rename vertices that pinch the horizon
*/
void qh_mark_dupridges(facetT *facetlist) {
  facetT *facet, *neighbor, **neighborp;
  int nummerge=0, skipdim;
  vertexT *apex;
  mergeT *merge, **mergep;
  

  trace4((qh ferr, "qh_mark_dupridges: identify duplicate ridges\n"));  
  FORALLfacet_(facetlist) {
    if (facet->dupridge) {
      FOREACHneighbor_(facet) {
        if (neighbor == qh_MERGEridge) {
	  facet->mergeridge= True;
	  continue;
	}
        if (neighbor->dupridge
	&& !qh_setin (neighbor->neighbors, facet)) { /* qh_MERGEridge */
	  qh_appendmergeset (facet, neighbor, MRGridge, NULL);
	  facet->mergeridge2= True;
	  facet->mergeridge= True;
	  nummerge++;
	}
      }
    }
  }
  if (!nummerge)