Преобразование Lat/Longs в координаты X/Y

Ответ 1

Вам нужно больше информации, чем просто одна пара lat/lon, чтобы это сделать.

На данном этапе предоставленная вами информация пропускает две вещи:

• Как большая область имеет покрытие вашего изображения (в терминах lat/lon)? Основываясь на том, что вы предоставили, я не знаю, показывает ли изображение область шириной в один метр или ширину километра.
• Какое место на вашем изображении ссылается на вашу координату (-37.803134, 145.132377)? Это один из углов? Где-то посередине?

Я также собираюсь предположить, что ваш образ выровнен на север/юг - например, он не имеет север, указывающий на верхний левый угол. Это будет усложнять ситуацию.

Самый простой способ - точно определить, какие координаты lat/lon соответствуют пикселю (0, 0) и пикселу (1017, 915). Затем вы можете определить пиксель, соответствующий заданной координате lat/lon, через interpolation.

Чтобы кратко описать этот процесс, представьте, что ваш (-37.803134, 145.132377) lat/lon соответствует вашему (0, 0) пикселю и что вы обнаружили, что ваш (1017, 915) пиксель соответствует лат /lon (-37.798917, 145.138535). Предполагая обычное соглашение с пикселем (0, 0) в нижнем левом углу, это означает, что север находится на изображении.

Затем, если вы заинтересованы в целевой координате (-37.801465, 145.134984), вы можете определить соответствующее количество пикселей вверх по изображению, которое оно выглядит следующим образом:

pixelY = ((targetLat - minLat) / (maxLat - minLat)) * (maxYPixel - minYPixel)
       = ((-37.801465 - -37.803134) / (-37.798917 - -37.803134)) * (915 - 0)
       = 362.138

То есть соответствующий пиксель составляет 362 пикселя от нижней части изображения. Затем вы можете сделать то же самое для горизонтального размещения пикселей, но вместо этого используйте долготы и X пикселей.

Часть ((targetLat - minLat) / (maxLat - minLat)) показывает, насколько вы находитесь между двумя опорными координатами, и давая 0, чтобы указать, что вы находитесь на первом, 1, чтобы указать второй, и числа между ними для указания местоположений между ними. Так, например, это произвело бы 0,25, чтобы указать, что вы находитесь на 25% пути к северу между двумя опорными координатами. Последний бит преобразует это в эквивалентные пиксели.

НТН!

EDIT Хорошо, основываясь на вашем комментарии, я могу быть немного более конкретным. Учитывая, что вы, кажется, используете верхний левый угол в качестве основной контрольной точки, я буду использовать следующие определения:

minLat = -37.803134
maxLat = -37.806232
MAP_HEIGHT = 916

Затем, если мы используем примерную координату (-37.804465, 145.134984), координата Y соответствующего пикселя относительно левого верхнего угла:

pixelY = ((targetLat - minLat) / (maxLat - minLat)) * (MAP_HEIGHT - 1)
       = ((-37.804465 - -37.803134) / (-37.806232 - -37.803134)) * 915
       = 393.11

Следовательно, соответствующий пиксель составляет 393 пикселя сверху. Я дам вам выработать горизонтальный эквивалент для себя - это в основном то же самое. ПРИМЕЧАНИЕ -1 с MAP_HEIGHT заключается в том, что если вы начинаете с нуля, максимальное число пикселей равно 915, а не 916.

РЕДАКТИРОВАТЬ: Одна вещь, которую я хотел бы воспользоваться возможностью, чтобы отметить, что это приближение. В действительности нет простой линейной зависимости между координатами широты и долготы и другими формами декартовых координат по ряду причин, включая проекции, которые используются при создании карт, и тот факт, что Земля не является идеальной сферой. В небольших районах это приближение достаточно близко, чтобы не иметь существенных различий, но в больших масштабах могут возникнуть расхождения. В зависимости от ваших потребностей, YMMV. (Моя благодарность за ура, ответ на который ниже напомнил мне, что это так).

Ответ 2

если вы ищете точное преобразование геодезических (lot, lan) в определенную декартовую координату (x, y метров от контрольной точки), вы можете сделать с моим фрагментом кода здесь, эта функция примет геодезическую координату в радиан и выводит результат в x, y

ввод:

• refLat, refLon: геодезическая координата который вы определили как 0,0 в картезианском координата (единица находится в радиусе)
• lat, lon: geodetic координаты, которые вы хотите вычислить его декартовую координату (единица находится в радиусе)
• xOffset, yOffset: результат в декартовой координаты x, y (единица измерения находится в метрах)

код:

#define GD_semiMajorAxis 6378137.000000
#define GD_TranMercB     6356752.314245
#define GD_geocentF      0.003352810664

void geodeticOffsetInv( double refLat, double refLon,
                        double lat,    double lon, 
                        double& xOffset, double& yOffset )
{
    double a = GD_semiMajorAxis;
    double b = GD_TranMercB;
    double f = GD_geocentF;

    double L     = lon-refLon
    double U1    = atan((1-f) * tan(refLat));
    double U2    = atan((1-f) * tan(lat));
    double sinU1 = sin(U1); 
    double cosU1 = cos(U1);
    double sinU2 = sin(U2);
    double cosU2 = cos(U2);

    double lambda = L;
    double lambdaP;
    double sinSigma;
    double sigma;
    double cosSigma;
    double cosSqAlpha;
    double cos2SigmaM;
    double sinLambda;
    double cosLambda;
    double sinAlpha;
    int iterLimit = 100;
    do {
        sinLambda = sin(lambda);
        cosLambda = cos(lambda);
        sinSigma = sqrt((cosU2*sinLambda) * (cosU2*sinLambda) + 
                        (cosU1*sinU2-sinU1*cosU2*cosLambda) * 
                        (cosU1*sinU2-sinU1*cosU2*cosLambda) );
        if (sinSigma==0)
        {
            xOffset = 0.0;
            yOffset = 0.0;
            return ;  // co-incident points
        }
        cosSigma    = sinU1*sinU2 + cosU1*cosU2*cosLambda;
        sigma       = atan2(sinSigma, cosSigma);
        sinAlpha    = cosU1 * cosU2 * sinLambda / sinSigma;
        cosSqAlpha  = 1 - sinAlpha*sinAlpha;
        cos2SigmaM  = cosSigma - 2*sinU1*sinU2/cosSqAlpha;
        if (cos2SigmaM != cos2SigmaM) //isNaN
        {
            cos2SigmaM = 0;  // equatorial line: cosSqAlpha=0 (§6)
        }
        double C = f/16*cosSqAlpha*(4+f*(4-3*cosSqAlpha));
        lambdaP = lambda;
        lambda = L + (1-C) * f * sinAlpha *
            (sigma + C*sinSigma*(cos2SigmaM+C*cosSigma*(-1+2*cos2SigmaM*cos2SigmaM)));
    } while (fabs(lambda-lambdaP) > 1e-12 && --iterLimit>0);

    if (iterLimit==0)
    {
        xOffset = 0.0;
        yOffset = 0.0;
        return;  // formula failed to converge
    }

    double uSq  = cosSqAlpha * (a*a - b*b) / (b*b);
    double A    = 1 + uSq/16384*(4096+uSq*(-768+uSq*(320-175*uSq)));
    double B    = uSq/1024 * (256+uSq*(-128+uSq*(74-47*uSq)));
    double deltaSigma = B*sinSigma*(cos2SigmaM+B/4*(cosSigma*(-1+2*cos2SigmaM*cos2SigmaM)-
        B/6*cos2SigmaM*(-3+4*sinSigma*sinSigma)*(-3+4*cos2SigmaM*cos2SigmaM)));
    double s = b*A*(sigma-deltaSigma);

    double bearing = atan2(cosU2*sinLambda,  cosU1*sinU2-sinU1*cosU2*cosLambda);
    xOffset = sin(bearing)*s;
    yOffset = cos(bearing)*s;
}

Ответ 3

Я не стал бы слишком беспокоиться о кривизне земли. Я раньше не использовал openstreetmap, но у меня был быстрый взгляд, и кажется, что они используют проекцию Меркатора.

Это просто означает, что они сплющили планету на вас на прямоугольник, делая X пропорциональным долготе и Y почти точно пропорционально Широта.

Итак, вы можете пойти дальше и использовать простые формулы Mac, и вы будете очень точны. Ваша Локатор будет намного меньше, чем пиксель, который стоит для небольшой карты, с которой вы имеете дело. Даже на карте размером с Викторию вы получите ошибку только 2-3%.

diverscuba23 указал, что вам нужно выбрать эллипсоид... openstreetmap использует WGS84, а также большинство современных карт. Однако будьте осторожны, что многие карты в Австралии используют более старый AGD66, который может отличаться на 100-200 метров или около того.

Ответ 4

double testClass::getX(double lon, int width)
{
    // width is map width
    double x = fmod((width*(180+lon)/360), (width +(width/2)));

    return x;
}

double testClass::getY(double lat, int height, int width)
{
    // height and width are map height and width
    double PI = 3.14159265359;
    double latRad = lat*PI/180;

    // get y value
    double mercN = log(tan((PI/4)+(latRad/2)));
    double y     = (height/2)-(width*mercN/(2*PI));
    return y;
}