Difference between revisions of "VTK/Examples/Python/MeshLabelImageColor"

From KitwarePublic
< VTK‎ | Examples‎ | Python
Jump to navigationJump to search
 
Line 4: Line 4:
 
==MeshLabelImageColor.py==
 
==MeshLabelImageColor.py==
 
<source lang="python">
 
<source lang="python">
 +
from __future__ import print_function
 +
 
import vtk
 
import vtk
+
 
 
input='labels.mhd'
 
input='labels.mhd'
+
 
 
# Prepare to read the file
 
# Prepare to read the file
+
 
 
readerVolume = vtk.vtkMetaImageReader()
 
readerVolume = vtk.vtkMetaImageReader()
 
readerVolume.SetFileName(input)
 
readerVolume.SetFileName(input)
 
readerVolume.Update()
 
readerVolume.Update()
+
 
+
 
 
# Extract the region of interest
 
# Extract the region of interest
 
voi = vtk.vtkExtractVOI()
 
voi = vtk.vtkExtractVOI()
Line 20: Line 22:
 
     voi.SetInput(readerVolume.GetOutput())
 
     voi.SetInput(readerVolume.GetOutput())
 
else:
 
else:
     voi.SetInputConnection(readerVolume.GetOutputPort())  
+
     voi.SetInputConnection(readerVolume.GetOutputPort())
 
#voi.SetVOI(0,517, 0,228, 0,392)
 
#voi.SetVOI(0,517, 0,228, 0,392)
 
voi.SetSampleRate(1,1,1)
 
voi.SetSampleRate(1,1,1)
Line 26: Line 28:
 
voi.Update()#necessary for GetScalarRange()
 
voi.Update()#necessary for GetScalarRange()
 
srange= voi.GetOutput().GetScalarRange()#needs Update() before!
 
srange= voi.GetOutput().GetScalarRange()#needs Update() before!
print "Range", srange
+
print("Range", srange)
+
 
+
 
 
##Prepare surface generation
 
##Prepare surface generation
 
contour = vtk.vtkDiscreteMarchingCubes() #for label images
 
contour = vtk.vtkDiscreteMarchingCubes() #for label images
Line 34: Line 36:
 
     contour.SetInput( voi.GetOutput() )
 
     contour.SetInput( voi.GetOutput() )
 
else:
 
else:
     contour.SetInputConnection( voi.GetOutputPort() )  
+
     contour.SetInputConnection( voi.GetOutputPort() )
 
#contour.ComputeNormalsOn()
 
#contour.ComputeNormalsOn()
+
 
+
 
 
#choose one label
 
#choose one label
 
index= 31
 
index= 31
+
 
print "Doing label", index
+
print("Doing label", index)
+
 
 
contour.SetValue(0, index)
 
contour.SetValue(0, index)
 
contour.Update() #needed for GetNumberOfPolys() !!!
 
contour.Update() #needed for GetNumberOfPolys() !!!
+
 
+
 
 
smoother= vtk.vtkWindowedSincPolyDataFilter()
 
smoother= vtk.vtkWindowedSincPolyDataFilter()
 
if vtk.VTK_MAJOR_VERSION <= 5:
 
if vtk.VTK_MAJOR_VERSION <= 5:
 
     smoother.SetInput(contour.GetOutput())
 
     smoother.SetInput(contour.GetOutput())
 
else:
 
else:
     smoother.SetInputConnection(contour.GetOutputPort())  
+
     smoother.SetInputConnection(contour.GetOutputPort())
 
smoother.SetNumberOfIterations(30) #this has little effect on the error!
 
smoother.SetNumberOfIterations(30) #this has little effect on the error!
 
#smoother.BoundarySmoothingOff()
 
#smoother.BoundarySmoothingOff()
Line 59: Line 61:
 
smoother.NonManifoldSmoothingOn()
 
smoother.NonManifoldSmoothingOn()
 
smoother.NormalizeCoordinatesOn()
 
smoother.NormalizeCoordinatesOn()
smoother.GenerateErrorScalarsOn()  
+
smoother.GenerateErrorScalarsOn()
 
#smoother.GenerateErrorVectorsOn()
 
#smoother.GenerateErrorVectorsOn()
 
smoother.Update()
 
smoother.Update()
+
 
 
smoothed_polys= smoother.GetOutput()
 
smoothed_polys= smoother.GetOutput()
 
smoother_error= smoothed_polys.GetPointData().GetScalars()
 
smoother_error= smoothed_polys.GetPointData().GetScalars()
+
 
 
##Find min and max z
 
##Find min and max z
 
se_range= smoother_error.GetRange()
 
se_range= smoother_error.GetRange()
print se_range
+
print(se_range)
 
minz= se_range[0] #min(smoother_error)
 
minz= se_range[0] #min(smoother_error)
 
maxz= se_range[1] #max(smoother_error)
 
maxz= se_range[1] #max(smoother_error)
 
if (maxz > 1):
 
if (maxz > 1):
     print "Big smoother error: min/max:", minz, maxz
+
     print("Big smoother error: min/max:", minz, maxz)
 
minz=  .3 #this way colours of different particles are comparable
 
minz=  .3 #this way colours of different particles are comparable
 
maxz= 1
 
maxz= 1
+
 
 
## Create the color map
 
## Create the color map
 
colorLookupTable= vtk.vtkLookupTable()
 
colorLookupTable= vtk.vtkLookupTable()
Line 84: Line 86:
 
#colorLookupTable.SetNumberOfColors(256) #256 default
 
#colorLookupTable.SetNumberOfColors(256) #256 default
 
colorLookupTable.Build()
 
colorLookupTable.Build()
+
 
 
##calc cell normal
 
##calc cell normal
 
triangleCellNormals= vtk.vtkPolyDataNormals()
 
triangleCellNormals= vtk.vtkPolyDataNormals()
Line 90: Line 92:
 
     triangleCellNormals.SetInput(smoothed_polys)
 
     triangleCellNormals.SetInput(smoothed_polys)
 
else:
 
else:
     triangleCellNormals.SetInputData(smoothed_polys)  
+
     triangleCellNormals.SetInputData(smoothed_polys)
 
triangleCellNormals.ComputeCellNormalsOn()
 
triangleCellNormals.ComputeCellNormalsOn()
 
triangleCellNormals.ComputePointNormalsOff()
 
triangleCellNormals.ComputePointNormalsOff()
Line 96: Line 98:
 
triangleCellNormals.AutoOrientNormalsOn()
 
triangleCellNormals.AutoOrientNormalsOn()
 
triangleCellNormals.Update() #creates vtkPolyData
 
triangleCellNormals.Update() #creates vtkPolyData
+
 
+
 
 
mapper= vtk.vtkPolyDataMapper()
 
mapper= vtk.vtkPolyDataMapper()
 
#mapper.SetInput(smoothed_polys) #this has no normals
 
#mapper.SetInput(smoothed_polys) #this has no normals
Line 103: Line 105:
 
     mapper.SetInput(triangleCellNormals.GetOutput()) #this is better for vis;-)
 
     mapper.SetInput(triangleCellNormals.GetOutput()) #this is better for vis;-)
 
else:
 
else:
     mapper.SetInputConnection(triangleCellNormals.GetOutputPort()) #this is better for vis;-)  
+
     mapper.SetInputConnection(triangleCellNormals.GetOutputPort()) #this is better for vis;-)
mapper.ScalarVisibilityOn()#show colour  
+
mapper.ScalarVisibilityOn()#show colour
 
mapper.SetScalarRange(minz, maxz)
 
mapper.SetScalarRange(minz, maxz)
 
#mapper.SetScalarModeToUseCellData() # contains the label eg. 31
 
#mapper.SetScalarModeToUseCellData() # contains the label eg. 31
 
mapper.SetScalarModeToUsePointData() #the smoother error relates to the verts
 
mapper.SetScalarModeToUsePointData() #the smoother error relates to the verts
 
mapper.SetLookupTable(colorLookupTable)
 
mapper.SetLookupTable(colorLookupTable)
+
 
+
 
 
# Take the isosurface data and create geometry
 
# Take the isosurface data and create geometry
+
 
 
actor = vtk.vtkLODActor()
 
actor = vtk.vtkLODActor()
+
 
 
actor.SetNumberOfCloudPoints( 100000 )
 
actor.SetNumberOfCloudPoints( 100000 )
 
actor.SetMapper(mapper)
 
actor.SetMapper(mapper)
+
 
+
 
 
# Create renderer
 
# Create renderer
+
 
 
ren = vtk.vtkRenderer()
 
ren = vtk.vtkRenderer()
ren.SetBackground( 0, 0, 0 )  
+
ren.SetBackground( 0, 0, 0 )
 
ren.AddActor(actor)
 
ren.AddActor(actor)
+
 
+
 
+
 
 
# Create a window for the renderer of size 250x250
 
# Create a window for the renderer of size 250x250
+
 
 
renWin = vtk.vtkRenderWindow()
 
renWin = vtk.vtkRenderWindow()
+
 
 
renWin.AddRenderer(ren)
 
renWin.AddRenderer(ren)
 
renWin.SetSize(250, 250)
 
renWin.SetSize(250, 250)
+
 
+
 
+
 
 
# Set an user interface interactor for the render window
 
# Set an user interface interactor for the render window
+
 
 
iren = vtk.vtkRenderWindowInteractor()
 
iren = vtk.vtkRenderWindowInteractor()
+
 
 
iren.SetRenderWindow(renWin)
 
iren.SetRenderWindow(renWin)
+
 
+
 
+
 
 
# Start the initialization and rendering
 
# Start the initialization and rendering
+
 
 
iren.Initialize()
 
iren.Initialize()
 
renWin.Render()
 
renWin.Render()
 
iren.Start()
 
iren.Start()
 +
 
</source>
 
</source>

Latest revision as of 02:45, 22 August 2015

This example takes a label image in Meta format and meshes a single label of it. It then smoothes the mesh and colours the vertices according to the displacement error introduced by the smoothing. The whole thing is then displayed. An example data can be found here Media:labels.tgz


MeshLabelImageColor.py

from __future__ import print_function

import vtk

input='labels.mhd'

# Prepare to read the file

readerVolume = vtk.vtkMetaImageReader()
readerVolume.SetFileName(input)
readerVolume.Update()


# Extract the region of interest
voi = vtk.vtkExtractVOI()
if vtk.VTK_MAJOR_VERSION <= 5:
    voi.SetInput(readerVolume.GetOutput())
else:
    voi.SetInputConnection(readerVolume.GetOutputPort())
#voi.SetVOI(0,517, 0,228, 0,392)
voi.SetSampleRate(1,1,1)
#voi.SetSampleRate(3,3,3)
voi.Update()#necessary for GetScalarRange()
srange= voi.GetOutput().GetScalarRange()#needs Update() before!
print("Range", srange)


##Prepare surface generation
contour = vtk.vtkDiscreteMarchingCubes() #for label images
if vtk.VTK_MAJOR_VERSION <= 5:
    contour.SetInput( voi.GetOutput() )
else:
    contour.SetInputConnection( voi.GetOutputPort() )
#contour.ComputeNormalsOn()


#choose one label
index= 31

print("Doing label", index)

contour.SetValue(0, index)
contour.Update() #needed for GetNumberOfPolys() !!!


smoother= vtk.vtkWindowedSincPolyDataFilter()
if vtk.VTK_MAJOR_VERSION <= 5:
    smoother.SetInput(contour.GetOutput())
else:
    smoother.SetInputConnection(contour.GetOutputPort())
smoother.SetNumberOfIterations(30) #this has little effect on the error!
#smoother.BoundarySmoothingOff()
#smoother.FeatureEdgeSmoothingOff()
#smoother.SetFeatureAngle(120.0)
#smoother.SetPassBand(.001)        #this increases the error a lot!
smoother.NonManifoldSmoothingOn()
smoother.NormalizeCoordinatesOn()
smoother.GenerateErrorScalarsOn()
#smoother.GenerateErrorVectorsOn()
smoother.Update()

smoothed_polys= smoother.GetOutput()
smoother_error= smoothed_polys.GetPointData().GetScalars()

##Find min and max z
se_range= smoother_error.GetRange()
print(se_range)
minz= se_range[0] #min(smoother_error)
maxz= se_range[1] #max(smoother_error)
if (maxz > 1):
    print("Big smoother error: min/max:", minz, maxz)
minz=  .3 #this way colours of different particles are comparable
maxz= 1

## Create the color map
colorLookupTable= vtk.vtkLookupTable()
colorLookupTable.SetTableRange(minz, maxz) #this does nothing, use mapper.SetScalarRange(minz, maxz)
colorLookupTable.SetHueRange(2/3.0, 1)
#colorLookupTable.SetSaturationRange(0, 0)
#colorLookupTable.SetValueRange(1, 0)
#colorLookupTable.SetNumberOfColors(256) #256 default
colorLookupTable.Build()

##calc cell normal
triangleCellNormals= vtk.vtkPolyDataNormals()
if vtk.VTK_MAJOR_VERSION <= 5:
    triangleCellNormals.SetInput(smoothed_polys)
else:
    triangleCellNormals.SetInputData(smoothed_polys)
triangleCellNormals.ComputeCellNormalsOn()
triangleCellNormals.ComputePointNormalsOff()
triangleCellNormals.ConsistencyOn()
triangleCellNormals.AutoOrientNormalsOn()
triangleCellNormals.Update() #creates vtkPolyData


mapper= vtk.vtkPolyDataMapper()
#mapper.SetInput(smoothed_polys) #this has no normals
if vtk.VTK_MAJOR_VERSION <= 5:
    mapper.SetInput(triangleCellNormals.GetOutput()) #this is better for vis;-)
else:
    mapper.SetInputConnection(triangleCellNormals.GetOutputPort()) #this is better for vis;-)
mapper.ScalarVisibilityOn()#show colour
mapper.SetScalarRange(minz, maxz)
#mapper.SetScalarModeToUseCellData() # contains the label eg. 31
mapper.SetScalarModeToUsePointData() #the smoother error relates to the verts
mapper.SetLookupTable(colorLookupTable)


# Take the isosurface data and create geometry

actor = vtk.vtkLODActor()

actor.SetNumberOfCloudPoints( 100000 )
actor.SetMapper(mapper)


# Create renderer

ren = vtk.vtkRenderer()
ren.SetBackground( 0, 0, 0 )
ren.AddActor(actor)



# Create a window for the renderer of size 250x250

renWin = vtk.vtkRenderWindow()

renWin.AddRenderer(ren)
renWin.SetSize(250, 250)



# Set an user interface interactor for the render window

iren = vtk.vtkRenderWindowInteractor()

iren.SetRenderWindow(renWin)



# Start the initialization and rendering

iren.Initialize()
renWin.Render()
iren.Start()