# SPDX-License-Identifier: Apache-2.0
# SPDX-FileCopyrightText: Copyright 2023-2024 TotalEnergies.
# SPDX-FileContributor: Martin Lemay, Paloma Martinez
import logging
import numpy as np
import numpy.typing as npt
from typing import ( Iterator, List, Sequence, Any, Union, Tuple )
from vtkmodules.util.numpy_support import ( numpy_to_vtk, vtk_to_numpy )
from vtkmodules.vtkCommonCore import vtkIdList, vtkPoints, reference, vtkDataArray, vtkLogger, vtkFloatArray
from vtkmodules.vtkCommonDataModel import ( vtkUnstructuredGrid, vtkMultiBlockDataSet, vtkPolyData, vtkDataSet,
vtkDataObject, vtkPlane, vtkCellTypes, vtkIncrementalOctreePointLocator,
VTK_TRIANGLE )
from vtkmodules.vtkFiltersCore import ( vtk3DLinearGridPlaneCutter, vtkPolyDataNormals, vtkPolyDataTangents )
from vtkmodules.vtkFiltersTexture import vtkTextureMapToPlane
from vtkmodules.vtkFiltersGeometry import vtkDataSetSurfaceFilter
from vtkmodules.vtkFiltersGeneral import vtkDataSetTriangleFilter
from geos.mesh.utils.multiblockHelpers import ( getBlockElementIndexesFlatten, getBlockFromFlatIndex )
from geos.utils.algebraFunctions import ( getAttributeMatrixFromVector, getAttributeVectorFromMatrix )
from geos.utils.geometryFunctions import ( getChangeOfBasisMatrix, CANONICAL_BASIS_3D )
from geos.utils.Logger import ( getLogger, Logger, VTKCaptureLog, RegexExceptionFilter )
from geos.utils.Errors import VTKError
__doc__ = """
Generic VTK utilities.
These methods include:
- extraction of a surface from a given elevation
- conversion from a list to vtkIdList
- conversion of vtk container into iterable
"""
[docs]
def isTriangulate( dataSet: vtkUnstructuredGrid ) -> bool:
"""Check if the mesh is fully triangulated.
Args:
dataSet (vtkUnstructuredGrid): The mesh to check
Returns:
bool: True if the mesh is triangulate only, False otherwise.
"""
cellTypes: vtkCellTypes = vtkCellTypes()
dataSet.GetCellTypes( cellTypes )
return all( cellTypes.GetCellType( cell ) == VTK_TRIANGLE for cell in range( cellTypes.GetNumberOfTypes() ) )
[docs]
def triangulateMesh(
dataSet: vtkUnstructuredGrid,
logger: Union[ Logger, None ] = None,
) -> vtkUnstructuredGrid:
"""Triangulate any type of dataset.
Args:
dataSet (vtkUnstructuredGrid): The mesh to triangulate
logger (Union[Logger, None], optional): A logger to manage the output messages.
Defaults to None, an internal logger is used.
Returns:
vtkUnstructuredGrid: Triangulated mesh
"""
vtkErrorLogger: Logger
if logger is None:
vtkErrorLogger = getLogger( "Triangulate Mesh vtkError Logger", True )
else:
vtkErrorLogger = logging.getLogger( f"{ logger.name } vtkError Logger" )
vtkErrorLogger.setLevel( logging.INFO )
vtkErrorLogger.addHandler( logger.handlers[ 0 ] )
vtkErrorLogger.propagate = False
vtkLogger.SetStderrVerbosity( vtkLogger.VERBOSITY_ERROR )
vtkErrorLogger.addFilter( RegexExceptionFilter() ) # will raise VTKError if captured VTK Error
with VTKCaptureLog() as capturedLog:
triangulateMeshFilter: vtkDataSetTriangleFilter = vtkDataSetTriangleFilter()
triangulateMeshFilter.SetInputData( dataSet )
triangulateMeshFilter.Update()
capturedLog.seek( 0 )
captured = capturedLog.read().decode()
if captured != "":
vtkErrorLogger.error( captured.strip() )
triangulatedMesh: vtkUnstructuredGrid = triangulateMeshFilter.GetOutput()
if triangulatedMesh.GetCellData() is None:
raise VTKError( "Something went wrong with VTK triangulation of the mesh." )
return triangulatedMesh
[docs]
def convertUnstructuredGridToPolyData(
dataSet: vtkUnstructuredGrid,
logger: Union[ Logger, None ] = None,
) -> vtkPolyData:
"""Convert vtkUnstructuredGrid to vtkPolyData.
..Warning:: Work only with triangulated mesh
Args:
dataSet (vtkUnstructuredGrid): Input mesh block
logger (Union[Logger, None], optional): A logger to manage the output messages.
Defaults to None, an internal logger is used.
Returns:
vtkPolyData: Extracted surface
"""
vtkErrorLogger: Logger
if logger is None:
vtkErrorLogger = getLogger( "Convert vtkUnstructuredGrid To vtkPolyData vtkError Logger", True )
else:
vtkErrorLogger = logging.getLogger( f"{ logger.name } vtkError Logger" )
vtkErrorLogger.setLevel( logging.INFO )
vtkErrorLogger.addHandler( logger.handlers[ 0 ] )
vtkErrorLogger.propagate = False
vtkLogger.SetStderrVerbosity( vtkLogger.VERBOSITY_ERROR )
vtkErrorLogger.addFilter( RegexExceptionFilter() ) # will raise VTKError if captured VTK Error
with VTKCaptureLog() as capturedLog:
extractSurfaceFilter: vtkDataSetSurfaceFilter = vtkDataSetSurfaceFilter()
extractSurfaceFilter.SetInputData( dataSet )
# fast mode should be used for rendering only
extractSurfaceFilter.FastModeOff()
# Delegation activated allow to accelerate the processing with unstructured mesh
# see https://vtk.org/doc/nightly/html/classvtkDataSetSurfaceFilter.html
extractSurfaceFilter.DelegationOn()
extractSurfaceFilter.Update()
capturedLog.seek( 0 )
captured = capturedLog.read().decode()
if captured != "":
vtkErrorLogger.error( captured.strip() )
extractedSurface: vtkPolyData = extractSurfaceFilter.GetOutput()
if extractedSurface is None:
raise VTKError( "Something went wrong with VTK convert vtu to vtp." )
return extractedSurface
[docs]
def toVtkIdList( data: List[ int ] ) -> vtkIdList:
"""Utility function transforming a list of ids into a vtkIdList.
Args:
data (list[int]): Id list
Returns:
result (vtkIdList): Vtk Id List corresponding to input data
"""
result = vtkIdList()
result.Allocate( len( data ) )
for d in data:
result.InsertNextId( d )
return result
[docs]
def vtkIter( vtkContainer: Union[ vtkIdList, vtkCellTypes ] ) -> Iterator[ Any ]:
"""Utility function transforming a vtk "container" into an iterable.
Args:
vtkContainer (vtkIdList | vtkCellTypes): A vtk container
Returns:
Iterator[ Any ]: The iterator
"""
if isinstance( vtkContainer, vtkIdList ):
for i in range( vtkContainer.GetNumberOfIds() ):
yield vtkContainer.GetId( i )
elif isinstance( vtkContainer, vtkCellTypes ):
for i in range( vtkContainer.GetNumberOfTypes() ):
yield vtkContainer.GetCellType( i )
[docs]
def getBounds(
input: Union[ vtkUnstructuredGrid,
vtkMultiBlockDataSet ] ) -> tuple[ float, float, float, float, float, float ]:
"""Get bounds of either single of composite data set.
Args:
input (Union[vtkUnstructuredGrid, vtkMultiBlockDataSet]): input mesh
Returns:
tuple[float, float, float, float, float, float]: tuple containing
bounds (xmin, xmax, ymin, ymax, zmin, zmax)
"""
if isinstance( input, vtkMultiBlockDataSet ):
return getMultiBlockBounds( input )
else:
return getMonoBlockBounds( input )
[docs]
def getMonoBlockBounds( input: vtkUnstructuredGrid, ) -> tuple[ float, float, float, float, float, float ]:
"""Get boundary box extrema coordinates for a vtkUnstructuredGrid.
Args:
input (vtkMultiBlockDataSet): input single block mesh
Returns:
tuple[float, float, float, float, float, float]: tuple containing
bounds (xmin, xmax, ymin, ymax, zmin, zmax)
"""
return input.GetBounds()
[docs]
def getMultiBlockBounds( input: vtkMultiBlockDataSet, ) -> tuple[ float, float, float, float, float, float ]:
"""Get boundary box extrema coordinates for a vtkMultiBlockDataSet.
Args:
input (vtkMultiBlockDataSet): input multiblock mesh
Returns:
tuple[float, float, float, float, float, float]: bounds.
"""
xmin, ymin, zmin = 3 * [ np.inf ]
xmax, ymax, zmax = 3 * [ -1.0 * np.inf ]
blockIndexes: list[ int ] = getBlockElementIndexesFlatten( input )
for blockIndex in blockIndexes:
block0: vtkDataObject = getBlockFromFlatIndex( input, blockIndex )
assert block0 is not None, "Mesh is undefined."
block: vtkDataSet = vtkDataSet.SafeDownCast( block0 )
bounds: tuple[ float, float, float, float, float, float ] = block.GetBounds()
xmin = bounds[ 0 ] if bounds[ 0 ] < xmin else xmin
xmax = bounds[ 1 ] if bounds[ 1 ] > xmax else xmax
ymin = bounds[ 2 ] if bounds[ 2 ] < ymin else ymin
ymax = bounds[ 3 ] if bounds[ 3 ] > ymax else ymax
zmin = bounds[ 4 ] if bounds[ 4 ] < zmin else zmin
zmax = bounds[ 5 ] if bounds[ 5 ] > zmax else zmax
return xmin, xmax, ymin, ymax, zmin, zmax
[docs]
def getBoundsFromPointCoords( cellPtsCoord: list[ npt.NDArray[ np.float64 ] ] ) -> Sequence[ float ]:
"""Compute bounding box coordinates of the list of points.
Args:
cellPtsCoord (list[npt.NDArray[np.float64]]): list of points
Returns:
Sequence[float]: bounding box coordinates (xmin, xmax, ymin, ymax, zmin, zmax)
"""
bounds: list[ float ] = [
np.inf,
-np.inf,
np.inf,
-np.inf,
np.inf,
-np.inf,
]
for ptsCoords in cellPtsCoord:
mins: npt.NDArray[ np.float64 ] = np.min( ptsCoords, axis=0 )
maxs: npt.NDArray[ np.float64 ] = np.max( ptsCoords, axis=0 )
for i in range( 3 ):
bounds[ 2 * i ] = float( min( bounds[ 2 * i ], mins[ i ] ) )
bounds[ 2 * i + 1 ] = float( max( bounds[ 2 * i + 1 ], maxs[ i ] ) )
return bounds
[docs]
def createSingleCellMesh( cellType: int, ptsCoord: npt.NDArray[ np.float64 ] ) -> vtkUnstructuredGrid:
"""Create a mesh that consists of a single cell.
Args:
cellType (int): cell type
ptsCoord (1DArray[np.float64]): cell point coordinates
Returns:
vtkUnstructuredGrid: output mesh
"""
nbPoints: int = ptsCoord.shape[ 0 ]
points: npt.NDArray[ np.float64 ] = np.vstack( ( ptsCoord, ) )
# Convert points to vtkPoints object
vtkpts: vtkPoints = vtkPoints()
vtkpts.SetData( numpy_to_vtk( points ) )
# create cells from point ids
cellsID: vtkIdList = vtkIdList()
for j in range( nbPoints ):
cellsID.InsertNextId( j )
# add cell to mesh
mesh: vtkUnstructuredGrid = vtkUnstructuredGrid()
mesh.SetPoints( vtkpts )
mesh.Allocate( 1 )
mesh.InsertNextCell( cellType, cellsID )
return mesh
[docs]
def createMultiCellMesh( cellTypes: list[ int ],
cellPtsCoord: list[ npt.NDArray[ np.float64 ] ],
sharePoints: bool = True ) -> vtkUnstructuredGrid:
"""Create a mesh that consists of multiple cells.
.. WARNING:: The mesh is not checked for conformity.
Args:
cellTypes (list[int]): Cell type
cellPtsCoord (list[1DArray[np.float64]]): List of cell point coordinates
sharePoints (bool): If True, cells share points, else a new point is created for each cell vertex
Returns:
vtkUnstructuredGrid: Output mesh
"""
assert len( cellPtsCoord ) == len( cellTypes ), "The lists of cell types of point coordinates must be of same size."
nbCells: int = len( cellPtsCoord )
mesh: vtkUnstructuredGrid = vtkUnstructuredGrid()
points: vtkPoints
cellVertexMapAll: list[ tuple[ int, ...] ]
points, cellVertexMapAll = createVertices( cellPtsCoord, sharePoints )
assert len( cellVertexMapAll ) == len(
cellTypes ), "The lists of cell types of cell point ids must be of same size."
mesh.SetPoints( points )
mesh.Allocate( nbCells )
# create mesh cells
for cellType, ptsId in zip( cellTypes, cellVertexMapAll, strict=True ):
# create cells from point ids
cellsID: vtkIdList = vtkIdList()
for ptId in ptsId:
cellsID.InsertNextId( ptId )
mesh.InsertNextCell( cellType, cellsID )
return mesh
[docs]
def createVertices( cellPtsCoord: list[ npt.NDArray[ np.float64 ] ],
shared: bool = True ) -> tuple[ vtkPoints, list[ tuple[ int, ...] ] ]:
"""Create vertices from cell point coordinates list.
Args:
cellPtsCoord (list[npt.NDArray[np.float64]]): list of cell point coordinates
shared (bool, optional): If True, collocated points are merged. Defaults to True.
Returns:
tuple[vtkPoints, list[tuple[int, ...]]]: tuple containing points and the
map of cell point ids
"""
# get point bounds
bounds: Sequence[ float ] = getBoundsFromPointCoords( cellPtsCoord )
points: vtkPoints = vtkPoints()
# use point locator to check for colocated points
pointsLocator = vtkIncrementalOctreePointLocator()
pointsLocator.InitPointInsertion( points, bounds )
cellVertexMapAll: list[ tuple[ int, ...] ] = []
ptId: reference = reference( 0 )
ptsCoords: npt.NDArray[ np.float64 ]
for ptsCoords in cellPtsCoord:
cellVertexMap: list[ int ] = []
pt: npt.NDArray[ np.float64 ] # 1DArray
for pt in ptsCoords:
if shared:
pointsLocator.InsertUniquePoint( pt.tolist(), ptId ) # type: ignore[arg-type]
else:
pointsLocator.InsertPointWithoutChecking( pt.tolist(), ptId, 1 ) # type: ignore[arg-type]
cellVertexMap += [ ptId.get() ] # type: ignore
cellVertexMapAll += [ tuple( cellVertexMap ) ]
return points, cellVertexMapAll
[docs]
def convertAttributeFromLocalToXYZForOneCell(
vector: npt.NDArray[ np.float64 ], localBasisVectors: tuple[ npt.NDArray[ np.float64 ], npt.NDArray[ np.float64 ],
npt.NDArray[ np.float64 ] ]
) -> npt.NDArray[ np.float64 ]:
"""Convert one cell attribute from local to XYZ basis.
.. Warning::
Vectors components are considered to be in GEOS output order such that \
v = ( M00, M11, M22, M12, M02, M01, M21, M20, M10 )
Args:
vector (npt.NDArray[np.float64]): The attribute expressed in local basis. The size of the vector must be 3, 9 or 6 (symmetrical case)
localBasisVectors (tuple[ npt.NDArray[np.float64], npt.NDArray[np.float64], npt.NDArray[np.float64]]): The local basis vectors expressed in the XYZ basis.
Returns:
vectorXYZ (npt.NDArray[np.float64]): The attribute expressed in XYZ basis.
"""
# Get components matrix from GEOS attribute vector
matrix3x3: npt.NDArray[ np.float64 ] = getAttributeMatrixFromVector( vector )
# Get transform matrix
transformMatrix: npt.NDArray[ np.float64 ] = getChangeOfBasisMatrix( localBasisVectors, CANONICAL_BASIS_3D )
# Apply transformation
arrayXYZ: npt.NDArray[ np.float64 ] = transformMatrix @ matrix3x3 @ transformMatrix.T
# Convert back to GEOS type attribute and return
return getAttributeVectorFromMatrix( arrayXYZ, vector.size )
[docs]
def getNormalVectors( surface: vtkPolyData ) -> npt.NDArray[ np.float64 ]:
"""Return the normal vectors of a surface mesh.
Args:
surface (vtkPolyData): The input surface.
Raises:
ValueError: No normal attribute found in the mesh. Use the computeNormals function beforehand.
Returns:
npt.NDArray[ np.float64 ]: The normal vectors of the input surface.
"""
normals: Union[ npt.NDArray[ np.float64 ],
None ] = surface.GetCellData().GetNormals() # type: ignore[no-untyped-call]
if normals is None:
raise ValueError( "No normal attribute found in the mesh. Use the computeNormals function beforehand." )
return vtk_to_numpy( normals )
[docs]
def getTangentsVectors( surface: vtkPolyData ) -> Tuple[ npt.NDArray[ np.float64 ], npt.NDArray[ np.float64 ] ]:
"""Return the tangential vectors of a surface.
Args:
surface (vtkPolyData): The input surface.
Raises:
ValueError: Tangent attribute not found in the mesh. Use computeTangents beforehand.
ValueError: Problem during the calculation of the second tangential vectors.
Returns:
Tuple[npt.NDArray[np.float64], npt.NDArray[np.float64]]: The tangents vectors of the input surface.
"""
# Get first tangential component
vtkTangents: Union[ vtkFloatArray, None ] = surface.GetCellData().GetTangents() # type: ignore[no-untyped-call]
tangents1: npt.NDArray[ np.float64 ]
try:
tangents1 = vtk_to_numpy( vtkTangents )
except AttributeError as err:
context: str = f"No tangential attribute found in the mesh. Use the computeTangents function beforehand.\n{ err }"
raise VTKError( context ) from err
else:
# Compute second tangential component
normals: npt.NDArray[ np.float64 ] = getNormalVectors( surface )
tangents2: npt.NDArray[ np.float64 ] = np.cross( normals, tangents1, axis=1 ).astype( np.float64 )
return ( tangents1, tangents2 )
[docs]
def getLocalBasisVectors(
surface: vtkPolyData,
logger: Union[ Logger, None ] = None,
) -> npt.NDArray[ np.float64 ]:
"""Return the local basis vectors for all cells of the input surface.
Args:
surface(vtkPolydata): The input surface.
logger (Union[Logger, None], optional): A logger to manage the output messages.
Defaults to None, an internal logger is used.
Returns:
npt.NDArray[np.float64]: Array with normal, tangential 1 and tangential 2 vectors.
"""
if logger is None:
logger = getLogger( "getLocalBasisVectors", True )
try:
normals: npt.NDArray[ np.float64 ] = getNormalVectors( surface )
surfaceWithNormals: vtkPolyData = surface
# ValueError raised if no normals found in the mesh
except ValueError:
# In that case, the normals are computed.
surfaceWithNormals = computeNormals( surface, logger )
normals = getNormalVectors( surfaceWithNormals )
# Tangents require normals to be present in the mesh
try:
tangents: Tuple[ npt.NDArray[ np.float64 ],
npt.NDArray[ np.float64 ] ] = getTangentsVectors( surfaceWithNormals )
# If no tangents is present in the mesh, they are computed on that surface
except VTKError:
surfaceWithTangents: vtkPolyData = computeTangents( surfaceWithNormals, logger )
tangents = getTangentsVectors( surfaceWithTangents )
return np.array( ( normals, *tangents ) )
[docs]
def computeNormals(
surface: vtkPolyData,
logger: Union[ Logger, None ] = None,
) -> vtkPolyData:
"""Compute and set the normals of a given surface.
Args:
surface (vtkPolyData): The input surface.
logger (Union[Logger, None], optional): A logger to manage the output messages.
Defaults to None, an internal logger is used.
Returns:
vtkPolyData: The surface with normal attribute.
"""
vtkErrorLogger: Logger
if logger is None:
vtkErrorLogger = getLogger( "Compute Surface Normals vtkError Logger", True )
else:
vtkErrorLogger = logging.getLogger( f"{ logger.name } vtkError Logger" )
vtkErrorLogger.setLevel( logging.INFO )
vtkErrorLogger.addHandler( logger.handlers[ 0 ] )
vtkErrorLogger.propagate = False
vtkLogger.SetStderrVerbosity( vtkLogger.VERBOSITY_ERROR )
vtkErrorLogger.addFilter( RegexExceptionFilter() ) # will raise VTKError if captured VTK Error
with VTKCaptureLog() as capturedLog:
normalFilter: vtkPolyDataNormals = vtkPolyDataNormals()
normalFilter.SetInputData( surface )
normalFilter.ComputeCellNormalsOn()
normalFilter.ComputePointNormalsOff()
normalFilter.Update()
capturedLog.seek( 0 )
captured = capturedLog.read().decode()
if captured != "":
vtkErrorLogger.error( captured.strip() )
outputSurface = normalFilter.GetOutput()
if outputSurface.GetCellData().GetNormals() is None:
raise VTKError( "Something went wrong with VTK calculation of the normals." )
return outputSurface
[docs]
def computeTangents(
triangulatedSurface: vtkPolyData,
logger: Union[ Logger, None ] = None,
) -> vtkPolyData:
"""Compute and set the tangents of a given surface.
.. Warning:: The computation of tangents requires a triangulated surface.
Args:
triangulatedSurface (vtkPolyData): The input surface. It should be triangulated beforehand.
logger (Union[Logger, None], optional): A logger to manage the output messages.
Defaults to None, an internal logger is used.
Returns:
vtkPolyData: The surface with tangent attribute
"""
# Need to compute texture coordinates required for tangent calculation
surface1: vtkPolyData = computeSurfaceTextureCoordinates( triangulatedSurface, logger )
# TODO: triangulate the surface before computation of the tangents if needed.
# Compute tangents
vtkErrorLogger: Logger
if logger is None:
vtkErrorLogger = getLogger( "Compute Surface Tangents vtkError Logger", True )
else:
vtkErrorLogger = logging.getLogger( f"{ logger.name } vtkError Logger" )
vtkErrorLogger.setLevel( logging.INFO )
vtkErrorLogger.addHandler( logger.handlers[ 0 ] )
vtkErrorLogger.propagate = False
vtkLogger.SetStderrVerbosity( vtkLogger.VERBOSITY_ERROR )
vtkErrorLogger.addFilter( RegexExceptionFilter() ) # will raise VTKError if captured VTK Error
with VTKCaptureLog() as capturedLog:
tangentFilter: vtkPolyDataTangents = vtkPolyDataTangents()
tangentFilter.SetInputData( surface1 )
tangentFilter.ComputeCellTangentsOn()
tangentFilter.ComputePointTangentsOff()
tangentFilter.Update()
surfaceOut: vtkPolyData = tangentFilter.GetOutput()
capturedLog.seek( 0 )
captured = capturedLog.read().decode()
if captured != "":
vtkErrorLogger.error( captured.strip() )
if surfaceOut is None:
raise VTKError( "Something went wrong in VTK calculation." )
# copy tangents attributes into filter input surface because surface attributes
# are not transferred to the output (bug that should be corrected by Kitware)
array: vtkDataArray = surfaceOut.GetCellData().GetTangents()
if array is None:
raise VTKError( "Attribute Tangents is not in the mesh." )
surface1.GetCellData().SetTangents( array )
surface1.GetCellData().Modified()
surface1.Modified()
return surface1
[docs]
def computeSurfaceTextureCoordinates(
surface: vtkPolyData,
logger: Union[ Logger, None ] = None,
) -> vtkPolyData:
"""Generate the 2D texture coordinates required for tangent computation. The dataset points are mapped onto a plane generated automatically.
Args:
surface (vtkPolyData): The input surface.
logger (Union[Logger, None], optional): A logger to manage the output messages.
Defaults to None, an internal logger is used.
Returns:
vtkPolyData: The input surface with generated texture map.
"""
# Need to compute texture coordinates required for tangent calculation
vtkErrorLogger: Logger
if logger is None:
vtkErrorLogger = getLogger( "Compute Surface Texture Coordinates vtkError Logger", True )
else:
vtkErrorLogger = logging.getLogger( f"{ logger.name } vtkError Logger" )
vtkErrorLogger.setLevel( logging.INFO )
vtkErrorLogger.addHandler( logger.handlers[ 0 ] )
vtkErrorLogger.propagate = False
vtkLogger.SetStderrVerbosity( vtkLogger.VERBOSITY_ERROR )
vtkErrorLogger.addFilter( RegexExceptionFilter() ) # will raise VTKError if captured VTK Error
with VTKCaptureLog() as capturedLog:
textureFilter: vtkTextureMapToPlane = vtkTextureMapToPlane()
textureFilter.SetInputData( surface )
textureFilter.AutomaticPlaneGenerationOn()
textureFilter.Update()
capturedLog.seek( 0 )
captured = capturedLog.read().decode()
if captured != "":
vtkErrorLogger.error( captured.strip() )
return textureFilter.GetOutput()