# GeoPoint distanceTo Classic List Threaded 4 messages Reply | Threaded
Open this post in threaded view
|

## GeoPoint distanceTo if getSRS() is projected, this method calculate the distance of (dx, dy, dz). but if getSRS() is geographic, this method only calculate the distance of (dx, dy). dz is ignored. is it right?
Reply | Threaded
Open this post in threaded view
|

## Re: GeoPoint distanceTo You are correct - that is what the code does. It should probably be ground distance and the projected case should probably be 2D only. I have opened a github ticket for it. Thanks for spotting. https://github.com/gwaldron/osgearth/issues/1668 Glenn Waldron / Pelican Mapping
Reply | Threaded
Open this post in threaded view
|

## Re: GeoPoint distanceTo In reply to this post by liugaijin //GeoData /** * Calculates the distance in meters from this geopoint to another. * dim = 3, 3d(spatial) distance; dim = 2, 2d(ground) distance */ double distanceTo(const GeoPoint& rhs, int dim = 3) const; //GeoData.cpp double GeoPoint::distanceTo(const GeoPoint& rhs, int dim) const {     if ( getSRS()->isProjected() && rhs.getSRS()->isProjected() )     {         if ( getSRS()->isEquivalentTo(rhs.getSRS()) )         {             osg::Vec3d vec = vec3d() - rhs.vec3d();             if (dim == 2)//ground distance                 vec.z() = 0;             return vec.length();         }         else         {             GeoPoint rhsT = rhs.transform(getSRS());             osg::Vec3d vec = vec3d() - rhsT.vec3d();            (dim == 2)//ground distance                 vec.z() = 0;             return vec.length();         }     }     else     {         // https://en.wikipedia.org/wiki/Geographical_distance#Ellipsoidal-surface_formulae        GeoPoint p1 = transform( getSRS()->getGeographicSRS() );         GeoPoint p2 = rhs.transform( p1.getSRS() );         double Re = getSRS()->getEllipsoid()->getRadiusEquator();         double Rp = getSRS()->getEllipsoid()->getRadiusPolar();         double F  = (Re-Rp)/Re; // flattening         double             lat1 = osg::DegreesToRadians(p1.y()),             lon1 = osg::DegreesToRadians(p1.x()),             lat2 = osg::DegreesToRadians(p2.y()),             lon2 = osg::DegreesToRadians(p2.x());         double B1 = atan( (1.0-F)*tan(lat1) );         double B2 = atan( (1.0-F)*tan(lat2) );         double P = (B1+B2)/2.0;         double Q = (B2-B1)/2.0;         double G = acos(sin(B1)*sin(B2) + cos(B1)*cos(B2)*cos(fabs(lon2-lon1)));         double             sinG = sin(G),             sinP = sin(P), sinQ = sin(Q),             cosP = cos(P), cosQ = cos(Q),             sinG2 = sin(G/2.0), cosG2 = cos(G/2.0);         double X = (G-sinG)*((sinP*sinP*cosQ*cosQ)/(cosG2*cosG2));         double Y = (G+sinG)*((cosP*cosP*sinQ*sinQ)/(sinG2*sinG2));         double dist = Re*(G-(F/2.0)*(X+Y));         // NaN could mean start/end points are the same         dist = osg::isNaN(dist)? 0.0 : dist;         if (dim == 3)//spatial distance         {             //altitudeMode of p1 and p2 must be ALTMODE_ABSOLUTE             dist = hypot(dist, p2.z() - p1.z());//refer to proj_lpz_dist in proj4         }         return dist;     } } This is my change.
Reply | Threaded
Open this post in threaded view
|

## Re: GeoPoint distanceTo Let's continue this discussion under the GitHub ticket - thanks Glenn Waldron / Pelican Mapping