Добавление местной координатной системы в GIS

Материал из GIS-Lab
Перейти к навигации Перейти к поиску
Эта страница является черновиком статьи.


Введение

Под местной системой координат (МСК) будет подразумеваться так называемая «городская» система, построенная независимо от государственной системы (ГСК) и включенная в неё заданием ключей перехода к СК-42 или СК-63. МСК крупных территорий, сравнимых с размерами субъектов Федерации, не являются предметом данной статьи, поскольку относятся к классическим картографическим проекциям.

Постановка задачи

Некоторые программы позволяют реализовать работу в МСК непосредственно. Так, в MapInfo любая проекция может быть дополнена аффинным преобразованием. ArcGIS в качестве МСК предлагает локальную проекцию: аналог ортометрической проекции на эллипсоиде, дополненный разворотом и масштабированием.

Другие программы, включая QGIS, работают только с «обыкновенными» проекциями, не допуская дополнительных геометрических преобразований. Задача статьи — показать, как сконструировать проекцию, позволяющую работать в МСК в таких средах, как QGIS или, скажем, бортовой софт приёмников GARMIN.

В качестве рабочей среды будем использовать командную строку UNIX. Это идеальный инструмент для экспериментирования, позволяющий непринуждённо сочетать PROJ.4, findkey и утилиты обработки текстовых потоков awk, pr.[1]

Подготовка данных

Тестовый пример основан на сгенерированных данных. В таблице X, Y — координаты пункта в государственной системе (ГСК), а именно в седьмой зоне СК-42, x, y — координаты в местной системе (МСК), p — вес пункта.

ID X Y x y p
1 7383477.64 6087377.60 1334.71 285.94 1.0
2 7382557.14 6081916.51 563.67 −5197.34 1.0
3 7386610.19 6088160.39 4444.27 1153.79 1.0
4 7381962.05 6090016.34 −252.07 2881.90 1.0

Подготовим файл исходных данных data0.dat:

1 7383477.64 6087377.60 1334.71   285.94 1.0
2 7382557.14 6081916.51  563.67 -5197.34 1.0
3 7386610.19 6088160.39 4444.27  1153.79 1.0
4 7381962.05 6090016.34 -252.07  2881.90 1.0

Убедимся в однородности данных обоих каталогов. Для этого вычислим параметры конформного преобразования:

$ findkey data0.dat key0.dat var0.dat

Сами по себе параметры нас не интересуют. Гораздо интереснее невязки, выведенные в файл var0.dat:

1 … -0.002 -0.001
2 … -0.017 0.013
3 … 0.031 0.017
4 … -0.012 -0.029

Малые значения невязок говорят о том, что можно начинать работать.

Пересчитаем координаты из ГСК в долготу/широту:

$ awk '{print $2, " ", $3}' data0.dat | proj -I -f "%.9f" +proj=tmerc +lat_0=0 +lon_0=39 +k=1 +x_0=7500000 +y_0=0 +ellps=krass > pt_longlat.dat

Команда awk извлекает из файла data0.dat вторую и третью колонки и передаёт их утилите proj, которая пересчитывает координаты из седьмой зоны СК-42 в долготу/широту и записывает результаты в файл pt_longlat.dat. Вот содержимое pt_longlat.dat:

37.183749701 54.896961118
37.171630976 54.847714659
37.232243035 54.904709276
37.159058387 54.920296878
  1. Команды примера воспроизводятся в среде MinGW, что входит в состав QGIS под MS Windows. Две особенности:
    • нужно заменить команду awk на gawk, или лучше в системе определить awk через gawk;
    • после некоторых команд (proj, pr) придётся добавить в пайп команду удаления лишних символов CR:
      | tr -d '\r'