1
0
Fork 0

some fixes for interpolating added points to the airport holes. fixes most gaps when allowing roads on airport landclass

add zero area triangle check when triangulating without extra points.
This commit is contained in:
Peter Sadrozinski 2015-02-08 06:52:39 -05:00
parent 314671886e
commit c5faeae357
10 changed files with 255 additions and 24 deletions

View file

@ -34,7 +34,7 @@
void TGConstruct::FixTJunctions( void ) {
int before, after;
std::vector<SGGeod> points;
std::vector<TGNode*> points;
tgRectangle bb;
// traverse each poly, and add intermediate nodes
@ -42,7 +42,7 @@ void TGConstruct::FixTJunctions( void ) {
for( unsigned int j = 0; j < polys_clipped.area_size(i); ++j ) {
tgPolygon current = polys_clipped.get_poly(i, j);
bb = current.GetBoundingBox();
nodes.get_geod_inside( bb.getMin(), bb.getMax(), points );
nodes.get_nodes_inside( bb.getMin(), bb.getMax(), points );
before = current.TotalNodes();
current = tgPolygon::AddColinearNodes( current, points );

View file

@ -6,6 +6,7 @@
#include "tg_accumulator.hxx"
#include "tg_contour.hxx"
#include "tg_polygon.hxx"
#include "tg_unique_tgnode.hxx"
tgContour tgContour::Snap( const tgContour& subject, double snap )
{
@ -621,6 +622,87 @@ static bool FindIntermediateNode( const SGGeod& start, const SGGeod& end,
return found_node;
}
static bool FindIntermediateNode( const SGGeod& start, const SGGeod& end,
const std::vector<TGNode*>& nodes, TGNode*& result,
double bbEpsilon, double errEpsilon )
{
bool found_node = false;
double m, m1, b, b1, y_err, x_err, y_err_min, x_err_min;
SGGeod p0 = start;
SGGeod p1 = end;
double xdist = fabs(p0.getLongitudeDeg() - p1.getLongitudeDeg());
double ydist = fabs(p0.getLatitudeDeg() - p1.getLatitudeDeg());
x_err_min = xdist + 1.0;
y_err_min = ydist + 1.0;
if ( xdist > ydist ) {
// sort these in a sensible order
SGGeod p_min, p_max;
if ( p0.getLongitudeDeg() < p1.getLongitudeDeg() ) {
p_min = p0;
p_max = p1;
} else {
p_min = p1;
p_max = p0;
}
m = (p_min.getLatitudeDeg() - p_max.getLatitudeDeg()) / (p_min.getLongitudeDeg() - p_max.getLongitudeDeg());
b = p_max.getLatitudeDeg() - m * p_max.getLongitudeDeg();
for ( int i = 0; i < (int)nodes.size(); ++i ) {
// cout << i << endl;
SGGeod current = nodes[i]->GetPosition();
if ( (current.getLongitudeDeg() > (p_min.getLongitudeDeg() + (bbEpsilon))) && (current.getLongitudeDeg() < (p_max.getLongitudeDeg() - (bbEpsilon))) ) {
y_err = fabs(current.getLatitudeDeg() - (m * current.getLongitudeDeg() + b));
if ( y_err < errEpsilon ) {
found_node = true;
if ( y_err < y_err_min ) {
result = nodes[i];
y_err_min = y_err;
}
}
}
}
} else {
// sort these in a sensible order
SGGeod p_min, p_max;
if ( p0.getLatitudeDeg() < p1.getLatitudeDeg() ) {
p_min = p0;
p_max = p1;
} else {
p_min = p1;
p_max = p0;
}
m1 = (p_min.getLongitudeDeg() - p_max.getLongitudeDeg()) / (p_min.getLatitudeDeg() - p_max.getLatitudeDeg());
b1 = p_max.getLongitudeDeg() - m1 * p_max.getLatitudeDeg();
for ( int i = 0; i < (int)nodes.size(); ++i ) {
SGGeod current = nodes[i]->GetPosition();
if ( (current.getLatitudeDeg() > (p_min.getLatitudeDeg() + (bbEpsilon))) && (current.getLatitudeDeg() < (p_max.getLatitudeDeg() - (bbEpsilon))) ) {
x_err = fabs(current.getLongitudeDeg() - (m1 * current.getLatitudeDeg() + b1));
if ( x_err < errEpsilon ) {
found_node = true;
if ( x_err < x_err_min ) {
result = nodes[i];
x_err_min = x_err;
}
}
}
}
}
return found_node;
}
static void AddIntermediateNodes( const SGGeod& p0, const SGGeod& p1, std::vector<SGGeod>& nodes, tgContour& result, double bbEpsilon, double errEpsilon )
{
SGGeod new_pt;
@ -639,6 +721,37 @@ static void AddIntermediateNodes( const SGGeod& p0, const SGGeod& p1, std::vecto
}
}
extern SGGeod InterpolateElevation( const SGGeod& dst_node, const SGGeod& start, const SGGeod& end );
static void AddIntermediateNodes( const SGGeod& p0, const SGGeod& p1, bool preserve3d, std::vector<TGNode*>& nodes, tgContour& result, double bbEpsilon, double errEpsilon )
{
TGNode* new_pt;
SGGeod new_geode;
SG_LOG(SG_GENERAL, SG_BULK, " " << p0 << " <==> " << p1 );
bool found_extra = FindIntermediateNode( p0, p1, nodes, new_pt, bbEpsilon, errEpsilon );
if ( found_extra ) {
if ( preserve3d ) {
// interpolate the new nodes elevation based on p0, p1
new_geode = InterpolateElevation( new_pt->GetPosition(), p0, p1 );
SG_LOG(SG_GENERAL, SG_ALERT, "INTERPOLATE ELVATION between " << p0 << " and " << p1 << " returned elvation " << new_geode.getElevationM() );
new_pt->SetElevation( new_geode.getElevationM() );
new_pt->SetFixedPosition( true );
}
AddIntermediateNodes( p0, new_pt->GetPosition(), preserve3d, nodes, result, bbEpsilon, errEpsilon );
result.AddNode( new_pt->GetPosition() );
SG_LOG(SG_GENERAL, SG_BULK, " adding = " << new_pt->GetPosition() );
AddIntermediateNodes( new_pt->GetPosition(), p1, preserve3d, nodes, result, bbEpsilon, errEpsilon );
}
}
tgContour tgContour::AddColinearNodes( const tgContour& subject, UniqueSGGeodSet& nodes )
{
SGGeod p0, p1;
@ -702,6 +815,37 @@ tgContour tgContour::AddColinearNodes( const tgContour& subject, std::vector<SGG
return result;
}
tgContour tgContour::AddColinearNodes( const tgContour& subject, bool preserve3d, std::vector<TGNode*>& nodes )
{
SGGeod p0, p1;
tgContour result;
for ( unsigned int n = 0; n < subject.GetSize()-1; n++ ) {
p0 = subject.GetNode( n );
p1 = subject.GetNode( n+1 );
// add start of segment
result.AddNode( p0 );
// add intermediate points
AddIntermediateNodes( p0, p1, preserve3d, nodes, result, SG_EPSILON*10, SG_EPSILON*4 );
}
p0 = subject.GetNode( subject.GetSize() - 1 );
p1 = subject.GetNode( 0 );
// add start of segment
result.AddNode( p0 );
// add intermediate points
AddIntermediateNodes( p0, p1, preserve3d, nodes, result, SG_EPSILON*10, SG_EPSILON*4 );
// maintain original hole flag setting
result.SetHole( subject.GetHole() );
return result;
}
// this is the opposite of FindColinearNodes - it takes a single SGGeode,
// and tries to find the line segment the point is colinear with
bool tgContour::FindColinearLine( const tgContour& subject, const SGGeod& node, SGGeod& start, SGGeod& end )

View file

@ -14,6 +14,8 @@
#include "clipper.hpp"
/* forward declarations */
class TGNode;
class tgPolygon;
typedef std::vector <tgPolygon> tgpolygon_list;
@ -112,6 +114,7 @@ public:
static bool IsInside( const tgContour& inside, const tgContour& outside );
static tgContour AddColinearNodes( const tgContour& subject, UniqueSGGeodSet& nodes );
static tgContour AddColinearNodes( const tgContour& subject, std::vector<SGGeod>& nodes );
static tgContour AddColinearNodes( const tgContour& subject, bool preserve3d, std::vector<TGNode*>& nodes );
static bool FindColinearLine( const tgContour& subject, const SGGeod& node, SGGeod& start, SGGeod& end );
// conversions

View file

@ -45,6 +45,19 @@ bool TGNodes::get_geod_inside( const SGGeod& min, const SGGeod& max, std::vector
return true;
}
bool TGNodes::get_nodes_inside( const SGGeod& min, const SGGeod& max, std::vector<TGNode*>& points ) const {
points.clear();
for ( unsigned int i = 0; i < tg_node_list.size(); i++ ) {
SGGeod const& pt = tg_node_list[i].GetPosition();
if ( IsAlmostWithin( pt, min, max ) ) {
points.push_back( &tg_node_list[i] );
}
}
return true;
}
bool TGNodes::get_geod_edge( const SGBucket& b, std::vector<SGGeod>& north, std::vector<SGGeod>& south, std::vector<SGGeod>& east, std::vector<SGGeod>& west ) const {
double north_compare = b.get_center_lat() + 0.5 * b.get_height();
double south_compare = b.get_center_lat() - 0.5 * b.get_height();
@ -100,7 +113,7 @@ void TGNodes::init_spacial_query( void )
// generate the tuple
Point pt( tg_node_list[i].GetPosition().getLongitudeDeg(), tg_node_list[i].GetPosition().getLatitudeDeg() );
double e( tg_node_list[i].GetPosition().getElevationM() );
Point_and_Elevation pande(pt, e);
Point_and_Elevation pande( pt, e, &tg_node_list[i] );
// and insert into tree
tg_kd_tree.insert( pande );
@ -142,6 +155,36 @@ bool TGNodes::get_geod_inside( const SGGeod& min, const SGGeod& max, std::vector
return true;
}
bool TGNodes::get_nodes_inside( const SGGeod& min, const SGGeod& max, std::vector<TGNode*>& points ) const {
points.clear();
// Have we generated the k-d tree?
if ( !kd_tree_valid ) {
SG_LOG(SG_GENERAL, SG_ALERT, "get_geod_inside called with invalid kdtree" );
exit(0);
return false;
}
// define an exact rectangulat range query (fuzziness=0)
Point ll( min.getLongitudeDeg() - fgPoint3_Epsilon, min.getLatitudeDeg() - fgPoint3_Epsilon );
Point ur( max.getLongitudeDeg() + fgPoint3_Epsilon, max.getLatitudeDeg() + fgPoint3_Epsilon );
Fuzzy_bb exact_bb(ll, ur);
// list of tuples as a result
std::list<Point_and_Elevation> result;
std::list<Point_and_Elevation>::iterator it;
// perform the query
tg_kd_tree.search(std::back_inserter( result ), exact_bb);
// and convert the tuples back into SGGeod
for ( it = result.begin(); it != result.end(); it++ ) {
points.push_back( boost::get<2>(*it) );
}
return true;
}
// This query finds all nodes along the tile borders (north, south, east and west)
bool TGNodes::get_geod_edge( const SGBucket& b, std::vector<SGGeod>& north, std::vector<SGGeod>& south, std::vector<SGGeod>& east, std::vector<SGGeod>& west ) const {
double north_compare = b.get_center_lat() + 0.5 * b.get_height();

View file

@ -15,9 +15,19 @@
#include <CGAL/Search_traits_adapter.h> /* Just use two dimensional lookup - elevation is a trait */
#include <boost/iterator/zip_iterator.hpp>
#include <simgear/compiler.h>
#include <simgear/bucket/newbucket.hxx>
#include <simgear/io/lowlevel.hxx>
#include "tg_unique_tgnode.hxx"
#define FG_PROXIMITY_EPSILON 0.000001
#define FG_COURSE_EPSILON 0.0001
typedef CGAL::Simple_cartesian<double> Kernel;
typedef Kernel::Point_2 Point;
typedef boost::tuple<Point, double> Point_and_Elevation;
typedef boost::tuple<Point, double, TGNode*> Point_and_Elevation;
//definition of the property map
#ifdef CGAL_OLD
@ -42,16 +52,6 @@ typedef CGAL::Fuzzy_iso_box<Traits> Fuzzy_bb;
typedef CGAL::Kd_tree<Traits> Tree;
#include <simgear/compiler.h>
#include <simgear/bucket/newbucket.hxx>
#include <simgear/io/lowlevel.hxx>
#include "tg_unique_tgnode.hxx"
#define FG_PROXIMITY_EPSILON 0.000001
#define FG_COURSE_EPSILON 0.0001
/* This class handles ALL of the nodes in a tile : 3d nodes in elevation data, 2d nodes generated from landclass, etc) */
class TGNodes {
public:
@ -123,6 +123,8 @@ public:
// Find all the nodes within a bounding box
bool get_geod_inside( const SGGeod& min, const SGGeod& max, std::vector<SGGeod>& points ) const;
bool get_nodes_inside( const SGGeod& min, const SGGeod& max, std::vector<TGNode*>& points ) const;
// Find a;; the nodes on the tile edges
bool get_geod_edge( const SGBucket& b, std::vector<SGGeod>& north, std::vector<SGGeod>& south, std::vector<SGGeod>& east, std::vector<SGGeod>& west ) const;

View file

@ -121,6 +121,7 @@ tgPolygon tgPolygon::AddColinearNodes( const tgPolygon& subject, std::vector<SGG
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for ( unsigned int c = 0; c < subject.Contours(); c++ ) {
result.AddContour( tgContour::AddColinearNodes( subject.GetContour(c), nodes ) );
@ -134,6 +135,22 @@ tgPolygon tgPolygon::AddColinearNodes( const tgPolygon& subject, UniqueSGGeodSet
return AddColinearNodes( subject, nodes.get_list() );
}
tgPolygon tgPolygon::AddColinearNodes( const tgPolygon& subject, std::vector<TGNode*>& nodes )
{
tgPolygon result;
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for ( unsigned int c = 0; c < subject.Contours(); c++ ) {
result.AddContour( tgContour::AddColinearNodes( subject.GetContour(c), subject.GetPreserve3D(), nodes ) );
}
return result;
}
// this is the opposite of FindColinearNodes - it takes a single SGGeode,
// and tries to find the line segment the point is colinear with
bool tgPolygon::FindColinearLine( const tgPolygon& subject, SGGeod& node, SGGeod& start, SGGeod& end )

View file

@ -86,6 +86,7 @@ void clipper_to_shapefile( ClipperLib::Paths polys, char* datasource );
// Forward Declaration:
class tgPolygon;
class tgChoppedPolygons;
class TGNode;
typedef std::vector <tgPolygon> tgpolygon_list;
typedef tgpolygon_list::iterator tgpolygon_list_iterator;
@ -413,6 +414,8 @@ public:
static tgPolygon AddColinearNodes( const tgPolygon& subject, std::vector<SGGeod>& nodes );
static bool FindColinearLine( const tgPolygon& subject, SGGeod& node, SGGeod& start, SGGeod& end );
static tgPolygon AddColinearNodes( const tgPolygon& subject, std::vector<TGNode*>& nodes );
// IO
void SaveToGzFile( gzFile& fp ) const;
void LoadFromGzFile( gzFile& fp );

View file

@ -9,6 +9,7 @@ tgPolygon tgPolygon::Snap( const tgPolygon& subject, double snap )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for (unsigned int c = 0; c < subject.Contours(); c++) {
result.AddContour( tgContour::Snap( subject.GetContour( c ), snap ) );
@ -24,6 +25,7 @@ tgPolygon tgPolygon::RemoveDups( const tgPolygon& subject )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for ( unsigned int c = 0; c < subject.Contours(); c++ ) {
result.AddContour( tgContour::RemoveDups( subject.GetContour( c ) ) );
@ -39,6 +41,7 @@ tgPolygon tgPolygon::RemoveBadContours( const tgPolygon& subject )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for ( unsigned int c = 0; c < subject.Contours(); c++ ) {
tgContour contour = subject.GetContour(c);
@ -59,6 +62,7 @@ tgPolygon tgPolygon::RemoveCycles( const tgPolygon& subject )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for ( unsigned int c = 0; c < subject.Contours(); c++ ) {
contours.clear();
@ -84,6 +88,7 @@ tgPolygon tgPolygon::SplitLongEdges( const tgPolygon& subject, double dist )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for ( unsigned c = 0; c < subject.Contours(); c++ )
{
@ -123,6 +128,7 @@ tgPolygon tgPolygon::StripHoles( const tgPolygon& subject )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
return result;
}
@ -148,6 +154,7 @@ tgPolygon tgPolygon::Simplify( const tgPolygon& subject )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
return result;
}
@ -160,6 +167,7 @@ tgPolygon tgPolygon::RemoveTinyContours( const tgPolygon& subject )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for ( unsigned int c = 0; c < subject.Contours(); c++ ) {
tgContour contour = subject.GetContour( c );
@ -183,6 +191,7 @@ tgPolygon tgPolygon::RemoveSpikes( const tgPolygon& subject )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
for ( unsigned int c = 0; c < subject.Contours(); c++ ) {
result.AddContour( tgContour::RemoveSpikes( subject.GetContour(c) ) );
@ -275,6 +284,7 @@ tgcontour_list tgPolygon::MergeSlivers( tgpolygon_list& polys, tgcontour_list& s
result.SetMaterial( polys[j].GetMaterial() );
result.SetTexParams( polys[j].GetTexParams() );
result.SetId( polys[j].GetId() );
result.SetPreserve3D( polys[j].GetPreserve3D() );
polys[j] = result;
done = true;
}

View file

@ -61,6 +61,8 @@ tgPolygon tgPolygon::Union( const tgPolygon& subject, tgPolygon& clip )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
return result;
}
@ -127,6 +129,8 @@ tgPolygon tgPolygon::Diff( const tgPolygon& subject, tgPolygon& clip )
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
return result;
}
@ -164,6 +168,8 @@ tgPolygon tgPolygon::Intersect( const tgPolygon& subject, const tgPolygon& clip
result.SetMaterial( subject.GetMaterial() );
result.SetTexParams( subject.GetTexParams() );
result.SetId( subject.GetId() );
result.SetPreserve3D( subject.GetPreserve3D() );
return result;
}

View file

@ -190,7 +190,10 @@ void tgPolygon::Tesselate()
SGGeod p1 = SGGeod::fromDeg( to_double(tri.vertex(1).x()), to_double(tri.vertex(1).y()) );
SGGeod p2 = SGGeod::fromDeg( to_double(tri.vertex(2).x()), to_double(tri.vertex(2).y()) );
/* Check for Zero Area before inserting */
if ( !SGGeod_isEqual2D( p0, p1 ) && !SGGeod_isEqual2D( p1, p2 ) && !SGGeod_isEqual2D( p0, p2 ) ) {
AddTriangle( p0, p1, p2 );
}
++count;
} else {