io.c

关于网格剖分的

    FOREACHvertex_(facet->vertices)
      fprintf (fp, " %d", qh_pointid (vertex->point));
    fprintf (fp, "\n");
    break;
  }
} /* printafacet */

/*-<a                             href="qh-c.htm#io"
  >-------------------------------</a><a name="printbegin">-</a>
  
  qh_printbegin(  )
    prints header for all output formats

  returns:
    checks for valid format
  
  notes:
    uses qh.visit_id for 3/4off
    changes qh.interior_point if printing centrums
    qh_countfacets clears facet->visitid for non-good facets
    
  see
    qh_printend() and qh_printafacet()
    
  design:
    count facets and related statistics
    print header for format
*/
void qh_printbegin (FILE *fp, int format, facetT *facetlist, setT *facets, boolT printall) {
  int numfacets, numsimplicial, numridges, totneighbors, numcoplanars;
  int i, num;
  facetT *facet, **facetp;
  vertexT *vertex, **vertexp;
  setT *vertices;
  pointT *point, **pointp, *pointtemp;

  qh printoutnum= 0;
  qh_countfacets (facetlist, facets, printall, &numfacets, &numsimplicial, 
      &totneighbors, &numridges, &numcoplanars);
  switch (format) {
  case qh_PRINTnone:
    break;
  case qh_PRINTarea:
    fprintf (fp, "%d\n", numfacets);
    break;
  case qh_PRINTcoplanars:
    fprintf (fp, "%d\n", numfacets);
    break;
  case qh_PRINTcentrums:
    if (qh CENTERtype == qh_ASnone)
      qh_clearcenters (qh_AScentrum);
    fprintf (fp, "%d\n%d\n", qh hull_dim, numfacets);
    break;
  case qh_PRINTfacets:
  case qh_PRINTfacets_xridge:
    if (facetlist)
      qh_printvertexlist (fp, "Vertices and facets:\n", facetlist, facets, printall);
    break;
  case qh_PRINTgeom: 
    if (qh hull_dim > 4)  /* qh_initqhull_globals also checks */
      goto LABELnoformat;
    if (qh VORONOI && qh hull_dim > 3)  /* PRINTdim == DROPdim == hull_dim-1 */
      goto LABELnoformat;
    if (qh hull_dim == 2 && (qh PRINTridges || qh DOintersections))
      fprintf (qh ferr, "qhull warning: output for ridges and intersections not implemented in 2-d\n");
    if (qh hull_dim == 4 && (qh PRINTinner || qh PRINTouter ||
			     (qh PRINTdim == 4 && qh PRINTcentrums)))
      fprintf (qh ferr, "qhull warning: output for outer/inner planes and centrums not implemented in 4-d\n");
    if (qh PRINTdim == 4 && (qh PRINTspheres))
      fprintf (qh ferr, "qhull warning: output for vertices not implemented in 4-d\n");
    if (qh PRINTdim == 4 && qh DOintersections && qh PRINTnoplanes)
      fprintf (qh ferr, "qhull warning: 'Gnh' generates no output in 4-d\n");
    if (qh PRINTdim == 2) {
      fprintf(fp, "{appearance {linewidth 3} LIST # %s | %s\n",
	      qh rbox_command, qh qhull_command);
    }else if (qh PRINTdim == 3) {
      fprintf(fp, "{appearance {+edge -evert linewidth 2} LIST # %s | %s\n",
	      qh rbox_command, qh qhull_command);
    }else if (qh PRINTdim == 4) {
      qh visit_id++;
      num= 0;
      FORALLfacet_(facetlist)    /* get number of ridges to be printed */
        qh_printend4geom (NULL, facet, &num, printall);
      FOREACHfacet_(facets)
        qh_printend4geom (NULL, facet, &num, printall);
      qh ridgeoutnum= num;
      qh printoutvar= 0;  /* counts number of ridges in output */
      fprintf (fp, "LIST # %s | %s\n", qh rbox_command, qh qhull_command);
    }
    if (qh PRINTdots) {
      qh printoutnum++;
      num= qh num_points + qh_setsize (qh other_points);
      if (qh DELAUNAY && qh ATinfinity)
	num--;
      if (qh PRINTdim == 4)
        fprintf (fp, "4VECT %d %d 1\n", num, num);
      else
	fprintf (fp, "VECT %d %d 1\n", num, num);
      for (i= num; i--; ) {
        if (i % 20 == 0)
          fprintf (fp, "\n");
	fprintf (fp, "1 ");
      }
      fprintf (fp, "# 1 point per line\n1 ");
      for (i= num-1; i--; ) {
        if (i % 20 == 0)
          fprintf (fp, "\n");
	fprintf (fp, "0 ");
      }
      fprintf (fp, "# 1 color for all\n");
      FORALLpoints {
        if (!qh DELAUNAY || !qh ATinfinity || qh_pointid(point) != qh num_points-1) {
	  if (qh PRINTdim == 4)
	    qh_printpoint (fp, NULL, point);
	  else
	    qh_printpoint3 (fp, point);
	}
      }
      FOREACHpoint_(qh other_points) {
	if (qh PRINTdim == 4)
	  qh_printpoint (fp, NULL, point);
	else
	  qh_printpoint3 (fp, point);
      }
      fprintf (fp, "0 1 1 1  # color of points\n");
    }
    if (qh PRINTdim == 4  && !qh PRINTnoplanes)
      /* 4dview loads up multiple 4OFF objects slowly */
      fprintf(fp, "4OFF %d %d 1\n", 3*qh ridgeoutnum, qh ridgeoutnum);
    qh PRINTcradius= 2 * qh DISTround;  /* include test DISTround */
    if (qh PREmerge) {
      maximize_(qh PRINTcradius, qh premerge_centrum + qh DISTround);
    }else if (qh POSTmerge)
      maximize_(qh PRINTcradius, qh postmerge_centrum + qh DISTround);
    qh PRINTradius= qh PRINTcradius;
    if (qh PRINTspheres + qh PRINTcoplanar)
      maximize_(qh PRINTradius, qh MAXabs_coord * qh_MINradius);
    if (qh premerge_cos < REALmax/2) {
      maximize_(qh PRINTradius, (1- qh premerge_cos) * qh MAXabs_coord);
    }else if (!qh PREmerge && qh POSTmerge && qh postmerge_cos < REALmax/2) {
      maximize_(qh PRINTradius, (1- qh postmerge_cos) * qh MAXabs_coord);
    }
    maximize_(qh PRINTradius, qh MINvisible); 
    if (qh JOGGLEmax < REALmax/2)
      qh PRINTradius += qh JOGGLEmax * sqrt (qh hull_dim);
    if (qh PRINTdim != 4 &&
	(qh PRINTcoplanar || qh PRINTspheres || qh PRINTcentrums)) {
      vertices= qh_facetvertices (facetlist, facets, printall);
      if (qh PRINTspheres && qh PRINTdim <= 3)
         qh_printspheres (fp, vertices, qh PRINTradius);
      if (qh PRINTcoplanar || qh PRINTcentrums) {
        qh firstcentrum= True;
        if (qh PRINTcoplanar&& !qh PRINTspheres) {
          FOREACHvertex_(vertices) 
            qh_printpointvect2 (fp, vertex->point, NULL,
				qh interior_point, qh PRINTradius);
	}
        FORALLfacet_(facetlist) {
	  if (!printall && qh_skipfacet(facet))
	    continue;
	  if (!facet->normal)
	    continue;
          if (qh PRINTcentrums && qh PRINTdim <= 3)
            qh_printcentrum (fp, facet, qh PRINTcradius);
	  if (!qh PRINTcoplanar)
	    continue;
          FOREACHpoint_(facet->coplanarset)
            qh_printpointvect2 (fp, point, facet->normal, NULL, qh PRINTradius);
          FOREACHpoint_(facet->outsideset)
            qh_printpointvect2 (fp, point, facet->normal, NULL, qh PRINTradius);
        }
        FOREACHfacet_(facets) {
	  if (!printall && qh_skipfacet(facet))
	    continue;
	  if (!facet->normal)
	    continue;
          if (qh PRINTcentrums && qh PRINTdim <= 3)
            qh_printcentrum (fp, facet, qh PRINTcradius);
	  if (!qh PRINTcoplanar)
	    continue;
          FOREACHpoint_(facet->coplanarset)
            qh_printpointvect2 (fp, point, facet->normal, NULL, qh PRINTradius);
          FOREACHpoint_(facet->outsideset)
            qh_printpointvect2 (fp, point, facet->normal, NULL, qh PRINTradius);
        }
      }
      qh_settempfree (&vertices);
    }
    qh visit_id++; /* for printing hyperplane intersections */
    break;
  case qh_PRINTids:
    fprintf (fp, "%d\n", numfacets);
    break;
  case qh_PRINTincidences:
    if (qh VORONOI)
      fprintf (qh ferr, "qhull warning: writing Delaunay.  Use 'p' or 'o' for Voronoi centers\n");
    qh printoutvar= qh vertex_id;  /* centrum id for non-simplicial facets */
    if (qh hull_dim <= 3)
      fprintf(fp, "%d\n", numfacets);
    else
      fprintf(fp, "%d\n", numsimplicial+numridges);
    break;
  case qh_PRINTinner:
  case qh_PRINTnormals:
  case qh_PRINTouter:
    if (qh CDDoutput)
      fprintf (fp, "%s | %s\nbegin\n    %d %d real\n", qh rbox_command, 
              qh qhull_command, numfacets, qh hull_dim+1);
    else
      fprintf (fp, "%d\n%d\n", qh hull_dim+1, numfacets);
    break;
  case qh_PRINTmathematica:  
    if (qh hull_dim > 3)  /* qh_initbuffers also checks */
      goto LABELnoformat;
    if (qh VORONOI)
      fprintf (qh ferr, "qhull warning: output is the Delaunay triangulation\n");
    fprintf(fp, "{\n");
    qh printoutvar= 0;   /* counts number of facets for notfirst */
    break;
  case qh_PRINTmerges:
    fprintf (fp, "%d\n", numfacets);
    break;
  case qh_PRINTpointintersect:
    fprintf (fp, "%d\n%d\n", qh hull_dim, numfacets);
    break;
  case qh_PRINTneighbors:
    fprintf (fp, "%d\n", numfacets);
    break;
  case qh_PRINToff:
  case qh_PRINTtriangles:
    if (qh VORONOI)
      goto LABELnoformat;
    num = qh hull_dim;
    if (format == qh_PRINToff || qh hull_dim == 2)
      fprintf (fp, "%d\n%d %d %d\n", num, 
        qh num_points+qh_setsize (qh other_points), numfacets, totneighbors/2);
    else { /* qh_PRINTtriangles */
      qh printoutvar= qh num_points+qh_setsize (qh other_points); /* first centrum */
      if (qh DELAUNAY)
        num--;  /* drop last dimension */
      fprintf (fp, "%d\n%d %d %d\n", num, qh printoutvar 
	+ numfacets - numsimplicial, numsimplicial + numridges, totneighbors/2);
    }
    FORALLpoints
      qh_printpointid (qh fout, NULL, num, point, -1);
    FOREACHpoint_(qh other_points)
      qh_printpointid (qh fout, NULL, num, point, -1);
    if (format == qh_PRINTtriangles && qh hull_dim > 2) {
      FORALLfacets {
	if (!facet->simplicial && facet->visitid)
          qh_printcenter (qh fout, format, NULL, facet);
      }
    }
    break;
  case qh_PRINTpointnearest:
    fprintf (fp, "%d\n", numcoplanars);
    break;
  case qh_PRINTpoints:
    if (!qh VORONOI)
      goto LABELnoformat;
    if (qh CDDoutput)
      fprintf (fp, "%s | %s\nbegin\n%d %d real\n", qh rbox_command,
             qh qhull_command, numfacets, qh hull_dim);
    else
      fprintf (fp, "%d\n%d\n", qh hull_dim-1, numfacets);
    break;
  case qh_PRINTvertices:
    fprintf (fp, "%d\n", numfacets);
    break;
  case qh_PRINTsummary:
  default:
  LABELnoformat:
    fprintf (qh ferr, "qhull internal error (qh_printbegin): can not use this format for dimension %d\n",
         qh hull_dim);
    qh_errexit (qh_ERRqhull, NULL, NULL);
  }
} /* printbegin */

/*-<a                             href="qh-c.htm#io"
  >-------------------------------</a><a name="printcenter">-</a>
  
  qh_printcenter( fp, string, facet )
    print facet->center as centrum or Voronoi center
    string may be NULL.  Don't include '%' codes.
    nop if qh CENTERtype neither CENTERvoronoi nor CENTERcentrum
    if upper envelope of Delaunay triangulation and point at-infinity
      prints qh_INFINITE instead;

  notes:
    defines facet->center if needed
    if format=PRINTgeom, adds a 0 if would otherwise be 2-d
*/
void qh_printcenter (FILE *fp, int format, char *string, facetT *facet) {
  int k, num;

  if (qh CENTERtype != qh_ASvoronoi && qh CENTERtype != qh_AScentrum)
    return;
  if (string)
    fprintf (fp, string, facet->id);
  if (qh CENTERtype == qh_ASvoronoi) {
    num= qh hull_dim-1;
    if (!facet->normal || !facet->upperdelaunay || !qh ATinfinity) {
      if (!facet->center)
        facet->center= qh_facetcenter (facet->vertices);
      for (k=0; k < num; k++)
        fprintf (fp, qh_REAL_1, facet->center[k]);
    }else {
      for (k=0; k < num; k++)
        fprintf (fp, qh_REAL_1, qh_INFINITE);
    }
  }else /* qh CENTERtype == qh_AScentrum */ {
    num= qh hull_dim;
    if (format == qh_PRINTtriangles && qh DELAUNAY) 
      num--;
    if (!facet->center) 
      facet->center= qh_getcentrum (facet);
    for (k=0; k < num; k++)
      fprintf (fp, qh_REAL_1, facet->center[k]);
  }
  if (format == qh_PRINTgeom && num == 2)
    fprintf (fp, " 0\n");
  else
    fprintf (fp, "\n");
} /* printcenter */

/*-<a                             href="qh-c.htm#io"
  >-------------------------------</a><a name="printcentrum">-</a>
  
  qh_printcentrum( fp, facet, radius )
    print centrum for a facet in OOGL format
    radius defines size of centrum
    2-d or 3-d only

  returns:
    defines facet->center if needed
*/
void qh_printcentrum (FILE *fp, facetT *facet, realT radius) {
  pointT *centrum, *projpt;
  boolT tempcentrum= False;
  realT xaxis[4], yaxis[4], normal[4], dist;
  realT green[3]={0, 1, 0};
  vertexT *apex;
  int k;
  
  if (qh CENTERtype == qh_AScentrum) {
    if (!facet->center)
      facet->center= qh_getcentrum (facet);
    centrum= facet->center;
  }else {
    centrum= qh_getcentrum (facet);
    tempcentrum= True;
  }
  fprintf (fp, "{appearance {-normal -edge normscale 0} ");
  if (qh firstcentrum) {
    qh firstcentrum= False;
    fprintf (fp, "{INST geom { define centrum CQUAD  # f%d\n\
-0.3 -0.3 0.0001     0 0 1 1\n\
 0.3 -0.3 0.0001     0 0 1 1\n\
 0.3  0.3 0.0001     0 0 1 1\n\
-0.3  0.3 0.0001     0 0 1 1 } transform { \n", facet->id);
  }else
    fprintf (fp, "{INST geom { : centrum } transform { # f%d\n", facet->id);
  apex= SETfirstt_(facet->vertices, vertexT);
  qh_distplane(apex->point, facet, &dist);
  projpt= qh_projectpoint(apex->point, facet, dist);
  for (k= qh hull_dim; k--; ) {
    xaxis[k]= projpt[k] - centrum[k];
    normal[k]= facet->normal[k];
  }
  if (qh hull_dim == 2) {
    xaxis[2]= 0;
    normal[2]= 0;
  }else if (qh hull_dim == 4) {
    qh_projectdim3 (xaxis, xaxis);
    qh_projectdim3 (normal, normal);
    qh_normalize2 (normal, qh PRINTdim, True, NULL, NULL);
  }
  qh_crossproduct (3, xaxis, normal, yaxis);
  fprintf (fp, "%8.4g %8.4g %8.4g 0\n", xaxis[0], xaxis[1], xaxis[2]);
  fprintf (fp, "%8.4g %8.4g %8.4g 0\n", yaxis[0], yaxis[1], yaxis[2]);
  fprintf (fp, "%8.4g %8.4g %8.4g 0\n", normal[0], normal[1], normal[2]);
  qh_printpoint3 (fp, centrum);
  fprintf (fp, "1 }}}\n"); 
  qh_memfree (projpt, qh normal_size);
  qh_printpointvect (fp, centrum, facet->normal, NULL, radius, green);
  if (tempcentrum)
    qh_memfree (centrum, qh normal_size);
} /* printcentrum */
  
/*-<a                             href="qh-c.htm#io"
  >-------------------------------</a><a name="printend">-</a>
  
  qh_printend( fp, format )
    prints trailer for all output formats