Как преобразовать геодезическое местоположение в позицию ECF, которая работает с моделью местности в цезии

Я пытаюсь поставить точку на вершине Эвереста в цезии. Моим наиболее вероятным кандидатом на вчерашний вечер был код, который я позаимствовал для преобразования геодезических данных в ecef (от PySatel.coord). После рассмотрения сегодня утром, это кажется правильным: 

a = 6378.137
b = 6356.7523142
esq = 6.69437999014 * 0.001
e1sq = 6.73949674228 * 0.001
f = 1 / 298.257223563


def geodetic2ecef(lat, lon, alt):
    """Convert geodetic coordinates to ECEF.

    Units are degrees and kilometers.
    """
    lat, lon = radians(lat), radians(lon)
    xi = sqrt(1 - esq * sin(lat))
    x = (a / xi + alt) * cos(lat) * cos(lon)
    y = (a / xi + alt) * cos(lat) * sin(lon)
    z = (a / xi * (1 - esq) + alt) * sin(lat)
    return x, y, z

Я вытащил lat / lon / alt для пика Mt. Эверест из Википедии. Я умножил координаты ECF, указанные в приведенном выше коде, на 1000 (м / км), прежде чем разместить объект в моем CZML. Я получаю местоположение ECF: [302995.41122130124, 5640733.98308375, 2981975.8695256836]. При использовании поставщика terrain по умолчанию (описанного в учебнике ) Эта точка значительно выше пика Mt. Эверест.

Вот соответствующий фрагмент CZML: 

{"position": 
  {"cartesian": [302995.41122130124, 5640733.98308375, 2981975.8695256836]}, 
 "id": "ellipsoid-1", 
 "ellipsoid": 
   {
     "radii": {"cartesian": [3545.5375159540376, 
                              164.44985193756034, 
                              164.62702908803794]}, 
     "material": {"solidColor": {"color": {"rgba": [0, 255, 0, 100]}}}
   }, 
 "orientation": {"unitQuaternion": [0.00014107125875577922, 
                                    -0.011462389405915903, 
                                    -0.010254110199791062, 
                                    -0.70702315200093502]}
}
2 3
geospatial coordinate-systems

2 ответа:

Здесь действует несколько факторов.

Во-первых, исходные данные, которые цезий использует для рельефа местности, могут иметь более низкую, чем ожидалось, высоту для вершины Эвереста. Мы используем набор данных CGIAR SRTM, поэтому этот пункт в их FAQ актуален:

Почему в некоторых горных регионах пики значительно ниже, чем они должны быть?

Как упоминалось ранее, многие исходные пустоты данных сосредоточены в горных районах и в заснеженных районах. Следовательно, множество вершин в высокогорных районах они фактически интерполируются. Без использования переменной высокого разрешения для интерполяции интерполяция не может определить, что пустота данных на самом деле является пиком, и стремится "сгладить" пик, что приводит к недооценке истинной высоты для этой области. Эта проблема в значительной степени решена в версии 4.
Они говорят, что он в значительной степени решен в версии v4, которую использует цезий, так что, надеюсь, этот первый фактор не является реальной проблемой.

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

В-третьих, Википедия дает высоту как высоту над средним уровнем моря (MSL). MSL-это сложная поверхность, с которой трудно работать математически, поэтому ваш geodetic2ecef не делает этого. Вместо этого он, как и цезий, предполагает, что высота относительно эллипсоида WGS84, который является гораздо более приятной поверхностью чтобы работать с ними.

У NGA есть веб-сайт, на котором можно найти высоту MSL над эллипсоидом WGS84, также известным как высота геоида: http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm96/intpt.html

Он сообщает, что для вершины Эвереста (27° 59' 17" N, 86° 55' 31" E) MSL составляет 28,73 метра ниже WGS84. Если вычесть это число из высоты пика, о котором сообщается в Википедии, вы должны, по крайней мере, приблизиться.

Эта страница имеет информация о программном вычислении высот геоида: http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm96/egm96.html

Я рекомендую интерполировать по 15-минутному файлу высоты геоида вместо вычисления высот из коэффициентов.

Пара других замечаний, не имеющих прямого отношения к вопросу:

  • цезий имеет код для преобразования LLA (мы называем его картографическим) в декартово. Смотрите Ellipsoid.cartographicToCartesian.
  • вы можете указать координаты в CZML в cartographicDegrees или cartographicRadians вместо картезианского, и тогда цезий сделает преобразование для вас автоматически. Однако вам все равно придется подстраиваться под геоид при указании высоты. Кроме того, не забывайте, что долгота-это первое.
6

Proj4js - порт знаменитой библиотеки proj4-может сделать эту работу за вас.

1