0

I'm working on plugin in QGIS which is creating 3d model from raster data. In that case, I need to create triangulation only in the vertical plane (no connecting vertical faces with neighbour points).

This is my code and example with the current effect.

    def delaunay_ogolny(self):
    import math
    X=[]
    Y=[]
    Z=[]
    for x in self.listqa:
        dada=x[0]
        X.append(dada)
        uu=x[1]
        Y.append(uu)
        dede=x[2]
        Z.append(dede)
    u = np.array(X)
    v = np.array(Y)
    x = u
    y = v
    z = np.array(Z)
    #Delaunay
    tri = Delaunay(np.array([u, v]).T,qhull_options="QJn")    #,furthest_site=True   #i put there QJn option
    #tri=Delaunay(self.listqa)
    self.faces = []
    self.points = []


    for vert in tri.simplices:
            self.faces.append([vert[0], vert[1], vert[2]])
    #print(self.faces)

    for i in range(x.shape[0]):
            self.points.append([x[i], y[i], z[i]])
    #print(self.points)
    return self.faces, self.points

After Delaunay function I create stl file from points and faces. Below is what I need but the algorithm only works properly when I got one straight line. I need it for all points.

enter image description here

EDIT1: @Ripi2

I tried your idea but it didn't work.

This is my code of iterating function where I iterate through all pixels and when the pixel have got value, my algorith add it to the list. After that I use this list in Delaunay triangulation.

This is the code:

    def iterator_shape(self):
    #selectedLayerIndex = self.dlg.comboBox.currentIndex()
    path=os.path.join(r'default\python\plugins\print_your_3d\TRASH\shape_to_raster.tif')
    layer = iface.addRasterLayer(sciezka, "qtaz")
    provider = layer.dataProvider()
    extent = provider.extent()
    rows = layer.height()
    cols = layer.width()
    block = provider.block(1, extent, cols, rows)
    xmin = extent.xMinimum()
    ymax = extent.yMaximum()
    xsize = layer.rasterUnitsPerPixelX()
    ysize = layer.rasterUnitsPerPixelY()
    k=1
    xinit = xmin + xsize / 2
    yinit = ymax - ysize / 2
    self.listqa = []
    for i in range(rows):
        for j in range(cols):
            x = xinit + j * xsize
            y = yinit
            k += 1
            if block.value(i, j)  == -3.4028234663852886e+38:
                continue
            elif block.value(i,j)>=self.dlg.doubleSpinBox.value():
                self.listqa.append([x, y, block.value(i, j)])
                self.listqa.append([x,y,block.value(i,j)-1])
                self.listqa.append([x,y+0.01,block.value(i,j)-1])   # @Ridi2 you mean that? 
        xinit = xmin + xsize / 2
        yinit -= ysize
    return self.listqa

This is example of raster which i would like to convert to 3d STL MODEL. 1 pixel have got X,Y and elevation value ( in code this is block.value)

enter image description here

Edit2: I need to make faces between that. This is better visualization made in ArcScene. enter image description here

Kaemer
  • 1
  • 1
  • This issue usually goes away when you have extra points in the lower part of the wall and a bit (1 mm) shifted towards the valley. – Ripi2 Apr 12 '20 at 17:23
  • Also be careful of the way you triangulate. Usually you just need a 2D triangulation (the Z is not needy for the process), and then points in the same vertical plane may mess the process (in the XY plane they are aligned or duplicated), A small bias left or right (that's the difficult part) of the vertical plane fixes issues. – Ripi2 Apr 12 '20 at 17:26
  • Thanks for response, I will check your advice. Soon I write how it go. – Kaemer Apr 12 '20 at 18:07
  • @Ripi2 i tried your idea and I edit my post. Could you check that? – Kaemer Apr 14 '20 at 20:55

0 Answers0