Comment puis-je trouver le lat / long qui est x km au nord de une donnée lat / long?

StackOverflow https://stackoverflow.com/questions/1125144

  •  13-09-2019
  •  | 
  •  

Question

J'ai un code C # qui génère google maps. Ce code se penche sur tous les points que je dois tracer sur la carte et fonctionne ensuite les limites d'un rectangle d'inclure ces points. Il passe ensuite cette limite à l'API Google Maps pour définir le niveau de zoom approprié pour afficher tous les points sur la carte.

Ce code fonctionne bien mais j'ai une nouvelle exigence.

L'un des points peut avoir une précision qui lui est associée. Si tel est le cas, alors je dessine un cercle autour du point du rayon défini sur la valeur de précision. Encore une fois cela fonctionne très bien mais mon contrôle des limites est maintenant ne fait pas ce que je veux faire. Je veux avoir la zone de délimitation comprennent le cercle complet.

Cela nécessite un algorithme pour prendre un point x et calculer le point y qui serait z mètres au nord de x et z aussi mètres au sud de x.

Est-ce que quelqu'un a cet algorithme, de préférence en C #. J'ai trouvé un algorithme générique mais je ne semblent pas avoir mis en œuvre correctement ce que les réponses que je am se sont 1000s de km à la dérive.

Ceci est l'exemple générique

Lat/lon given radial and distance

A point {lat,lon} is a distance d out on the tc radial from point 1 if:

     lat=asin(sin(lat1)*cos(d)+cos(lat1)*sin(d)*cos(tc))
     IF (cos(lat)=0)
        lon=lon1      // endpoint a pole
     ELSE
        lon=mod(lon1-asin(sin(tc)*sin(d)/cos(lat))+pi,2*pi)-pi
     ENDIF

Et voici ma traduction C #.

  // Extend a Point North/South by the specified distance
    public static Point ExtendPoint(Point _pt, int _distance, int _bearing )
    {
        Decimal lat = 0.0;
        Decimal lng = 0.0;

        lat = Math.Asin(Math.Sin(_pt.Lat) * Math.Cos(_distance) + Math.Cos(_pt.Lat) * 
            Math.Sin(_distance) * Math.Cos(_bearing));

         if (Math.Cos(lat) == 0)
         {
            lng = _pt.Lng;      // endpoint a pole
         }
         else 
         {
             lng = (
                 (_pt.Lng - Math.Asin(Math.Sin(_bearing) * Math.Sin(_distance) / Math.Cos(lat)) 
                 + Math.PI) % (2 * Math.PI)) - Math.PI;
         }

         ret = new Point(lat,lng);
         return ret;
    }

J'appelle cette fonction avec un palier de 0 à calculer la nouvelle position nord et une valeur de 180 pour calculer la nouvelle position du sud.

Quelqu'un peut-il soit voir ce que je l'ai fait mal ou peut-être fournir un algorithme de travail connu?

Était-ce utile?

La solution

Si vous avez une latitude et la longitude, vous pouvez calculer la latitude et la longitude correctes d'un changement x km en latitude comme ceci:

new-lat = ((old-km-north + x-km-change)/40,075) * 360)
           ^ is the ratio of the                  ^ times the ratio of the circle
           of the earth the change                by 360 to get the total ratio 
           covers.                                covered in degrees.

La même chose peut appliquer à la longitude. Si vous avez la distance totale, plus le changement que vous pouvez calculer les degrés totaux de la même façon.

new-long = ((old-km-east + x-km-change)/40,075) * 360)
           ^ is the ratio of the                  ^ times the ratio of the circle
           of the earth the change                by 360 to get the total ratio 
           covers.                                covered in degrees.

Encore une fois, ces calculs doivent travailler, mais je suis en cours d'exécution de l'intuition pure, mais la logique ne semble vrai.

Edit: Comme le souligne Skizz 40,075 doit être ajusté à la circonférence de la terre à une latitude donnée en utilisant 2.pi.r.cos (LAT) ou 40074.cos (lat)

Autres conseils

J'ai un morceau très similaire de code. Il me a obtenu des résultats très proches par rapport à une autre mise en œuvre.

Je pense que le problème avec le vôtre est que vous utilisez la « distance » que la distance linéaire en mètres au lieu de la distance angulaire en radians.

/// <summary>
/// Calculates the end-point from a given source at a given range (meters) and bearing (degrees).
/// This methods uses simple geometry equations to calculate the end-point.
/// </summary>
/// <param name="source">Point of origin</param>
/// <param name="range">Range in meters</param>
/// <param name="bearing">Bearing in degrees</param>
/// <returns>End-point from the source given the desired range and bearing.</returns>
public static LatLonAlt CalculateDerivedPosition(LatLonAlt source, double range, double bearing)
{
    double latA = source.Latitude * UnitConstants.DegreesToRadians;
    double lonA = source.Longitude * UnitConstants.DegreesToRadians;
    double angularDistance = range / GeospatialConstants.EarthRadius;
    double trueCourse = bearing * UnitConstants.DegreesToRadians;

    double lat = Math.Asin(
        Math.Sin(latA) * Math.Cos(angularDistance) + 
        Math.Cos(latA) * Math.Sin(angularDistance) * Math.Cos(trueCourse));

    double dlon = Math.Atan2(
        Math.Sin(trueCourse) * Math.Sin(angularDistance) * Math.Cos(latA), 
        Math.Cos(angularDistance) - Math.Sin(latA) * Math.Sin(lat));

    double lon = ((lonA + dlon + Math.PI) % UnitConstants.TwoPi) - Math.PI;

    return new LatLonAlt(
        lat * UnitConstants.RadiansToDegrees, 
        lon * UnitConstants.RadiansToDegrees, 
        source.Altitude);
}

public const double EarthRadius = 6378137.0;   //  WGS-84 ellipsoid parameters

et LatLonAlt est en degrés / mètres (conversion a lieu en interne). Ajuster au besoin.

Je suppose que vous pouvez comprendre ce que la valeur de UnitConstants.DegreesToRadians est:)

Pour les paresseux, (comme moi;)) une solution de copier-coller, la version de Erich Mirabal avec des changements mineurs:

using System.Device.Location; // add reference to System.Device.dll
public static class GeoUtils
{
    /// <summary>
    /// Calculates the end-point from a given source at a given range (meters) and bearing (degrees).
    /// This methods uses simple geometry equations to calculate the end-point.
    /// </summary>
    /// <param name="source">Point of origin</param>
    /// <param name="range">Range in meters</param>
    /// <param name="bearing">Bearing in degrees</param>
    /// <returns>End-point from the source given the desired range and bearing.</returns>
    public static GeoCoordinate CalculateDerivedPosition(this GeoCoordinate source, double range, double bearing)
    {
        var latA = source.Latitude * DegreesToRadians;
        var lonA = source.Longitude * DegreesToRadians;
        var angularDistance = range / EarthRadius;
        var trueCourse = bearing * DegreesToRadians;

        var lat = Math.Asin(
            Math.Sin(latA) * Math.Cos(angularDistance) +
            Math.Cos(latA) * Math.Sin(angularDistance) * Math.Cos(trueCourse));

        var dlon = Math.Atan2(
            Math.Sin(trueCourse) * Math.Sin(angularDistance) * Math.Cos(latA),
            Math.Cos(angularDistance) - Math.Sin(latA) * Math.Sin(lat));

        var lon = ((lonA + dlon + Math.PI) % (Math.PI*2)) - Math.PI;

        return new GeoCoordinate(
            lat * RadiansToDegrees,
            lon * RadiansToDegrees,
            source.Altitude);
    }

    private const double DegreesToRadians = Math.PI/180.0;
    private const double RadiansToDegrees = 180.0/ Math.PI;
    private const double EarthRadius = 6378137.0;
}

Utilisation:

[TestClass]
public class CalculateDerivedPositionUnitTest
{
    [TestMethod]
    public void OneDegreeSquareAtEquator()
    {
        var center = new GeoCoordinate(0, 0);
        var radius = 111320;
        var southBound = center.CalculateDerivedPosition(radius, -180);
        var westBound = center.CalculateDerivedPosition(radius, -90);
        var eastBound = center.CalculateDerivedPosition(radius, 90);
        var northBound = center.CalculateDerivedPosition(radius, 0);

        Console.Write($"leftBottom: {southBound.Latitude} , {westBound.Longitude} rightTop: {northBound.Latitude} , {eastBound.Longitude}");
    }
}

Je ne sais pas si je manque quelque chose ici, mais je pense que la question pourrait être reformulée comme «J'ai un lat / point de longitude, et je veux trouver le point x mètres au nord et au sud de x mètres ce point. "

Si c'est la question, alors vous n'avez pas besoin de trouver une nouvelle longitude (qui rend les choses plus simples), vous avez juste besoin d'une nouvelle latitude. Un degré de latitude est à 60 miles nautiques à long partout sur la Terre, et un mile nautique est 1.852 mètres. Ainsi, pour de nouvelles latitudes x mètres au nord et au sud:

north_lat = lat + x / (1852 * 60)
north_lat = min(north_lat, 90)

south_lat = lat - x / (1852 * 60)
south_lat = max(south_lat, -90)

Ce n'est pas tout à fait exact parce que la Terre est pas une sphère parfaite avec exactement 60 miles nautiques entre chaque degré de latitude. Cependant, les autres réponses supposent que les lignes de latitude sont à égale distance, donc je suppose que vous ne vous inquiétez pas à ce sujet. Si vous êtes intéressé à combien erreur qui pourrait introduire, il y a une belle table sur Wikipédia qui montre « la distance de la surface par 1 ° changement de latitude » pour différentes latitudes à ce lien:

http://en.wikipedia.org/wiki/Latitude#Degree_length

Il y a des problèmes avec les deux équations sur Ed William site plutôt impressionnant de ... mais je ne les ai pas analyser de voir pourquoi.

Une troisième équation je trouve ici semble donner des résultats appropriés .

est ici le cas de test en php ... la troisième équation est correcte, les deux premières donnent des valeurs incorrectes d'une manière extravagante pour la longitude.

<?php
            $lon1 = -108.553412; $lat1 = 35.467155; $linDistance = .5; $bearing = 170;
            $lon1 = deg2rad($lon1); $lat1 = deg2rad($lat1);
            $distance = $linDistance/6371;  // convert dist to angular distance in radians
            $bearing = deg2rad($bearing);

            echo "lon1: " . rad2deg($lon1) . " lat1: " . rad2deg($lat1) . "<br>\n";

// doesn't work
            $lat2 = asin(sin($lat1) * cos($distance) + cos($lat1) * sin($distance) * cos($bearing) );
            $dlon = atan2(sin($bearing) * sin($distance) * cos($lat1), cos($distance) - sin($lat1) * sin($lat2));
            $lon2 = (($lon1 - $dlon + M_PI) % (2 * M_PI)) - M_PI;  // normalise to -180...+180

            echo "lon2: " . rad2deg($lon2) . " lat2: " . rad2deg($lat2) . "<br>\n";

// same results as above
            $lat3 = asin( (sin($lat1) * cos($distance)) + (cos($lat1) * sin($distance) * cos($bearing)));
            $lon3 = (($lon1 - (asin(sin($bearing) * sin($distance) / cos($lat3))) + M_PI) % (2 * M_PI)) - M_PI;

            echo "lon3: " . rad2deg($lon3) . " lat3: " . rad2deg($lat3) . "<br>\n";

// gives correct answer... go figure
            $lat4 = asin(sin($lat1) * cos($linDistance/6371) + cos($lat1) * sin($linDistance/6371) * cos($bearing) );
            $lon4 = $lon1 + atan2( (sin($bearing) * sin($linDistance/6371) * cos($lat1) ), (cos($linDistance/6371) - sin($lat1) * sin($lat2)));

            echo "lon4: " . rad2deg($lon4) . " lat4: " . rad2deg($lat4) . "<br>\n";
?>

Remarque J'ai reçu par e-mail de l'auteur (Ed Williams) des deux premières équations:

  

De mes "notes de mise en œuvre":

     

Note sur la fonction mod. Cela semble être mis en œuvre différemment   langues différentes, avec des conventions différentes quant à savoir si le signe de la   résultat suit le signe du diviseur ou du dividende. (Nous voulons que le signe   suivre le diviseur ou euclidien. Le fmod de C et de% de Java ne fonctionnent pas.)   Dans ce document, Mod (y, x) est le reste de la division par y x et toujours   se situe dans la plage de 0 <= mod      

Si vous avez une fonction de plancher (int dans Excel), qui renvoie étage (x) =   « Le plus grand nombre entier inférieur ou égal à x » par exemple, étage (-2,3) = - 3 et   étage (2,3) = 2

mod(y,x) = y - x*floor(y/x)
  

Les éléments suivants doivent travailler en l'absence d'un plancher Fonction- quel que soit   si "int" ou tronque tours vers le bas:

mod=y - x * int(y/x)
if ( mod < 0) mod = mod + x
  

php est comme fmod en C et le fait "mal" pour mes besoins.

Il est plus précis si vous devez d'abord reprojeter à UTM, puis vérifiez la distance.

Hope this helps

scroll top