More bug fixes for airport generation.

src/Airports/GenAirports/build.cxx

@@ -48,6 +48,7 @@
 #include <Geometry/trinodes.hxx>
 #include <Polygon/index.hxx>
 #include <Polygon/polygon.hxx>
+#include <Polygon/split.hxx>
 #include <Triangulate/trieles.hxx>
 
 #include "convex_hull.hxx"
@@ -141,18 +142,18 @@ void add_intermediate_nodes( int contour, const Point3D& start,
     Point3D p0 = start;
     Point3D p1 = end;
 
-    cout << "  add_intermediate_nodes()" << endl;
+    // cout << "  add_intermediate_nodes()" << endl;
     printf("   %.7f %.7f %.7f <=> %.7f %.7f %.7f\n",
 	   p0.x(), p0.y(), p0.z(), p1.x(), p1.y(), p1.z() );
 
     double xdist = fabs(p0.x() - p1.x());
     double ydist = fabs(p0.y() - p1.y());
-    cout << "xdist = " << xdist << "  ydist = " << ydist << endl;
+    // cout << "xdist = " << xdist << "  ydist = " << ydist << endl;
     x_err_min = xdist + 1.0;
     y_err_min = ydist + 1.0;
 
     if ( xdist > ydist ) {
-	cout << "use y = mx + b" << endl;
+	// cout << "use y = mx + b" << endl;
 
 	// sort these in a sensible order
 	Point3D p_min, p_max;
@@ -168,28 +169,28 @@ void add_intermediate_nodes( int contour, const Point3D& start,
 	b = p_max.y() - m * p_max.x();
 
 	// if ( temp ) {
-	cout << "m = " << m << " b = " << b << endl;
+	// cout << "m = " << m << " b = " << b << endl;
 	// }
 
 	current = nodes.begin();
 	last = nodes.end();
 	counter = 0;
 	for ( ; current != last; ++current ) {
-	    cout << counter << endl;
+	    // cout << counter << endl;
 
 	    if ( (current->x() > (p_min.x() + FG_EPSILON)) 
 		 && (current->x() < (p_max.x() - FG_EPSILON)) ) {
 
-		printf("found a potential candidate %.7f %.7f %.7f\n",
+		// printf("found a potential candidate %.7f %.7f %.7f\n",
-		       current->x(), current->y(), current->z() );
+		//        current->x(), current->y(), current->z() );
 
 		y_err = fabs(current->y() - (m * current->x() + b));
-		cout << "y_err = " << y_err << endl;
+		// cout << "y_err = " << y_err << endl;
 
 		if ( y_err < tgAirportEpsilon ) {
-		    cout << "FOUND EXTRA SEGMENT NODE (Y)" << endl;
+		    // cout << "FOUND EXTRA SEGMENT NODE (Y)" << endl;
-		    cout << p_min << " < " << *current << " < "
+		    // cout << p_min << " < " << *current << " < "
-		         << p_max << endl;
+		    //      << p_max << endl;
 		    found_extra = true;
 		    if ( y_err < y_err_min ) {
 			extra_index = counter;
@@ -200,7 +201,7 @@ void add_intermediate_nodes( int contour, const Point3D& start,
 	    ++counter;
 	}
     } else {
-	cout << "use x = m1 * y + b1" << endl;
+	// cout << "use x = m1 * y + b1" << endl;
 
 	// sort these in a sensible order
 	Point3D p_min, p_max;
@@ -217,7 +218,7 @@ void add_intermediate_nodes( int contour, const Point3D& start,
 
 	// bool temp = true;
 	// if ( temp ) {
-	cout << "  m1 = " << m1 << " b1 = " << b1 << endl;
+	// cout << "  m1 = " << m1 << " b1 = " << b1 << endl;
 	// }
 
 	// cout << "  should = 0 = " << fabs(p_min.x() - (m1 * p_min.y() + b1)) << endl;;
@@ -230,20 +231,20 @@ void add_intermediate_nodes( int contour, const Point3D& start,
 	    if ( (current->y() > (p_min.y() + FG_EPSILON)) 
 		 && (current->y() < (p_max.y() - FG_EPSILON)) ) {
 		
-		printf("found a potential candidate %.7f %.7f %.7f\n",
+		// printf("found a potential candidate %.7f %.7f %.7f\n",
-		       current->x(), current->y(), current->z() );
+		//        current->x(), current->y(), current->z() );
 
 		x_err = fabs(current->x() - (m1 * current->y() + b1));
-		cout << "x_err = " << x_err << endl;
+		// cout << "x_err = " << x_err << endl;
 
 		// if ( temp ) {
-		cout << "  (" << counter << ") x_err = " << x_err << endl;
+		// cout << "  (" << counter << ") x_err = " << x_err << endl;
 		// }
 
 		if ( x_err < tgAirportEpsilon ) {
-		    cout << "FOUND EXTRA SEGMENT NODE (X)" << endl;
+		    // cout << "FOUND EXTRA SEGMENT NODE (X)" << endl;
-		    cout << p_min << " < " << *current << " < "
+		    // cout << p_min << " < " << *current << " < "
-		         << p_max << endl;
+		    //      << p_max << endl;
 		    found_extra = true;
 		    if ( x_err < x_err_min ) {
 			extra_index = counter;
@@ -257,8 +258,8 @@ void add_intermediate_nodes( int contour, const Point3D& start,
 
     if ( found_extra ) {
 	// recurse with two sub segments
-	cout << "dividing " << p0 << " " << nodes[extra_index]
+	// cout << "dividing " << p0 << " " << nodes[extra_index]
-	     << " " << p1 << endl;
+	//      << " " << p1 << endl;
 	add_intermediate_nodes( contour, p0, nodes[extra_index], tmp_nodes, 
 				result );
 
@@ -600,12 +601,14 @@ void build_airport( string airport_raw, string_list& runways_raw,
 	cout << "result_b = " << result_b.contours() << endl;
 	accum = polygon_union( runway_b, accum );
 
+#if 0
 	// after clip, but before removing T intersections
 	char tmpa[256], tmpb[256];
 	sprintf( tmpa, "a%d", i );
 	sprintf( tmpb, "b%d", i );
 	write_polygon( result_a, tmpa );
 	write_polygon( result_b, tmpb );
+#endif
 
 	// print runway points
 	cout << "clipped runway pts (a)" << endl;
@@ -643,7 +646,7 @@ void build_airport( string airport_raw, string_list& runways_raw,
     // generate convex hull
     FGPolygon hull = convex_hull(apt_pts);
     FGPolygon base_poly = polygon_diff( hull, accum );
-    write_polygon( base_poly, "base-raw" );
+    // write_polygon( base_poly, "base-raw" );
 
     // add segments to polygons to remove any possible "T"
     // intersections
@@ -697,19 +700,21 @@ void build_airport( string airport_raw, string_list& runways_raw,
 	cout << "add nodes/remove dups runway = " << k << endl;
 	rwy_polys[k] = add_nodes_to_poly( rwy_polys[k], tmp_nodes );
 
+#if 0
 	char tmp[256];
 	sprintf( tmp, "r%d", k );
 	write_polygon( rwy_polys[k], tmp );
+#endif
 
         rwy_polys[k] = remove_dups( rwy_polys[k] );
         rwy_polys[k] = remove_bad_contours( rwy_polys[k] );
     }
     cout << "add nodes/remove dups base " << endl;
     base_poly = add_nodes_to_poly( base_poly, tmp_nodes );
-    write_polygon( base_poly, "base-add" );
+    // write_polygon( base_poly, "base-add" );
     base_poly = remove_dups( base_poly );
     base_poly = remove_bad_contours( base_poly );
-    write_polygon( base_poly, "base-final" );
+    // write_polygon( base_poly, "base-final" );
 
     // tesselate the polygons and prepair them for final output
 
@@ -880,6 +885,7 @@ void build_airport( string airport_raw, string_list& runways_raw,
     write_index( objpath, b, name );
 
     string holepath = root + "/AirportArea";
-    long int poly_index = poly_index_next();
+    // long int poly_index = poly_index_next();
-    write_boundary( holepath, b, hull, poly_index );
+    // write_boundary( holepath, b, hull, poly_index );
+    split_polygon( holepath, HoleArea, hull );
 }

src/Lib/Geometry/poly_support.cxx

@@ -360,9 +360,10 @@ void polygon_tesselate( const FGPolygon &p,
     // edge list (e), and a triangle neighbor list (n).
     // no new points on boundary (Y), no internal segment
     // splitting (YY), no quality refinement (q)
+    // Quite (Q)
 
     string tri_options;
-    tri_options = "pzYYen";
+    tri_options = "pzYYenQ";
     cout << "Triangulation with options = " << tri_options << endl;
 
     triangulate( (char *)tri_options.c_str(), &in, &out, &vorout );
@@ -431,6 +432,13 @@ void polygon_tesselate( const FGPolygon &p,
 // will modify the points_inside list for your polygon.
 
 FGPolygon polygon_tesselate_alt( FGPolygon &p ) {
+    FGPolygon result;
+    result.erase();
+
+    // Bail right away if polygon is empty
+    if ( p.contours() == 0 ) {
+	return result;
+    }
 
     // 1.  Robustly find a point inside each contour that is not
     //     inside any other contour
@@ -447,8 +455,6 @@ FGPolygon polygon_tesselate_alt( FGPolygon &p ) {
 
     // 3.  Convert the tesselated output to a list of tringles.
     //     basically a polygon with a contour for every triangle
-    FGPolygon result;
-    result.erase();
     for ( int i = 0; i < (int)trieles.size(); ++i ) {
 	FGTriEle t = trieles[i];
 	Point3D p1 = nodes[ t.get_n1() ];
@@ -614,9 +620,10 @@ static void contour_tesselate( FGContourNode *node, const FGPolygon &p,
     // edge list (e), and a triangle neighbor list (n).
     // no new points on boundary (Y), no internal segment
     // splitting (YY), no quality refinement (q)
+    // Quite (Q)
 
     string tri_options;
-    tri_options = "pzYYen";
+    tri_options = "pzYYenQ";
     cout << "Triangulation with options = " << tri_options << endl;
 
     triangulate( (char *)tri_options.c_str(), &in, &out, &vorout );
@@ -716,10 +723,12 @@ static Point3D point_inside_contour( FGContourNode *node, const FGPolygon &p ) {
 
     // build list of hole points
     for ( int i = 0; i < node->get_num_kids(); ++i ) {
-	contour_num = node->get_kid(i)->get_contour_num();
+        if ( node->get_kid( i ) != NULL ) {
-	hole_pts.push_back( p.get_point_inside( contour_num ) );
+	    contour_num = node->get_kid(i)->get_contour_num();
-	point_list contour = p.get_contour( contour_num );
+	    hole_pts.push_back( p.get_point_inside( contour_num ) );
-	hole_polys.add_contour( contour, 1 );
+	    point_list contour = p.get_contour( contour_num );
+	    hole_polys.add_contour( contour, 1 );
+        }
     }
 
     triele_list elelist;
@@ -828,12 +837,16 @@ static void build_contour_tree( FGContourNode *node,
     for ( int i = 0; i < node->get_num_kids(); ++i ) {
 	for ( int j = 0; j < node->get_num_kids(); ++j ) {
 	    if ( i != j ) {
-		if ( p.is_inside( i, j ) ) {
+		if ( (node->get_kid(i) != NULL)&&(node->get_kid(j) != NULL) ) {
-		    // need to remove contour j from the kid list
+		    if ( p.is_inside( i, j ) ) {
-		    avail[ node->get_kid( i ) -> get_contour_num() ] = 1;
+		        // need to remove contour j from the kid list
-		    node->remove_kid( i );
+		        avail[ node->get_kid( i ) -> get_contour_num() ] = 1;
-		    cout << "removing kid " << i << " which is inside of kid "
+		        node->remove_kid( i );
-			 << j << endl;
+		        cout << "removing kid " << i 
+                             << " which is inside of kid " << j << endl;
+		    }
+		} else {
+		    // one of these kids is already NULL, skip
 		}
 	    } else {
 		// doesn't make sense to check if a contour is inside itself

src/Lib/TriangleJRS/tri_support.c

@@ -34,7 +34,7 @@ void write_tri_data( struct triangulateio *out ) {
     fprintf(node, "%d 2 %d 0\n", 
 	    out->numberofpoints, out->numberofpointattributes);
     for (i = 0; i < out->numberofpoints; ++i) {
-	fprintf(node, "%d %.6f %.6f %.2f\n", 
+	fprintf(node, "%d %.8f %.8f %.2f\n", 
 		i, out->pointlist[2*i], out->pointlist[2*i + 1], 0.0);
     }
     fclose(node);
@@ -47,7 +47,7 @@ void write_tri_data( struct triangulateio *out ) {
 	    fprintf(ele, "%d ", out->trianglelist[i * out->numberofcorners + j]);
         }
         for (j = 0; j < out->numberoftriangleattributes; ++j) {
-	    fprintf(ele, "%.6f ", 
+	    fprintf(ele, "%.8f ", 
 		    out->triangleattributelist[i 
 					      * out->numberoftriangleattributes
 					      + j]
@@ -67,12 +67,12 @@ void write_tri_data( struct triangulateio *out ) {
     }
     fprintf(fp, "%d\n", out->numberofholes);
     for (i = 0; i < out->numberofholes; ++i) {
-	fprintf(fp, "%d %.6f %.6f\n", 
+	fprintf(fp, "%d %.8f %.8f\n", 
 		i, out->holelist[2*i], out->holelist[2*i + 1]);
     }
     fprintf(fp, "%d\n", out->numberofregions);
     for (i = 0; i < out->numberofregions; ++i) {
-	fprintf(fp, "%d %.6f %.6f %.6f\n", 
+	fprintf(fp, "%d %.8f %.8f %.8f\n", 
 		i, out->regionlist[4*i], out->regionlist[4*i + 1],
 		out->regionlist[4*i + 2]);
     }