QMapControl и яндекс карты

[ipkis]

Захотелось отображать яндекс карты в QMapControl, для этого создал yandexMapAdapter

Код: Выделить всё

    YandexMapAdapter::YandexMapAdapter(   )
    : TileMapAdapter("vec01.maps.yandex.net", "/tiles?l=map&v=3.340.0&x=%2&y=%3&z=%1&lang=ru_RU", 256, 0, 21)
{}

Все прекрасно отображается, однако координаты по Y имеют какое-то смещение, один и тот же дом имеет разные координаты при использовании имеющегося googleMapAdapter и моего yandexMapAdapter. Скорее всего нужно переопределить эти функции:

Код: Выделить всё

 QPoint TileMapAdapter::coordinateToDisplay(const QPointF& coordinate) const
{
    qreal x = (coordinate.x()+180) * (mNumberOfTiles*mTileSize)/360.; // coord to pixel!
    qreal y = (1-(log(tan(PI/4+deg_rad(coordinate.y())/2)) /PI)) /2  * (mNumberOfTiles*mTileSize);

    return QPoint(int(x), int(y));
}

QPointF TileMapAdapter::displayToCoordinate(const QPoint& point) const
{
    qreal longitude = (point.x()*(360/(mNumberOfTiles*mTileSize)))-180;
    qreal latitude = rad_deg(atan(sinh((1-point.y()*(2/(mNumberOfTiles*mTileSize)))*PI)));

    return QPointF(longitude, latitude);

}

Но какие формулы нужно использовать для яндекс карт? может кто-то уже это делал, если не трудно подскажите куда копать или дайте ссылку на формулы.

[freeExec]

Яндекс основан на эллипсоиде, а не на сфере как у гугла.

[ipkis]

так, получается долгота совпадает, подскажите формулу для вычисления широты

[trir]

в ГОСТ'е посмотри

[Игорь Белов]

В общепринятых координатах зависимость y(lat):

Код: Выделить всё

y = a * (atanh(sin(lat)) - e * atanh(e * sin(lat)))

Здесь параметры a, e — соответственно экваториальный радиус эллипсоида и его эксцентриситет.
Радиус a далее сокращается при вычислении тайловых/пиксельных координат.
Эксцентриситет же разный:

• сфера Гугла и др.: e = 0
• эллипсоид Яндекса: e = 0.08181919084262149433480245175385667

При e ≠ 0 обратная функция lat(y) в лоб не считается. Либо итеративный поиск корня, либо разложение в тригонометрический ряд с последующим обращением. С первым, надеюсь, сами справитесь. Второе выложу, если ну прям очень надо.

[ipkis]

с первым справлюсь, а второе ну очень надо Спасибо за помощь!

[Игорь Белов]

Второе: тригонометрический ряд

Вычислим широту на сфере χ (хи):

Код: Выделить всё

chi = asin(tanh(y/a))

Широта на сфере χ

chi.png (642 байт) 11503 просмотра

Широта на эллипсоиде φ в виде тригонометрического ряда по χ:

Широта на эллипсоиде φ

phi.png (9.45 КБ) 11503 просмотра

Здесь e — эксцентриситет. Квадрат эксцентриситета для эллипсоида WGS 84:
e² = 0.006694379990141316996137233540044784

[stout]

Вот формула Боуринга для перехода от изометрической широты к геодезической

isometric2Geodic.png (47.73 КБ) 11443 просмотра

А вот график ошибки этой формулы в секундах дуги (φ(true) - φ(approx))*20624.806247

isometric2GeodicERROR.png (11.18 КБ) 11443 просмотра

Вертикальная шкала в единицах 10 знака после запятой.

[stout]

Скорее всего для Игорь Белов, это элементарщина, но мне кажется надо добавить, что сумму тригонометрического ряда удобно находить по алгоритму Кленшоу

Clenshaw formula (sin).png (36.11 КБ) 11370 просмотров

[Игорь Белов]

t = 2 cos 2χ ?

[stout]

t = 2 cos 2χ ?

Да, конечно! Вы правы.
В качестве частичного искупления оставлю здесь ссылку на неплохую статью (на мой взгляд)
SOME APPLICATIONS OF CLENSHAW'S RECURRENCE FORMULA IN MAP PROJECTIONS
в качестве добавки в этой статье описан алгоритм proj_etmerc.c из PROJ

[ipkis]

Я выложил решение тут https://ru.stackoverflow.com/questions/ ... 417#972417