/** * gts_graph_partition_print_stats: * @partition: a list of @GtsGraph representing a partition. * @fp: a file pointer. * * Writes to @fp a summary of the properties of @partition. */ void gts_graph_partition_print_stats (GSList * partition, FILE * fp) { GtsRange weight; GSList * i; g_return_if_fail (partition != NULL); g_return_if_fail (fp != NULL); gts_range_init (&weight); i = partition; while (i) { gts_range_add_value (&weight, gts_graph_weight (i->data)); i = i->next; } gts_range_update (&weight); fprintf (fp, "# parts: %d\n" "# edge cuts: %5d edge cuts weight: %5g\n" "# weight: ", g_slist_length (partition), gts_graph_partition_edges_cut (partition), gts_graph_partition_edges_cut_weight (partition)); gts_range_print (&weight, fp); fputc ('\n', fp); }
/** * gts_bb_tree_triangle_distance: * @tree: a bounding box tree. * @t: a #GtsTriangle. * @distance: a #GtsBBoxDistFunc. * @delta: spatial scale of the sampling to be used. * @range: a #GtsRange to be filled with the results. * * Given a triangle @t, points are sampled regularly on its surface * using @delta as increment. The distance from each of these points * to the closest object of @tree is computed using @distance and the * gts_bb_tree_point_distance() function. The fields of @range are * filled with the number of points sampled, the minimum, average and * maximum value and the standard deviation. */ void gts_bb_tree_triangle_distance (GNode * tree, GtsTriangle * t, GtsBBoxDistFunc distance, gdouble delta, GtsRange * range) { GtsPoint * p1, * p2, * p3, * p; GtsVector p1p2, p1p3; gdouble l1, t1, dt1; guint i, n1; g_return_if_fail (tree != NULL); g_return_if_fail (t != NULL); g_return_if_fail (distance != NULL); g_return_if_fail (delta > 0.); g_return_if_fail (range != NULL); gts_triangle_vertices (t, (GtsVertex **) &p1, (GtsVertex **) &p2, (GtsVertex **) &p3); gts_vector_init (p1p2, p1, p2); gts_vector_init (p1p3, p1, p3); gts_range_init (range); p = GTS_POINT (gts_object_new (GTS_OBJECT_CLASS (gts_point_class ()))); l1 = sqrt (gts_vector_scalar (p1p2, p1p2)); n1 = l1/delta + 1; dt1 = 1.0/(gdouble) n1; t1 = 0.0; for (i = 0; i <= n1; i++, t1 += dt1) { gdouble t2 = 1. - t1; gdouble x = t2*p1p3[0]; gdouble y = t2*p1p3[1]; gdouble z = t2*p1p3[2]; gdouble l2 = sqrt (x*x + y*y + z*z); guint j, n2 = (guint) (l2/delta + 1); gdouble dt2 = t2/(gdouble) n2; x = t2*p1->x + t1*p2->x; y = t2*p1->y + t1*p2->y; z = t2*p1->z + t1*p2->z; t2 = 0.0; for (j = 0; j <= n2; j++, t2 += dt2) { p->x = x + t2*p1p3[0]; p->y = y + t2*p1p3[1]; p->z = z + t2*p1p3[2]; gts_range_add_value (range, gts_bb_tree_point_distance (tree, p, distance, NULL)); } } gts_object_destroy (GTS_OBJECT (p)); gts_range_update (range); }
static void angle_stats (GtsEdge * e, GtsRange * angle) { GSList * i; GtsTriangle * t1 = NULL, * t2 = NULL; i = e->triangles; while (i) { if (GTS_IS_FACE (i->data)) { if (!t1) t1 = i->data; else if (!t2) t2 = i->data; else return; } i = i->next; } if (!t1 || !t2) return; gts_range_add_value (angle, fabs (gts_triangles_angle (t1, t2))); }
/** * gts_bb_tree_segment_distance: * @tree: a bounding box tree. * @s: a #GtsSegment. * @distance: a #GtsBBoxDistFunc. * @delta: spatial scale of the sampling to be used. * @range: a #GtsRange to be filled with the results. * * Given a segment @s, points are sampled regularly on its length * using @delta as increment. The distance from each of these points * to the closest object of @tree is computed using @distance and the * gts_bb_tree_point_distance() function. The fields of @range are * filled with the number of points sampled, the minimum, average and * maximum value and the standard deviation. */ void gts_bb_tree_segment_distance (GNode * tree, GtsSegment * s, gdouble (*distance) (GtsPoint *, gpointer), gdouble delta, GtsRange * range) { GtsPoint * p1, * p2, * p; GtsVector p1p2; gdouble l, t, dt; guint i, n; g_return_if_fail (tree != NULL); g_return_if_fail (s != NULL); g_return_if_fail (distance != NULL); g_return_if_fail (delta > 0.); g_return_if_fail (range != NULL); p1 = GTS_POINT (s->v1); p2 = GTS_POINT (s->v2); gts_vector_init (p1p2, p1, p2); gts_range_init (range); p = GTS_POINT (gts_object_new (GTS_OBJECT_CLASS (gts_point_class()))); l = sqrt (gts_vector_scalar (p1p2, p1p2)); n = (guint) (l/delta + 1); dt = 1.0/(gdouble) n; t = 0.0; for (i = 0; i <= n; i++, t += dt) { p->x = p1->x + t*p1p2[0]; p->y = p1->y + t*p1p2[1]; p->z = p1->z + t*p1p2[2]; gts_range_add_value (range, gts_bb_tree_point_distance (tree, p, distance, NULL)); } gts_object_destroy (GTS_OBJECT (p)); gts_range_update (range); }
/* stripe - Turns the input surface into triangle strips and outputs a Geomview representation of the result. */ int main (int argc, char * argv[]) { GtsSurface * s; GSList * strips = NULL, * i; gboolean verbose = FALSE; int c = 0; GtsFile * fp; /* parse options using getopt */ while (c != EOF) { #ifdef HAVE_GETOPT_LONG static struct option long_options[] = { {"help", no_argument, NULL, 'h'}, {"verbose", no_argument, NULL, 'v'} }; int option_index = 0; switch ((c = getopt_long (argc, argv, "hv", long_options, &option_index))) { #else /* not HAVE_GETOPT_LONG */ switch ((c = getopt (argc, argv, "hv"))) { #endif /* not HAVE_GETOPT_LONG */ case 'v': /* verbose */ verbose = TRUE; break; case 'h': /* help */ fprintf (stderr, "Usage: stripe [OPTION] < FILE\n" "Turns the input surface into triangle strips and outputs a\n" "Geomview representation of the result.\n" "\n" " -v --verbose print statistics about the surface and strips\n" " -h --help display this help and exit\n" "\n" "Report bugs to %s\n", GTS_MAINTAINER); return 0; /* success */ break; case '?': /* wrong options */ fprintf (stderr, "Try `stripe --help' for more information.\n"); return 1; /* failure */ } } /* read surface in */ s = gts_surface_new (gts_surface_class (), gts_face_class (), gts_edge_class (), gts_vertex_class ()); fp = gts_file_new (stdin); if (gts_surface_read (s, fp)) { fputs ("stripe: file on standard input is not a valid GTS file\n", stderr); fprintf (stderr, "stdin:%d:%d: %s\n", fp->line, fp->pos, fp->error); return 1; /* failure */ } gts_file_destroy (fp); if (verbose) gts_surface_print_stats (s, stderr); strips = gts_surface_strip (s); /* if verbose on print stats */ if (verbose) { GtsRange l; gts_range_init (&l); i = strips; while (i) { gts_range_add_value (&l, g_slist_length (i->data)); i = i->next; } gts_range_update (&l); fprintf (stderr, "# Strips: %d\n# length : ", l.n); gts_range_print (&l, stderr); fputc ('\n', stderr); } puts ("LIST {\n"); i = strips; while (i) { GList * j = i->data; GtsTriangle * oldt = NULL; GtsColor c; c.r = rand ()/(gdouble) RAND_MAX; c.g = rand ()/(gdouble) RAND_MAX; c.b = rand ()/(gdouble) RAND_MAX; while (j) { GtsTriangle * t = j->data; GtsPoint * p1 = GTS_POINT (GTS_SEGMENT (t->e1)->v1), * p2 = GTS_POINT (GTS_SEGMENT (t->e1)->v2), * p3 = GTS_POINT (gts_triangle_vertex (t)); printf ("OFF 3 1 3\n%g %g %g\n%g %g %g\n%g %g %g\n3 0 1 2 %g %g %g\n", p1->x, p1->y, p1->z, p2->x, p2->y, p2->z, p3->x, p3->y, p3->z, c.r, c.g, c.b); if (oldt) { GtsSegment * cs = GTS_SEGMENT (gts_triangles_common_edge (t, oldt)); GtsPoint * op1 = GTS_POINT (GTS_SEGMENT (oldt->e1)->v1), * op2 = GTS_POINT (GTS_SEGMENT (oldt->e1)->v2), * op3 = GTS_POINT (gts_triangle_vertex (oldt)); printf ("VECT 1 3 0 3 0 %g %g %g %g %g %g %g %g %g\n", (op1->x + op2->x + op3->x)/3., (op1->y + op2->y + op3->y)/3., (op1->z + op2->z + op3->z)/3., (GTS_POINT (cs->v1)->x + GTS_POINT (cs->v2)->x)/2., (GTS_POINT (cs->v1)->y + GTS_POINT (cs->v2)->y)/2., (GTS_POINT (cs->v1)->z + GTS_POINT (cs->v2)->z)/2., (p1->x + p2->x + p3->x)/3., (p1->y + p2->y + p3->y)/3., (p1->z + p2->z + p3->z)/3.); } oldt = t; j = j->next; } i = i->next; } puts ("}\n"); return 0; /* success */ }
static void wave_run (GfsSimulation * sim) { GfsDomain * domain = GFS_DOMAIN (sim); GfsWave * wave = GFS_WAVE (sim); SolidFluxParams par; par.div = gfs_variable_from_name (domain->variables, "P"); g_assert (par.div); par.p = &sim->advection_params; par.fv = gfs_temporary_variable (domain); gfs_simulation_refine (sim); gfs_simulation_init (sim); while (sim->time.t < sim->time.end && sim->time.i < sim->time.iend) { gdouble tstart = gfs_clock_elapsed (domain->timer); gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gfs_event_do, sim); /* get global timestep */ gfs_domain_face_traverse (domain, FTT_XYZ, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttFaceTraverseFunc) gfs_face_reset_normal_velocity, NULL); gfs_simulation_set_timestep (sim); gdouble dt = sim->advection_params.dt; gdouble g = sim->physical_params.g/sim->physical_params.L; gdouble tnext = sim->tnext; /* spatial advection */ guint ik, ith; for (ik = 0; ik < wave->nk; ik++) { FttVector cg; group_velocity (ik, 0, &cg, wave->ntheta, g); gfs_domain_face_traverse (domain, FTT_XYZ, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttFaceTraverseFunc) set_group_velocity, &cg); if (wave->alpha_s > 0.) { /* stability criterion for GSE diffusion */ gdouble cfl = sim->advection_params.cfl; sim->advection_params.cfl = MIN (cfl, 2./(4.*wave->alpha_s*M_PI/wave->ntheta)); /* fixme: this should be: sim->advection_params.cfl = MIN (cfl, sqrt(3.)/(wave->alpha_s*2.*M_PI/wave->ntheta)); */ gfs_simulation_set_timestep (sim); sim->advection_params.cfl = cfl; } else gfs_simulation_set_timestep (sim); /* subcycling */ guint n = rint (dt/sim->advection_params.dt); g_assert (fabs (sim->time.t + sim->advection_params.dt*n - tnext) < 1e-12); while (n--) { for (ith = 0; ith < wave->ntheta; ith++) { FttVector cg; group_velocity (ik, ith, &cg, wave->ntheta, g); gfs_domain_face_traverse (domain, FTT_XYZ, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttFaceTraverseFunc) set_group_velocity, &cg); GfsVariable * t = GFS_WAVE (sim)->F[ik][ith]; sim->advection_params.v = t; gfs_domain_traverse_leaves (domain, (FttCellTraverseFunc) solid_flux, &par); gfs_tracer_advection_diffusion (domain, &sim->advection_params, NULL); sim->advection_params.fv = par.fv; gfs_domain_traverse_merged (domain, (GfsMergedTraverseFunc) gfs_advection_update, &sim->advection_params); if (wave->alpha_s > 0.) gse_alleviation_diffusion (domain, t, &cg, sim->advection_params.dt); gfs_domain_bc (domain, FTT_TRAVERSE_LEAFS, -1, t); gfs_domain_cell_traverse (domain, FTT_POST_ORDER, FTT_TRAVERSE_NON_LEAFS, -1, (FttCellTraverseFunc) t->fine_coarse, t); } gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) redo_some_events, sim); gfs_simulation_adapt (sim); } } sim->advection_params.dt = dt; /* source terms */ if (wave->source) (* wave->source) (wave); sim->time.t = sim->tnext = tnext; sim->time.i++; gts_range_add_value (&domain->timestep, gfs_clock_elapsed (domain->timer) - tstart); gts_range_update (&domain->timestep); gts_range_add_value (&domain->size, gfs_domain_size (domain, FTT_TRAVERSE_LEAFS, -1)); gts_range_update (&domain->size); } gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gfs_event_do, sim); gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gts_object_destroy, NULL); gts_object_destroy (GTS_OBJECT (par.fv)); }
static void gfs_skew_symmetric_run (GfsSimulation * sim) { GfsVariable * p, * res = NULL, * gmac[FTT_DIMENSION]; GfsDomain * domain; GSList * i; domain = GFS_DOMAIN (sim); p = gfs_variable_from_name (domain->variables, "P"); g_assert (p); FttComponent c; for (c = 0; c < FTT_DIMENSION; c++) gmac[c] = gfs_temporary_variable (domain); gfs_variable_set_vector (gmac, FTT_DIMENSION); gfs_simulation_refine (sim); gfs_simulation_init (sim); i = domain->variables; while (i) { if (GFS_IS_VARIABLE_RESIDUAL (i->data)) res = i->data; i = i->next; } gfs_simulation_set_timestep (sim); GfsVariable ** u = gfs_domain_velocity (domain); GfsVariable ** velfaces = GFS_SKEW_SYMMETRIC(sim)->velfaces; GfsVariable ** velold = GFS_SKEW_SYMMETRIC(sim)->velold; FaceData fd = { velfaces, velold, u, p, &sim->advection_params.dt, GFS_SKEW_SYMMETRIC(sim)->beta}; if (sim->time.i == 0) { gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttCellTraverseFunc) reset_unold, &fd); gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttCellTraverseFunc) get_face_values, &fd); gfs_mac_projection (domain, &sim->projection_params, sim->advection_params.dt/2., p, sim->physical_params.alpha, gmac, NULL); gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttCellTraverseFunc) get_velfaces, &fd); gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttCellTraverseFunc) initialize_unold, &fd); } while (sim->time.t < sim->time.end && sim->time.i < sim->time.iend) { gdouble tstart = gfs_clock_elapsed (domain->timer); gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gfs_event_do, sim); gfs_skew_symmetric_momentum (sim, &fd, gmac); gfs_mac_projection (domain, &sim->projection_params, sim->advection_params.dt/2., p, sim->physical_params.alpha, gmac, NULL); gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttCellTraverseFunc) correct_face_velocity, NULL); gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gfs_event_half_do, sim); gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttCellTraverseFunc) get_velfaces, &fd); gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttCellTraverseFunc) get_cell_values, &fd); gfs_domain_cell_traverse (domain, FTT_POST_ORDER, FTT_TRAVERSE_NON_LEAFS, -1, (FttCellTraverseFunc) gfs_cell_coarse_init, domain); gfs_simulation_adapt (sim); sim->time.t = sim->tnext; sim->time.i++; gfs_simulation_set_timestep (sim); gfs_advance_tracers (sim, sim->advection_params.dt); gts_range_add_value (&domain->timestep, gfs_clock_elapsed (domain->timer) - tstart); gts_range_update (&domain->timestep); gts_range_add_value (&domain->size, gfs_domain_size (domain, FTT_TRAVERSE_LEAFS, -1)); gts_range_update (&domain->size); } gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gfs_event_do, sim); gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gts_object_destroy, NULL); for (c = 0; c < FTT_DIMENSION; c++) gts_object_destroy (GTS_OBJECT (gmac[c])); }
static void gfs_porous_run (GfsSimulation * sim) { GfsVariable * p, * pmac, * res = NULL, * g[FTT_DIMENSION], * gmac[FTT_DIMENSION]; GfsVariable ** gc = sim->advection_params.gc ? g : NULL; GfsDomain * domain; GfsPorous *por; GSList * i; domain = GFS_DOMAIN (sim); por = GFS_POROUS (sim); p = gfs_variable_from_name (domain->variables, "P"); g_assert (p); pmac = gfs_variable_from_name (domain->variables, "Pmac"); g_assert (pmac); FttComponent c; for (c = 0; c < FTT_DIMENSION; c++) { gmac[c] = gfs_temporary_variable (domain); if (sim->advection_params.gc) g[c] = gfs_temporary_variable (domain); else g[c] = gmac[c]; } gfs_variable_set_vector (gmac, FTT_DIMENSION); gfs_variable_set_vector (g, FTT_DIMENSION); gfs_simulation_refine (sim); gfs_simulation_init (sim); i = domain->variables; while (i) { if (GFS_IS_VARIABLE_RESIDUAL (i->data)) res = i->data; i = i->next; } gfs_simulation_set_timestep (sim); if (sim->time.i == 0) { /*inserted changes inside this function*/ gfs_approximate_projection_por (domain, por, &sim->approx_projection_params, sim->advection_params.dt, p, sim->physical_params.alpha, res, g, NULL); gfs_simulation_set_timestep (sim); gfs_advance_tracers (sim, sim->advection_params.dt/2.); } else if (sim->advection_params.gc) gfs_update_gradients_por (domain, por, p, sim->physical_params.alpha, g); while (sim->time.t < sim->time.end && sim->time.i < sim->time.iend) { gdouble tstart = gfs_clock_elapsed (domain->timer); gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gfs_event_do, sim); /*inserted changes */ gfs_pre_projection (domain, por, FTT_DIMENSION); if (sim->advection_params.linear) { /* linearised advection */ gfs_domain_face_traverse (domain, FTT_XYZ, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttFaceTraverseFunc) gfs_face_reset_normal_velocity, NULL); gfs_domain_face_traverse (domain, FTT_XYZ, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, -1, (FttFaceTraverseFunc) gfs_face_interpolated_normal_velocity, sim->u0); } else gfs_predicted_face_velocities (domain, FTT_DIMENSION, &sim->advection_params); gfs_variables_swap (p, pmac); gfs_mac_projection_por (domain, por, &sim->projection_params, sim->advection_params.dt/2., p, sim->physical_params.alpha, gmac, NULL); gfs_variables_swap (p, pmac); gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gfs_event_half_do, sim); gfs_centered_velocity_advection_diffusion (domain, FTT_DIMENSION, &sim->advection_params, gmac, sim->time.i > 0 || !gc ? gc : gmac, sim->physical_params.alpha); if (gc) { gfs_source_darcy_implicit (domain, sim->advection_params.dt); gfs_correct_centered_velocities (domain, FTT_DIMENSION, sim->time.i > 0 ? gc : gmac, -sim->advection_params.dt); /*inserted changes*/ gfs_post_projection (domain, por, FTT_DIMENSION); } else if (gfs_has_source_coriolis (domain)) { gfs_correct_centered_velocities (domain, FTT_DIMENSION, gmac, sim->advection_params.dt); gfs_source_darcy_implicit (domain, sim->advection_params.dt); gfs_correct_centered_velocities (domain, FTT_DIMENSION, gmac, -sim->advection_params.dt); /*inserted changes*/ gfs_post_projection (domain, por, FTT_DIMENSION); } gfs_domain_cell_traverse (domain, FTT_POST_ORDER, FTT_TRAVERSE_NON_LEAFS, -1, (FttCellTraverseFunc) gfs_cell_coarse_init, domain); gfs_simulation_adapt (sim); /*inserted changes */ gfs_approximate_projection_por (domain, por, &sim->approx_projection_params, sim->advection_params.dt, p, sim->physical_params.alpha, res, g, NULL); /*inserted changes */ sim->time.t = sim->tnext; sim->time.i++; gfs_simulation_set_timestep (sim); gfs_advance_tracers (sim, sim->advection_params.dt); gts_range_add_value (&domain->timestep, gfs_clock_elapsed (domain->timer) - tstart); gts_range_update (&domain->timestep); gts_range_add_value (&domain->size, gfs_domain_size (domain, FTT_TRAVERSE_LEAFS, -1)); gts_range_update (&domain->size); } gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gfs_event_do, sim); gts_container_foreach (GTS_CONTAINER (sim->events), (GtsFunc) gts_object_destroy, NULL); for (c = 0; c < FTT_DIMENSION; c++) { gts_object_destroy (GTS_OBJECT (gmac[c])); if (sim->advection_params.gc) gts_object_destroy (GTS_OBJECT (g[c])); } }
int main (int argc, char * argv[]) { GtsSurface * s; GtsFile * fp; GtsFace * first = NULL; int c = 0; gboolean verbose = FALSE; if (!setlocale (LC_ALL, "POSIX")) g_warning ("cannot set locale to POSIX"); colormap = colormap_red_blue (); /* default */ /* parse options using getopt */ while (c != EOF) { #ifdef HAVE_GETOPT_LONG static struct option long_options[] = { {"cmap", required_argument, NULL, 'c'}, {"help", no_argument, NULL, 'h'}, {"verbose", no_argument, NULL, 'v'}, { NULL } }; int option_index = 0; switch ((c = getopt_long (argc, argv, "hvc:", long_options, &option_index))) { #else /* not HAVE_GETOPT_LONG */ switch ((c = getopt (argc, argv, "hvc:"))) { #endif /* not HAVE_GETOPT_LONG */ case 'c': { /* cmap */ FILE * fptr = fopen (optarg, "rt"); if (!fptr) { fprintf (stderr, "traverse: cannot open colormap file `%s'.\n", optarg); return 1; } colormap = colormap_read (fptr); fclose (fptr); break; } case 'v': /* verbose */ verbose = TRUE; break; case 'h': /* help */ fprintf (stderr, "Usage: traverse [OPTION] < file.gts > file.oogl\n" "Output an OOGL (geomview) surface colored according to the (graph) distance\n" "from a random face to the others\n" "\n" " -c FILE --cmap=FILE load FILE as colormap\n" " -v --verbose print statistics about the surface\n" " -h --help display this help and exit\n" "\n" "Reports bugs to %s\n", GTS_MAINTAINER); return 0; /* success */ break; case '?': /* wrong options */ fprintf (stderr, "Try `traverse --help' for more information.\n"); return 1; /* failure */ } } s = gts_surface_new (gts_surface_class (), GTS_FACE_CLASS (depth_face_class ()), gts_edge_class (), gts_vertex_class ()); fp = gts_file_new (stdin); if (gts_surface_read (s, fp)) { fputs ("traverse: file on standard input is not a valid GTS file\n", stderr); fprintf (stderr, "stdin:%d:%d: %s\n", fp->line, fp->pos, fp->error); return 1; /* failure */ } if (verbose) gts_surface_print_stats (s, stderr); gts_surface_foreach_face (s, (GtsFunc) pick_first_face, &first); gts_range_init (&depth_range); if (first) { GtsSurfaceTraverse * t = gts_surface_traverse_new (s, first); GtsFace * f; guint level; while ((f = gts_surface_traverse_next (t, &level))) { DEPTH_FACE (f)->depth = level; gts_range_add_value (&depth_range, level); } gts_surface_traverse_destroy (t); } gts_range_update (&depth_range); if (verbose) { fputs ("distance: ", stderr); gts_range_print (&depth_range, stderr); fputc ('\n', stderr); } gts_surface_write_oogl (s, stdout); return 0; }