qhull.c

关于网格剖分的

      isoutside= True;
  }

  if (isoutside) {
    if (!bestfacet->outsideset 
    || !qh_setlast (bestfacet->outsideset)) {
      qh_setappend(&(bestfacet->outsideset), point);
      if (!bestfacet->newfacet) {
        qh_removefacet (bestfacet);  /* make sure it's after qh facet_next */
        qh_appendfacet (bestfacet);
      }
#if !qh_COMPUTEfurthest
      bestfacet->furthestdist= bestdist;
#endif
    }else {
#if qh_COMPUTEfurthest
      zinc_(Zcomputefurthest);
      qh_distplane (oldfurthest, bestfacet, &dist);
      if (dist < bestdist) 
	qh_setappend(&(bestfacet->outsideset), point);
      else
	qh_setappend2ndlast(&(bestfacet->outsideset), point);
#else
      if (bestfacet->furthestdist < bestdist) {
	qh_setappend(&(bestfacet->outsideset), point);
	bestfacet->furthestdist= bestdist;
      }else
	qh_setappend2ndlast(&(bestfacet->outsideset), point);
#endif
    }
    qh num_outside++;
    trace4((qh ferr, "qh_partitionpoint: point p%d is outside facet f%d\n",
	  qh_pointid(point), bestfacet->id));
  }else if (bestdist >= -qh MAXcoplanar) {
    zzinc_(Zcoplanarpart);
    if (qh DELAUNAY)
      qh_precision ("nearly incident point");
    if (qh KEEPcoplanar + qh KEEPnearinside || bestdist > qh max_outside) 
      qh_partitioncoplanar (point, bestfacet, &bestdist);
  }else if (qh KEEPnearinside && bestdist > -qh NEARinside) {
    zinc_(Zpartnear);
    qh_partitioncoplanar (point, bestfacet, &bestdist);
  }else {
    zinc_(Zpartinside);
    trace4((qh ferr, "qh_partitionpoint: point p%d is inside all facets, closest to f%d dist %2.2g\n",
	  qh_pointid(point), bestfacet->id, bestdist));
    if (qh KEEPinside)	  
      qh_partitioncoplanar (point, bestfacet, &bestdist);
  }
} /* partitionpoint */

/*-<a                             href="qh-c.htm#qhull"
  >-------------------------------</a><a name="partitionvisible">-</a>
  
  qh_partitionvisible( allpoints, numoutside )
    partitions points in visible facets to qh.newfacet_list
    qh.visible_list= visible facets
    for visible facets
      1st neighbor (if any) points to a horizon facet or a new facet
    if allpoints (not used),
      repartitions coplanar points

  returns:
    updates outside sets and coplanar sets of qh.newfacet_list
    updates qh.num_outside (count of outside points)
  
  notes:
    qh.findbest_notsharp should be clear (extra work if set)

  design:
    for all visible facets with outside set or coplanar set
      select a newfacet for visible facet
      if outside set
        partition outside set into new facets
      if coplanar set and keeping coplanar/near-inside/inside points
        if allpoints
          partition coplanar set into new facets, may be assigned outside
        else
          partition coplanar set into coplanar sets of new facets
    for each deleted vertex
      if allpoints
        partition vertex into new facets, may be assigned outside
      else
        partition vertex into coplanar sets of new facets
*/
void qh_partitionvisible(/*visible_list*/ boolT allpoints, int *numoutside) {
  facetT *visible, *newfacet;
  pointT *point, **pointp;
  int coplanar=0, size;
  unsigned count;
  vertexT *vertex, **vertexp;
  
  if (qh ONLYmax)
    maximize_(qh MINoutside, qh max_vertex);
  *numoutside= 0;
  FORALLvisible_facets {
    if (!visible->outsideset && !visible->coplanarset)
      continue;
    newfacet= visible->f.replace;
    count= 0;
    while (newfacet && newfacet->visible) {
      newfacet= newfacet->f.replace;
      if (count++ > qh facet_id)
	qh_infiniteloop (visible);
    }
    if (!newfacet)
      newfacet= qh newfacet_list;
    if (visible->outsideset) {
      size= qh_setsize (visible->outsideset);
      *numoutside += size;
      qh num_outside -= size;
      FOREACHpoint_(visible->outsideset) 
        qh_partitionpoint (point, newfacet);
    }
    if (visible->coplanarset && (qh KEEPcoplanar + qh KEEPinside + qh KEEPnearinside)) {
      size= qh_setsize (visible->coplanarset);
      coplanar += size;
      FOREACHpoint_(visible->coplanarset) {
        if (allpoints)
          qh_partitionpoint (point, newfacet);
        else
          qh_partitioncoplanar (point, newfacet, NULL);
      }
    }
  }
  FOREACHvertex_(qh del_vertices) {
    if (vertex->point) {
      if (allpoints)
        qh_partitionpoint (vertex->point, qh newfacet_list);
      else
        qh_partitioncoplanar (vertex->point, qh newfacet_list, NULL);
    }
  }
  trace1((qh ferr,"qh_partitionvisible: partitioned %d points from outsidesets and %d points from coplanarsets\n", *numoutside, coplanar));
} /* partitionvisible */



/*-<a                             href="qh-c.htm#qhull"
  >-------------------------------</a><a name="precision">-</a>
  
  qh_precision( reason )
    restart on precision errors if not merging and if 'QJn'
*/
void qh_precision (char *reason) {

  if (qh ALLOWrestart && !qh PREmerge && !qh MERGEexact) {
    if (qh JOGGLEmax < REALmax/2) {
      trace0((qh ferr, "qh_precision: qhull restart because of %s\n", reason));
      longjmp(qh restartexit, qh_ERRprec);
    }
  }
} /* qh_precision */

/*-<a                             href="qh-c.htm#qhull"
  >-------------------------------</a><a name="printsummary">-</a>
  
  qh_printsummary( fp )
    prints summary to fp

  notes:
    not in io.c so that user_eg.c can prevent io.c from loading

  design:
    determine number of points, vertices, and coplanar points
    print summary
*/
void qh_printsummary(FILE *fp) {
  realT ratio, outerplane, innerplane;
  float cpu;
  int size, id, nummerged, numvertices, numcoplanars= 0;
  facetT *facet;
  char *s;

  size= qh num_points + qh_setsize (qh other_points);
  numvertices= qh num_vertices - qh_setsize (qh del_vertices);
  id= qh_pointid (qh GOODpointp);
  FORALLfacets {
    if (facet->coplanarset)
      numcoplanars += qh_setsize( facet->coplanarset);
  }
  if (id >=0 && qh STOPcone-1 != id && -qh STOPpoint-1 != id)
    size--;
  if (qh STOPcone || qh STOPpoint)
      fprintf (fp, "\nAt a premature exit due to 'TVn', 'TCn', or precision error.");
  if (qh VORONOI) {
    if (qh UPPERdelaunay)
      fprintf (fp, "\n\
Furthest-site Voronoi vertices by the convex hull of %d points in %d-d:\n\n", size, qh hull_dim);
    else
      fprintf (fp, "\n\
Voronoi diagram by the convex hull of %d points in %d-d:\n\n", size, qh hull_dim);
    fprintf(fp, "  Number of Voronoi regions%s: %d\n",
              qh ATinfinity ? " and at-infinity" : "", numvertices);
    if (numcoplanars) 
      fprintf(fp, "  Number of nearly incident points: %d\n", numcoplanars); 
    else if (size > numvertices) 
      fprintf(fp, "  Total number of nearly incident points: %d\n", size - numvertices); 
    fprintf(fp, "  Number of Voronoi vertices: %d\n", qh num_good);
    fprintf(fp, "  Number of facets in hull: %d\n", qh num_facets - qh num_visible);
  }else if (qh DELAUNAY) {
    if (qh UPPERdelaunay)
      fprintf (fp, "\n\
Furthest-site Delaunay triangulation by the convex hull of %d points in %d-d:\n\n", size, qh hull_dim);

    else
      fprintf (fp, "\n\
Delaunay triangulation by the convex hull of %d points in %d-d:\n\n", size, qh hull_dim);
    fprintf(fp, "  Number of input sites%s: %d\n", 
              qh ATinfinity ? " and at-infinity" : "", numvertices);
    if (numcoplanars) 
      fprintf(fp, "  Number of nearly incident points: %d\n", numcoplanars); 
    else if (size > numvertices) 
      fprintf(fp, "  Total number of nearly incident points: %d\n", size - numvertices); 
    fprintf(fp, "  Number of Delaunay facets: %d\n", qh num_good);
    fprintf(fp, "  Number of facets in hull: %d\n", qh num_facets - qh num_visible);
  }else if (qh HALFspace) {
    fprintf (fp, "\n\
Halfspace intersection by the convex hull of %d points in %d-d:\n\n", size, qh hull_dim);
    fprintf(fp, "  Number of nonredundant halfspaces: %d\n", numvertices);
    if (numcoplanars) {
      if (qh KEEPinside && qh KEEPcoplanar)
      	s= "similar and redundant";
      else if (qh KEEPinside)
        s= "redundant";
      else
        s= "similar"; 
      fprintf(fp, "  Number of %s halfspaces: %d\n", s, numcoplanars);
    } 
    fprintf(fp, "  Number of intersection points: %d\n", qh num_facets - qh num_visible);
    if (qh num_good)
      fprintf(fp, "  Number of 'good' intersections: %d\n", qh num_good);
  }else {
    fprintf (fp, "\n\
Convex hull of %d points in %d-d:\n\n", size, qh hull_dim);
    fprintf(fp, "  Number of vertices: %d\n", numvertices);
    if (numcoplanars) {
      if (qh KEEPinside && qh KEEPcoplanar)
      	s= "coplanar and interior";
      else if (qh KEEPinside)
        s= "interior";
      else
        s= "coplanar"; 
      fprintf(fp, "  Number of %s points: %d\n", s, numcoplanars);
    } 
    fprintf(fp, "  Number of facets: %d\n", qh num_facets - qh num_visible);
    if (qh num_good)
      fprintf(fp, "  Number of 'good' facets: %d\n", qh num_good);
  }
  fprintf(fp, "\nStatistics for: %s | %s", 
                      qh rbox_command, qh qhull_command);
  if (qh ROTATErandom != INT_MIN)
    fprintf(fp, " QR%d\n\n", qh ROTATErandom);
  else
    fprintf(fp, "\n\n");
  fprintf(fp, "  Number of points processed: %d\n", zzval_(Zprocessed));
  fprintf(fp, "  Number of hyperplanes created: %d\n", zzval_(Zsetplane));
  fprintf(fp, "  Number of distance tests for qhull: %d\n", zzval_(Zpartition)+
      zzval_(Zpartitionall)+zzval_(Znumvisibility)+zzval_(Zpartcoplanar));
#if 0  /* NOTE: must print before printstatistics() */
  {realT stddev, ave;
  fprintf(fp, "  average new facet balance: %2.2g\n",
	  wval_(Wnewbalance)/zval_(Zprocessed));
  stddev= qh_stddev (zval_(Zprocessed), wval_(Wnewbalance), 
                                 wval_(Wnewbalance2), &ave);
  fprintf(fp, "  new facet standard deviation: %2.2g\n", stddev);
  fprintf(fp, "  average partition balance: %2.2g\n",
	  wval_(Wpbalance)/zval_(Zpbalance));
  stddev= qh_stddev (zval_(Zpbalance), wval_(Wpbalance), 
                                 wval_(Wpbalance2), &ave);
  fprintf(fp, "  partition standard deviation: %2.2g\n", stddev);
  }
#endif
  nummerged= zzval_(Ztotmerge) - zzval_(Zcyclehorizon) + zzval_(Zcyclefacettot);
  if (nummerged) {
    fprintf(fp,"  Number of merged facets: %d\n", nummerged);
    fprintf(fp,"  Number of distance tests for merging: %d\n",zzval_(Zbestdist)+
          zzval_(Zcentrumtests)+zzval_(Zdistconvex)+zzval_(Zdistcheck)+
          zzval_(Zdistzero));
  }
  if (!qh RANDOMoutside && qh QHULLfinished) {
    cpu= qh hulltime;
    cpu /= qh_SECticks;
    wval_(Wcpu)= cpu;
    fprintf (fp, "  CPU seconds to compute hull (after input): %2.4g\n", cpu);
  }
  if (qh RERUN) {
    if (!qh PREmerge && !qh MERGEexact)
      fprintf(fp, "  Percentage of runs with precision errors: %4.1f\n",
	   zzval_(Zretry)*100.0/qh build_cnt);  /* careful of order */
  }else if (qh JOGGLEmax < REALmax/2) {
    if (zzval_(Zretry))
      fprintf(fp, "  After %d retries, input joggled by: %2.2g\n",
         zzval_(Zretry), qh JOGGLEmax);
    else
      fprintf(fp, "  Input joggled by: %2.2g\n", qh JOGGLEmax);
  }
  if (qh totarea != 0.0) 
    fprintf(fp, "  %s facet area:   %2.8g\n", 
	    zzval_(Ztotmerge) ? "Approximate" : "Total", qh totarea);
  if (qh totvol != 0.0) 
    fprintf(fp, "  %s volume:       %2.8g\n", 
	    zzval_(Ztotmerge) ? "Approximate" : "Total", qh totvol);
  if (qh MERGING) {
    qh_outerinner (NULL, &outerplane, &innerplane);
    if (outerplane > 2 * qh DISTround) {
      fprintf(fp, "  Maximum distance of %spoint above facet: %2.2g", 
	    (qh QHULLfinished ? "" : "merged "), outerplane);
      ratio= outerplane/(qh ONEmerge + qh DISTround);
      if (ratio > 0.05 && qh ONEmerge > qh MINoutside && qh JOGGLEmax > REALmax/2)
        fprintf (fp, " (%.1fx)\n", ratio);
      else
        fprintf (fp, "\n");
    }
    if (innerplane < -2 * qh DISTround) {
      fprintf(fp, "  Maximum distance of %svertex below facet: %2.2g", 
	    (qh QHULLfinished ? "" : "merged "), innerplane);
      ratio= -innerplane/(qh ONEmerge+qh DISTround);
      if (ratio > 0.05 && qh JOGGLEmax > REALmax/2)
        fprintf (fp, " (%.1fx)\n", ratio);
      else
        fprintf (fp, "\n");
    }
  }
  fprintf(fp, "\n");
} /* printsummary */