Add runway shoulders support.

This was harder than expected as the shoulders were often partially
covered by runway polygons. This was caused by inaccuracies introduced
by the calculation of dimensions from the runway center.

Now we create shoulders if set in apt.dat. If no shoulders are set,
a 1m broad stripe is added to concrete and asphalt runways to make
the blending runway-airport surface more visually appealing.

src/Airports/GenAirports850/airport.cxx

@@ -438,6 +438,15 @@ void Airport::BuildBtg(const string& root, const string_list& elev_src )
         SG_LOG(SG_GENERAL, SG_DEBUG, "no pavements");
     }
 
+    // Build runway shoulders here
+    for (i=0; i<runways.size(); i++ )
+    {
+        if ( runways[i]->GetsShoulder() )
+        {
+            runways[i]->BuildShoulder( altitude, &rwy_polys, &rwy_tps, &accum );
+        }
+    }
+
     // build the base and clearing if there's a boundary
     if (boundary)
     {

src/Airports/GenAirports850/runway.hxx

@@ -39,7 +39,18 @@ public:
         return ( Point3D( (rwy.lon[0]+rwy.lon[1])/2.0f, (rwy.lat[0]+rwy.lat[1])/2.0f, 0.0f) );
     }
 
+    bool GetsShoulder()
+    {
+        if (rwy.surface < 3 || rwy.shoulder > 0)
+            return true;
+            else return false;
+    }
+
     int BuildBtg( float alt_m, superpoly_list* rwy_polys, texparams_list* texparams, superpoly_list* rwy_lights, ClipPolyType* accum, TGPolygon* apt_base, TGPolygon* apt_clearing );
+    void BuildShoulder( float alt_m,
+                                superpoly_list *rwy_polys,
+                                texparams_list *texparams,
+                                ClipPolyType *accum );
     
 private:
     struct TGRunway {
@@ -84,7 +95,6 @@ private:
         return ( gen_wgs84_area( Point3D(GetStart()), Point3D(GetEnd()), rwy.length + 2.0*length_extend, displ1, displ2, rwy.width + 2.0*width_extend, rwy.heading, alt_m, false) );
     }
 
-
     void gen_rw_designation( const std::string& material,
                              TGPolygon poly, double heading, string rwname,
                              double &start_pct, double &end_pct,

src/Airports/GenAirports850/rwy_gen.cxx

@@ -22,6 +22,7 @@
 #include <simgear/compiler.h>
 #include <simgear/constants.h>
 #include <simgear/debug/logstream.hxx>
+#include <simgear/math/sg_geodesy.hxx>
 #include "runway.hxx"
 
 #include <stdlib.h>
@@ -332,3 +333,78 @@ for ( int rwhalf=0; rwhalf<2; ++rwhalf ){
 
 }
 }
+
+void Runway::BuildShoulder( float alt_m,
+                            superpoly_list *rwy_polys,
+                            texparams_list *texparams,
+                            ClipPolyType *accum )
+{
+    string shoulder_surface = "";
+    double shoulder_width;
+    if (rwy.shoulder > 0){  // Add a shoulder to the runway
+        shoulder_width = rwy.width * 0.15;
+        if (shoulder_width > 8.0){
+            shoulder_width = 8.0;
+        }
+        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
+        shoulder_width = 1;
+        if (rwy.surface == 1){
+            shoulder_surface = "pa_shoulder_f";
+        } else if (rwy.surface == 2){
+            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 );
+
+    double lat, lon,r;
+    // Create both shoulder sides
+    for (int i=0; i<2; ++i){
+        double step;
+        double lat, lon,r;
+
+        if (i == 0){
+            step= (rwy.width + shoulder_width)*0.5;
+        } else if (i == 1) {
+            step= -(rwy.width + shoulder_width)*0.5;
+        }
+        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 shoulder1 = Point3D( lon, lat, 0.0f );
+
+        geo_direct_wgs_84 ( alt_m, rwy.lat[1], rwy.lon[1], left_hdg,
+                            step, &lat, &lon, &r );
+
+        Point3D shoulder2 = Point3D( lon, lat, 0.0f );
+
+        TGPolygon shoulder = gen_wgs84_area( shoulder1, shoulder2, 0.0, 0.0, 0.0, shoulder_width, rwy.heading, alt_m, false);
+
+        TGSuperPoly sp;
+        TGTexParams tp;
+        TGPolygon clipped = tgPolygonDiff( shoulder, *accum );
+        TGPolygon split = tgPolygonSplitLongEdges( clipped, 400.0 );
+
+        sp.erase();
+        sp.set_poly( split );
+        sp.set_material( shoulder_surface );
+        rwy_polys->push_back( sp );
+
+        *accum = tgPolygonUnion( shoulder, *accum );
+        tp = TGTexParams( shoulder.get_pt(0,0), shoulder_width , rwy.length + 2, rwy.heading );
+        if (i == 1){
+            tp.set_maxu(0);
+            tp.set_minu(1);
+        }
+        texparams->push_back( tp );
+    }
+}