Определение линейных размеров полигонов в ArcGis

ArcGIS 8.x,9.x,10.x (Arcview, ArcEditor, Arcinfo).
Борис
Интересующийся
Сообщения: 36
Зарегистрирован: 22 ноя 2005, 12:32
Репутация: 0

Определение линейных размеров полигонов в ArcGis

Сообщение Борис » 08 окт 2007, 11:43

Добрый день!

Есть полигональная шейп-тема.
Помогите, пожалуйста, определить размер полигонов - то есть минимальный (ширину) и максимальный (длину) размеры. Дело осложняется тем, что полигоны - не прямоугольной формы, а всякие изрезанные и пр.

AndreyL
Завсегдатай
Сообщения: 484
Зарегистрирован: 17 авг 2006, 14:04
Репутация: 0
Откуда: Новосибирск

Сообщение AndreyL » 14 ноя 2007, 17:07

Если тема еще интересует, то можно попробовать вписать каждый полигон в прямоугольник минимальной площади, повернутый на определенный угол. Угол найти минимизируя площадь прямоугольника.

Борис
Интересующийся
Сообщения: 36
Зарегистрирован: 22 ноя 2005, 12:32
Репутация: 0

Сообщение Борис » 17 ноя 2007, 23:55

Да, тема еще интересует (очень:).
AndreyL, поясните, пожалуйста, подробнее, как это сделать.
Лучше всего - в ArcGis, но есть и ArcView, и MapInfo.

AndreyL
Завсегдатай
Сообщения: 484
Зарегистрирован: 17 авг 2006, 14:04
Репутация: 0
Откуда: Новосибирск

Сообщение AndreyL » 18 ноя 2007, 16:40

Есть подозрение, что придется задействовать математические пакеты (или хотя бы Excel) – решение, похоже, численное. Идея то проста: поворот группы точек (всех вертексов полигона) вокруг центра координат на заданный угол сложности не представляет (для пущей страсти можно сначала параллельно перенести ближе к центру). Далее посчитать минимумы и максимумы по Х и по Y и площадь (или периметр, если хотите) тоже не сложно. Теперь осталось подобрать нужный угол поворота, а вот это задачка численная.

Борис
Интересующийся
Сообщения: 36
Зарегистрирован: 22 ноя 2005, 12:32
Репутация: 0

Сообщение Борис » 20 ноя 2007, 12:07

AndreyL, а не могли бы Вы бодробнее (конкретнее) описать последовательность действий. С общим ходом рассуждений согласен, но как все это вополтить... :roll: ???

AndreyL
Завсегдатай
Сообщения: 484
Зарегистрирован: 17 авг 2006, 14:04
Репутация: 0
Откуда: Новосибирск

Сообщение AndreyL » 20 ноя 2007, 13:29

Вот так примерно можно на Вольфрамовской Математике (точки, естественно, для примера):
points = Table[{Random[], Random[]*1.5}, {i, 25}];

ClearAll[x]
mat = {{Cos[x], Sin[x]}, {-Sin[x], Cos[x]}};

V = (Max[Transpose[points.mat][[1]]] -
Min[Transpose[points.mat][[1]]])*(Max[
Transpose[points.mat][[2]]] - Min[Transpose[points.mat][[2]]]) //
Simplify;
rez = FindMinimum[V, {x, 0, Pi}]

rect = {{Min[Transpose[points.mat][[1]]],
Min[Transpose[points.mat][[2]]]}, {Min[Transpose[points.mat][[1]]],
Max[Transpose[points.mat][[2]]]}, {Max[Transpose[points.mat][[1]]],
Max[Transpose[points.mat][[2]]]}, {Max[Transpose[points.mat][[1]]],
Min[Transpose[points.mat][[2]]]}, {Min[Transpose[points.mat][[1]]],
Min[Transpose[points.mat][[2]]]}} /. rez[[2]];
mat2 = {{Cos[-x], Sin[-x]}, {-Sin[-x], Cos[-x]}};
rrect = rect.mat2 /. rez[[2]];

gr1 = ListPlot[points, PlotStyle -> PointSize[0.02],
DisplayFunction -> Identity];
gr2 = ListPlot[rrect, PlotJoined -> True, DisplayFunction -> Identity];
Show[gr1, gr2, AspectRatio -> Automatic, DisplayFunction -> $DisplayFunction];

Борис
Интересующийся
Сообщения: 36
Зарегистрирован: 22 ноя 2005, 12:32
Репутация: 0

Сообщение Борис » 21 ноя 2007, 22:19

Спасибо за ответ :) , разбираюсь :? , а не могли бы Вы пояснить на конкретном примере. Последовательность действий такая: 1)есть полигон; 2)узлы полигона превращены в точечный шейп; 3)определены координаты всех узлов; 4)получены,например, следующие аттрибуты:

polygon vertice _____POINT_X__________ POINT_Y
1_____ 1_____ 545649,517515_____ 6023352,726300
1_____ 2_____ 721286,473581_____ 5968844,705450
1_____ 3_____ 642552,665689_____ 5955520,522580
1_____ 4_____ 577143,040672_____ 5980957,598980
1_____ 5_____ 522635,019824_____ 6001549,517960
1_____ 6_____ 545649,517515_____ 6023352,726300

как быть дальше? как на автомате прочесть их Математикой (ведь это только один полигон (колонка "polygon", а на самом деле их много). Принцип построения атрибут.таблицы простой: колонка номера полигона, далее номер узла в полигоне, далее координаты.

AndreyL
Завсегдатай
Сообщения: 484
Зарегистрирован: 17 авг 2006, 14:04
Репутация: 0
Откуда: Новосибирск

Сообщение AndreyL » 22 ноя 2007, 13:18

Конвертора DBase-таблиц в Математику я не знаю, но можно конвертировать таблицу в текстовый файл – такое Математика нормально читает. А потом в самой Математике разбить таблицу по признаку polygon:

dat = Import[FileName, "Table"];
list = Transpose[dat][[1]];
UnVal = {};
n = Length
  • ;
    Do[If[Position[UnVal, list[]] == {}, UnVal = Insert[UnVal, list[], -1]], {i, n}]

    PointsTable = Table[{}, {i, Length[UnVal]}];
    Do[k = Position[UnVal, list[]][[1, 1]];
    PointsTable[[k]] = Insert[PointsTable[[k]], {dat[[i, 3]], dat[[i, 4]]}, -1], {i, n}]
    PointsTable


    А дальше по вкусу – можно создать табличку с вертексами всех прямоугольников, табличку с длиной, шириной, площадью, периметром и углом поворота каждого прямоугольника и т.д.

geologic
Гуру
Сообщения: 852
Зарегистрирован: 15 сен 2005, 13:19
Репутация: 5
Откуда: москва
Контактная информация:

Сообщение geologic » 22 ноя 2007, 18:53

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

AndreyL
Завсегдатай
Сообщения: 484
Зарегистрирован: 17 авг 2006, 14:04
Репутация: 0
Откуда: Новосибирск

Сообщение AndreyL » 23 ноя 2007, 09:25

geologic писал(а):Подобрать угол поворота наверно будет сложно, первое что на ум приходит - тупо поворачивать с заданной дискретностью и проверять, макс-не макс.
Да упаси меня Господи от таких подвигов. Тогда уж проще выбросить ГИС вместе с компьютером, и решать задачу на бумажной карте, вооружившись рейсшиной и калькулятором. Мы же четыре поста назад соорудили код автоматического поиска угла. Только ГИС, как и любая информационная система, на серьезные подвиги неспособна, поэтому задействуется математическое программное обеспечение, хотя алгоритм оптимизации можно и на VBA написать (если задача не единичная и делать кучу конвертов ленивей, чем писать код).

geologic
Гуру
Сообщения: 852
Зарегистрирован: 15 сен 2005, 13:19
Репутация: 5
Откуда: москва
Контактная информация:

Сообщение geologic » 23 ноя 2007, 12:26

Осенило :) Определить длинную ось можно, найдя максимальное расстояние в группе точек. Это можно сделать простым запросом SQL в реляционной БД, к таблице с колонками X, Y, IDA (индекс "бывшего" полигона, признак группы). Группа, кстати, не обязательно должна быть бывшим полигоном, любое облако годится.

Определить азимут длинной оси можно, но механизм сложнее. Вычислить вторую ось как перпендикуляр... Здесь мысль тормозит пока :(

AndreyL
Завсегдатай
Сообщения: 484
Зарегистрирован: 17 авг 2006, 14:04
Репутация: 0
Откуда: Новосибирск

Сообщение AndreyL » 23 ноя 2007, 14:43

geologic писал(а):Осенило :) Определить длинную ось можно, найдя максимальное расстояние в группе точек.
В принципе Вы правы, если нужно найти прямоугольник с максимальной длиной. Можно поставить задачу найти прямоугольник с минимальной шириной. Можно искать прямоугольник с минимальной длиной или максимальной шириной (типа минимакса и максимина). Мы пока ищем прямоугольник с минимальной площадью. Это вопрос – что мы должны понимать под длиной и шириной. И кроме Бориса на него никто не ответит.

Борис! Я тут пока накидал нечто, проводящее расчеты после загрузки и разбора таблицы, надеюсь, разберетесь и подправите по Вашим требованиям:
mm = Length[PointsTable]
vertTab = Table[{}, {i, mm}];
rectTab = Table[{}, {i, mm}];

Do[
points = PointsTable[[mmi]];
rot = Transpose[points.mat];
V = (Max[rot[[1]]] - Min[rot[[1]]])*(Max[rot[[2]]] - Min[rot[[2]]]) //
Simplify;
rez = FindMinimum[V, {x, 0, Pi}];

rot = Transpose[points.mat] /. rez[[2]];
rect = {{Min[rot[[1]]], Min[rot[[2]]]}, {Min[rot[[1]]],
Max[rot[[2]]]}, {Max[rot[[1]]], Max[rot[[2]]]}, {Max[rot[[1]]],
Min[rot[[2]]]}};

rrect = rect.mat2 /. rez[[2]];

perim = 2*(Max[rot[[1]]] - Min[rot[[1]]] + Max[rot[[2]]] - Min[rot[[2]]]);
midPoint = {Sum[rrect[[i, 1]], {i, Length[rrect]}]/Length[rrect],
Sum[rrect[[i, 2]], {i, Length[rrect]}]/Length[rrect]};
angle =
If[(Max[rot[[2]]] - Min[rot[[2]]]) > (Max[rot[[1]]] - Min[rot[[1]]]),
x/Pi*180, x/Pi*180 + 90] /. rez[[2]];
angle = If[angle < 0, angle + 360, angle];
vertTab[[mmi]] =
Table[{UnVal[[mmi]], rrect[[i, 1]], rrect[[i, 2]]}, {i, Length[rrect]}];
rectTab[[mmi]] = {UnVal[[mmi]], midPoint[[1]], midPoint[[2]], rez[[1]],
perim, angle},
{mmi, mm}]


Кстати, этот код, похоже, универсальный в отношении того, что мы понимаем под длиной и шириной – все будет зависеть от целевой функции V. Например, для того решения, что предложил geologic целевая функция будет:
V = -(Max[rot[[1]]] - Min[rot[[1]]]) // Simplify;

geologic
Гуру
Сообщения: 852
Зарегистрирован: 15 сен 2005, 13:19
Репутация: 5
Откуда: москва
Контактная информация:

Сообщение geologic » 23 ноя 2007, 15:02

Ну, выборка одно дело, расчёть - другое. Задача ставилась рассчитать для всех длину и ширину. С длиной разобрались, что это просто - у вас свой подход сложился, у меня свой. А вот насчет ширины...

Я бы прежде поставил вопрос "что такое ширина", соответственно. В том смысле, что самому интересно услышать. задача и правда часто стоит, и не только для замкнутых контуров, измерить ширину, мол.

Борис
Интересующийся
Сообщения: 36
Зарегистрирован: 22 ноя 2005, 12:32
Репутация: 0

Сообщение Борис » 24 ноя 2007, 23:53

Всем спасибо, погружаюсь в мир вольфрама 8)
а про ширину --- да, ширина - это (в бытовом плане ( :? ) ) то, что перпенидикулярно длине. Задачка, думаю, и впрямь может быть довольно широкоупотребимая, так что пытаюсь сам и призываю всех желающих создать скрипт (VBA или с коннектом на ту же Математику)
:wink:

AndreyL
Завсегдатай
Сообщения: 484
Зарегистрирован: 17 авг 2006, 14:04
Репутация: 0
Откуда: Новосибирск

Сообщение AndreyL » 26 ноя 2007, 14:06

Борис писал(а):пытаюсь сам и призываю всех желающих создать скрипт (VBA или с коннектом на ту же Математику)
Коннект на Математику не нужен, если уж писать скрипт, то писать. Математика используется, когда задача единичная – решил и забыл. Или для исследования задачи, чтобы потом перенести это в тот же скрипт.
Если сделаете оболочку, то готов помочь с поиском оптимального значения, поскольку задачу я исследовал пока отвечал Вам.

Ответить

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