Трансформация MapInfo.prj в WKT и proj

[dab]

Для данных в формате Mapinfo (TAB) с системой координат

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

"МСК-31 зона 2", 8, 9999, 3, 23.57, -140.95, -79.8, 0, -0.35, -0.79, -0.22, 0, 7, 38.48333333333, 0, 1, 2250000, -5212900.56

утилита gdalsrsinfo предлагает

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

PROJCS["unnamed",
    GEOGCS["unnamed",
        DATUM["MIF 9999,3,23.57,-140.95,-79.8,0,-0.35,-0.79,-0.22,0",
            SPHEROID["Krassovsky",6378245,298.3]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",38.48333333333],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",2250000],
    PARAMETER["false_northing",-5212900.56],
    UNIT["Meter",1.0]]

А как на самом деле должны выглядеть эти параметры в формате OGC WKT?

[SergeyRyzhkov]

В смысле?
Датум по идее должен быть нормально описан

DATUM["МСК-31 зона 2",
SPHEROID["Krassowsky 1940",6378245,298.3,
AUTHORITY["EPSG","7024"]],
TOWGS84[23.57,-140.95,-79.8,0,-0.35,-0.79,-0.22]
],

[dab]

В смысле?
Датум по идее должен быть нормально описан ...

В том смысле, что после прочтения статьи непонятно как в формате OGC WKT представить систему с пользовательским датумом из Mapinfo (с 9999).
А утилита gdalsrsinfo выдаёт некорректный результат.
Вот я и хочу просто увидеть корректное описание в формате OGC WKT для

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

"МСК-31 зона 2", 8, 9999, 3, 23.57, -140.95, -79.8, 0, -0.35, -0.79, -0.22, 0, 7, 38.48333333333, 0, 1, 2250000, -5212900.56

[SergeyRyzhkov]

В статье же написано вроде, цитирую:
" ... Вам может встретиться описание СК, где в качестве номера датума будет стоять "9999" или подобный номер, которого нет в таблице ..."
Вопрос к статье? Или как представить датум в WKT? Или к gdalsrsinfo ?
Выше я привел как должно быть в WKT корректно описан датум, на основании записи СК в формате MapInfo

[dab]

... Вопрос к статье? Или как представить датум в WKT? Или к gdalsrsinfo ?

Вопрос: как представить в формате OGC WKT систему координат, имеющуюся в указанном ниже формате Mapinfo

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

"МСК-31 зона 2", 8, 9999, 3, 23.57, -140.95, -79.8, 0, -0.35, -0.79, -0.22, 0, 7, 38.48333333333, 0, 1, 2250000, -5212900.56

[SergeyRyzhkov]

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

PROJCS["МСК-31 зона 2",
    GEOGCS["MI Custum",
        DATUM["MI Custum",
                      SPHEROID["Krassowsky 1940",6378245,298.3,AUTHORITY["EPSG","7024"]],
                      TOWGS84[23.57,-140.95,-79.8,0,-0.35,-0.79,-0.22]
        ],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",38.48333333333],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",2250000],
    PARAMETER["false_northing",-5212900.56],
    UNIT["Meter",1.0]]

[dab]

Сергей, спасибо за ответы.
После некоторых манипуляций в QGIS, я остановился всё-таки на другом описании

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

PROJCS["unnamed",
    GEOGCS["Krassovsky, 1942",
        DATUM["unknown",
                SPHEROID["krass",6378245,298.3],
                TOWGS84[23.57,-140.95,-79.8,0,0.35,0.79,-0.22]
            ],
            PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",38.48333333333],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",2250000],
    PARAMETER["false_northing",-5212900.56],
    UNIT["Meter",1]]

Такое описание обеспечивает одинаковый результат по сравнению с пересчётом и в Mapinfo и в QGIS с использованием PROJ

[SergeyRyzhkov]

Вы правы, для QGIS надо так описать
TOWGS84[23.57,-140.95,-79.8,0,0.35,0.79,-0.22]

А для MapInfo, если бы она понимала WKT
TOWGS84[23.57,-140.95,-79.8,0,-0.35,-0.79,-0.22]

И дело не в правильности параметров или записи, а в том какой метод использует ПО для трансформации датумов (при 7-и параметрах ессно):

Position Vector или Coordinate Frame Rotation (пишу по памяти могу ошибиться)

[ericsson]

См. http://proj.maptools.org/gen_parms.html секцию "Datum transformation to WGS84".
Собственно, на сколько я понимаю, в QGIS используется http://epsg.io/9606-method а в Mapinfo http://epsg.io/9607-method

[Boris]

А давайте посмотрим, что думает Mapinfo по этому поводу:

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

Print CoordSysStringToWKT$("CoordSys Earth Projection 8, 9999, 3, 23.57, -140.95, -79.8, 0, -0.35, -0.79, -0.22, 0, 7, 38.48333333333, 0, 1, 2250000, -5212900.56")

выдает:

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

PPROJCS["_MI_0",
    GEOGCS["unnamed",
        DATUM["unnamed",
            SPHEROID["World Geodetic System of 1984, GEM 10C",6378137,298.257223563]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433]
    ],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",38.48333333333],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",2250000],
    PARAMETER["false_northing",-5212900.56],
    UNIT["METER",1]
]

а для стандартной проекции 1001 эллипсоид указан верно, но +TOWGS=, все равно не осилил:

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

Print CoordSysStringToWKT$("CoordSys Earth Projection 8, 1001, 7, 38.48333333333, 0, 1, 2250000, -5212900.56")

получаем:

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

PROJCS["_MI_1",
    GEOGCS["Pulkovo 1942",
        DATUM["Pulkovo_1942",
            SPHEROID["Krassovsky - 1940/1948",6378245,298.300000000038,
                AUTHORITY["EPSG","7024"]],
            AUTHORITY["EPSG","6284"]
        ],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4284"]
    ],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",38.48333333333],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",2250000],
    PARAMETER["false_northing",-5212900.56],
    UNIT["METER",1]
]

[Донецков]

Вы правы, для QGIS надо так описать
TOWGS84[23.57,-140.95,-79.8,0,0.35,0.79,-0.22]
А для MapInfo, если бы она понимала WKT
TOWGS84[23.57,-140.95,-79.8,0,-0.35,-0.79,-0.22]

тогда вопрос, что нужно указывать при работе непосредственно с библиотекой proj.dll? (модуль который работает с этой библиотекой)
по моему должно быть: TOWGS84[23.57,-140.95,-79.8,0,0.35,0.79,-0.22]
проверял на известных данных, если добавить "-", то ошибки увеличиваются...
хотя встречается и такое:

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

+proj=tmerc +lat_0=0 +lon_0=39 +k=1 +x_0=7500000 +y_0=0 +ellps=krass +towgs84=23.57,-140.95,-79.8,0,-0.35,-0.79,-0.22 +units=m +no_defs 

например тут

[SergeyRyzhkov]

PROJ.4 для трансформации датумов использует метод Position Vector (EPSG 9606)

MapInfo для трансформации датумов использует метод Coordinate frame rotation (EPSG 9607)

По ссылке выше явно указан метод http://epsg.io/9607-method , поэтому знаки с минусом

[Донецков]

Точно... указано, что

Method: Coordinate Frame Rotation (geog2D domain)

, и странно, что для всех приведенных параметров для РФ используется только он..., а например для Азербайджана есть и

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

Method: Position Vector transformation (geog2D domain)

[surh_rus]

Вы правы, для QGIS надо так описать
TOWGS84[23.57,-140.95,-79.8,0,0.35,0.79,-0.22]

А для MapInfo, если бы она понимала WKT
TOWGS84[23.57,-140.95,-79.8,0,-0.35,-0.79,-0.22]

И дело не в правильности параметров или записи, а в том какой метод использует ПО для трансформации датумов (при 7-и параметрах ессно):

В файле res4 указано
Зона 1, 1.5 градусная
+proj=tmerc +ellps=krass +towgs84=24,-123,-94,0.02,-0.25,-0.13,1.1 +units=m +lon_0=66.08333333 +lat_0=0 +k_0=1 +x_0=1500000 +y_0=-6000000
Вношу проекцию в QGIS. Все таки, с + или - эти параметры должны быть?

[trir]

исходя из:

PROJ.4 для трансформации датумов использует метод Position Vector (EPSG 9606)

MapInfo для трансформации датумов использует метод Coordinate frame rotation (EPSG 9607)

коэффициенты вращения нужно уможать на -1 при переходе от MapInfo к PROJ.4

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

res = Ellipsoid.toProj & " +towgs84=" & X & "," & Y & "," & Z & "," & (-1 * RotationX) & "," & (-1 * RotationY) & "," & (-1 * RotationZ) & "," & Scale_ppm