Заданная точка (широта, долгота), расстояние и азимут, как получить новую широту и долготу

Я нашел кусок кода в Интернете. Он рассчитывает минимальный ограничивающий прямоугольник по заданной точке широты и расстояния.

private static void GetlatLon(double LAT, double LON, double distance, double angle, out double newLon, out double newLat)
        {     
            double dx = distance * 1000 * Math.Sin(angle * Math.PI / 180.0);
            double dy = distance * 1000 * Math.Cos(angle * Math.PI / 180.0);
            double ec = 6356725 + 21412 * (90.0 - LAT) / 90.0;
            double ed = ec * Math.Cos(LAT * Math.PI / 180);
            newLon = (dx / ed + LON * Math.PI / 180.0) * 180.0 / Math.PI;
            newLat = (dy / ec + LAT * Math.PI / 180.0) * 180.0 / Math.PI;

            }

    public static void GetRectRange(double centorlatitude, double centorLogitude, double distance,
                                  out double maxLatitude, out double minLatitude, out double maxLongitude, out double minLongitude)
{       
            GetlatLon(centorlatitude, centorLogitude, distance, 0, out temp, out maxLatitude);
            GetlatLon(centorlatitude, centorLogitude, distance, 180, out temp, out minLatitude);
            GetlatLon(centorlatitude, centorLogitude, distance, 90, out minLongitude, out temp);
            GetlatLon(centorlatitude, centorLogitude, distance, 270, out maxLongitude, out temp);
}



double ec = 6356725 + 21412 * (90.0 - LAT) / 90.0; //why?
double ed = ec * Math.Cos(LAT * Math.PI / 180);    // why?
dx / ed                                            //why?
dy / ec                                            //why?

6378137 - радиус экватора, 6356725 - полярный радиус, 21412 =6378137 -6356725. по ссылке я знаю немного о значениях. Но эти четыре строки, я не знаю почему. Не могли бы вы помочь, чтобы дать больше информации? Не могли бы вы помочь, дайте мне знать вывод формулы?

По ссылке в разделе "Точка назначения с учетом расстояния и направления от начальной точки" он дает другую формулу для получения результата. Каков вывод формулы?

По этой ссылке я знаю происхождение формулы Haversine, это очень информативно. Я не думаю, что формула в разделе "Пункт назначения с учетом расстояния и пеленга от начальной точки" является просто простым изменением Haversine.

Большое спасибо!

2 ответа

Решение

Это яркий пример того, почему комментирование вашего кода делает его более читабельным и понятным. Математически вы смотрите на следующее:

double ec = 6356725 + 21412 * (90,0 - LAT) / 90,0; //Зачем?

Это мера эксцентриситета для учета экваториальной выпуклости в некотором роде. 21412 Это, как вы знаете, разница в радиусе земли между экватором и полюсом. 6356725 это полярный радиус. (90.0 - LAT) / 90.0 является 1 на экваторе и 0 на полюсе. Формула просто оценивает, сколько выпуклости присутствует на любой данной широте.

double ed = ec * Math.Cos (LAT * Math.PI / 180); // Зачем?

(LAT * Math.PI / 180) это перевод широты из градусов в радианы. cos (0) = 1 а также cos(1) = 0Таким образом, на экваторе вы применяете полную величину эксцентриситета, в то время как на полюсе вы не применяете ничего. Аналогично предыдущей строке.

dx / ed // почему?

dy / ec // почему?

Выше, кажется, дробные добавления к расстоянию в обоих x а также y направления, относящиеся к выпуклости в любой широте / долготе, используемые в newLonnewLat расчет, чтобы прибыть на новое место.

Я не исследовал найденный вами фрагмент кода, но математически это именно то, что происходит. Надеюсь, это приведет вас в правильном направлении.


Пример Хаверсин в C

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

double m2ft (double l) {            /* convert meters to feet       */
    return l/(1200.0/3937.0);
}

double ft2smi (double l) {          /* convert feet to statute miles*/
    return l/5280.0;
}

double km2smi (double l) {          /* convert km to statute mi.    */
    return ft2smi(m2ft( l * 1000.0 ));
}

static const double deg2rad = 0.017453292519943295769236907684886;
static const double earth_rad_m = 6372797.560856;

typedef struct pointd {
    double lat;
    double lon;
} pointd;

/* Computes the arc, in radian, between two WGS-84 positions.
   The result is equal to Distance(from,to)/earth_rad_m
      = 2*asin(sqrt(h(d/earth_rad_m )))
   where:
      d is the distance in meters between 'from' and 'to' positions.
      h is the haversine function: h(x)=sin²(x/2)
   The haversine formula gives:
      h(d/R) = h(from.lat-to.lat)+h(from.lon-to.lon)+cos(from.lat)*cos(to.lat)
   http://en.wikipedia.org/wiki/Law_of_haversines
 */
double arcradians (const pointd *from, const pointd *to)
{
    double latitudeArc  = (from-> lat - to-> lat) * deg2rad;
    double longitudeArc = (from-> lon - to-> lon) * deg2rad;

    double latitudeH = sin (latitudeArc * 0.5);
    latitudeH *= latitudeH;

    double lontitudeH = sin (longitudeArc * 0.5);
    lontitudeH *= lontitudeH;

    double tmp = cos (from-> lat * deg2rad) * cos (to-> lat * deg2rad);

    return 2.0 * asin (sqrt (latitudeH + tmp*lontitudeH));
}

/* Computes the distance, in meters, between two WGS-84 positions.
   The result is equal to earth_rad_m*ArcInRadians(from,to)
 */
double dist_m (const pointd *from, const pointd *to) {
    return earth_rad_m * arcradians (from, to);
}

int main (int argc, char **argv) {

    if (argc < 5 ) {
        fprintf (stderr, "Error: insufficient input, usage: %s (lat,lon) (lat,lon)\n", argv[0]);
        return 1;
    }

    pointd points[2];

    points[0].lat = strtod (argv[1], NULL);
    points[0].lon = strtod (argv[2], NULL);

    points[1].lat = strtod (argv[3], NULL);
    points[1].lon = strtod (argv[4], NULL);

    printf ("\nThe distance in meters from 1 to 2 (smi): %lf\n\n", km2smi (dist_m (&points[0], &points[1])/1000.0) );

    return 0;
}

/* Results/Example.
    ./bin/gce 31.77 -94.61 31.44 -94.698
    The distance in miles from Nacogdoches to Lufkin, Texas (smi): 23.387997 miles
 */

Я предполагаю, что 6356725 как-то связан с радиусом Земли. Проверьте этот ответ, а также взгляните на формулу Haversine.

Другие вопросы по тегам