Я работаю с набором данных граничных файлов FSA, предоставленных Статистическим управлением Канады. Данные и полные сведения доступны здесь, набор 2016 г.
Чтобы данные были совместимы с остальной частью того, над чем я работаю, мне нужно преобразовать из предоставленного формата в широты / долготы.
Подробности предоставленного формата:
- Проекция: конформный конус Ламберта
- Ложное восточное направление: 6200000.000000
- Ложное северное направление: 3000000.000000
- Центральный меридиан: -91.866667
- Стандартная параллель 1: 49.000000
- Стандартная параллель 2: 77.000000
- Широта происхождения: 63.390675
- Линейная единица: метр (1.000000)
- Базовая: Североамериканская 1983 (NAD83)
- Первичный меридиан: Гринвич
- Угловая единица: градус
- Сфероид: GRS 1980
Я использую rgeos :: gCentroid для получения центроида каждого FSA, поэтому в решении не нужно принимать во внимание объекты sf или sp.
Все FSA, представляющие наибольший интерес, - вседовольно далеко от центра широты / долготы 63,390674, -91,8666667, поэтому приближения, подходящие для коротких расстояний, не работают для меня.
Из решения, детализированного здесь Я адаптировал следующий код:
# Read Shape File
fsa <- sf::st_read('~/lfsa000b16a_e.shp')
# Compute Centroid per FSA and store in dataframe
fsa_centroids = gCentroid(as(fsa, 'Spatial'), byid = TRUE, id = fsa$CFSAUID)
fsa_df <- data.frame(fsa_centroids)
fsa_df['fsa'] = rownames(fsa_df)
# Earth Radius
R = 6371007.1810
# Central Lat/Lon
lat0 = 63.390675
lon0 = -91.866667
# Remove false Easting/Northing
fsa_df$x = fsa_df$x - 6200000
fsa_df$y = fsa_df$y - 3000000
# Sphere-approximation calculations
fsa_df$lat = asin(sin(lat0*pi/180)*cos(fsa_df$y/R) + cos(lat0*pi/180)*sin(fsa_df$y/R)*cos(0)) * 180/pi
fsa_df$lon = lon0 + atan2(sin(pi/2)*sin(fsa_df$x/R*cos(fsa_df$lat*pi/180)), cos(fsa_df$x/R)-sin(fsa_df$lat*pi/180)^2) * 180/pi
(извиненияс такими бессмысленными терминами, как cos(0), я хотел сохранить полную формулу и я был неуверен, какой подшипник был 0 рад, поэтому я возился)
Моя проблема двоякая: это решение основано на совершенной сфереПредположение и данные использует сфероид GRS 1980.Кроме того, поскольку мы пересекаем сфероид, имеет значение порядок, в котором вы идете на восток / запад и север / юг (хотя оба порядка имеют очень высокие ошибки в приближенном решении).
Есть ли более простой способ или библиотекиЯ могу использовать это преобразование?