merge.c

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

    return;
  trace4((qh ferr, "qh_mark_dupridges: rebuild normals for same cycles\n"));  
  apex= SETfirst_(facetlist->vertices);
  skipdim= qh_mindiff (apex->point, qh interior_point, qh hull_dim);
  FORALLfacet_(facetlist) {            /* gets rid of qh_MERGEridge */
    if (facet->mergeridge && !facet->mergeridge2)   
      qh_makeridges (facet);
  }
  FOREACHmerge_(qh facet_mergeset) {   /* restore the missing neighbors */
    if (merge->type == MRGridge) {
      qh_setappend (&merge->facet2->neighbors, merge->facet1);
      qh_makeridges (merge->facet1);   /* and the missing ridges */
    }
  }
  trace1((qh ferr, "qh_mark_dupridges: found %d duplicated ridges\n", 
                nummerge));
} /* mark_dupridges */

/*-------------------------------------------
-maydropneighbor -- drop neighbor relationship if no ridge between facet and neighbor
  bumps qh visit_id
returns:
  appends degenerate facets to facet_mergeset
  won't cause redundant facets since vertex inclusion is the same
  may drop vertex and neighbor if no ridge
*/
void qh_maydropneighbor (facetT *facet) {
  ridgeT *ridge, **ridgep;
  realT angledegen= qh_ANGLEdegen;
  facetT *neighbor, **neighborp;

  qh visit_id++;
  trace4((qh ferr, "qh_maydropneighbor: test f%d for no ridges to a neighbor\n",
	  facet->id));
  FOREACHridge_(facet->ridges) {
    ridge->top->visitid= qh visit_id;
    ridge->bottom->visitid= qh visit_id;
  }
  FOREACHneighbor_(facet) {
    if (neighbor->visitid != qh visit_id) {
      trace0((qh ferr, "qh_maydropneighbor: facets f%d and f%d are no longer neighbors during p%d\n",
	    facet->id, neighbor->id, qh furthest_id));
      zinc_(Zdropneighbor);
      qh_setdel (facet->neighbors, neighbor);
      neighborp--;  /* repeat, deleted a neighbor */
      qh_setdel (neighbor->neighbors, facet);
      if (qh_setsize (neighbor->neighbors) < qh hull_dim) {
        zinc_(Zdropdegen);
        qh_appendmergeset (neighbor, neighbor, MRGdegen, &angledegen);
        trace2((qh ferr, "qh_maydropneighbors: f%d is degenerate.\n", neighbor->id));
      }
    }
  }
  if (qh_setsize (facet->neighbors) < qh hull_dim) {
    zinc_(Zdropdegen);
    qh_appendmergeset (facet, facet, MRGdegen, &angledegen);
    trace2((qh ferr, "qh_maydropneighbors: f%d is degenerate.\n", facet->id));
  }
} /* maydropneighbor */


/*---------------------------------------
-merge_degenredundant- merge all degenerate and redundant facets
  merges in qh degen_mergeset from qh_degen_redundant_neighbors()
returns:
  number of merges performed
  resets facet->degenerate/redundant
  if deleted (visible) facet has no neighbors, sets ->f.replace to NULL
notes:
  redundant merges happen before degenerate ones
  merging and renaming vertices can result in degen/redundant facets
*/
int qh_merge_degenredundant () {
  int size;
  mergeT *merge;
  facetT *bestneighbor, *facet1, *facet2;
  realT dist, mindist, maxdist;
  vertexT *vertex, **vertexp;
  int nummerges= 0;
  mergeType mergetype;
  boolT wasmerge;

  while ((merge= (mergeT*)qh_setdellast (qh degen_mergeset))) {
    facet1= merge->facet1;
    facet2= merge->facet2;
    mergetype= merge->type;
    qh_memfree (merge, sizeof(mergeT));
    if (facet1->visible)
      continue;
    facet1->degenerate= False; 
    facet1->redundant= False; 
    if (qh TRACEmerge-1 == zzval_(Ztotmerge))
      qhmem.IStracing= qh IStracing= qh TRACElevel;
    wasmerge= True;
    if (mergetype == MRGredundant) {
      zinc_(Zneighbor);
      while (facet2->visible) {
        if (!facet2->f.replace) {
          fprintf (qh ferr, "qhull internal error (qh_merge_degenredunant): f%d redundant but f%d has no replacement\n",
	       facet1->id, facet2->id);
          qh_errexit2 (qh_ERRqhull, facet1, facet2);
        }
        facet2= facet2->f.replace;
      }
      if (facet1 == facet2) {
	qh_degen_redundant_facet (facet1); /* in case of others */
	continue;
      }
      trace2((qh ferr, "qh_merge_degenredundant: facet f%d is contained in f%d, will merge\n",
	    facet1->id, facet2->id));
      qh_mergefacet(facet1, facet2, NULL, NULL, !qh_MERGEapex);
      /* merge distance is already accounted for */
      nummerges++;
    }else {  /* mergetype == MRGdegen, other merges may have fixed */
      if (!(size= qh_setsize (facet1->neighbors))) {
        zinc_(Zdelfacetdup);
        trace2((qh ferr, "qh_merge_degenredundant: facet f%d has no neighbors.  Deleted\n", facet1->id));
        qh_willdelete (facet1, NULL);
        FOREACHvertex_(facet1->vertices) {
  	  qh_setdel (vertex->neighbors, facet1);
	  if (!SETfirst_(vertex->neighbors)) {
	    zinc_(Zdegenvertex);
	    trace2((qh ferr, "qh_merge_degenredundant: deleted v%d because f%d has no neighbors\n",
         	 vertex->id, facet1->id));
	    vertex->deleted= True;
	    qh_setappend (&qh del_vertices, vertex);
	  }
        }
        nummerges++;
      }else if (size < qh hull_dim) {
        bestneighbor= qh_findbestneighbor(facet1, &dist, &mindist, &maxdist);
        trace2((qh ferr, "qh_merge_degenredundant: facet f%d has %d neighbors, merge into f%d dist %2.2g\n",
	      facet1->id, size, bestneighbor->id, dist));
        qh_mergefacet(facet1, bestneighbor, &mindist, &maxdist, !qh_MERGEapex);
        nummerges++;
        if (qh PRINTstatistics) {
	  zinc_(Zdegen);
	  wadd_(Wdegentot, dist);
	  wmax_(Wdegenmax, dist);
        }
      }	/* else, another merge fixed the degeneracy and redundancy tested */
    }
  }
  return nummerges;
} /* merge_degenredundant */

/*------------------------------------------------------------------
-merge_nonconvex- merge nonconvex ridge
*/
void qh_merge_nonconvex (facetT *facet1, facetT *facet2, mergeType mergetype) {
  facetT *bestfacet, *bestneighbor, *neighbor;
  realT dist, dist2, mindist, mindist2, maxdist, maxdist2;

  if (qh TRACEmerge-1 == zzval_(Ztotmerge))
    qhmem.IStracing= qh IStracing= qh TRACElevel;
  trace3((qh ferr, "qh_merge_nonconvex: merge #%d for f%d and f%d type %d\n",
      zzval_(Ztotmerge) + 1, facet1->id, facet2->id, mergetype));
  /* concave or coplanar */
  if (!facet1->newfacet) {
    bestfacet= facet2;   /* avoid merging old facet if new is ok */
    facet2= facet1;
    facet1= bestfacet;
  }else
    bestfacet= facet1;
  bestneighbor= qh_findbestneighbor(bestfacet, &dist, &mindist, &maxdist);
  neighbor= qh_findbestneighbor(facet2, &dist2, &mindist2, &maxdist2);
  if (dist < dist2) {
    qh_mergefacet(bestfacet, bestneighbor, &mindist, &maxdist, !qh_MERGEapex);
  }else if (qh AVOIDold && !facet2->newfacet
  && ((mindist >= -qh MAXcoplanar && maxdist <= qh max_outside)
       || dist * 1.5 < dist2)) {
    zinc_(Zavoidold);
    wadd_(Wavoidoldtot, dist);
    wmax_(Wavoidoldmax, dist);
    trace2((qh ferr, "qh_merge_nonconvex: avoid merging old facet f%d dist %2.2g.  Use f%d dist %2.2g instead\n",
           facet2->id, dist2, facet1->id, dist2));
    qh_mergefacet(bestfacet, bestneighbor, &mindist, &maxdist, !qh_MERGEapex);
  }else {
    qh_mergefacet(facet2, neighbor, &mindist2, &maxdist2, !qh_MERGEapex);
    dist= dist2;
  }
  if (qh PRINTstatistics) {
    if (mergetype == MRGanglecoplanar) {
      zinc_(Zacoplanar);
      wadd_(Wacoplanartot, dist);
      wmax_(Wacoplanarmax, dist);
    }else if (mergetype == MRGconcave) {
      zinc_(Zconcave);
      wadd_(Wconcavetot, dist);
      wmax_(Wconcavemax, dist);
    }else { /* MRGcoplanar */
      zinc_(Zcoplanar);
      wadd_(Wcoplanartot, dist);
      wmax_(Wcoplanarmax, dist);
    }
  }
} /* merge_nonconvex */

/*-------------------------------------------------
-mergecycle- merge a ->f.samecycle into newfacet with ->normal
  samecycle facets are simplicial from apex
  newfacet is non-simplicial horizon
  initializes vertex neighbors on first merge
returns:
  samecycle deleted (on visible_list)
  newfacet at end of facet_list
  deleted vertices on qh del_vertices
notes: 
  similar to mergefacet
*/
void qh_mergecycle (facetT *samecycle, facetT *newfacet) {
  int traceonce= False, tracerestore= 0;
  vertexT *apex;
#ifndef qh_NOtrace
  facetT *same;
#endif

  if (!qh VERTEXneighbors)
    qh_vertexneighbors();
  zzinc_(Ztotmerge);
  if (qh POSTmerging && zzval_(Ztotmerge) > qh mergereport + qh REPORTfreq2)
    qh_tracemerging();
#ifndef qh_NOtrace
  if (qh TRACEmerge == zzval_(Ztotmerge))
    qhmem.IStracing= qh IStracing= qh TRACElevel;
  trace2((qh ferr, "qh_mergecycle: merge #%d for facets from cycle f%d into coplanar horizon f%d\n", 
        zzval_(Ztotmerge), samecycle->id, newfacet->id));
  if (newfacet == qh tracefacet) {
    tracerestore= qh IStracing;
    qh IStracing= 4;
    fprintf (qh ferr, "qh_mergecycle: ========= trace merge %d of samecycle %d into trace f%d, furthest is p%d\n",
	       zzval_(Ztotmerge), samecycle->id, newfacet->id,  qh furthest_id);
    traceonce= True;
  }
  if (qh IStracing >=4) {
    fprintf (qh ferr, "  same cycle:");
    FORALLsame_cycle_(samecycle)
      fprintf(qh ferr, " f%d", same->id);
    fprintf (qh ferr, "\n");
  }
  if (qh IStracing >=4)
    qh_errprint ("MERGING CYCLE", samecycle, newfacet, NULL, NULL);
#endif /* !qh_NOtrace */
  apex= SETfirst_(samecycle->vertices);
  qh_makeridges (newfacet);
  qh_mergecycle_neighbors (samecycle, newfacet);
  qh_mergecycle_ridges (samecycle, newfacet);
  qh_mergecycle_vneighbors (samecycle, newfacet);
  if (SETfirst_(newfacet->vertices) != apex) 
    qh_setaddnth (&newfacet->vertices, 0, apex);  /* apex has last id */
  if (!newfacet->newfacet)
    qh_newvertices (newfacet->vertices);
  qh_mergecycle_facets (samecycle, newfacet);
  qh_tracemerge (samecycle, newfacet);
  /* check for degen_redundant_neighbors after qh_forcedmerges() */
  if (traceonce) {
    fprintf (qh ferr, "qh_mergecycle: end of trace facet\n");
    qh IStracing= tracerestore;
  }
} /* mergecycle */

/*-------------------------------------------------
-mergecycle_all- merge ->samecycle's of coplanar facets into horizon
  don't merge facets with ->mergeridge (these already have ->normal)
  facets are simplicial from apex
  ->cycledone == False
returns:
  all newfacets merged into coplanar horizon facets
    deleted vertices on qh del_vertices
    sets wasmerge
*/
void qh_mergecycle_all (facetT *facetlist, boolT *wasmerge) {
  facetT *facet, *same, *prev, *horizon;
  facetT *samecycle= NULL, *nextfacet, *nextsame;
  int cycles=0, total=0, facets, nummerge;

  trace2((qh ferr, "qh_mergecycle_all: begin\n"));
  for (facet= facetlist; facet && (nextfacet= facet->next); facet= nextfacet) {
    if (facet->normal)
      continue;
    if (!facet->mergehorizon) {
      fprintf (qh ferr, "qh_mergecycle_all: f%d without normal\n", facet->id);
      qh_errexit (qh_ERRqhull, facet, NULL);
    }
    horizon= SETfirst_(facet->neighbors);
    if (facet->f.samecycle == facet) {
      zinc_(Zonehorizon);  
      /* merge distance done in qh_findhorizon */
      qh_mergefacet (facet, horizon, NULL, NULL, qh_MERGEapex);
    }else {
      samecycle= facet;
      facets= 0;
      prev= facet;
      for (same= facet->f.samecycle; same;  /* FORALLsame_cycle_(facet) */
	   same= (same == facet ? NULL :nextsame)) { /* ends at facet */
	nextsame= same->f.samecycle;
        if (same->cycledone || same->visible)
          qh_infiniteloop (same);
        same->cycledone= True;
        if (same->normal) { 
          prev->f.samecycle= same->f.samecycle; /* unlink ->mergeridge */
	  same->f.samecycle= NULL;
        }else {
          prev= same;
	  facets++;
	}
      }
      while (nextfacet && nextfacet->cycledone)  /* will delete samecycle */
	nextfacet= nextfacet->next;
      horizon->f.newcycle= NULL;
      qh_mergecycle (samecycle, horizon);
      nummerge= horizon->nummerge + facets;
      if (nummerge > qh_MAXnummerge) 
      	horizon->nummerge= qh_MAXnummerge;
      else
        horizon->nummerge= nummerge;
      zzinc_(Zcyclehorizon);
      total += facets;
      zzadd_(Zcyclefacettot, facets);
      zmax_(Zcyclefacetmax, facets);
    }
    cycles++;
  }
  if (cycles)
    *wasmerge= True;
  trace1((qh ferr, "qh_mergecycle_all: merged %d same cycles or facets into coplanar horizons\n", cycles));
} /* mergecycle_all */

/*-------------------------------------------------
-mergecycle_facets- finish merge of samecycle into newfacet
returns:
  samecycle prepended to visible_list for later deletion and partitioning
     same->f.replace == newfacet
  newfacet moved to end of qh facet_list
    makes newfacet a newfacet (get's facet1->id if it was old)
    sets newfacet->newmerge set
    clears newfacet->center (unless merging into a large facet)
    clears newfacet->tested and ridge->tested for facet1
  adds neighboring facets to facet_mergeset if redundant or degenerate
*/
void qh_mergecycle_facets (facetT *samecycle, facetT *newfacet) {
  facetT *same, *next;
  
  trace4((qh ferr, "qh_mergecycle_facets: make newfacet new and samecycle deleted\n"));  
  qh_removefacet(newfacet);  /* append as a newfacet to end of qh facet_list */
  qh_appendfacet(newfacet);
  newfacet->newfacet= True;
  newfacet->simplicial= False;
  newfacet->newmerge= True;
  
  for (same= samecycle->f.samecycle; same; same= (same == samecycle ?  NULL : next)) {
    next= same->f.samecycle;  /* reused by willdelete */
    qh_willdelete (same, newfacet);
  }
  if (newfacet->center 
      && qh_setsize (newfacet->vertices) <= qh hull_dim + qh_MAXnewcentrum) {
    qh_memfree (newfacet->center, qh normal_size);
    newfacet->center= NULL;
  }
  trace3((qh ferr, "qh_mergecycle_facets: merged facets from cycle f%d into f%d\n", 
             samecycle->id, newfacet->id));
} /* mergecycle_facets */

/*-------------------------------------------------
-mergecycle_neighbors- add neighbors for ->f.samecycle to newfacet
  newfacet not in samecycle
  usually, facets are new, simplicial facets without internal ridges 
     not so if horizon facet is coplanar to two different samecycles
  newfacet has ridges
returns:
  newfacet with updated neighbors and vice-versa
  all neighbors of newfacet marked with qh visit_id
    cycle marked with qh visit_id-1
  ridge updated for simplicial neighbors of samecycle with a ridge
notes:
  similar to mergeneighbors
*/
void qh_mergecycle_neighbors(facetT *samecycle, facetT *newfacet) {
  facetT *same, *neighbor, **neighborp;
  int delneighbors= 0, newneighbors= 0;
  int samevisitid;