Помогите с привязкой карты по трём точкам

[flex]

Здравствуйте..

У меня следующий вопрос. Необходимо произвести привязку карты по трём точкам.

Точки имеют формат,

X,Y (экранные координаты) + широта, долгота

Линии широты параллельны друг другу, линии долготы также. Известно, что отсканированная карта может быть повёрнута, на какой-то угол.

Подскажите, как можно выяснить, повёрнута ли карта?

Слышал, что по трём точкам, используя аффинное преобразование, можно вычислить коэффициенты привязки, но как найти наклон карты, я так и не выяснил.

Зная наклон, сама привязка не вызывает труда.

До того, как использовать полиномиальное преобразование, хотелось бы сначала развернуть карту, чтобы линии широты и долготы шли параллельно границам экрана.

Правильно ли я понимаю проблему, или карту разворачивать (координатные оси) и вовсе не нужно?

Заранее благодарен.

[Максим Дубинин]

Если у вас есть 3 точки - используйте их, аффинное преобразование включает поворот, отдельно его вычислять не надо.

[flex]

Если у вас есть 3 точки - используйте их, аффинное преобразование включает поворот, отдельно его вычислять не надо.

Спасибо за ответ.

Полином какого порядка лучше брать (точек может быть и больше), если есть экранная координата и географическое положение?

Если я правильно понимаю, то для полиномиального преобразования действительна эта формула:

x0 = a0 + a1x + a2y

y0 = b0 + b1x + b2y

Если в эту формулу вместо коэффициентов a подставить экранные координаты по оси X, а в b по оси Y, то куда пойдут широта и долгота?

Если у меня есть три точки,

X = 10, Y = 10, Longitude = -150, Latitude = 80
X = 200, Y = 200, Lon = -170, Lat = -80
X = 1000, Y = 100, Lon = 140, Lat = 50

то как будет записано уравнение?

Далее по заданным X, Y скажем X = 50, Y = 50 должны быть получены соответсвующие значения широты и долготы.

[Максим Дубинин]

Полином 1-го порядка = аффинное преобразование.

Полином какого порядка лучше брать (точек может быть и больше), если есть экранная координата и географическое положение?

Вопрос странновато сформулирован, "экранная" как вы ее называете координата и ее реальное положение есть всегда, неважно какое преобразование вы используете, полиномиальное или аффинное.

То какой полином брать зависит от:
1. что вам нужно получить в резлультате (то есть насколько сильно допустимо деформирование растра), смысловое ограничение
2. сколько опорных точек в наличие, практическое ограничение

Вы читали вот это?
http://gis-lab.info/qa/polynom.html

[flex]

Полином 1-го порядка = аффинное преобразование.

Вопрос странновато сформулирован, "экранная" как вы ее называете координата и ее реальное положение есть всегда, неважно какое преобразование вы используете, полиномиальное или аффинное.

То какой полином брать зависит от:
1. что вам нужно получить в резлультате (то есть насколько сильно допустимо деформирование растра), смысловое ограничение
2. сколько опорных точек в наличие, практическое ограничение

Вы читали вот это?
http://gis-lab.info/qa/polynom.html

Да, я прочитал эту статью. Не совсем понятно, как связать положение на карте (пиксел) c широтой и долготой.

Есть выборка точек, скажем:

X1 = 10, Y1 = 15 Lon = -150, Lat = -80
Y2 = 20, Y2 = 23 Lon = -170, Lat = 60
Y3 = 100, Y3 = 100 Lon = 140, Lat = 50

И по ним нужно взять полином. Но как это связать с формулой в статье, не совсем понятно. Т.е. далее, если у меня

X4 = 80, Y4 = 70 какие будут широта и долгота?

Отвечая на вопросы:

1. Получить в результате широту, долготу.
Т.е. мышь движеться по карте и передаёт экранные координаты, которые должны быть переведены в широту и долготу. Карта может быть повёрнута. Но известно, что линии широты и долготы не деформированы.
2. Опорных точек может быть "бесконечно" много, но для этого я хотел брать среднее арифметическое.

[flex]

http://gis-lab.info/qa/polynom.html

Кстати не совсем ясно из статьи, как вычисляется коэффициет k.

Т.е. по логике а должно быть a0, a1, a2 - для полинома первого порядка.

Но если на первом шаге итерации

i = 0, j = 0

k = 0

то на втором шаге

i = 1, j = 0

k = (i*i+j)/2 + j = 1/2, т.е. никак 1

на третьем шаге

i = 1, j = 1

k = (1+1)/2+1 = 2

Почему индекс второго шага не правильный?

0.5 != 1

Или я что-то не так делаю?

[Максим Дубинин]

Есть выборка точек, скажем:

X1 = 10, Y1 = 15 Lon = -150, Lat = -80
Y2 = 20, Y2 = 23 Lon = -170, Lat = 60
Y3 = 100, Y3 = 100 Lon = 140, Lat = 50

Надо решить систему уравнений и найти все 6 коэффициентов:
x1' = a*x1 + b*y1 + c
y1' = d*x1 + e*y1 + f
x2' = a*x2 + b*y2 + c
y2' = d*x2 + e*y2 + f
x3' = a*x3 + b*y3 + c
y3' = d*x3 + e*y3 + f

Отвечая на вопросы:

Это были не вопросы, это были соображения от чего зависит выбор полинома.

[Максим Дубинин]

Почему индекс второго шага не правильный?

Похоже на мелкий глюк унаследованный из эрдасовской документации, вообще обычно aij индекс используют для бинарных полиномов, вот например: http://mathworld.wolfram.com/BivariatePolynomial.html

[flex]

Почему индекс второго шага не правильный?

Похоже на мелкий глюк унаследованный из эрдасовской документации, вообще обычно aij индекс используют для бинарных полиномов, вот например: http://mathworld.wolfram.com/BivariatePolynomial.html

Спасибо за ответ.

Вы не видели кстати, есть ли программная реализация полиномов различных порядков.

Т.к. если я хочу использовать скажем полином пятого порядка, то не так то просто это запрограммировать.

[Максим Дубинин]

5 полином врядли пригоден для привязки материалов, слишком большие искажения.

По поводу алгоритмов, я бы для начала посмотрел GRASS или GDAL, они с исходниками на C и тот и тот - GPL, надо только найти где именно прописана трансформация.

[Максим Дубинин]

на случай если не нашли

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

grass-6.2.1\imagery\i.rectify\crs.c

необходимые процедуры:
exactdet - точек ровно столько сколько нужно
calcls - точек больше чем нужно

единственно не уверен, годятся ли они для любого порядка, это надо внимательно разбираться

[flex]

на случай если не нашли

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

grass-6.2.1\imagery\i.rectify\crs.c

необходимые процедуры:
exactdet - точек ровно столько сколько нужно
calcls - точек больше чем нужно

единственно не уверен, годятся ли они для любого порядка, это надо внимательно разбираться

Огромное спасибо. Я посмотрю.

Ещё один вопрос. Если у меня вот такой полином

a0 + a1*x + a2*y + a3*x^2 + a4*x*y + a5*y^2 = 15

a0' + a1' * x + a2' * y + a3' * x^2 + a4' *x*y + a5'*y^2 = 25

и коэффициенты я нашёл. Как найти теперь x и y.

Просто если я по координатам x, y могу выдавать широту и долготу, то как сделать обратную операцию. Т.е. как в моем случае, широта долгота известны и равны 15, 25, все коэффициенты известны, но x и y нужно найти.

Еслиб полином был первой степени, то вот такое уравнение

a0 + a1*x + a2*y = 15
a0' + a1'*x + a2'*y = 25

я конечно решить смогу.

[Максим Дубинин]

Боюсь четко ответить на вопрос я не смогу, математику забыл, внесу "пищу" для размышления. Просто решать систему уравнений - не очень эффективно, как я понимаю, допустим преобразование аффинное и вы нашли коэффициенты a0,a1,a2,b0,b1,b2

X = a0+a1x+a2y
Y = b0+b1x+b2y

однако, если вы хотите делать обратный пересчет, то есть, зная X и Y вычислить x и y, то коэффициенты будут ДРУГИМИ!

Пример прямого преобразования:
X = 50 + 2x + 0y
Y = 50 + 0x + 2y

А вот обратное будет уже:
x = -25 + 0.5X + 0Y
y = -25 + 0X + 0.5Y

Это не мешает, конечно уравнения в лоб как вы предлагаете. Соответственно с полиномами второй степени, надо применять тот же подход, матричную алгебру, перевычислять коэффициенты и заново решать уравнения с новыми коэффициентами. Как именно это сделать сразу не скажу, надо порыться в книжках. Для аффинного преобразования формулы для вычисления обратного преобразования можно найти здесь: http://www.posc.org/Epicentre.2_2/DataM ... _cs35.html
Там же есть и для полиномиального, но не очень понятные, надо разбираться.

[flex]

Боюсь четко ответить на вопрос я не смогу, математику забыл, внесу "пищу" для размышления. Просто решать систему уравнений - не очень эффективно, как я понимаю, допустим преобразование аффинное и вы нашли коэффициенты a0,a1,a2,b0,b1,b2

X = a0+a1x+a2y
Y = b0+b1x+b2y

однако, если вы хотите делать обратный пересчет, то есть, зная X и Y вычислить x и y, то коэффициенты будут ДРУГИМИ!

Пример прямого преобразования:
X = 50 + 2x + 0y
Y = 50 + 0x + 2y

А вот обратное будет уже:
x = -25 + 0.5X + 0Y
y = -25 + 0X + 0.5Y

Это не мешает, конечно уравнения в лоб как вы предлагаете. Соответственно с полиномами второй степени, надо применять тот же подход, матричную алгебру, перевычислять коэффициенты и заново решать уравнения с новыми коэффициентами. Как именно это сделать сразу не скажу, надо порыться в книжках. Для аффинного преобразования формулы для вычисления обратного преобразования можно найти здесь: http://www.posc.org/Epicentre.2_2/DataM ... _cs35.html
Там же есть и для полиномиального, но не очень понятные, надо разбираться.

Спасибо за ответ.

Если честно, я не смог найти хороших различных источников, где бы про это можно было бы почитать. Поэтому пытаюсь делать в лоб. Спасибо за ссылку. Я её посмотрел, буду с этим разбираться. Если есть какие-то интересные книги на примете, вообще про преобразование, полиномы или что-то связаное с ГИС и калибровкой, то я был бы очень признателен за любую информацию.

Ещё одно решение в лоб.

У нас есть X, Y и широта, долгота.

Вычисляем сначала коэффициенты для X, Y -> широта долгота, и коэффициенты для пары широта, долгота -> X, Y. Т.е. две независимых друг от друга системы. У нас же есть табличные значения.

Проблема - нет зависимости между переменными, и при незначительных отклонениях (неправильно задана одна из точек), может быть так, что
если X = 10, Y = 20 -> широта = 20 и долгота = 10
То пытаясь выяснить X, Y по широте, долготе получим
широта = 20 и долгота = 10 -> X = 123, Y = 830
Что конечно не приемлимо. Поэтому если решать, как было выше указано, т.е. для первого порядка и сохранить зависимость, т.е. выражать x, y через широту долготу и обратно,

x = -25 + 0.5X + 0Y
y = -25 + 0X + 0.5Y

получим

X = 10, Y = 20 -> широта = 20 и долгота = 10
широта = 20 и долгота = 10 -> X = 9.9, Y = 20.01

Прошу на цифры не обращать внимания, это просто пример. Но в реальности так примерно и получилось.

[Максим Дубинин]

Честно говоря, последние ваши соображения из последнего сообщения я не очень понял.

Если честно, я не смог найти хороших различных источников, где бы про это можно было бы почитать.

Их нет, я думаю, для этого и был создан этот сайт, создать такие источники В связи с этим, прочитайте, пожалуйста, все-таки личное сообщение.

По сути: потратив недюжинное количество времени и нервных клеток, я таки вывел алгоритм вычисления коэффициентов для полиномов 1-го порядка (аффинное преобразование) и 2-го порядков, для порогового количества точек, причем как прямое, так и обратное. Все это - матричная алгебра, я проверил в R, результат прекрасно совпал в тем набором коэф-тов, который вычисляет например ERDAS.

Для того, чтобы реализовать этот подход программно, вам надо найти алгоритмы вычисления двух штук: инвертирования матрицы и умножения матриц. Если это проблем для вас не составит, то все остальное - детали. Если готовы в это окунуться - пишите на почту.