Форумы GIS-Lab.info

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

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

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

Серьёзные люди, имея полные данные по исходным данным и измерениям в старой сети и в новой, используют коллокацию для совместного переуравнивания. По результатам получается поле деформаций вместе с оценками погрешностей. В нашем (то есть в Вашем) случае можно только симулировать эту деятельность.

Создайте файл разностей долгот и широт.
Постройте грид разностей долгот b1 и грид разностей широт b2. Предпочтительно использовать для этого функции Грина для сплайнов на сфере.
Добавьте нулевые гриды погрешностей той же размерности b3 и b4.
Создайте виртуальный четырёхканальный растр:
Код: Выделить всё
```
gdalbuildvrt -separate germany.vrt b1.nc b2.nc b3.nc b4.nc
```

Конвертируйте его в формат 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).

Во-первых: спасибо за ответ!Есть ещё вопросы...

Игорь Белов писал(а): ↑20 фев 2020, 20:42Постройте грид разностей долгот b1 и грид разностей широт b2. Предпочтительно использовать для этого функции Грина для сплайнов на сфере.

Это можно сделать штатными инструментами QGIS? Нигде не могу найти такие модули.

Игорь Белов писал(а): ↑20 фев 2020, 20:42Добавьте нулевые гриды погрешностей той же размерности b3 и b4.

Т.е. просто создаю растр с нулевыми значениями пикселей в тех же размерах пикселей и в том же экстенте?

aspirant-1 писал(а): ↑02 мар 2020, 20:40Это можно сделать штатными инструментами QGIS?

Нет. Рекомендую модуль greenspline из пакета GMT.

aspirant-1 писал(а): ↑02 мар 2020, 20:40просто создаю растр с нулевыми значениями пикселей в тех же размерах пикселей и в том же экстенте?

Да.

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

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

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

После интерполяции результат явно не верный, и нет пространственной привязки. Либо я не правильно входной файл подготовил, либо код неверный.

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 (40.9 КБ) 9357 просмотров

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

Игорь Белов писал(а): ↑03 мар 2020, 07:48Рекомендую модуль greenspline

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

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

там ровно то, чего Вы просили, RTFM

gamm писал(а): ↑03 мар 2020, 17:04там ровно то, чего Вы просили, RTFM

Значит формат входных данных правильный? Ок, поэкспериментирую. Как я понял, там надо добавить ограничение на область интерполяции.

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

-Rxmin/xmax[/ymin/ymax]

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)".

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

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

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

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.

5P. S. Скорее всего, сеть неравномерная, где-то густо, где-то пусто. Перед использованием greenspline или surface не забывайте прогнать данные через модуль blockmean

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

Игорь Белов писал(а): ↑03 мар 2020, 18:51На выходе грид в замечательном формате netCDF.

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

2) 2D интерполяция идет в градусах, но на сфере? Из документации непонятно. В общем исходники нужно читать, как всегда. Там правда, брат

aspirant-1 писал(а): ↑03 мар 2020, 19:49что означает вот этот кусочек кода -I5m/3m ?

Шаг грида в угловых минутах 5′ по X и 3′ по Y:

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

-I5m/3m

Форумы GIS-Lab.info

NTv2

NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2

Re: NTv2