Bug 550 / 454 work (not yet fully fixed)

Adding a spherical projection mode which will work at the poles and
+180/-180 line correctly. Not yet enabled since it has other issues to
be addressed.

diff --git a/src/GUI/MapWidget.cxx b/src/GUI/MapWidget.cxx
index b063bc79f8d8b0a6cf17e5800f430aa0dc12b782..393a1c4e6fc078f34dc112136cf454e0c4505134 100644 (file)
--- a/src/GUI/MapWidget.cxx
+++ b/src/GUI/MapWidget.cxx
@@ -387,7 +387,7 @@ MapWidget::MapWidget(int x, int y, int maxX, int maxY) :
   _width = maxX - x;
   _height = maxY - y;
   _hasPanned = false;
-  _orthoAzimuthProject = false;
+  _projection = PROJECTION_SAMSON_FLAMSTEED;
   _magneticHeadings = false;
   
   MapData::setFont(legendFont);
@@ -1621,26 +1621,45 @@ SGVec2d MapWidget::project(const SGGeod& geod) const
   SGVec2d p;
   double r = earth_radius_lat(geod.getLatitudeRad());
   
-  if (_orthoAzimuthProject) {
-    // http://mathworld.wolfram.com/OrthographicProjection.html
-    double cosTheta = cos(geod.getLatitudeRad());
-    double sinDLambda = sin(geod.getLongitudeRad() - _projectionCenter.getLongitudeRad());
-    double cosDLambda = cos(geod.getLongitudeRad() - _projectionCenter.getLongitudeRad());
-    double sinTheta1 = sin(_projectionCenter.getLatitudeRad());
-    double sinTheta = sin(geod.getLatitudeRad());
-    double cosTheta1 = cos(_projectionCenter.getLatitudeRad());
-    
-    p = SGVec2d(cosTheta * sinDLambda,
-                (cosTheta1 * sinTheta) - (sinTheta1 * cosTheta * cosDLambda)) * r * currentScale();
+    switch (_projection) {
+    case PROJECTION_SAMSON_FLAMSTEED:
+    {
+        // Sanson-Flamsteed projection, relative to the projection center
+        double lonDiff = geod.getLongitudeRad() - _projectionCenter.getLongitudeRad(),
+        latDiff = geod.getLatitudeRad() - _projectionCenter.getLatitudeRad();
+        
+        p = SGVec2d(cos(geod.getLatitudeRad()) * lonDiff, latDiff) * r * currentScale();
+        break;
+    }
+            
+    case PROJECTION_ORTHO_AZIMUTH:
+    {
+        // http://mathworld.wolfram.com/OrthographicProjection.html
+        double cosTheta = cos(geod.getLatitudeRad());
+        double sinDLambda = sin(geod.getLongitudeRad() - _projectionCenter.getLongitudeRad());
+        double cosDLambda = cos(geod.getLongitudeRad() - _projectionCenter.getLongitudeRad());
+        double sinTheta1 = sin(_projectionCenter.getLatitudeRad());
+        double sinTheta = sin(geod.getLatitudeRad());
+        double cosTheta1 = cos(_projectionCenter.getLatitudeRad());
+        
+        p = SGVec2d(cosTheta * sinDLambda,
+                    (cosTheta1 * sinTheta) - (sinTheta1 * cosTheta * cosDLambda)) * r * currentScale();
+        break;
+    }
+            
+    case PROJECTION_SPHERICAL:
+    {
+        SGVec3d cartCenter = SGVec3d::fromGeod(_projectionCenter);
+        SGVec3d cartPt = SGVec3d::fromGeod(geod) - cartCenter;
+        
+        // rotate relative to projection center
+        SGQuatd orient = SGQuatd::fromLonLat(_projectionCenter);
+        cartPt = orient.rotateBack(cartPt);
+        return SGVec2d(cartPt.y(), cartPt.x()) * currentScale();
+        break;
+    }
+    } // of projection mode switch
     
-  } else {
-    // Sanson-Flamsteed projection, relative to the projection center
-    double lonDiff = geod.getLongitudeRad() - _projectionCenter.getLongitudeRad(),
-      latDiff = geod.getLatitudeRad() - _projectionCenter.getLatitudeRad();
-
-    p = SGVec2d(cos(geod.getLatitudeRad()) * lonDiff, latDiff) * r * currentScale();
-      
-  }
   
 // rotate as necessary
   double cost = cos(_upHeading * SG_DEGREES_TO_RADIANS),
@@ -1658,24 +1677,43 @@ SGGeod MapWidget::unproject(const SGVec2d& p) const
   SGVec2d ur(cost * p.x() - sint * p.y(),
              sint * p.x() + cost * p.y());
 
-  double r = earth_radius_lat(_projectionCenter.getLatitudeRad());
-  SGVec2d unscaled = ur * (1.0 / (currentScale() * r));
-
-  if (_orthoAzimuthProject) {
-      double phi = length(p);
-      double c = asin(phi);
-      double sinTheta1 = sin(_projectionCenter.getLatitudeRad());
-      double cosTheta1 = cos(_projectionCenter.getLatitudeRad());
-      
-      double lat = asin(cos(c) * sinTheta1 + ((unscaled.y() * sin(c) * cosTheta1) / phi));
-      double lon = _projectionCenter.getLongitudeRad() + 
+  
+
+    switch (_projection) {
+    case PROJECTION_SAMSON_FLAMSTEED:
+    {
+        double r = earth_radius_lat(_projectionCenter.getLatitudeRad());
+        SGVec2d unscaled = ur * (1.0 / (currentScale() * r));
+        double lat = unscaled.y() + _projectionCenter.getLatitudeRad();
+        double lon = (unscaled.x() / cos(lat)) + _projectionCenter.getLongitudeRad();
+        return SGGeod::fromRad(lon, lat);
+    }
+        
+    case PROJECTION_ORTHO_AZIMUTH:
+    {
+        double r = earth_radius_lat(_projectionCenter.getLatitudeRad());
+        SGVec2d unscaled = ur * (1.0 / (currentScale() * r));
+        
+        double phi = length(p);
+        double c = asin(phi);
+        double sinTheta1 = sin(_projectionCenter.getLatitudeRad());
+        double cosTheta1 = cos(_projectionCenter.getLatitudeRad());
+        
+        double lat = asin(cos(c) * sinTheta1 + ((unscaled.y() * sin(c) * cosTheta1) / phi));
+        double lon = _projectionCenter.getLongitudeRad() +
         atan((unscaled.x()* sin(c)) / (phi  * cosTheta1 * cos(c) - unscaled.y() * sinTheta1 * sin(c)));
-      return SGGeod::fromRad(lon, lat);
-  } else {
-      double lat = unscaled.y() + _projectionCenter.getLatitudeRad();
-      double lon = (unscaled.x() / cos(lat)) + _projectionCenter.getLongitudeRad();
-      return SGGeod::fromRad(lon, lat);
-  }
+        return SGGeod::fromRad(lon, lat);
+    }
+        
+    case PROJECTION_SPHERICAL:
+    {
+        SGVec2d unscaled = ur * (1.0 / currentScale());
+        SGQuatd orient = SGQuatd::fromLonLat(_projectionCenter);
+        SGVec3d cartCenter = SGVec3d::fromGeod(_projectionCenter);
+        SGVec3d cartPt = orient.rotate(SGVec3d(unscaled.x(), unscaled.y(), 0.0));
+        return SGGeod::fromCart(cartPt + cartCenter);
+    }
+    } // of projection mode switch
 }
 
 double MapWidget::currentScale() const

diff --git a/src/GUI/MapWidget.hxx b/src/GUI/MapWidget.hxx
index 877ab7d53a907bbab8828ff6cb54ac19706eec28..0a826f761b61884c79c0db3b2a2460ef0d173e55 100644 (file)
--- a/src/GUI/MapWidget.hxx
+++ b/src/GUI/MapWidget.hxx
@@ -33,6 +33,12 @@ public:
     
   void setProperty(SGPropertyNode_ptr prop);
 private:
+    enum Projection {
+        PROJECTION_SAMSON_FLAMSTEED,
+        PROJECTION_ORTHO_AZIMUTH,
+        PROJECTION_SPHERICAL
+    };
+    
   int zoom() const;
   
   void handlePan(int x, int y);
@@ -103,7 +109,7 @@ private:
   bool _hasPanned;
   
   SGGeod _projectionCenter;
-  bool _orthoAzimuthProject;
+  Projection _projection;
   SGGeod _aircraft;
   SGGeod _clickGeod;
   SGVec2d _hitLocation;