vasya @quas read-only
Пользователь
2 июня 2014 в 16:06

Разработка → Геодезия: по полю на «питоне» из песочницы

Qt*, Python*
image

Доброго времени, Хабр!
Немного истории: в ходе учебы программированию, искал я себе реальную задачу, да такую чтобы с пользой. Нашел. Увидел как знакомый геодезист, на работе, считает объем земельного участка. Очень долго и нудно…

Геодезический расчет объемов:
При возведении жилых сооружений, высокотехнологичных помещений, автомобильных и железных дорог, а так – же в целях определения объемов строительных материалов и подсчета объема земляных работ, требуется помощь геодезистов. Они “отстреливают” территорию, разбивая всю площадь на так называемую геосетку, далее полученые точки из прибора выгружаются в autoCAD и высчитывают объем всей территории. Ниже пример геосетки:

image

Каждая точка имеет свои координаты (x,y,z-высота) относительно балтийского моря. По порядку берутся 4 точки как показано на предыдущем изображении и вычисляется исходя из координат этих точек объем участка земли, так вычисляются все объемы получившихся на сетке фигур, плюсуются и в итоге мы имеем объем всей территории.

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

Несколько дней я был заморочен на алгоритме поиска соседей, в результате вышло не совсем то чего я ждал, алгоритм работал как хотел в связи с тем, что сетка не ровная как в теории, а кривая как ходил помощник геодезиста по полю и ставил отметки, затем я понял что и AutoCAD очень просто сортирует по времени создания точки. Далее немного включившись в AutoCAD я заметил интересный порядок выгрузки точек в XML файл, и нашел возможность назначения каждой точки ID, теперь мы назначили каждой точке имя с левого угла до правого края и так-же понимаясь на ряд выше про именовали все точки (это куда быстрее, чем два дня считать объем территории площадью 500х700 метров). Ниже пример такого готового файла с именованными по порядку точками:

image

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

Плюсую все эти объемы и получаю за 0.2 секунды результат к которому он идет пару дней, проверяли на трех проектах пока данные сходятся, программой даже более точнее получается объем. Продолжаем ее тестировать. Теперь на днях я решил прикрутить этот код к юзабельному интерфейсу, с PyQt4, мой первый опыт работы с написанием графики, ниже я опубликую код который относится только к GUI.

from PyQt4 import QtGui
import sys
import cvgLeicaXmlReader
import cvgMath

class myWindow(QtGui.QMainWindow):

    def __init__(self, parent=None):
        QtGui.QWidget.__init__(self, parent)
        self.setWindowTitle('CVG2014')
        self.setFixedSize(350, 350)
        self.setWindowIcon(QtGui.QIcon('static/Icon.png'))
        self.setStyleSheet("QMainWindow {background-image: url(static/background.png);}")
        self.directory = ''
        self.file = ''

        self.labelFilename = QtGui.QLabel('Select .XML file with points', self)
        self.labelFilename.setFixedWidth(300)
        self.labelFilename.setFixedHeight(25)
        self.labelFilename.move(10, 5)
        self.labelFilename.setStyleSheet("QLabel { background-color: white; \
                                           border: 1px solid grey; \
                                           color: grey;}")

        self.SB_WidthOfXAxis = QtGui.QSpinBox(self)
        self.SB_WidthOfXAxis.move(10, 35)
        self.SB_WidthOfXAxis.setFixedWidth(50)
        self.SB_WidthOfXAxis.setMaximum(9999)

        self.labelPointsWidth = QtGui.QLabel('Length of points on X axis', self)
        self.labelPointsWidth.setFixedWidth(250)
        self.labelPointsWidth.setFixedHeight(25)
        self.labelPointsWidth.move(65, 47.5)

        self.SB_HeightAboveSeaLevel = QtGui.QDoubleSpinBox(self)
        self.SB_HeightAboveSeaLevel.move(10, 75)
        self.SB_HeightAboveSeaLevel.setFixedWidth(50)
        self.SB_HeightAboveSeaLevel.setRange(-9999.99, 9999.99)

        self.labelPointsSeaLevel = QtGui.QLabel('Height above sea level', self)
        self.labelPointsSeaLevel.setFixedWidth(250)
        self.labelPointsSeaLevel.setFixedHeight(25)
        self.labelPointsSeaLevel.move(65, 87)

        self.buttonOpenFile = QtGui.QPushButton('...', self)
        self.buttonOpenFile.setFixedWidth(30)
        self.buttonOpenFile.setFixedHeight(27)
        self.buttonOpenFile.move(311, 4)
        self.buttonOpenFile.clicked.connect(self.getXmlFile)

        self.buttonGetVolume = QtGui.QPushButton('RUN', self)
        self.buttonGetVolume.setFixedWidth(52)
        self.buttonGetVolume.setFixedHeight(35)
        self.buttonGetVolume.move(9, 115)
        self.buttonGetVolume.clicked.connect(self.getVolume)

        self.showVolume = QtGui.QLabel('Get Vol', self)
        self.showVolume.setFixedWidth(140)
        self.showVolume.setFixedHeight(32)
        self.showVolume.move(72, 117)
        self.showVolume.setStyleSheet("QLabel { background-color: white; \
                                           border: 1px solid grey; \
                                           color: grey;}")


    def getVolume(self):
        xml_file = self.getFileName()
        points = cvgLeicaXmlReader.getPointsFromXmlFile(xml_file)
        if(xml_file and points):


            # length_of_points = len(points)
            QUANTITY_POINTS_AT_X_AXIS = self.SB_WidthOfXAxis.value()
            STATIC_HEIGHT = self.SB_HeightAboveSeaLevel.value()

            rows = cvgLeicaXmlReader.getRowsFromPoints(points, QUANTITY_POINTS_AT_X_AXIS)
            quads = cvgLeicaXmlReader.getAllQuads(rows)

            volumes = []
            if (STATIC_HEIGHT or QUANTITY_POINTS_AT_X_AXIS) != 0:
                for quadrangle in quads:
                    Quadrangle_type = cvgMath.getTypeQuadrangle(quadrangle)
                    v = cvgMath.getVolumeQuadrangle(quadrangle, Quadrangle_type, STATIC_HEIGHT)
                    volumes.append(v)
            else:
                volumes = 0

            volumes = 0 if STATIC_HEIGHT == 0 else (round(sum(volumes), 3))
            result = '-'+str(volumes) if STATIC_HEIGHT < 0 else str(volumes)
            result = '0' if result == '-0' else result
            self.showVolume.setStyleSheet("QLabel { background-color: white; \
                                           border: 1px solid grey; \
                                           color: grey;}")
            self.showVolume.setText(result)
        else:
            self.showVolume.setText('Select the correct file!')
            self.showVolume.setStyleSheet("QLabel { background-color: white; \
                                           border: 1px solid grey; \
                                           color: red; \
                                           font-weight: bold}")
    
    def getFileName(self):
        return self.file

    def getXmlFile(self):
        sender = self.sender()
        path = QtGui.QFileDialog.getOpenFileName(sender, 'Open Xml file with points', self.directory, 'XML *.xml')
        fileName = path[path.rfind('/')+1:]
        self.directory = path[:path.rfind('/')]
        if(len(path) > 54):
            start = len(path)-54
            pathSlice = path[start:]
            pathSlice = pathSlice[pathSlice.find('/'):]
            pathSlice = '..'+pathSlice
        else:
            pathSlice = path
        self.labelFilename.setText(pathSlice)
        print(len(path))
        self.file = path


if __name__ == '__main__':
    app = QtGui.QApplication(sys.argv)
    win = myWindow()
    win.show()
    app.exec_()


С остальными модулями кому интересно, можно ознакомиться на гитхабе
Так же я загрузил .rar архив со скомпилированным кодом в exe для windows.

image

Что требуется от нас: открываем XML файл с точками (два файла для примера лежат так-же в на гите (2 корректных и один битый, для теста)), указываем количество точек по ширине поля, ниже указываем высоту относительно уровня моря, и ждем кнопку RUN, получаем объем в поле с право от кнопки запуска.

Если честно я не знаю, можно ли применять этот инструмент на производстве. Но свое желание я удовлетворил. Всем спасибо за внимание.

Файл с 65 точками (ширина 13 точек)
Файл с 78 точками (ширина 13 точек)
Модуль чтения XML файла
Модуль с математическими вычислениями
Модуль визуального оформления
Архив со скомпилированным кодом в Exe
vasya @quas
карма
1,0
рейтинг 0,0
Пользователь
Реклама помогает поддерживать и развивать наши сервисы

Подробнее
Реклама

Самое читаемое Разработка

Комментарии (12)

  • –2
    Для вычисления площади хитрых фигур метод Монте-Карло должен отлично подойти.
    • +2
      Монте-Карло имеет смысл для многомерных интергралов, т.к. в задачах высоких размерностей измельчение шага интегрированиея в 10 раз в каждом измерении приведет к росту объема расчетов в 10^n раз, где n — размерность пространства.
      Для таких простых задач из треугольничков, считать лучше обычными квадратурами.
  • 0
    Я так и не понял какой именно объем считается. каждого стобика, ограниченного четырёхугольником, от уровня моря?
    Другими словами, объем всего куска земли находящегося под данной поверхностью и над уровнем моря?

    Думаю, проверка на квадрат, прямоугольник — это лишняя и бессмысленная оптимизация, учитывая что формулы Пифагора и Герона в принципе одинаковой сложности, а так же количество прямоугольников до смешного минимально для современных компьютеров.

    Что бы не мучиться с поиском соседей, можно скормить этот массив точек алгоритму триангуляции Делоне, потом просто отбросить все треугольники со сшиком длинными сторонами (которые могут появиться по краям участка, т.к. триангуляция всегда выпуклая фигура)
    • 0
      И, мне кажется, ссылка на термин «геосетка» немного о другом.
      • 0
        С геосеткой все в порядке, именно ее предварительно используют на поле, и затем по ней на пересечениях отстреливают точки. А в нашем, конкретном, случае ей не пользуются, но в итоге после измерений получается подобие именно такой сетки. Я не силен в геодезии, извините.
    • 0
      Да, именно, объем всего куска земли, и относительно уровня моря, под либо над, но сперва расчитывается объем как Вы сказали каждого столбика, видимо для экономии расходного матерьяла, а если визуализировать все эти столбики в пространстве, вся картина будет похожа на «эквалайзер» (я это вижу так), а в дальнейшем объемы всех столбцов плюсуются.

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

      С треангуляцией Делоне, когда мучался с алгоритмом поиска соседей, столкнулся, но на первый взгляд она меня испугала, не разобрался я, и как раз в тот период я нашел более простой выход для себя (именование точек), и пошел по более легкому пути, который сохранил мне время, но думаю я еще вернусь к Делоне.
  • +1
    А почему объем на скриншоте отрицательный? Не сводите геодезистов с ума :)
    • 0
      Отрицательное значение в объеме (которое допустимо, если об этом речь) лишь показатель того, что область находиться ниже уровня моря. Изначально для вычислений геодезисту дается статическая высота, которую я и использую в расчетах в своем коде, и соответственно она может быть и отрицательной, отсюда и отрицательный объем на выходе.
    • 0
      может это Израиль?
      • 0
        это Казахстан.
      • 0
        Тогда уже Австралия, где вода закручивается в другую сторону :)
  • 0
    Теперь на днях я решил прикрутить этот код к юзабельному интерфейсу, с PyQt4, мой первый опыт работы с написанием графики…
    self.labelFilename.move(10, 5)


    Изучите как работает QLayout в Qt, который позволит избежать ручного размещения компонентов по координатам и добавит опциональную «резиновость» интерфейсу.
    И заодно можно посмотреть QDesigner, чтобы понять как работает каждый из видов Layout, да и в целом я советую разрабатывать UI там, чтобы скрыть сотни строк шаблонного кода.

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