[DEM]

Собсно по мотивам темы
Перевод формата STL в твердотельную геометрию (Rhino, Salome и другие программы)
Решено сделать программу, т.к. та сетка которая получалась в Salome, делалась довольно таки долго и с непредсказуемыми последствиями.
Файлы в Salome загружались очень долго и т.д.
В качестве программы в которой будет создаваться модель, принят FreeCad соответственно программирование будет на Python.
Сразу оговорюсь, Python принят в качестве языка программирования т.к. потом к нему собираемся через SWIG подключить библиотеку SCAD API на C++.
На данный момент столкнулся с задачкой которую в лоб решить не могу.
У меня есть список отрезков вот такого вида

Код:

[Выделить все]

1я плоскость ['6/41', '6/40', '4/35', '30/31', '5/36', '1/34', '4/30', '5/35', '32/33', '31/32', '33/34', '1/41', '37/38', '38/39', '36/37', '39/40']
2я плоскость
['19/86', '19/87', '102/103', '15/71', '15/70', '95/96', '85/86', '4/104', '3/29', '61/62', '82/83', '84/85', '83/84', '92/93', '65/66', '30/31', '3/101', '90/91', '72/73', '75/76', '71/72', '28/29', '22/99', '32/33', '54/55', '74/75', '24/99', '24/98', '22/100', '21/97', '79/80', '21/100', '73/74', '25/26', '76/77', '2/67', '10/101', '10/102', '59/60', '16/77', '16/78', '1/34', '13/56', '26/27', '1/54', '4/30', '18/82', '18/81', '14/96', '78/79', '103/104', '87/88', '57/58', '58/59', '2/25', '27/28', '88/89', '93/94', '11/66', '94/95', '69/70', '13/57', '55/56', '23/97', '17/80', '17/81', '23/98', '62/63', '20/92', '20/91', '68/69', '14/68', '60/61', '64/65', '11/67', '89/90', '33/34', '12/63', '12/64', '31/32']
3я плоскость
['10/52', '7/42', '10/53', '44/45', '43/44', '42/43', '45/46', '7/53', '48/49', '50/51', '8/47', '8/46', '9/48', '9/47', '51/52', '49/50']

Собсно надо получить списки замкнутых контуров с обходом точек по часовой стрелке....
Типа [19/87, 87/88, 88/89 (и т.д.)]

Координаты точек я так понимаю особо то и не потребуются....
Либо алгоритм нужен, либо код на Paython 2.7
Ну иль просто подсказка в какую сторону копать.....

[VVA]

Сообщение от DEM Собсно надо получить списки замкнутых контуров с обходом точек по часовой стрелке

В качестве идеи Convex Hull с использованием Graham Scan Algorithm



Вроде как Phyton есть готовые решения
graham scan algorithm python

[DEM]

Координаты вообще не хочется использовать....
Тут аналитически можно определить.
Использовать что нибудь типа Граф....

[trir]

1. Строим граф
2. Получаем список вершин, у которых >2 связей
3. Для всех рёбер из этих вершин находим кратчайшее расстояние между вершинами (не считая само ребро конечно)
4. ???
5. PROFIT!

[DEM]

Сообщение от trir 1. Строим граф 2. Получаем список вершин, у которых >2 связей 3. Для всех рёбер из этих вершин находим кратчайшее расстояние между вершинами (не считая само ребро конечно) 4. ??? 5. PROFIT!

В виде кода можно?

[Brandashmыg]

1. А почему 2.7?
2. В Shapely есть возможность создать полигон по точкам и потом получать от него всё что угодно, в том числе координаты наружного контура.

Код:

[Выделить все]

from shapely.geometry import Polygon
polygon = Polygon([(0, 0), (1, 1), (1, 0)])
list(polygon.exterior.coords) #собственно координаты контура

А ещё проще каждый полигон при импорте сделать полигоном, а потом их просто объединить (a.union(b))

Почему лучше так, а не повторением алгоритма обхода точек самостоятельно - библиотеки в питоне обычно быстрее, т.к. "узкие" участки переписаны на С++.

[DEM]

Brandashmыg
1. Потому как Freecad на 2.7 а для построения и редактирования схем я его хочу использовать.
2. Да тут как бы можно и без библиотеки обойтись, ведь задача чисто аналитическая.

[Sleekka]

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

Нужно использовать именно массив отрезков или полигонов.

[Brandashmыg]

Почему аналитическая? Выше пример привели - всё через координаты точек. Сами отрезки не нужны. Но тогда может какую-либо область "срезать". Поэтому и проще идти через объединение полигонов - тогда и отверстия и "выгрызы" останутся.

[DEM]

Ну вот я граф нарисовал, и там все довольно простг получается...

[Brandashmыg]

По количеству рёбер? У 62 узла и у 310 одинаковое количество отрезков.

[Sleekka]

Сообщение от DEM Ну вот я граф нарисовал, и там все довольно простг получается...

Граф рисуется по точкам, но не по отрезкам, следовательно достоверность совпадения графа с контуром не 100%-я вопрос лишь когда и на каких данных вылезут различия.

[DEM]

Brandashmыg
Так просто чтобы сделать как ты показал, придется координаты преобразовывать в 2d....
Sleekka
http://pythontutor.ru/lessons/graphs/
http://pythontutor.ru/lessons/dfs/
Реализация графа в ребрах представлена...

[Brandashmыg]

DEM
Shapely поддерживает 3д
The Polygon constructor takes two positional parameters. The first is an ordered sequence of (x, y[, z])

[DEM]

Сообщение от Brandashmыg По количеству рёбер? У 62 узла и у 310 одинаковое количество отрезков.

Чего-то не понял тебя....
Представленно 3 плоскости....
В них различные контуры представлены.
Часть узлов уже удалена, т.к. они не нужны.
Это те узлы которые находятся внутри контура....

----- добавлено через ~9 мин. -----
Зы. А по части Shapely ты по моему не прав....
Там все таки надо по порядку задавать грани...
Но там вроде как есть создание полигона по точкам...
Да и в чистом пайтоне я встречал алгоритм для 2d

[Brandashmыg]

Сообщение от DEM Там все таки надо по порядку задавать грани...

Каждый КЭ СРАЗУ после чтения загоняешь в полигон, потом делаешь для всех union и получаешь обобщённый полигон плиты.
А для чего это всё? Опалубку вытаскивать?

-----------------
Пример:
https://github.com/kuvbur/asf-reinfo...aster/parse.py
класс asf функция get_element

[DEM]

Brandashmыg
Да хочу получить геометрическую модель а не сеточную.
Чтобы потом проще править было.
Сетку потом буду делать в gmsh, или перекидывать в dxf и в Робот загонять.
Вообще-то идея создание своего препроцессора и постпроцессора.

Потом можно будет и армирование просматривать, в том же Freecade и сразу распечатку делать.
Много чего задумано....

----- добавлено через ~8 ч. -----

Сообщение от Brandashmыg 1. А почему 2.7? Код: [Выделить все] from shapely.geometry import Polygon polygon = Polygon([(0, 0), (1, 1), (1, 0)]) list(polygon.exterior.coords) #собственно координаты контура

Так сперва то надо создать полигон....
Надо определить последовательность точек...

[Brandashmыg]

Для четырёх узлового КЭ порядок точек будет не важен. Посмотри код на гитхабе - там просто список скармливаешь ему и всё.

[trir]

у меня проблемы с п.3 - не могу добиться нужного результата...

[Brandashmыg]

Сообщение от trir у меня проблемы с п.3 - не могу добиться нужного результата...

Сообщение от trir 3. Для всех рёбер из этих вершин находим кратчайшее расстояние между вершинами (не считая само ребро конечно)

А почему самое короткое ребро обязательно внешнее?

[trir]

кратчайшее расстояние между вершинами

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

[DEM]

Brandashmыg
Ты гнешь свою линию и всё...
Я говорю, что у меня гораздо быстрее сейчас все отрабатывается и безо всяких библиотек..
Сейчас делаю поиск пути по аналогу с графами, все это будет отрабатываться гораздо быстрее чем графические библиотеки.
У меня же не четырех угольники, у меня контур уже найден, его надо отсортировать по верштнам.
trir
В п. 3 должно получится 3 полигона....
Я сейчас взял за аналог работу с графами.
Обхожу смежные вершины.
Возможны и вот такие варианты(см рисунок)
Тогда надо создать 2 полигона....

[Brandashmыg]

Сообщение от DEM Brandashmыg Ты гнешь свою линию и всё...

Гну) Буду очень признателен, если покажете итоговый код для нахождения контура с помощью графов.

[DEM]

Brandashmыg
Думаю сегодня получу....

[trir]

В п. 3 должно получится 3 полигона....

это просто у тебя subdivision surface, а в общем случае их может быть сколько угодно, а может и не быть

[DEM]

trir
Не замкнутых не будет.....
По крайней мере пока планирую так...


----- добавлено через ~2 ч. -----
У меня почему то с цифрами не срабатывает....
Ищу почему же.....

Код:

[Выделить все]

def find_all_paths(graph, start, end, path=[]):
    path = path + [start]
    if start == end:
        return [path]
    if not graph.has_key(start):
        return []

    paths = []
    for node in graph[start]:
        if node not in path:
            newpaths = find_all_paths(graph, node, end, path)
            for newpath in newpaths:
                paths.append(newpath)
    return paths


def find_path(graph, start, end, path=[]):
    path = path + [start]
    if start == end:
        return path
    if not graph.has_key(start):
        return None
    for node in graph[start]:
        if node not in path:
            newpath = find_path(graph, node, end, path)
            if newpath: return newpath
    return None
graf2={}
rr0=['6/41', '6/40', '4/35', '30/31', '5/36', '1/34', '4/30', '5/35', '32/33', '31/32', '33/34', '1/41', '37/38', '38/39', '36/37', '39/40']
rr=['19/86', '19/87', '102/103', '15/71', '15/70', '95/96', '85/86', '4/104', '3/29', '61/62', '82/83', '84/85', '83/84', '92/93', '65/66', '30/31', '3/101', '90/91', '72/73', '75/76', '71/72', '28/29', '22/99', '32/33', '54/55', '74/75', '24/99', '24/98', '22/100', '21/97', '79/80', '21/100', '73/74', '25/26', '76/77', '2/67', '10/101', '10/102', '59/60', '16/77', '16/78', '1/34', '13/56', '26/27', '1/54', '4/30', '18/82', '18/81', '14/96', '78/79', '103/104', '87/88', '57/58', '58/59', '2/25', '27/28', '88/89', '93/94', '11/66', '94/95', '69/70', '13/57', '55/56', '23/97', '17/80', '17/81', '23/98', '62/63', '20/92', '20/91', '68/69', '14/68', '60/61', '64/65', '11/67', '89/90', '33/34', '12/63', '12/64', '31/32']
for ff in rr:
    nnn = ff.split("/")
    if graf2.get(nnn[0])==None:
        graf2[nnn[0]]=[nnn[1]]
    else:
        graf2[nnn[0]] = graf2[nnn[0]]+[nnn[1]]
    if graf2.get(nnn[1]) == None:
        graf2[nnn[1]] = [nnn[0]]
    else:
        graf2[nnn[1]] = graf2[nnn[1]] + [nnn[0]]

print graf2
print find_all_paths(graf2, '24', '99')

[Brandashmыg]

Мой вариант. line = rezult[0].exterior - наружная граница.

Код:

[Выделить все]

import numpy as np
import shapely.geometry as sh
import shapely.ops as sho

class msh(object):
    def __init__(self, filename):
        self.__read_asf(filename)
    def __parse_string(self, string, start_string, end_string):
        out = string.split(start_string)[1].split(end_string)[0]
        out = out.strip().split("\n")[1:]
        return out
    def __get_index_node(self, node):
        """ 
        Индекс точки в массиве data по её номеру
        """
        index = int(np.where(self.nodes[:,[0]] == node)[0][0])
        return index
    def __get_index_element(self, element):
        """
        Индекс элемента в массиве data по его номеру
        """
        key = self.n_elements[element]
        index = int(np.where(self.elements[key][:,[0]] == element)[0][0])
        return key, index 
    def __get_coordinates_node(self,node):
        """
        Координаты точки по её номеру
        """
        index = self.__get_index_node(node)
        coordinates = self.nodes[index,1:]
        return coordinates
    def get_node(self,node):
        """
        Получить точку в виде объекта shapely
        """
        coordinates = self.__get_coordinates_node(node)
        return sh.Point(coordinates) 
    def get_list_element_by_type(self, type_element):
        element_list = []
        for element in self.elements[type_element]:
            element_list.append(element[0])    
        return element_list
    def get_element(self,element):
        """
        Полигон в виде объекта shapely по его номеру
        """
        key, index = self.__get_index_element(element)
        nodes = self.elements[key][index][3:]
        for i, node in enumerate(nodes):
            if i == 0:
                coordinates = self.__get_coordinates_node(node)
            else:
                coordinates = np.vstack((coordinates, self.__get_coordinates_node(node)))        
        if key == 1:
            elements = sh.LineString(coordinates)
        elif key == 2:
            elements = sh.Polygon(coordinates)
        elif key == 3:
            elements = sh.Polygon(coordinates)
        return elements
    def axis_surf(self, polygon):
        x = polygon.exterior.coords[0][0]
        y = polygon.exterior.coords[0][1]
        z = polygon.exterior.coords[0][2]
        xx = 1
        yy = 1
        zz = 1
        for coord in polygon.exterior.coords[1:]:
            if x != coord[0]:
                xx = 0
            if y != coord[1]:
                yy = 0
            if z != coord[2]:
                yy = 0
        if xx == 0:
            axis = 0
        elif yy == 0:
            axis = 1
        elif zz == 0:
            axis = 2
        return axis
    def __read_asf(self,filename):
        """
        Читает файл
        """
        msh = open(filename, 'r').read().split('$EndMeshFormat')[1]
        #Получаем два массива - с точками и элементами
        nodes_str = self.__parse_string(msh, '$Nodes', '$EndNodes')
        elements_str = self.__parse_string(msh, '$Elements', '$EndElements')
        
        #Создаём массив с точками
        self.nodes = np.fromstring(nodes_str[0], dtype=float, sep=' ')
        for node in nodes_str[1:]:
            a = np.fromstring(node, dtype=float, sep=' ')
            self.nodes =  np.vstack((self.nodes, a))
        #Теперь - элементы
        #Для каждого типа элемента пропишем количество вершин в словарь и рассортируем элементы по номеру
        self.elements = {1:2, 2:3, 3:4, 4:4, 5:8}
        #Итак далее, в примере только типы 1 и 3, на них пока и остановимся.
        self.n_elements = {} #словарь для быстрого поиска типа элемента по его номеру
        for i, element in enumerate(elements_str):
            a = np.fromstring(element, dtype=float, sep=' ')
            number_element = int(a[0]) 
            elm_type = int(a[1])
            self.n_elements[number_element] = elm_type
            if type(self.elements[elm_type]) == int:
               self.elements[elm_type] = a
            else:
                self.elements[elm_type] = np.vstack((self.elements[elm_type], a))
          
               
filename = "G://0555.msh"
test = msh(filename)
poly_list = test.get_list_element_by_type(3)
poly_x = []
for polygon in poly_list:
    poly = test.get_element(polygon)
    if test.axis_surf(poly) == 0:
        if poly.is_valid:
            poly_x.append(poly)
        else:
            print("err")
rezult = sho.cascaded_union(poly_x)
line = rezult[0].exterior

Причём решение честное - я не использовал элементы 1-го типа, предполагая, что они идёт по контуру. Только 3-го типа.

[DEM]

А на 20 тысячах элементов не зависнет?

[Brandashmыg]

Давай файл - проверим)

[DEM]

Brandashmыg
Так в СКАДе же сейчас любой файл можно сохранить в msh...

[DEM]

Сообщение от Brandashmыg Давай файл - проверим)

Ну вот, надо все контуры создать....
PS. Опять проблемы с shapely не получается поставить на Анаконду2 и на FreeCad.