Is it possible to create Delaunay triangulation only in vertical (Z dimension)?

255 views Asked by At

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

0

There are 0 answers