1
0
Fork 0

Lot's pf success on the quest to a full file parse.

Just three airports crashed, and no hangs.  One of the crashes is still NZSP
I don't know if I'll fix this soon.  It's the south pole problem with stuff
crossing the IDL.
- Added dome debug support for generating shapefiles from TGPolygon.
  This allows directly creating shapefiles instead of chopping a poly, then
  having to run poly2ogr on it.  The API is flexible, allowing multiple
  datafiles, and layers open simutaneously, and the features in the shapefiles
  have a description field that can be set.
- Added a restart-id= which allows a bash script to check the exit value of genapts,
  and automatically restart on the airport after the one that crashed.
  I've included an example script I have used.
- Stability fixes:
  1) actually, this is both stability, and graphical fix.  Some airport linear
     features and polygon contours were using bezier nodes for straight lines.
     When the control point lies just outside the line, the parser would attempt to
     draw a 'hook'  unfortunately, some of my calcs return NAN, and sometimes you just
     get crazy looking polys to attempt to clip - not good.
     Fix is to check cubic and quadratic nodes to see if they should really be linear.
     Substitute linear segments where apropriate.
  2) Along the same line as 1, There are some really short bezier segments, that were
     getting subdivided into the hard coded 8 segments.  Now, if the approximately linear
     distance of the curve is less than 4 meters, the curve is broken into 1/2 meter segments.
     Likewise, if the curve is greater than 800 meters, the curve is broken into 100 meter
     segments.
  3) Sometimes, during clipping, degenerate ontours are created.  Although these are cleaned
     before triangulation, there were still issues with generating the accumulation poly with
     the union operator, causing a crash.  I now simplify, and reduce degeneracy on all polys
     before clipping against the accumulator.
- Extra special bonus - successful completeion of an airport is printed in green :)
This commit is contained in:
Peter Sadrozinski 2012-02-14 18:08:11 -05:00 committed by Christian Schmitt
parent a40b72b15b
commit 071d98b3bc
18 changed files with 845 additions and 340 deletions

View file

@ -24,6 +24,7 @@ target_link_libraries(genapts
TriangleJRS
${SIMGEAR_CORE_LIBRARIES}
${SIMGEAR_CORE_LIBRARY_DEPENDENCIES}
${GDAL_LIBRARY}
${GPC_LIBRARY}
${NEWMAT_LIBRARY})

View file

@ -25,6 +25,7 @@ target_link_libraries(genapts850
Polygon Geometry
Array Optimize Output poly2tri
TriangleJRS
${GDAL_LIBRARY}
${SIMGEAR_CORE_LIBRARIES}
${SIMGEAR_CORE_LIBRARY_DEPENDENCIES}
${GPC_LIBRARY}

View file

@ -29,6 +29,12 @@
#include "elevations.hxx"
#include <stdio.h>
string SGLOG_GREEN = "\033[0;32m";
string SGLOG_NORMAL = "\033[0m";
Airport::Airport( int c, char* def)
{
int numParams;
@ -426,6 +432,14 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
}
}
#if 0
if (icao == "SSOK")
{
SG_LOG(SG_GENERAL, SG_ALERT, "Injecting error at icao " << icao );
exit(-1);
}
#endif
// Starting to clip the polys
gettimeofday(&build_start, NULL);
@ -435,7 +449,17 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
for ( unsigned int i=0; i<features.size(); i++ )
{
SG_LOG(SG_GENERAL, SG_INFO, "Build Feature Poly " << i << " of " << features.size() << " : " << features[i]->GetDescription() );
features[i]->BuildBtg( altitude, &line_polys, &line_tps, &line_accum, &rwy_lights );
#if 0
if ( i == 22 ) {
features[i]->BuildBtg( altitude, &line_polys, &line_tps, &line_accum, &rwy_lights, true );
} else {
//features[i]->BuildBtg( altitude, &line_polys, &line_tps, &line_accum, &rwy_lights, false );
}
#else
features[i]->BuildBtg( altitude, &line_polys, &line_tps, &line_accum, &rwy_lights, false );
#endif
}
}
else
@ -526,6 +550,7 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
gettimeofday(&log_time, NULL);
SG_LOG( SG_GENERAL, SG_ALERT, "Finished building pavements for " << icao << " at " << ctime(&log_time.tv_sec) );
#if 0
// Build runway shoulders here
for ( unsigned int i=0; i<runways.size(); i++ )
{
@ -533,9 +558,10 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
if ( runways[i]->GetsShoulder() )
{
runways[i]->BuildShoulder( altitude, &rwy_polys, &rwy_tps, &accum );
runways[i]->BuildShoulder( altitude, &rwy_polys, &rwy_tps, &accum, &apt_base, &apt_clearing );
}
}
#endif
// build the base and clearing if there's a boundary
if (boundary)
@ -671,6 +697,7 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
poly = remove_spikes( poly );
poly = remove_dups( poly );
poly = remove_bad_contours( poly );
poly = remove_small_contours( poly );
rwy_polys[k].set_poly( poly );
}
@ -688,6 +715,7 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
poly = remove_spikes( poly );
poly = remove_dups( poly );
poly = remove_bad_contours( poly );
poly = remove_small_contours( poly );
pvmt_polys[k].set_poly( poly );
}
@ -701,10 +729,12 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
poly = remove_cycles( poly );
poly = remove_dups( poly );
poly = remove_bad_contours( poly );
poly = tgPolygonSimplify( poly );
poly = remove_tiny_contours( poly );
poly = remove_spikes( poly );
poly = remove_dups( poly );
poly = remove_bad_contours( poly );
poly = remove_tiny_contours( poly );
line_polys[k].set_poly( poly );
}
@ -728,6 +758,7 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
base_poly = remove_spikes( base_poly );
base_poly = remove_dups( base_poly );
base_poly = remove_bad_contours( base_poly );
base_poly = remove_small_contours( base_poly );
gettimeofday(&cleanup_end, NULL);
timersub(&cleanup_end, &cleanup_start, &cleanup_time);
@ -778,7 +809,7 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
TGPolygon poly = pvmt_polys[i].get_poly();
#if 0
if ( i == 1 ) {
if ( i == 0 ) {
SG_LOG(SG_GENERAL, SG_INFO, "Problem poly: " << poly );
tgChopNormalPolygon( "/home/pete", "Base", poly, false );
@ -814,10 +845,10 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
TGPolygon poly = line_polys[i].get_poly();
#if 0
if ( i == 627 ) {
if ( i == 282 ) {
SG_LOG(SG_GENERAL, SG_INFO, "Problem poly: " << poly );
tgChopNormalPolygon( "/home/pete", "Base", poly, false );
tgChopNormalPolygon( "/home/pete/", "Base", poly, false );
verbose_triangulation = true;
} else {
verbose_triangulation = false;
@ -839,7 +870,7 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
{
SG_LOG(SG_GENERAL, SG_INFO, "Problem poly: " << base_poly );
tgChopNormalPolygon( "/home/pete", "Base", base_poly, false );
tgChopNormalPolygon( "/home/pete/", "Base", base_poly, false );
verbose_triangulation = true;
}
#endif
@ -1099,7 +1130,7 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
// Now build the fitted airport surface ...
// calculation min/max coordinates of airport area
SG_LOG(SG_GENERAL, SG_INFO, " calculation min/max coordinates of airport area");
SG_LOG(SG_GENERAL, SG_DEBUG, " calculation min/max coordinates of airport area");
Point3D min_deg(9999.0, 9999.0, 0), max_deg(-9999.0, -9999.0, 0);
for ( unsigned int j = 0; j < nodes.get_node_list().size(); ++j )
@ -1127,7 +1158,7 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
}
}
SG_LOG(SG_GENERAL, SG_INFO, "Before extending for lights: min = " << min_deg << " max = " << max_deg );
SG_LOG(SG_GENERAL, SG_DEBUG, "Before extending for lights: min = " << min_deg << " max = " << max_deg );
// extend the min/max coordinates of airport area to cover all
// lights as well
@ -1159,8 +1190,6 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
}
}
SG_LOG(SG_GENERAL, SG_DEBUG, " done " );
// Extend the area a bit so we don't have wierd things on the edges
double dlon = max_deg.lon() - min_deg.lon();
double dlat = max_deg.lat() - min_deg.lat();
@ -1571,4 +1600,6 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
// write_boundary( holepath, b, hull, poly_index );
tgChopNormalPolygon( holepath, "Hole", divided_base, true );
tgChopNormalPolygon( holepath, "Airport", apt_clearing, false );
SG_LOG( SG_GENERAL, SG_ALERT, SGLOG_GREEN << "\nSUCCESS generating " << icao << SGLOG_NORMAL << "\n" );
}

View file

@ -58,6 +58,11 @@ public:
}
}
int NumFeatures( void )
{
return features.size();
}
void SetBoundary( ClosedPoly* bndry )
{
boundary = bndry;

View file

@ -9,6 +9,12 @@
#include <simgear/debug/logstream.hxx>
// TEMP...
inline double LinearDistance( Point3D p0, Point3D p1 )
{
return sqrt( pow( fabs( p1.x()-p0.x() ), 2 ) +
pow( fabs( p1.y()-p0.y() ), 2 ) );
}
inline Point3D CalculateLinearLocation( Point3D p0, Point3D p1, double t )
{
Point3D result;
@ -21,6 +27,11 @@ inline Point3D CalculateLinearLocation( Point3D p0, Point3D p1, double t )
return result;
}
inline double QuadraticDistance( Point3D p0, Point3D cp, Point3D p1 )
{
return LinearDistance( p0, cp ) + LinearDistance( cp, p1 );
}
inline Point3D CalculateQuadraticLocation( Point3D p0, Point3D cp, Point3D p1, double t )
{
Point3D result;
@ -34,7 +45,11 @@ inline Point3D CalculateQuadraticLocation( Point3D p0, Point3D cp, Point3D p1, d
return result;
}
// TODO: Should be in a math library
inline double CubicDistance( Point3D p0, Point3D cp0, Point3D cp1, Point3D p1 )
{
return LinearDistance( p0, cp0 ) + LinearDistance( cp0, cp1 ) + LinearDistance( cp1, p1 );
}
inline Point3D CalculateCubicLocation( Point3D p0, Point3D cp0, Point3D cp1, Point3D p1, double t )
{
Point3D result;

View file

@ -100,6 +100,7 @@ void ClosedPoly::AddNode( BezNode* node )
}
}
#if 0
void ClosedPoly::CreateConvexHull( void )
{
TGPolygon convexHull;
@ -122,6 +123,7 @@ void ClosedPoly::CreateConvexHull( void )
SG_LOG(SG_GENERAL, SG_ALERT, "Boundary size too small: " << boundary->size() << ". Ignoring..." );
}
}
#endif
void ClosedPoly::CloseCurContour()
{
@ -132,7 +134,7 @@ void ClosedPoly::CloseCurContour()
if (cur_feature)
{
SG_LOG(SG_GENERAL, SG_DEBUG, "We still have an active linear feature - add the first node to close it");
cur_feature->Finish(true);
cur_feature->Finish(true, features.size() );
features.push_back(cur_feature);
cur_feature = NULL;
@ -145,7 +147,7 @@ void ClosedPoly::CloseCurContour()
boundary = cur_contour;
// generate the convex hull from the bezcontour node locations
CreateConvexHull();
// CreateConvexHull();
cur_contour = NULL;
}
@ -158,18 +160,18 @@ void ClosedPoly::CloseCurContour()
void ClosedPoly::ConvertContour( BezContour* src, point_list *dst )
{
BezNode* prevNode;
BezNode* curNode;
BezNode* nextNode;
Point3D prevLoc;
Point3D curLoc;
Point3D nextLoc;
Point3D cp1;
Point3D cp2;
int curve_type = CURVE_LINEAR;
unsigned int i;
double total_dist;
double meter_dist = 1.0f/96560.64f;
int num_segs = BEZIER_DETAIL;
SG_LOG(SG_GENERAL, SG_DEBUG, "Creating a contour with " << src->size() << " nodes");
@ -177,23 +179,11 @@ void ClosedPoly::ConvertContour( BezContour* src, point_list *dst )
dst->empty();
// iterate through each bezier node in the contour
for (i=0; i <= src->size()-1; i++)
for (unsigned int i = 0; i <= src->size()-1; i++)
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\nHandling Node " << i << "\n\n");
if (i == 0)
{
// set prev node to last in the contour, as all contours must be closed
prevNode = src->at( src->size()-1 );
}
else
{
// otherwise, it's just the previous index
prevNode = src->at( i-1 );
}
curNode = src->at(i);
if (i < src->size() - 1)
{
nextNode = src->at(i+1);
@ -204,37 +194,6 @@ void ClosedPoly::ConvertContour( BezContour* src, point_list *dst )
nextNode = src->at(0);
}
// determine the type of curve from prev (just to get correct prev location)
// once we start drawing the curve from cur to next, we can just remember the prev loc
if (prevNode->HasNextCp())
{
// curve from prev is cubic or quadratic
if(curNode->HasPrevCp())
{
// curve from prev is cubic : calculate the last location on the curve
prevLoc = CalculateCubicLocation( prevNode->GetLoc(), prevNode->GetNextCp(), curNode->GetPrevCp(), curNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (BEZIER_DETAIL-1) );
}
else
{
// curve from prev is quadratic : use prev node next cp
prevLoc = CalculateQuadraticLocation( prevNode->GetLoc(), prevNode->GetNextCp(), curNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (BEZIER_DETAIL-1) );
}
}
else
{
// curve from prev is quadratic or linear
if( curNode->HasPrevCp() )
{
// curve from prev is quadratic : calculate the last location on the curve
prevLoc = CalculateQuadraticLocation( prevNode->GetLoc(), curNode->GetPrevCp(), curNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (BEZIER_DETAIL-1) );
}
else
{
// curve from prev is linear : just use prev node location
prevLoc = prevNode->GetLoc();
}
}
// now determine how we will iterate from current node to next node
if( curNode->HasNextCp() )
{
@ -245,12 +204,14 @@ void ClosedPoly::ConvertContour( BezContour* src, point_list *dst )
curve_type = CURVE_CUBIC;
cp1 = curNode->GetNextCp();
cp2 = nextNode->GetPrevCp();
total_dist = CubicDistance( curNode->GetLoc(), cp1, cp2, nextNode->GetLoc() );
}
else
{
// curve is quadratic using current nodes cp as the cp
curve_type = CURVE_QUADRATIC;
cp1 = curNode->GetNextCp();
total_dist = QuadraticDistance( curNode->GetLoc(), cp1, nextNode->GetLoc() );
}
}
else
@ -261,28 +222,94 @@ void ClosedPoly::ConvertContour( BezContour* src, point_list *dst )
// curve is quadratic using next nodes cp as the cp
curve_type = CURVE_QUADRATIC;
cp1 = nextNode->GetPrevCp();
total_dist = QuadraticDistance( curNode->GetLoc(), cp1, nextNode->GetLoc() );
}
else
{
// curve is linear
curve_type = CURVE_LINEAR;
total_dist = LinearDistance( curNode->GetLoc(), nextNode->GetLoc() );
}
}
#if 0
double num_meters = total_dist / meter_dist;
// If total distance is < 4 meters, then we need to modify num Segments so that each segment >= 0.5 meter
if (num_meters < 4.0f)
{
num_segs = ((int)num_meters + 1) * 2;
}
else if (num_meters > 800.0f)
{
num_segs = num_meters / 100.0f + 1;
}
else
{
num_segs = 8;
}
#endif
double num_meters = total_dist / meter_dist;
if (num_meters < 4.0f)
{
if (curve_type != CURVE_LINEAR)
{
// If total distance is < 4 meters, then we need to modify num Segments so that each segment >= 1/2 meter
num_segs = ((int)num_meters + 1) * 2;
SG_LOG(SG_GENERAL, SG_DEBUG, "Segment from (" << curNode->GetLoc().x() << "," << curNode->GetLoc().y() << ") to (" << nextNode->GetLoc().x() << "," << nextNode->GetLoc().y() << ")" );
SG_LOG(SG_GENERAL, SG_DEBUG, " Distance is " << num_meters << " ( < 4.0) so num_segs is " << num_segs );
}
else
{
num_segs = 1;
}
}
else if (num_meters > 800.0f)
{
// If total distance is > 800 meters, then we need to modify num Segments so that each segment <= 100 meters
num_segs = num_meters / 100.0f + 1;
SG_LOG(SG_GENERAL, SG_DEBUG, "Segment from (" << curNode->GetLoc().x() << "," << curNode->GetLoc().y() << ") to (" << nextNode->GetLoc().x() << "," << nextNode->GetLoc().y() << ")" );
SG_LOG(SG_GENERAL, SG_DEBUG, " Distance is " << num_meters << " ( > 100.0) so num_segs is " << num_segs );
}
else
{
if (curve_type != CURVE_LINEAR)
{
num_segs = 8;
// num_segs = 16;
SG_LOG(SG_GENERAL, SG_DEBUG, "Segment from (" << curNode->GetLoc().x() << "," << curNode->GetLoc().y() << ") to (" << nextNode->GetLoc().x() << "," << nextNode->GetLoc().y() << ")" );
SG_LOG(SG_GENERAL, SG_DEBUG, " Distance is " << num_meters << " (OK) so num_segs is " << num_segs );
}
else
{
// make sure linear segments don't got over 100m
//num_segs = 1;
num_segs = num_meters / 100.0f + 1;
}
}
// if only one segment, revert to linear
if (num_segs == 1)
{
curve_type = CURVE_LINEAR;
}
// initialize current location
curLoc = curNode->GetLoc();
if (curve_type != CURVE_LINEAR)
{
for (int p=0; p<BEZIER_DETAIL; p++)
for (int p=0; p<num_segs; p++)
{
// calculate next location
if (curve_type == CURVE_QUADRATIC)
{
nextLoc = CalculateQuadraticLocation( curNode->GetLoc(), cp1, nextNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (p+1) );
nextLoc = CalculateQuadraticLocation( curNode->GetLoc(), cp1, nextNode->GetLoc(), (1.0f/num_segs) * (p+1) );
}
else
{
nextLoc = CalculateCubicLocation( curNode->GetLoc(), cp1, cp2, nextNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (p+1) );
nextLoc = CalculateCubicLocation( curNode->GetLoc(), cp1, cp2, nextNode->GetLoc(), (1.0f/num_segs) * (p+1) );
}
// add the pavement vertex
@ -298,8 +325,38 @@ void ClosedPoly::ConvertContour( BezContour* src, point_list *dst )
SG_LOG(SG_GENERAL, SG_DEBUG, " add bezier node (type " << curve_type << ") at (" << curLoc.x() << "," << curLoc.y() << ")");
}
// now set set cur location for the next iteration
curLoc = nextLoc;
}
}
else
{
// nextLoc = nextNode->GetLoc();
// // just add the one vertex - linear
// dst->push_back( curLoc );
// SG_LOG(SG_GENERAL, SG_DEBUG, "adding Linear Anchor node at (" << curLoc.x() << "," << curLoc.y() << ")");
if (num_segs > 1)
{
for (int p=0; p<num_segs; p++)
{
// calculate next location
nextLoc = CalculateLinearLocation( curNode->GetLoc(), nextNode->GetLoc(), (1.0f/num_segs) * (p+1) );
// add the feature vertex
dst->push_back( curLoc );
if (p==0)
{
SG_LOG(SG_GENERAL, SG_DEBUG, "adding Linear anchor node at (" << curLoc.x() << "," << curLoc.y() << ")");
}
else
{
SG_LOG(SG_GENERAL, SG_DEBUG, " add linear node at (" << curLoc.x() << "," << curLoc.y() << ")");
}
// now set set prev and cur locations for the next iteration
prevLoc = curLoc;
curLoc = nextLoc;
}
}
@ -307,13 +364,18 @@ void ClosedPoly::ConvertContour( BezContour* src, point_list *dst )
{
nextLoc = nextNode->GetLoc();
// just add the one vertex - linear
// just add the one vertex - dist is small
dst->push_back( curLoc );
SG_LOG(SG_GENERAL, SG_DEBUG, "adding Linear Anchor node at (" << curLoc.x() << "," << curLoc.y() << ")");
curLoc = nextLoc;
}
}
}
}
#if 0
void ExpandPoint( Point3D *prev, Point3D *cur, Point3D *next, double expand_by, double *heading, double *offset )
{
double offset_dir;
@ -352,7 +414,9 @@ void ExpandPoint( Point3D *prev, Point3D *cur, Point3D *next, double expand_by,
SG_LOG(SG_GENERAL, SG_DEBUG, "heading is " << *heading << " distance is " << *offset );
}
#endif
#if 0
void ClosedPoly::ExpandContour( point_list& src, TGPolygon& dst, double dist )
{
point_list expanded_boundary;
@ -443,6 +507,7 @@ void ClosedPoly::ExpandContour( point_list& src, TGPolygon& dst, double dist )
dst.add_contour( expanded_boundary, 9 );
}
#endif
// finish the poly - convert to TGPolygon, and tesselate
void ClosedPoly::Finish()
@ -465,7 +530,7 @@ void ClosedPoly::Finish()
// and add it to the geometry
pre_tess.add_contour( dst_contour, 0 );
// The convert the hole contours
// Then convert the hole contours
for (unsigned int i=0; i<holes.size(); i++)
{
dst_contour.clear();
@ -535,24 +600,9 @@ int ClosedPoly::BuildBtg( float alt_m, superpoly_list* rwy_polys, texparams_list
SG_LOG(SG_GENERAL, SG_DEBUG, "BuildBtg: original poly has " << pre_tess.contours() << " contours");
// do this before clipping and generating the base
pre_tess = remove_bad_contours( pre_tess );
pre_tess = remove_dups( pre_tess );
pre_tess = tgPolygonSimplify( pre_tess );
pre_tess = reduce_degeneracy( pre_tess );
//for (int c=0; c<pre_tess.contours(); c++)
//{
// for (int pt=0; pt<pre_tess.contour_size(c); pt++)
// {
// SG_LOG(SG_GENERAL, SG_DEBUG, "BuildBtg: contour " << c << " pt " << pt << ": (" << pre_tess.get_pt(c, pt).x() << "," << pre_tess.get_pt(c, pt).y() << ")" );
// }
//}
//SG_LOG(SG_GENERAL, SG_DEBUG, "BuildBtg: original poly has " << pre_tess.contours() << " contours");
//for (int i=0; i<pre_tess.contours(); i++)
//{
// SG_LOG(SG_GENERAL, SG_DEBUG, "BuildBtg: original countour " << i << " has " << pre_tess.contour_size(i) << " points" );
//}
// grow pretess by a little bit
//pre_tess = tgPolygonExpand( pre_tess, 0.05); // 5cm
@ -566,9 +616,9 @@ int ClosedPoly::BuildBtg( float alt_m, superpoly_list* rwy_polys, texparams_list
}
#if 1
TGPolygon clipped = tgPolygonDiffClipper( pre_tess, *accum );
#else
TGPolygon clipped = tgPolygonDiff( pre_tess, *accum );
#else
TGPolygon clipped = tgPolygonDiffClipper( pre_tess, *accum );
#endif
SG_LOG(SG_GENERAL, SG_DEBUG, "BuildBtg: clipped poly has " << clipped.contours() << " contours");
for (int i=0; i<clipped.contours(); i++)
@ -587,9 +637,9 @@ int ClosedPoly::BuildBtg( float alt_m, superpoly_list* rwy_polys, texparams_list
rwy_polys->push_back( sp );
SG_LOG(SG_GENERAL, SG_DEBUG, "clipped = " << clipped.contours());
#if 1
*accum = tgPolygonUnionClipper( pre_tess, *accum );
#else
*accum = tgPolygonUnion( pre_tess, *accum );
#else
*accum = tgPolygonUnionClipper( pre_tess, *accum );
#endif
tp = TGTexParams( pre_tess.get_pt(0,0), 5.0, 5.0, texture_heading );
texparams->push_back( tp );

View file

@ -6,26 +6,33 @@
#include <Geometry/poly_support.hxx>
// for debugging clipping errors
#include <Polygon/chop.hxx>
#include "beznode.hxx"
#include "linearfeature.hxx"
#include "math.h"
void LinearFeature::ConvertContour( BezContour* src, bool closed )
{
BezNode* prevNode;
BezNode* curNode;
BezNode* nextNode;
Point3D prevLoc;
Point3D curLoc;
Point3D nextLoc;
Point3D cp1;
Point3D cp2;
int curve_type = CURVE_LINEAR;
double total_dist, linear_dist;
double meter_dist = 1.0f/96560.64f;
double theta1, theta2;
int num_segs = BEZIER_DETAIL;
Marking* cur_mark = NULL;
Lighting* cur_light = NULL;
SG_LOG(SG_GENERAL, SG_DEBUG, " LinearFeature::ConvertContour - Creating a contour with " << src->size() << " nodes");
// clear anything in the point list
@ -34,19 +41,6 @@ void LinearFeature::ConvertContour( BezContour* src, bool closed )
// iterate through each bezier node in the contour
for (unsigned int i=0; i <= src->size()-1; i++)
{
SG_LOG(SG_GENERAL, SG_DEBUG, " LinearFeature::ConvertContour: Handling Node " << i << "\n\n");
if (i == 0)
{
// set prev node to last in the contour, as all contours must be closed
prevNode = src->at( src->size()-1 );
}
else
{
// otherwise, it's just the previous index
prevNode = src->at( i-1 );
}
curNode = src->at(i);
if (i < src->size() - 1)
@ -132,38 +126,6 @@ void LinearFeature::ConvertContour( BezContour* src, bool closed )
}
////////////////////////////////////////////////////////////////////////////////////
// determine the type of curve from prev (just to get correct prev location)
// once we start drawing the curve from cur to next, we can just remember the prev loc
if (prevNode->HasNextCp())
{
// curve from prev is cubic or quadratic
if(curNode->HasPrevCp())
{
// curve from prev is cubic : calculate the last location on the curve
prevLoc = CalculateCubicLocation( prevNode->GetLoc(), prevNode->GetNextCp(), curNode->GetPrevCp(), curNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (BEZIER_DETAIL-1) );
}
else
{
// curve from prev is quadratic : use prev node next cp
prevLoc = CalculateQuadraticLocation( prevNode->GetLoc(), prevNode->GetNextCp(), curNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (BEZIER_DETAIL-1) );
}
}
else
{
// curve from prev is quadratic or linear
if( curNode->HasPrevCp() )
{
// curve from prev is quadratic : calculate the last location on the curve
prevLoc = CalculateQuadraticLocation( prevNode->GetLoc(), curNode->GetPrevCp(), curNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (BEZIER_DETAIL-1) );
}
else
{
// curve from prev is linear : just use prev node location
prevLoc = prevNode->GetLoc();
}
}
// now determine how we will iterate from current node to next node
if( curNode->HasNextCp() )
{
@ -174,12 +136,17 @@ void LinearFeature::ConvertContour( BezContour* src, bool closed )
curve_type = CURVE_CUBIC;
cp1 = curNode->GetNextCp();
cp2 = nextNode->GetPrevCp();
total_dist = CubicDistance( curNode->GetLoc(), cp1, cp2, nextNode->GetLoc() );
theta1 = SGMiscd::rad2deg( CalculateTheta( curNode->GetLoc(), cp1, nextNode->GetLoc()) );
theta2 = SGMiscd::rad2deg( CalculateTheta( curNode->GetLoc(), cp2, nextNode->GetLoc()) );
}
else
{
// curve is quadratic using current nodes cp as the cp
curve_type = CURVE_QUADRATIC;
cp1 = curNode->GetNextCp();
total_dist = QuadraticDistance( curNode->GetLoc(), cp1, nextNode->GetLoc() );
theta1 = SGMiscd::rad2deg( CalculateTheta( curNode->GetLoc(), cp1, nextNode->GetLoc()) );
}
}
else
@ -190,28 +157,123 @@ void LinearFeature::ConvertContour( BezContour* src, bool closed )
// curve is quadratic using next nodes cp as the cp
curve_type = CURVE_QUADRATIC;
cp1 = nextNode->GetPrevCp();
total_dist = QuadraticDistance( curNode->GetLoc(), cp1, nextNode->GetLoc() );
theta1 = SGMiscd::rad2deg( CalculateTheta( curNode->GetLoc(), cp1, nextNode->GetLoc()) );
}
else
{
// curve is linear
curve_type = CURVE_LINEAR;
total_dist = LinearDistance( curNode->GetLoc(), nextNode->GetLoc() );
}
}
// One more test - some people are using bezier curves to draw straight lines - this can cause a bit of havoc...
// Sometimes, the control point lies just beyond the final point. We try to make a 'hook' at the end, which makes some really bad polys
// Just convert the entire segment to linear
// this can be detected in quadratic curves (current issue in LFKJ) when the contol point lies within the line generated from point 1 to point 2
// theat close to 180 at the control point to the cur node and next node
if ( curve_type == CURVE_QUADRATIC )
{
if ( (abs(theta1 - 180.0) < 5.0 ) || (abs(theta1) < 5.0 ) || (isnan(theta1)) )
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\nLinearFeature: Quadtratic curve with cp in line : convert to linear: " << description << ": theta is " << theta1 );
curve_type = CURVE_LINEAR;
}
else
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\nLinearFeature: Quadtratic curve withOUT cp in line : keep quadtratic: " << description << ": theta is " << theta1 );
}
}
if ( curve_type == CURVE_CUBIC )
{
if ( (abs(theta1 - 180.0) < 5.0 ) || (abs(theta1) < 5.0 ) || (isnan(theta1)) )
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\nLinearFeature: Cubic curve with cp1 in line : " << description << ": theta is " << theta1 );
if ( (abs(theta2 - 180.0) < 5.0 ) || (abs(theta2) < 5.0 ) || (isnan(theta2)) )
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\n and cp2 in line : " << description << ": theta is " << theta2 << " CONVERTING TO LINEAR" );
curve_type = CURVE_LINEAR;
}
else
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\n BUT cp2 NOT in line : " << description << ": theta is " << theta2 );
}
}
else
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\nLinearFeature: Cubic curve withOUT cp1 in line : keep quadtratic: " << description << ": theta is " << theta1 );
if ( (abs(theta2 - 180.0) < 5.0 ) || (abs(theta2) < 5.0 ) || (isnan(theta2)) )
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\n BUT cp2 IS in line : " << description << ": theta is " << theta2 );
}
else
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\n AND cp2 NOT in line : " << description << ": theta is " << theta2 );
}
}
}
double num_meters = total_dist / meter_dist;
if (num_meters < 4.0f)
{
if (curve_type != CURVE_LINEAR)
{
// If total distance is < 4 meters, then we need to modify num Segments so that each segment >= 1/2 meter
num_segs = ((int)num_meters + 1) * 2;
SG_LOG(SG_GENERAL, SG_DEBUG, "Segment from (" << curNode->GetLoc().x() << "," << curNode->GetLoc().y() << ") to (" << nextNode->GetLoc().x() << "," << nextNode->GetLoc().y() << ")" );
SG_LOG(SG_GENERAL, SG_DEBUG, " Distance is " << num_meters << " ( < 4.0) so num_segs is " << num_segs );
}
else
{
num_segs = 1;
}
}
else if (num_meters > 800.0f)
{
// If total distance is > 800 meters, then we need to modify num Segments so that each segment <= 100 meters
num_segs = num_meters / 100.0f + 1;
SG_LOG(SG_GENERAL, SG_DEBUG, "Segment from (" << curNode->GetLoc().x() << "," << curNode->GetLoc().y() << ") to (" << nextNode->GetLoc().x() << "," << nextNode->GetLoc().y() << ")" );
SG_LOG(SG_GENERAL, SG_DEBUG, " Distance is " << num_meters << " ( > 100.0) so num_segs is " << num_segs );
}
else
{
if (curve_type != CURVE_LINEAR)
{
num_segs = 8;
SG_LOG(SG_GENERAL, SG_DEBUG, "Segment from (" << curNode->GetLoc().x() << "," << curNode->GetLoc().y() << ") to (" << nextNode->GetLoc().x() << "," << nextNode->GetLoc().y() << ")" );
SG_LOG(SG_GENERAL, SG_DEBUG, " Distance is " << num_meters << " (OK) so num_segs is " << num_segs );
}
else
{
num_segs = 1;
}
}
// if only one segment, revert to linear
if (num_segs == 1)
{
curve_type = CURVE_LINEAR;
}
// initialize current location
curLoc = curNode->GetLoc();
if (curve_type != CURVE_LINEAR)
{
for (int p=0; p<BEZIER_DETAIL; p++)
for (int p=0; p<num_segs; p++)
{
// calculate next location
if (curve_type == CURVE_QUADRATIC)
{
nextLoc = CalculateQuadraticLocation( curNode->GetLoc(), cp1, nextNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (p+1) );
nextLoc = CalculateQuadraticLocation( curNode->GetLoc(), cp1, nextNode->GetLoc(), (1.0f/num_segs) * (p+1) );
}
else
{
nextLoc = CalculateCubicLocation( curNode->GetLoc(), cp1, cp2, nextNode->GetLoc(), (1.0f/BEZIER_DETAIL) * (p+1) );
nextLoc = CalculateCubicLocation( curNode->GetLoc(), cp1, cp2, nextNode->GetLoc(), (1.0f/num_segs) * (p+1) );
}
// add the feature vertex
@ -227,26 +289,16 @@ void LinearFeature::ConvertContour( BezContour* src, bool closed )
}
// now set set prev and cur locations for the next iteration
prevLoc = curLoc;
curLoc = nextLoc;
}
}
else
{
// For linear features, sometime long linear lines confuse the tesselator. Add intermediate nodes to keep the rectangles from
// getting too long.
double az1 = 0.0f;
double az2 = 0.0f;
double dist = 0.0f;
// calculate linear distance to determine how many segments we want
Point3D destLoc = nextNode->GetLoc();
geo_inverse_wgs_84( curLoc.y(), curLoc.x(), destLoc.y(), destLoc.x(), &az1, &az2, &dist);
if (dist > 100.0)
if (num_segs > 1)
{
int num_segs = (dist / 100.0f) + 1;
for (int p=0; p<num_segs; p++)
{
// calculate next location
@ -265,10 +317,7 @@ void LinearFeature::ConvertContour( BezContour* src, bool closed )
}
// now set set prev and cur locations for the next iteration
prevLoc = curLoc;
curLoc = nextLoc;
SG_LOG(SG_GENERAL, SG_DEBUG, "Set prevLoc = (" << prevLoc.x() << "," << prevLoc.y() << ") and curLoc = (" << curLoc.x() << "," << curLoc.y() << ")" );
}
}
else
@ -280,16 +329,11 @@ void LinearFeature::ConvertContour( BezContour* src, bool closed )
SG_LOG(SG_GENERAL, SG_DEBUG, "adding Linear Anchor node at (" << curLoc.x() << "," << curLoc.y() << ")");
prevLoc = curLoc;
curLoc = nextLoc;
SG_LOG(SG_GENERAL, SG_DEBUG, "Set prevLoc = (" << prevLoc.x() << "," << prevLoc.y() << ") and curLoc = (" << curLoc.x() << "," << curLoc.y() << ")" );
}
}
}
// TEST TEST TEST : This should do it
#if 1
if (closed)
{
SG_LOG(SG_GENERAL, SG_DEBUG, "Closed COntour : adding last node at (" << curLoc.x() << "," << curLoc.y() << ")");
@ -297,8 +341,6 @@ void LinearFeature::ConvertContour( BezContour* src, bool closed )
// need to add the markings for last segment
points.push_back( curLoc );
}
#endif
// TEST TEST TEST
// check for marking that goes all the way to the end...
if (cur_mark)
@ -402,6 +444,20 @@ Point3D LinearFeature::OffsetPointMiddle( Point3D *prev, Point3D *cur, Point3D *
// straight line blows up math - dist should be exactly as given
dist = offset_by;
}
else if ( isnan(theta) ) {
SG_LOG(SG_GENERAL, SG_DEBUG, "\nLinearFeature: (theta is NAN) " << description );
// find the direction to the next point
geo_inverse_wgs_84( cur->y(), cur->x(), next->y(), next->x(), &next_dir, &az2, &dist);
offset_dir = next_dir - 90.0;
while (offset_dir < 0.0)
{
offset_dir += 360.0;
}
// straight line blows up math - dist should be exactly as given
dist = offset_by;
}
else
{
SG_LOG(SG_GENERAL, SG_DEBUG, "\nLinearFeature: (theta NOT close to 180) " << description << ": theta is " << theta );
@ -498,7 +554,9 @@ Point3D midpoint( Point3D p0, Point3D p1 )
return Point3D( (p0.x() + p1.x()) / 2, (p0.y() + p1.y()) / 2, (p0.z() + p1.z()) / 2 );
}
int LinearFeature::Finish( bool closed )
#define DEBUG_LF (0)
int LinearFeature::Finish( bool closed, unsigned int idx )
{
TGPolygon poly;
TGPolygon normals_poly;
@ -516,6 +574,15 @@ int LinearFeature::Finish( bool closed )
double light_delta = 0;
double pt_x = 0, pt_y = 0;
#if DEBUG_LF
void* ds_id;
void* l_id;
// Create a datasource for each linear feature
char ds_name[128];
sprintf(ds_name, "./lf_debug/%04d_%s", idx, description.c_str());
ds_id = tgShapefileOpenDatasource( ds_name );
#endif
// create the inner and outer boundaries to generate polys
// this generates 2 point lists for the contours, and remembers
@ -644,6 +711,13 @@ int LinearFeature::Finish( bool closed )
exit(1);
}
#if DEBUG_LF
// Create a new layer in the datasource for each Mark
char layer_name[128];
sprintf( layer_name, "%04d_%s", i, material.c_str() );
l_id = tgShapefileOpenLayer( ds_id, layer_name );
#endif
last_end_v = 0.0f;
for (unsigned int j = marks[i]->start_idx; j <= marks[i]->end_idx; j++)
{
@ -681,6 +755,12 @@ int LinearFeature::Finish( bool closed )
poly.add_node( 0, cur_outer );
poly.add_node( 0, cur_inner );
#if DEBUG_LF
char feature_name[128];
sprintf( feature_name, "%04d", j);
tgShapefileCreateFeature( ds_id, l_id, poly, feature_name );
#endif
sp.erase();
sp.set_poly( poly );
sp.set_material( material );
@ -699,6 +779,11 @@ int LinearFeature::Finish( bool closed )
}
}
#if DEBUG_LF
// Close the datasource
tgShapefileCloseDatasource( ds_id );
#endif
// now generate the supoerpoly list for lights with constant distance between lights (depending on feature type)
for (unsigned int i=0; i<lights.size(); i++)
{
@ -831,15 +916,23 @@ int LinearFeature::Finish( bool closed )
return 1;
}
int LinearFeature::BuildBtg(float alt_m, superpoly_list* line_polys, texparams_list* line_tps, ClipPolyType* line_accum, superpoly_list* lights )
int LinearFeature::BuildBtg(float alt_m, superpoly_list* line_polys, texparams_list* line_tps, ClipPolyType* line_accum, superpoly_list* lights, bool debug )
{
TGPolygon poly;
TGPolygon clipped;
if (debug) {
sglog().setLogLevels( SG_GENERAL, SG_BULK );
} else {
sglog().setLogLevels( SG_GENERAL, SG_INFO );
}
SG_LOG(SG_GENERAL, SG_DEBUG, "\nLinearFeature::BuildBtg: " << description);
for ( unsigned int i = 0; i < marking_polys.size(); i++)
{
poly = marking_polys[i].get_poly();
poly = tgPolygonSimplify( poly );
poly = remove_tiny_contours( poly );
SG_LOG(SG_GENERAL, SG_DEBUG, "LinearFeature::BuildBtg: clipping poly " << i << " of " << marking_polys.size() );
#if 1
@ -848,9 +941,43 @@ int LinearFeature::BuildBtg(float alt_m, superpoly_list* line_polys, texparams_l
clipped = tgPolygonDiffClipper( poly, *line_accum );
#endif
// clean the poly before union with accum
clipped = reduce_degeneracy( clipped );
marking_polys[i].set_poly( clipped );
line_polys->push_back( marking_polys[i] );
#if LF_DEBUG
if ( (debug) && ( i == 78 ) ) {
void* ds_id;
void* l_id;
SG_LOG(SG_GENERAL, SG_INFO, "Problem poly: " << poly );
char ds_name[128];
sprintf(ds_name, "./lf_debug/problem");
ds_id = tgShapefileOpenDatasource( ds_name );
char layer_name[128];
sprintf( layer_name, "problem");
l_id = tgShapefileOpenLayer( ds_id, layer_name );
char feature_name[128];
sprintf( feature_name, "prob");
tgShapefileCreateFeature( ds_id, l_id, poly, feature_name );
sprintf( layer_name, "accum");
l_id = tgShapefileOpenLayer( ds_id, layer_name );
sprintf( feature_name, "accum");
tgShapefileCreateFeature( ds_id, l_id, *line_accum, feature_name );
tgShapefileCloseDatasource( ds_id );
}
#endif
SG_LOG(SG_GENERAL, SG_DEBUG, "LinearFeature::BuildBtg: union poly " << i << " of " << marking_polys.size() );
#if 1
*line_accum = tgPolygonUnion( poly, *line_accum );
#else

View file

@ -91,8 +91,8 @@ public:
contour.push_back( b );
}
int Finish( bool closed );
int BuildBtg( float alt_m, superpoly_list* line_polys, texparams_list* line_tps, ClipPolyType* line_accum, superpoly_list* lights );
int Finish( bool closed, unsigned int idx );
int BuildBtg( float alt_m, superpoly_list* line_polys, texparams_list* line_tps, ClipPolyType* line_accum, superpoly_list* lights, bool debug );
int BuildBtg( float alt_m, superpoly_list* line_polys, texparams_list* line_tps, Polygons* line_accum, superpoly_list* lights );
private:

View file

@ -25,6 +25,7 @@
#include <Polygon/index.hxx>
#include <Geometry/util.hxx>
#include <Geometry/poly_support.hxx>
#include "beznode.hxx"
#include "closedpoly.hxx"
@ -50,7 +51,7 @@ using namespace std;
static void usage( int argc, char **argv ) {
SG_LOG(SG_GENERAL, SG_ALERT,
"Usage " << argv[0] << " --input=<apt_file> "
<< "--work=<work_dir> [ --start-id=abcd ] [ --nudge=n ] "
<< "--work=<work_dir> [ --start-id=abcd ] [ --restart-id=abcd ] [ --nudge=n ] "
<< "[--min-lon=<deg>] [--max-lon=<deg>] [--min-lat=<deg>] [--max-lat=<deg>] "
<< "[--clear-dem-path] [--dem-path=<path>] [--max-slope=<decimal>] "
<< "[ --airport=abcd ] [--tile=<tile>] [--chunk=<chunk>] [--verbose] [--help]");
@ -90,6 +91,8 @@ static void help( int argc, char **argv, const string_list& elev_src ) {
cout << "a valid airport code eg. --airport-id=KORD, or a starting airport can be specified using --start-id=abcd \n";
cout << "where once again abcd is a valid airport code. In this case, all airports in the file subsequent to the \n";
cout << "start-id are done. This is convienient when re-starting after a previous error. \n";
cout << "If you want to restart with the airport after a problam icao, use --restart-id=abcd, as this works the same as\n";
cout << " with the exception that the airport abcd is skipped \n";
cout << "\nAn input area may be specified by lat and lon extent using min and max lat and lon. \n";
cout << "Alternatively, you may specify a chunk (10 x 10 degrees) or tile (1 x 1 degree) using a string \n";
cout << "such as eg. w080n40, e000s27. \n";
@ -134,7 +137,9 @@ int main(int argc, char **argv)
string work_dir = "";
string input_file = "";
string start_id = "";
string restart_id = "";
string airport_id = "";
string last_apt_file = "./last_apt.txt";
int arg_pos;
for (arg_pos = 1; arg_pos < argc; arg_pos++)
@ -156,10 +161,18 @@ int main(int argc, char **argv)
{
start_id = arg.substr(11);
}
else if ( arg.find("--restart-id=") == 0 )
{
restart_id = arg.substr(13);
}
else if ( arg.find("--nudge=") == 0 )
{
nudge = atoi( arg.substr(8).c_str() );
}
else if ( arg.find("--last_apt_file=") == 0 )
{
last_apt_file = arg.substr(16);
}
else if ( arg.find("--min-lon=") == 0 )
{
min_lon = atof( arg.substr(10).c_str() );
@ -270,6 +283,9 @@ int main(int argc, char **argv)
string counter_file = airportareadir+"/poly_counter";
poly_index_init( counter_file );
// Initialize shapefile support (for debugging)
tgShapefileInit();
sg_gzifstream in( input_file );
if ( !in.is_open() )
{
@ -289,7 +305,7 @@ int main(int argc, char **argv)
SG_LOG(SG_GENERAL, SG_INFO, "Finished Adding airport - now parse");
// and start the parser
parser->Parse();
parser->Parse( last_apt_file );
}
else if ( start_id != "" )
{
@ -302,7 +318,23 @@ int main(int argc, char **argv)
parser->AddAirports( position, min_lat, min_lon, max_lat, max_lon );
// parse all the airports that were found
parser->Parse();
parser->Parse( last_apt_file );
}
else if ( restart_id != "" )
{
SG_LOG(SG_GENERAL, SG_INFO, "move forward airport after " << restart_id );
// scroll forward in datafile
position = parser->FindAirport( restart_id );
// add all remaining airports within boundary
parser->AddAirports( position, min_lat, min_lon, max_lat, max_lon );
// but remove the restart id - it's broken
parser->RemoveAirport( restart_id );
// parse all the airports that were found
parser->Parse( last_apt_file );
}
else
{
@ -310,12 +342,13 @@ int main(int argc, char **argv)
parser->AddAirports( 0, min_lat, min_lon, max_lat, max_lon );
// and parser them
parser->Parse();
parser->Parse( last_apt_file );
}
delete parser;
SG_LOG(SG_GENERAL, SG_INFO, "Done");
exit(0);
return 0;
}

View file

@ -306,7 +306,18 @@ void Parser::AddAirports( long start_pos, float min_lat, float min_lon, float ma
}
}
void Parser::Parse()
void Parser::RemoveAirport( string icao )
{
for (unsigned int i = 0; i < airport_icaos.size(); i++ ) {
if (airport_icaos[i] == icao) {
parse_positions.erase(parse_positions.begin()+i);
airport_icaos.erase(airport_icaos.begin()+i);
break;
}
}
}
void Parser::Parse( string last_apt_file )
{
char tmp[2048];
struct timeval parse_start;
@ -335,6 +346,11 @@ void Parser::Parse()
SG_LOG( SG_GENERAL, SG_ALERT, "Parsing airport " << airport_icaos[i] << " at " << parse_positions[i] << " start time " << ctime(&parse_start.tv_sec) );
in.seekg(parse_positions[i], ios::beg);
// save the airport we are working on
char command[256];
sprintf( command, "echo %s > %s", airport_icaos[i].c_str(), last_apt_file.c_str() );
system( command );
while ( !in.eof() && (cur_state != STATE_DONE ) )
{
in.getline(tmp, 2048);
@ -760,7 +776,7 @@ int Parser::ParseLine(char* line)
}
if (cur_airport)
{
cur_feat->Finish( true );
cur_feat->Finish( true, cur_airport->NumFeatures() );
cur_airport->AddFeature( cur_feat );
}
cur_feat = NULL;
@ -794,7 +810,7 @@ int Parser::ParseLine(char* line)
}
if (cur_airport)
{
cur_feat->Finish( false );
cur_feat->Finish( false, cur_airport->NumFeatures() );
cur_airport->AddFeature( cur_feat );
}
}

View file

@ -88,7 +88,8 @@ public:
long FindAirport( string icao );
void AddAirport( string icao );
void AddAirports( long start_pos, float min_lat, float min_lon, float max_lat, float max_lon );
void Parse( void );
void RemoveAirport( string icao );
void Parse( string last_apt_file );
private:
bool IsAirportDefinition( char* line, string icao );

View file

@ -163,6 +163,8 @@ int Runway::BuildBtg( float alt_m, superpoly_list* rwy_polys, texparams_list* te
if (apt_base)
{
// If we have shoulders, we need to grow even further...
// generate area around runways
base = gen_runway_area_w_extend( 0.0, 20.0, -rwy.overrun[0], -rwy.overrun[1], 20.0);

View file

@ -50,7 +50,9 @@ public:
void BuildShoulder( float alt_m,
superpoly_list *rwy_polys,
texparams_list *texparams,
ClipPolyType *accum );
ClipPolyType *accum,
TGPolygon* apt_base,
TGPolygon* apt_clearing );
private:
struct TGRunway {
@ -152,6 +154,7 @@ private:
superpoly_list gen_ssalx( float alt_m, const string& kind, bool recip );
superpoly_list gen_malsx( float alt_m, const string& kind, bool recip );
};
typedef std::vector <Runway *> RunwayList;

View file

@ -23,6 +23,8 @@
#include <simgear/constants.h>
#include <simgear/debug/logstream.hxx>
#include <simgear/math/sg_geodesy.hxx>
#include "beznode.hxx"
#include "runway.hxx"
#include <stdlib.h>
@ -153,8 +155,7 @@ void Runway::gen_rwy( double alt_m,
TGPolygon runway = gen_runway_w_mid( alt_m, 0, 0 );
TGPolygon runway_half;
for ( int rwhalf=0; rwhalf<2; ++rwhalf ){
for ( int rwhalf=0; rwhalf<2; ++rwhalf ){
if (rwhalf == 0) {
@ -195,7 +196,6 @@ for ( int rwhalf=0; rwhalf<2; ++rwhalf ){
double heading = 0.0;
string rwname;
//
// Displaced threshold if it exists
//
@ -292,8 +292,6 @@ for ( int rwhalf=0; rwhalf<2; ++rwhalf ){
gen_rw_designation( material, runway_half, heading,
rwname, start1_pct, end1_pct, rwy_polys, texparams, accum );
if (rwy.marking[rwhalf] > 1){
// Generate remaining markings depending on type of runway
gen_rw_marking( runway_half,
@ -329,33 +327,34 @@ for ( int rwhalf=0; rwhalf<2; ++rwhalf ){
gen_runway_overrun( runway_half, rwhalf,
material,
rwy_polys, texparams, accum );
}
}
}
void Runway::BuildShoulder( float alt_m,
superpoly_list *rwy_polys,
texparams_list *texparams,
ClipPolyType *accum )
ClipPolyType *accum,
TGPolygon* apt_base,
TGPolygon* apt_clearing )
{
TGPolygon base, safe_base;
string shoulder_surface = "";
double shoulder_width = 0.0f;
if (rwy.shoulder > 0){ // Add a shoulder to the runway
//shoulder_width = rwy.width * 0.15;
//if (shoulder_width > 11.0){
// shoulder_width = 11.0;
//}
if (rwy.shoulder > 0){
// Add a shoulder to the runway
shoulder_width = 11.0f;
if (rwy.shoulder == 1){
shoulder_surface = "pa_shoulder";
} else if (rwy.shoulder == 2){
shoulder_surface = "pc_shoulder";
} else SG_LOG(SG_GENERAL, SG_ALERT, "Unknown shoulder surface code = " << rwy.shoulder );
} else if (rwy.shoulder == 0){ // We add a fake shoulder if the runway has an asphalt or concrete surface
} else {
SG_LOG(SG_GENERAL, SG_ALERT, "Unknown shoulder surface code = " << rwy.shoulder );
}
} else {
// We add a fake shoulder if the runway has an asphalt or concrete surface
shoulder_width = 1.0;
if (rwy.surface == 1){
shoulder_surface = "pa_shoulder_f";
@ -363,48 +362,104 @@ void Runway::BuildShoulder( float alt_m,
shoulder_surface = "pc_shoulder_f";
}
}
SG_LOG(SG_GENERAL, SG_DEBUG, "Shoulder width = " << shoulder_width );
SG_LOG(SG_GENERAL, SG_DEBUG, "Shoulder surface is: " << shoulder_surface );
if (shoulder_width > 0.0f) {
// we need to break these shoulders up into managable pieces, as we're asked to triangulate
// 3-4 km long by 1m wide strips - jrs-can't handle it.
double max_dist = (double)shoulder_width * 25.0f;
int numSegs = (rwy.length / max_dist) + 1;
double dist = rwy.length / (double)numSegs;
// Create both shoulder sides
for (int i=0; i<2; ++i){
double step;
double lat = 0.0f, lon = 0.0f, r;
/* nudge the shoulders so the really long lines overlap the runway a bit */
/* If the are 'equal' there's a good chance roundoff error can create a */
/* REALY thin long polygon, which causes a segfault */
if (i == 0){
step= (rwy.width + shoulder_width) * 0.5;
} else if (i == 1) {
step= -(rwy.width + shoulder_width) * 0.5;
if (max_dist > 100.0f) {
max_dist = 100.0f;
}
double left_hdg = rwy.heading - 90.0;
if ( left_hdg < 0 ) { left_hdg += 360.0; }
geo_direct_wgs_84 ( alt_m, rwy.lat[0], rwy.lon[0], left_hdg, step, &lat, &lon, &r );
Point3D ref = Point3D( lon, lat, 0.0f );
for (int j=0; j<numSegs; j++)
{
geo_direct_wgs_84 ( alt_m, ref.y(), ref.x(), rwy.heading, (j*dist), &lat, &lon, &r );
TGPolygon shoulderSegment = gen_wgs84_rect( lat, lon, rwy.heading, dist, shoulder_width+1.0f );
int numSegs = (rwy.length / max_dist) + 1;
double dist = rwy.length / numSegs;
TGPolygon poly;
TGSuperPoly sp;
TGTexParams tp;
double lat = 0.0f;
double lon = 0.0f;
double r = 0.0f;
Point3D inner_start, inner_end;
Point3D outer_start, outer_end;
Point3D curInnerLoc, nextInnerLoc;
Point3D curOuterLoc, nextOuterLoc;
// Create two paralell lines from start position to end position, and interpolate in between
// many airports line the taxiway directly to the corner of the runway. This can create problems,
// so extend the shoulders 0.5 meters past each end of the runway
for (int i=0; i<2; i++) {
double rev_hdg = rwy.heading - 180.0;
if ( rev_hdg < 0 ) { rev_hdg += 360.0; }
if (i == 0) {
// left side
double left_hdg = rwy.heading - 90.0;
if ( left_hdg < 0 ) { left_hdg += 360.0; }
geo_direct_wgs_84 ( 0, rwy.lat[0], rwy.lon[0], left_hdg, rwy.width*.5, &lat, &lon, &r );
geo_direct_wgs_84 ( 0, lat, lon, rev_hdg, 0.5f, &lat, &lon, &r );
inner_start = Point3D( lon, lat, 0.0f );
geo_direct_wgs_84 ( 0, lat, lon, rwy.heading, rwy.length+1.0f, &lat, &lon, &r );
inner_end = Point3D( lon, lat, 0.0f );
geo_direct_wgs_84 ( 0, rwy.lat[0], rwy.lon[0], left_hdg, rwy.width*.5 + shoulder_width, &lat, &lon, &r );
geo_direct_wgs_84 ( 0, lat, lon, rev_hdg, 0.5f, &lat, &lon, &r );
outer_start = Point3D( lon, lat, 0.0f );
geo_direct_wgs_84 ( 0, lat, lon, rwy.heading, rwy.length+1.0f, &lat, &lon, &r );
outer_end = Point3D( lon, lat, 0.0f );
} else {
// right side
double right_hdg = rwy.heading + 90.0;
if ( right_hdg > 360 ) { right_hdg -= 360.0; }
geo_direct_wgs_84 ( 0, rwy.lat[0], rwy.lon[0], right_hdg, rwy.width*.5, &lat, &lon, &r );
geo_direct_wgs_84 ( 0, lat, lon, rev_hdg, 0.5f, &lat, &lon, &r );
inner_start = Point3D( lon, lat, 0.0f );
geo_direct_wgs_84 ( 0, lat, lon, rwy.heading, rwy.length+1.0f, &lat, &lon, &r );
inner_end = Point3D( lon, lat, 0.0f );
geo_direct_wgs_84 ( 0, rwy.lat[0], rwy.lon[0], right_hdg, rwy.width*.5 + shoulder_width, &lat, &lon, &r );
geo_direct_wgs_84 ( 0, lat, lon, rev_hdg, 0.5f, &lat, &lon, &r );
outer_start = Point3D( lon, lat, 0.0f );
geo_direct_wgs_84 ( 0, lat, lon, rwy.heading, rwy.length+1.0f, &lat, &lon, &r );
outer_end = Point3D( lon, lat, 0.0f );
}
curInnerLoc = inner_start;
curOuterLoc = outer_start;
for (int p=0; p<numSegs; p++)
{
// calculate next locations
nextInnerLoc = CalculateLinearLocation( inner_start, inner_end, (1.0f/numSegs) * (p+1) );
nextOuterLoc = CalculateLinearLocation( outer_start, outer_end, (1.0f/numSegs) * (p+1) );
// Generate a poly
poly.erase();
if (i == 0 ) {
poly.add_node( 0, curInnerLoc );
poly.add_node( 0, nextInnerLoc );
poly.add_node( 0, nextOuterLoc );
poly.add_node( 0, curOuterLoc );
} else {
poly.add_node( 0, curOuterLoc );
poly.add_node( 0, nextOuterLoc );
poly.add_node( 0, nextInnerLoc );
poly.add_node( 0, curInnerLoc );
}
#if 1
TGPolygon clipped = tgPolygonDiff( shoulderSegment, *accum );
TGPolygon clipped = tgPolygonDiff( poly, *accum );
#else
TGPolygon clipped = tgPolygonDiffClipper( shoulderSegment, *accum );
TGPolygon clipped = tgPolygonDiffClipper( poly, *accum );
#endif
sp.erase();
sp.set_poly( clipped );
@ -412,17 +467,33 @@ void Runway::BuildShoulder( float alt_m,
rwy_polys->push_back( sp );
#if 1
*accum = tgPolygonUnion( shoulderSegment, *accum );
*accum = tgPolygonUnion( poly, *accum );
#else
*accum = tgPolygonUnionClipper( shoulderSegment, *accum );
*accum = tgPolygonUnionClipper( poly, *accum );
#endif
tp = TGTexParams( shoulderSegment.get_pt(0,0), -shoulder_width, dist, rwy.heading );
if (i == 0){
tp = TGTexParams( poly.get_pt(0,0), shoulder_width, dist, rwy.heading );
tp.set_maxv(dist);
// reverse u direction for right side
if ( i == 1 ) {
tp.set_maxu(0);
tp.set_minu(1);
}
texparams->push_back( tp );
// Add to base / safe base
base = tgPolygonExpand( poly, 20.0f );
safe_base = tgPolygonExpand( poly, 50.0f );
// add this to the airport clearing
*apt_clearing = tgPolygonUnion(safe_base, *apt_clearing);
// and add the clearing to the base
*apt_base = tgPolygonUnion( base, *apt_base );
// now set cur locations for the next iteration
curInnerLoc = nextInnerLoc;
curOuterLoc = nextOuterLoc;
}
}
}

View file

@ -14,6 +14,7 @@ target_link_libraries(fgfs-construct
Polygon Geometry
Array Optimize Output landcover poly2tri
TriangleJRS
${GDAL_LIBRARY}
${SIMGEAR_CORE_LIBRARIES}
${SIMGEAR_CORE_LIBRARY_DEPENDENCIES}
${GPC_LIBRARY})

View file

@ -1,4 +1,6 @@
if (GDAL_FOUND)
include_directories(${GDAL_INCLUDE_DIR})
endif (GDAL_FOUND)
add_library(Geometry STATIC
contour_tree.cxx

View file

@ -1452,6 +1452,7 @@ TGPolygon remove_bad_contours( const TGPolygon &poly ) {
// remove any too small contours
TGPolygon remove_tiny_contours( const TGPolygon &poly ) {
double min_area = SG_EPSILON*SG_EPSILON;
TGPolygon result;
result.erase();
@ -1460,8 +1461,34 @@ TGPolygon remove_tiny_contours( const TGPolygon &poly ) {
double area=poly.area_contour(i);
if (-SG_EPSILON*SG_EPSILON<area && area<SG_EPSILON*SG_EPSILON) {
// cout << "tossing a bad contour " << i << " (too small)" << endl;
if (-min_area<area && area<min_area) {
SG_LOG(SG_GENERAL, SG_INFO, "remove_tiny_contours " << i << " area is " << area << ": removing");
continue;
} else {
SG_LOG(SG_GENERAL, SG_DEBUG, "remove_tiny_contours NO - " << i << " area is " << area << " requirement is " << min_area);
}
/* keeping the contour */
int flag = poly.get_hole_flag( i );
result.add_contour( contour, flag );
}
return result;
}
// remove any too small contours
TGPolygon remove_small_contours( const TGPolygon &poly ) {
double min_area = 200*SG_EPSILON*SG_EPSILON;
TGPolygon result;
result.erase();
for ( int i = 0; i < poly.contours(); ++i ) {
point_list contour = poly.get_contour( i );
double area=poly.area_contour(i);
if (-min_area<area && area<min_area) {
SG_LOG(SG_GENERAL, SG_INFO, "remove_small_contours " << i << " area is " << area << ": removing");
continue;
}
@ -1473,4 +1500,113 @@ TGPolygon remove_tiny_contours( const TGPolygon &poly ) {
return result;
}
// seperate
#include <ogrsf_frmts.h>
const char* format_name="ESRI Shapefile";
void tgShapefileInit( void )
{
OGRRegisterAll();
}
void* tgShapefileOpenDatasource( const char* datasource_name )
{
OGRDataSource *datasource;
OGRSFDriver *ogrdriver;
ogrdriver = OGRSFDriverRegistrar::GetRegistrar()->GetDriverByName(format_name);
if (!ogrdriver) {
SG_LOG(SG_GENERAL, SG_ALERT, "Unknown datasource format driver: " << format_name);
exit(1);
}
datasource = ogrdriver->CreateDataSource(datasource_name, NULL);
if (!datasource) {
SG_LOG(SG_GENERAL, SG_ALERT, "Unable to create datasource: " << datasource_name);
exit(1);
}
return (void*)datasource;
}
void* tgShapefileOpenLayer( void* ds_id, const char* layer_name ) {
OGRDataSource* datasource = (OGRDataSource *)ds_id;
OGRLayer* layer;
OGRSpatialReference srs;
srs.SetWellKnownGeogCS("WGS84");
layer = datasource->CreateLayer( layer_name, &srs, wkbPolygon25D, NULL);
if (!layer) {
SG_LOG(SG_GENERAL, SG_ALERT, "Creation of layer '" << layer_name << "' failed");
return NULL;
}
OGRFieldDefn descriptionField("ID", OFTString);
descriptionField.SetWidth(128);
if( layer->CreateField( &descriptionField ) != OGRERR_NONE ) {
SG_LOG(SG_GENERAL, SG_ALERT, "Creation of field 'Description' failed");
}
return (void*)layer;
}
void tgShapefileCreateFeature( void* ds_id, void* l_id, const TGPolygon &poly, const char* description )
{
OGRDataSource* datasource = (OGRDataSource *)ds_id;
OGRLayer* layer = (OGRLayer*)l_id;
OGRPolygon* polygon = new OGRPolygon();
for ( int i = 0; i < poly.contours(); i++ ) {
bool skip_ring=false;
point_list contour = poly.get_contour( i );
if (contour.size()<3) {
SG_LOG(SG_GENERAL, SG_DEBUG, "Polygon with less than 3 points");
skip_ring=true;
}
// FIXME: Current we ignore the hole-flag and instead assume
// that the first ring is not a hole and the rest
// are holes
OGRLinearRing *ring=new OGRLinearRing();
for (unsigned int pt = 0; pt < contour.size(); pt++) {
OGRPoint *point=new OGRPoint();
point->setX( contour[pt].x() );
point->setY( contour[pt].y() );
point->setZ( 0.0 );
ring->addPoint(point);
}
ring->closeRings();
if (!skip_ring) {
polygon->addRingDirectly(ring);
}
OGRFeature* feature = NULL;
feature = new OGRFeature( layer->GetLayerDefn() );
feature->SetField("ID", description);
feature->SetGeometry(polygon);
if( layer->CreateFeature( feature ) != OGRERR_NONE )
{
SG_LOG(SG_GENERAL, SG_ALERT, "Failed to create feature in shapefile");
}
OGRFeature::DestroyFeature(feature);
}
}
void tgShapefileCloseLayer( void* l_id )
{
OGRLayer* layer = (OGRLayer *)l_id;
//OGRLayer::DestroyLayer( layer );
}
void tgShapefileCloseDatasource( void* ds_id )
{
OGRDataSource* datasource = (OGRDataSource *)ds_id;
OGRDataSource::DestroyDataSource( datasource );
}

View file

@ -107,6 +107,16 @@ TGPolygon remove_bad_contours( const TGPolygon &poly );
// remove any too small contours
TGPolygon remove_tiny_contours( const TGPolygon &poly );
TGPolygon remove_small_contours( const TGPolygon &poly );
// Write Polygons to Shapefile Support
void tgShapefileInit( void );
void* tgShapefileOpenDatasource( const char* datasource_name );
void* tgShapefileOpenLayer( void* ds_id, const char* layer_name );
void tgShapefileCreateFeature( void* ds_id, void* l_id, const TGPolygon &poly, const char* feature_name );
void tgShapefileCloseLayer( void* l_id );
void tgShapefileCloseDatasource( void* ds_id );
#endif // _POLY_SUPPORT_HXX