LineOnMesh
VTKExamples/Python/DataManipulation/LineOnMesh
Description¶
Create a terrain with regularly spaced points and smooth the terrain (see SmoothMeshgrid. Then, a vtkCellLocator is used to find the intersection between vertical lines and the surface (using the IntersectWithLine method). Each intersection is saved as x, y, z coordinates and used to plot a spline on the surface.
- Contributed by Michka Popoff
Code¶
LineOnMesh.py
#!/usr/bin/env python import vtk import random import numpy # Make a 32 x 32 grid size = 32 # Define z values for the topography topography = numpy.zeros([size, size]) for i in range(size): for j in range(size): topography[i][j] = random.randrange(0, 5) # Define points, triangles and colors colors = vtk.vtkUnsignedCharArray() colors.SetNumberOfComponents(3) points = vtk.vtkPoints() triangles = vtk.vtkCellArray() # Build the meshgrid manually count = 0 for i in range(size-1): for j in range(size-1): z1 = topography[i][j] z2 = topography[i][j+1] z3 = topography[i+1][j] # Triangle 1 points.InsertNextPoint(i, j, z1) points.InsertNextPoint(i, (j+1), z2) points.InsertNextPoint((i+1), j, z3) triangle = vtk.vtkTriangle() triangle.GetPointIds().SetId(0, count) triangle.GetPointIds().SetId(1, count + 1) triangle.GetPointIds().SetId(2, count + 2) triangles.InsertNextCell(triangle) z1 = topography[i][j+1] z2 = topography[i+1][j+1] z3 = topography[i+1][j] # Triangle 2 points.InsertNextPoint(i, (j+1), z1) points.InsertNextPoint((i+1), (j+1), z2) points.InsertNextPoint((i+1), j, z3) triangle = vtk.vtkTriangle() triangle.GetPointIds().SetId(0, count + 3) triangle.GetPointIds().SetId(1, count + 4) triangle.GetPointIds().SetId(2, count + 5) count += 6 triangles.InsertNextCell(triangle) # Add some color r = [int(i/float(size)*255),int(j/float(size)*255),0] colors.InsertNextTypedTuple(r) colors.InsertNextTypedTuple(r) colors.InsertNextTypedTuple(r) colors.InsertNextTypedTuple(r) colors.InsertNextTypedTuple(r) colors.InsertNextTypedTuple(r) # Create a polydata object trianglePolyData = vtk.vtkPolyData() # Add the geometry and topology to the polydata trianglePolyData.SetPoints(points) trianglePolyData.GetPointData().SetScalars(colors) trianglePolyData.SetPolys(triangles) # Clean the polydata so that the edges are shared ! cleanPolyData = vtk.vtkCleanPolyData() cleanPolyData.SetInputData(trianglePolyData) # Use a filter to smooth the data (will add triangles and smooth) smooth_loop = vtk.vtkLoopSubdivisionFilter() smooth_loop.SetNumberOfSubdivisions(3) smooth_loop.SetInputConnection(cleanPolyData.GetOutputPort()) # Create a mapper and actor for smoothed dataset mapper = vtk.vtkPolyDataMapper() mapper.SetInputConnection(smooth_loop.GetOutputPort()) actor_loop = vtk.vtkActor() actor_loop.SetMapper(mapper) actor_loop.GetProperty().SetInterpolationToFlat() # Update the pipeline so that vtkCellLocator finds cells ! smooth_loop.Update() # Define a cellLocator to be able to compute intersections between lines # and the surface locator = vtk.vtkCellLocator() locator.SetDataSet(smooth_loop.GetOutput()) locator.BuildLocator() maxloop = 1000 dist = 20.0/maxloop tolerance = 0.001 # Make a list of points. Each point is the intersection of a vertical line # defined by p1 and p2 and the surface. points = vtk.vtkPoints() for i in range(maxloop): p1 = [2+i*dist, 16, -1] p2 = [2+i*dist, 16, 6] # Outputs (we need only pos which is the x, y, z position # of the intersection) t = vtk.mutable(0) pos = [0.0, 0.0, 0.0] pcoords = [0.0, 0.0, 0.0] subId = vtk.mutable(0) locator.IntersectWithLine(p1, p2, tolerance, t, pos, pcoords, subId) # Add a slight offset in z pos[2] += 0.01 # Add the x, y, z position of the intersection points.InsertNextPoint(pos) # Create a spline and add the points spline = vtk.vtkParametricSpline() spline.SetPoints(points) functionSource = vtk.vtkParametricFunctionSource() functionSource.SetUResolution(maxloop) functionSource.SetParametricFunction(spline) # Map the spline mapper = vtk.vtkPolyDataMapper() mapper.SetInputConnection(functionSource.GetOutputPort()) # Define the line actor actor = vtk.vtkActor() actor.SetMapper(mapper) actor.GetProperty().SetColor([1.0, 0.0, 0.0]) actor.GetProperty().SetLineWidth(3) # Visualize renderer = vtk.vtkRenderer() renderWindow = vtk.vtkRenderWindow() renderWindow.AddRenderer(renderer) renderWindowInteractor = vtk.vtkRenderWindowInteractor() renderWindowInteractor.SetRenderWindow(renderWindow) # Add actors and render renderer.AddActor(actor) renderer.AddActor(actor_loop) renderer.SetBackground(1, 1, 1) # Background color white renderWindow.SetSize(800, 800) renderWindow.Render() renderWindowInteractor.Start()