Вытащить данные с ЕМГИС Тольятти (новый вопрос)

ArcGIS 8.x,9.x,10.x (Arcview, ArcEditor, Arcinfo).
Ответить
Аватара пользователя
romanshuvalov
Новоприбывший
Сообщения: 8
Зарегистрирован: 13 авг 2016, 10:20
Репутация: 3

Вытащить данные с ЕМГИС Тольятти (новый вопрос)

Сообщение romanshuvalov » 06 апр 2018, 18:39

Актуальный вопрос см. ниже

---

Добрый день.

Данные с ЕМГИС Тольятти теперь можно использовать в OSM (подробности).

Для удобного использования OSM-мапперами нужно загрузить подложку в JOSM. Но изображения хранятся в совершенно непонятном виде с неизвестной проекцией.

Вот что я вижу на сервере:
http://tech.mfc63.ru/arcgis/rest/servic ... MapServer/

Обратите внимание на Spatial Reference:

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

PROJCS["Local Conditional",GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Local"],PARAMETER["False_Easting",27478.878],PARAMETER["False_Northing",19463.17],PARAMETER["Scale_Factor",1.0],PARAMETER["Azimuth",-1.272472],PARAMETER["Longitude_Of_Center",49.3933174],PARAMETER["Latitude_Of_Center",53.4992761],UNIT["Meter",1.0]] 
Внизу есть кнопочка "Export map", но в той форме невозможно указать нужную мне проекцию, она всегда сбрасывается на "свою" локальную. Есть мнение, что данные вообще не привязаны к географическим координатам, этого мы не знаем.

Я не без помощи товарищей с телеграм-канала @RUOSM написал скрипт, который высчитывает bbox для каждого тайла, отправляет post-запрос и вытаскивает картинку. Однако, во-первых, немного не угадал с проекцией, а во-вторых, что более важно, система координат имеет ПОВОРОТ. Поэтому даже если я внесу правильные координаты bbox-а, он мне все равно выдаст повернутое изображение.

Что подскажете?

Вариант "подогнать все вручную" оставляем только на крайний случай.
Последний раз редактировалось romanshuvalov 09 апр 2018, 16:08, всего редактировалось 2 раза.

Monstria
Активный участник
Сообщения: 133
Зарегистрирован: 17 май 2011, 06:22
Репутация: 50
Откуда: Нижний Новгород

Re: Вытащить данные с ЕМГИС Тольятти

Сообщение Monstria » 06 апр 2018, 19:18

Если знаком с NET - вот мои функции для работы с OSM тайлами

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

Imports System.Drawing

Public Class Projection

    Public Shared Function LatLonToTile(lat As Double, lon As Double, zoom As Integer) As Point
        Dim x As Single = ((lon + 180.0) / 360.0 * (1 << zoom))
        Dim y As Single = ((1.0 - Math.Log(Math.Tan(lat * Math.PI / 180.0) + 1.0 / Math.Cos(lat * Math.PI / 180.0)) / Math.PI) / 2.0 * (1 << zoom))
        Return New Point(x, y)
    End Function
    Public Shared Function LatLonToTile(lat As Double, lon As Double, zoom As Integer, fixed As Boolean) As Point
        Dim x As Single = ((lon + 180.0) / 360.0 * (1 << zoom))
        Dim y As Single = ((1.0 - Math.Log(Math.Tan(lat * Math.PI / 180.0) + 1.0 / Math.Cos(lat * Math.PI / 180.0)) / Math.PI) / 2.0 * (1 << zoom))
        Return New Point(Math.Floor(x), Math.Floor(y))
    End Function


    ''' <summary>
    ''' Расчет точки координат в пикселях
    ''' </summary>
    ''' <param name="lat"></param>
    ''' <param name="lon"></param>
    ''' <param name="zoom"></param>
    ''' <returns></returns>
    ''' <remarks></remarks>
    Shared Function LatLonToPixelF(lat As Double, lon As Double, zoom As Double) As PointF
        Dim p As New PointF
        Dim centrPoint As Double = Math.Pow(2, zoom + 7)
        Dim totalPixel As Double = 2 * centrPoint
        Dim pixelsPerLngDegree As Double = totalPixel / 360
        Dim pixelsPerLngRadian As Double = totalPixel / (2 * Math.PI)
        Dim siny As Double = Math.Min(Math.Max(Math.Sin(lat * (Math.PI / 180)), -0.9999), 0.9999)
        p.X = Math.Round(centrPoint + lon * pixelsPerLngDegree)
        p.Y = Math.Round(centrPoint - 0.5 * Math.Log((1 + siny) / (1 - siny)) * pixelsPerLngRadian)
        Return p
    End Function

    ''' <summary>
    ''' Пиксель карты в координаты Lat,Lon
    ''' </summary>
    ''' <param name="px"></param>
    ''' <param name="py"></param>
    ''' <param name="zoom"></param>
    ''' <returns>Lat,Lon</returns>
    ''' <remarks></remarks>
    Shared Function PixelsToLatLon(px As Double, py As Double, zoom As Integer) As PointF
        Dim r() As Double = PixelsToMeters(px, py, zoom)
        r = MetersToLatLon(r(0), r(1))
        Dim res As New PointF(r(0), r(1))
        Return res
    End Function

    Shared Function MetersToLatLon(x As Double, y As Double) As Double()
        Dim lon As Double = (x / originShift) * 180
        Dim lat As Double = (y / originShift) * 180
        lat = 180 / Math.PI * (2 * Math.Atan(Math.Exp(lat * Math.PI / 180)) - Math.PI / 2)
        Return {lat, lon}
    End Function


    Shared Function PixelsToMeters(px As Double, py As Double, zoom As Integer) As Double()
        Dim res = Resolution(zoom)
        Dim x As Double = px * res - originShift
        Dim y As Double = Math.Abs(py * res - originShift)
        Return {x, y}
    End Function

    Shared Function Resolution(zoom As Integer) As Double
        Return initialResolution / Math.Pow(2, zoom)
    End Function

    Const tileSize = 256.0
    Const initialResolution = 2 * Math.PI * 6378137 / tileSize
    Const originShift = 2 * Math.PI * 6378137 / 2

End Class

Monstria
Активный участник
Сообщения: 133
Зарегистрирован: 17 май 2011, 06:22
Репутация: 50
Откуда: Нижний Новгород

Re: Вытащить данные с ЕМГИС Тольятти

Сообщение Monstria » 06 апр 2018, 19:25

Добавлю... для получения тайлов - строка запроса:

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

String.Format("https://{0}.tile.openstreetmap.org/{1}/{2}/{3}.png", serverType(severId), zoom, x, y)
Где serverType:

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

Private Shared serverType As New List(Of Char) From {"a", "b", "c"}

Аватара пользователя
romanshuvalov
Новоприбывший
Сообщения: 8
Зарегистрирован: 13 авг 2016, 10:20
Репутация: 3

Re: Вытащить данные с ЕМГИС Тольятти

Сообщение romanshuvalov » 06 апр 2018, 19:46

Monstria писал(а):
06 апр 2018, 19:18
Если знаком с NET - вот мои функции для работы с OSM тайлами
Я точно то же самое написал на Си. В текущей задаче оно не помогает (вернее, я пользовался этой утилитой, чтобы получить координты тайла, но дальше начинаются проблемы - см. первый пост).

P.S. Аркгиса у меня нет, работать с ним не умею и вообще не знаю, бывает ли он для линукса. Надеялся обойтись найденной формой загрузки и всякими gdal-утилитами.

Аватара пользователя
romanshuvalov
Новоприбывший
Сообщения: 8
Зарегистрирован: 13 авг 2016, 10:20
Репутация: 3

Re: Вытащить данные с ЕМГИС Тольятти

Сообщение romanshuvalov » 07 апр 2018, 00:46

Проблему решил.

Решение сложное и далеко не идеальное: вытаскивал фрагменты карты, из JSON-ответа вытаскивал реальные координаты bbox-а (они не совпадали с запрашиваемыми), затем через gdal_translate и gdalwarp получал нужный тайл 1024х1024, затем через тот же gdalwarp нарезал его на 16 256х256 уже готовых тайлов. Повторять в цикле до посинения.

Методы сглаживания gdalwarp давали дичайшие артефакты, пришлось обходиться без фильтра (nearest). Качество пострадало, но с таким подходом по-другому никак.

Аватара пользователя
romanshuvalov
Новоприбывший
Сообщения: 8
Зарегистрирован: 13 авг 2016, 10:20
Репутация: 3

Re: Вытащить данные с ЕМГИС Тольятти (готово)

Сообщение romanshuvalov » 09 апр 2018, 16:07

Возможно ли заставить сделать полный вывод геометрии, без "more..."?

http://tech.mfc63.ru/arcgis/rest/servic ... lse&f=html

Аватара пользователя
Denis Rykov
Гуру
Сообщения: 3376
Зарегистрирован: 11 апр 2008, 21:09
Репутация: 529
Ваше звание: Author
Контактная информация:

Re: Вытащить данные с ЕМГИС Тольятти (новый вопрос)

Сообщение Denis Rykov » 09 апр 2018, 16:13

Запросить в JSON?
Spatial is now, more than ever, just another column- The Geometry Column.

Аватара пользователя
romanshuvalov
Новоприбывший
Сообщения: 8
Зарегистрирован: 13 авг 2016, 10:20
Репутация: 3

Re: Вытащить данные с ЕМГИС Тольятти (новый вопрос)

Сообщение romanshuvalov » 11 апр 2018, 01:42

Хм, точно. Туплю. Спасибо.

Ответить

Вернуться в «ArcGIS»

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

Сейчас этот форум просматривают: нет зарегистрированных пользователей и 5 гостей