NTv2

Системы координат, проекции, преобразования, привязка
aspirant-1
Активный участник
Сообщения: 127
Зарегистрирован: 04 май 2012, 21:23
Репутация: 12
Откуда: Пенза

NTv2

Сообщение aspirant-1 » 20 фев 2020, 19:17

Здравствуйте! Вопрос по созданию файлов .gsb для трансформации по NTv2 (извините, может некорректно выразился).

Как создать файл трансформации? У меня есть данные по координатным парам в разных СК (около 100 пар пунктов ГГС). Хотелось бы создать сетку, а не использовать пресловутые 7 параметров.

Есть ли какое-либо, позволяющее быстро создать такой грид?

Аватара пользователя
Игорь Белов
Гуру
Сообщения: 2227
Зарегистрирован: 04 янв 2011, 22:00
Репутация: 1500
Откуда: Казань

Re: NTv2

Сообщение Игорь Белов » 20 фев 2020, 20:42

Серьёзные люди, имея полные данные по исходным данным и измерениям в старой сети и в новой, используют коллокацию для совместного переуравнивания. По результатам получается поле деформаций вместе с оценками погрешностей. В нашем (то есть в Вашем) случае можно только симулировать эту деятельность.
  1. Создайте файл разностей долгот и широт.
  2. Постройте грид разностей долгот b1 и грид разностей широт b2. Предпочтительно использовать для этого функции Грина для сплайнов на сфере.
  3. Добавьте нулевые гриды погрешностей той же размерности b3 и b4.
  4. Создайте виртуальный четырёхканальный растр:

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

    gdalbuildvrt -separate germany.vrt b1.nc b2.nc b3.nc b4.nc
  5. Конвертируйте его в формат NTv2:

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

    gdal_translate -of NTv2 germany.vrt germany.gsb -mo VERSION=NTv2.0 -mo SYSTEM_F=DHDN90 -mo SYSTEM_T=ETRS89 -mo MAJOR_F=6377397.155 -mo MINOR_F=6356078.963 -mo MAJOR_T=6378137 -mo MINOR_T=6356752.314 -mo SUB_NAME=DHDN90 -mo CREATED=06-11-09 -mo UPDATED=06-11-09
Пример имитирует создание известной трансформации DHDN to ETRS89 (8).
The purpose of computing is insight, not numbers

aspirant-1
Активный участник
Сообщения: 127
Зарегистрирован: 04 май 2012, 21:23
Репутация: 12
Откуда: Пенза

Re: NTv2

Сообщение aspirant-1 » 02 мар 2020, 20:40

Во-первых: спасибо за ответ!Есть ещё вопросы...
Игорь Белов писал(а):
20 фев 2020, 20:42
Постройте грид разностей долгот b1 и грид разностей широт b2. Предпочтительно использовать для этого функции Грина для сплайнов на сфере.
Это можно сделать штатными инструментами QGIS? Нигде не могу найти такие модули.
Игорь Белов писал(а):
20 фев 2020, 20:42
Добавьте нулевые гриды погрешностей той же размерности b3 и b4.
Т.е. просто создаю растр с нулевыми значениями пикселей в тех же размерах пикселей и в том же экстенте?

Аватара пользователя
Игорь Белов
Гуру
Сообщения: 2227
Зарегистрирован: 04 янв 2011, 22:00
Репутация: 1500
Откуда: Казань

Re: NTv2

Сообщение Игорь Белов » 03 мар 2020, 07:48

aspirant-1 писал(а):
02 мар 2020, 20:40
Это можно сделать штатными инструментами QGIS?
Нет. Рекомендую модуль greenspline из пакета GMT.
aspirant-1 писал(а):
02 мар 2020, 20:40
просто создаю растр с нулевыми значениями пикселей в тех же размерах пикселей и в том же экстенте?
Да.
The purpose of computing is insight, not numbers

aspirant-1
Активный участник
Сообщения: 127
Зарегистрирован: 04 май 2012, 21:23
Репутация: 12
Откуда: Пенза

Re: NTv2

Сообщение aspirant-1 » 03 мар 2020, 15:42

Скачал, установил. Все запускается. Создал файл вот такого вида

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

44.84362995277778 53.23528129166667 0.05757
44.86153118055556 53.12470514722222 0.05443
45.03488886111111 53.080609952777785 0.05417
45.29735655277778 53.23538368888889 0.05509
44.93713768333333 53.244917497222225 0.05810
"Долгота широта разница широт"

Вводил код из примеров, естественно названия файлов другие.

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

greenspline -V -Rg -Sp -D3 -I1 -GP1994.nc mag_obs_1990.txt
После интерполяции результат явно не верный, и нет пространственной привязки. Либо я не правильно входной файл подготовил, либо код неверный.

gamm
Гуру
Сообщения: 4044
Зарегистрирован: 15 окт 2010, 08:33
Репутация: 1049
Ваше звание: программист
Откуда: Казань

Re: NTv2

Сообщение gamm » 03 мар 2020, 16:30

aspirant-1 писал(а):
03 мар 2020, 15:42
После интерполяции результат явно не верный
Не очень понятно, об чем вы. В результате сетка, и отметки по осям (на всю сферу), все в hdf. Вы чего-то другого ожидали? я вашу табличку из сообщения отправил интерполировать ...

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

  group name       otype dclass       dim
0     /  lat H5I_DATASET  FLOAT       181
1     /  lon H5I_DATASET  FLOAT       361
2     /    z H5I_DATASET  FLOAT 361 x 181
grid.jpg
grid.jpg (40.9 КБ) 7166 просмотров

aspirant-1
Активный участник
Сообщения: 127
Зарегистрирован: 04 май 2012, 21:23
Репутация: 12
Откуда: Пенза

Re: NTv2

Сообщение aspirant-1 » 03 мар 2020, 16:44

Эмм...А это интерполировались именно значения разностей широт? Которые 0.05757...? Просто я тоже получил такую же картинку. Но, как мне показалась, в QGIS она загрузилась без координат. А оказывается результат интерполяции 5 точек почти на мир? Правильно ли я понимаю подписи и шкалы осей?

gamm
Гуру
Сообщения: 4044
Зарегистрирован: 15 окт 2010, 08:33
Репутация: 1049
Ваше звание: программист
Откуда: Казань

Re: NTv2

Сообщение gamm » 03 мар 2020, 17:02

Игорь Белов писал(а):
03 мар 2020, 07:48
Рекомендую модуль greenspline
Игорь, а форматом вывода можно управлять? Он то HDF5 выводит, то черте-что (пробовал область поменять, в первых байтах получил CDF вместо HDF, но читалка read.cdffile.list его не распознает)

gamm
Гуру
Сообщения: 4044
Зарегистрирован: 15 окт 2010, 08:33
Репутация: 1049
Ваше звание: программист
Откуда: Казань

Re: NTv2

Сообщение gamm » 03 мар 2020, 17:04

aspirant-1 писал(а):
03 мар 2020, 16:44
А оказывается результат интерполяции 5 точек почти на мир
там ровно то, чего Вы просили, RTFM

aspirant-1
Активный участник
Сообщения: 127
Зарегистрирован: 04 май 2012, 21:23
Репутация: 12
Откуда: Пенза

Re: NTv2

Сообщение aspirant-1 » 03 мар 2020, 17:10

gamm писал(а):
03 мар 2020, 17:04
там ровно то, чего Вы просили, RTFM
Значит формат входных данных правильный? Ок, поэкспериментирую. Как я понял, там надо добавить ограничение на область интерполяции.

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

-Rxmin/xmax[/ymin/ymax]

gamm
Гуру
Сообщения: 4044
Зарегистрирован: 15 окт 2010, 08:33
Репутация: 1049
Ваше звание: программист
Откуда: Казань

Re: NTv2

Сообщение gamm » 03 мар 2020, 17:43

aspirant-1 писал(а):
03 мар 2020, 17:10
Как я понял, там надо добавить ограничение на область интерполяции.
Ага, поскольку "The two shorthands -Rg and -Rd stand for global domain (0/360 and -180/+180 in longitude respectively, with -90/+90 in latitude)".

я попробовал менять область, тогда формат вывода поменялся зачем-то, я его не смог прочитать ...

Аватара пользователя
Игорь Белов
Гуру
Сообщения: 2227
Зарегистрирован: 04 янв 2011, 22:00
Репутация: 1500
Откуда: Казань

Re: NTv2

Сообщение Игорь Белов » 03 мар 2020, 18:51

Сделайте так:

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

gmt greenspline dlat.txt -R44.8/45.35/53/53.25 -I5m/3m -St0.2 -D4 -Gb2.nc
На выходе грид в замечательном формате netCDF. GDAL его знает, это главное.

P. S. Скорее всего, сеть неравномерная, где-то густо, где-то пусто. Перед использованием greenspline или surface не забывайте прогнать данные через модуль blockmean.
Вложения
gmt_map.png
gmt_map.png (1.42 КБ) 7108 просмотров
The purpose of computing is insight, not numbers

aspirant-1
Активный участник
Сообщения: 127
Зарегистрирован: 04 май 2012, 21:23
Репутация: 12
Откуда: Пенза

Re: NTv2

Сообщение aspirant-1 » 03 мар 2020, 19:49

5P. S. Скорее всего, сеть неравномерная, где-то густо, где-то пусто. Перед использованием greenspline или surface не забывайте прогнать данные через модуль blockmean
Да, Вы правы, сеть действительно неоднородна. Это связано с известными причинами - уничтоженные ГГС и т.п., сами представляете. Больше спасибо за консультацию. А что означает вот этот кусочек кода -I5m/3m ? Размер ячейки? NetCDF без проблем отображается в QGIS.

gamm
Гуру
Сообщения: 4044
Зарегистрирован: 15 окт 2010, 08:33
Репутация: 1049
Ваше звание: программист
Откуда: Казань

Re: NTv2

Сообщение gamm » 03 мар 2020, 20:06

Игорь Белов писал(а):
03 мар 2020, 18:51
На выходе грид в замечательном формате netCDF.
1) вопрос не в этом, вопрос в том, как управляется формат вывода, и на кой его менять. Единственное, что я нашел в доках - если задать сетку для интерполяции самому (-Nnet_file), то можно получить выхлоп в ASCII.

2) 2D интерполяция идет в градусах, но на сфере? Из документации непонятно. В общем исходники нужно читать, как всегда. Там правда, брат :mrgreen:
Последний раз редактировалось gamm 03 мар 2020, 20:41, всего редактировалось 1 раз.

Аватара пользователя
Игорь Белов
Гуру
Сообщения: 2227
Зарегистрирован: 04 янв 2011, 22:00
Репутация: 1500
Откуда: Казань

Re: NTv2

Сообщение Игорь Белов » 03 мар 2020, 20:11

aspirant-1 писал(а):
03 мар 2020, 19:49
что означает вот этот кусочек кода -I5m/3m ?
Шаг грида в угловых минутах 5′ по X и 3′ по Y:
The purpose of computing is insight, not numbers

Ответить

Вернуться в «Координаты и привязка»

Кто сейчас на конференции

Сейчас этот форум просматривают: Bing [Bot] и 9 гостей