def __init__(self):

        pypes.pypeScript.__init__(self)

        self.Surface = None
        self.DeformedSurface = None
        self.SourcePoints = vtk.vtkPoints()
        self.TargetPoints = vtk.vtkPoints()
        self.DisplacementNorms = vtk.vtkDoubleArray()
        self.Displacements = vtk.vtkDoubleArray()
        self.Displacements.SetNumberOfComponents(3)
        self.SourceSpheres = vtk.vtkPolyData()
        self.TargetSpheres = vtk.vtkPolyData()
        self.SourceSpheres.SetPoints(self.SourcePoints)
        self.TargetSpheres.SetPoints(self.TargetPoints)
        self.SourceSpheres.GetPointData().SetScalars(self.DisplacementNorms)
        self.SourceSpheres.GetPointData().SetVectors(self.Displacements)
        self.vmtkRenderer = None
        self.OwnRenderer = 0
        self.DisplayDeformed = False
        self.SurfaceMapper = None
        self.Opacity = 1.0
        self.SourceSpheresActor = None
        self.TargetSpheresActor = None

        self.SetScriptName("vmtkthinplatesplinedeformation")
        self.SetInputMembers(
            [
                ["Surface", "i", "vtkPolyData", 1, "", "the input surface", "vmtksurfacereader"],
                ["Opacity", "opacity", "float", 1, "(0.0,1.0)", "object opacities in the scene"],
                ["vmtkRenderer", "renderer", "vmtkRenderer", 1, "", "external renderer"],
            ]
        )
        self.SetOutputMembers([["DeformedSurface", "o", "vtkPolyData", 1, "", "", "vmtksurfacewriter"]])

    def __init__(self):

        pypes.pypeScript.__init__(self)
        
        self.Surface = None
        self.ResolutionArrayName = 'ResolutionArray'
        self.RBFType = 'biharmonic'
        self.Spheres = vtk.vtkPolyData()
        self.SphereIds = vtk.vtkIdList()
        self.vmtkRenderer = None
        self.OwnRenderer = 0
        self.DisplayArray = False
        self.SurfaceMapper = None
        self.CurrentSphereId = -1        
        self.SphereWidget = None
        self.Opacity = 1.
        self.SpheresActor = None
        self.ScalarBarActor = None
        self.InteractionMode = 0
        self.ExamineSurface = None
        self.ExamineSpheres = vtk.vtkPolyData()
        self.ExamineSpheresActor = None
        self.ExamineText = None
        
        self.SetScriptName('vtksurfaceresolution')
        self.SetInputMembers([
            ['Surface','i','vtkPolyData',1,'','the input surface','vmtksurfacereader'],
            ['ResolutionArrayName','resolutionarray','str',1,'','array storing the desired edge length'],
            ['RBFType','rbftype','str',1,'["thinplatespline","biharmonic","triharmonic"]','the type of RBF interpolation'],
            ['Opacity','opacity','float',1,'(0.0,1.0)','object opacities in the scene'],
            ['vmtkRenderer','renderer','vmtkRenderer',1,'','external renderer']
            ])
        self.SetOutputMembers([
            ['Surface','o','vtkPolyData',1,'','','vmtksurfacewriter']
            ])

def CreateCoords_versore(o, r):
    """ Ritorna una lista di attori contenenti i il sistema di coordinate:
    o = origine
    r = versore"""
    points = []
    Lines=[]
    Polygon = vtk.vtkPolyData()
    ac=[]
    
    points = vtk.vtkPoints()
            
    points.SetNumberOfPoints(4)

    points.SetPoint(0, self.midPoint)
    points.SetPoint(1, [self.FrenetBinormalArray[0]+self.midPoint[0], self.FrenetBinormalArray[1]+self.midPoint[1], self.FrenetBinormalArray[2]+self.midPoint[2]])
    points.SetPoint(2, [self.FrenetNormalArray[0]+self.midPoint[0], self.FrenetNormalArray[1]+self.midPoint[1], self.FrenetNormalArray[2]+self.midPoint[2]])
    points.SetPoint(3, [self.FrenetTangentArray[0]+self.midPoint[0], self.FrenetTangentArray[1]+self.midPoint[1], self.FrenetTangentArray[2]+self.midPoint[2]])
     
    points.SetPoint(0, o)
    points.SetPoint(1, [o[0]+r[0], o[1]         , o[2]])
    points.SetPoint(2, [o[0]        , o[1]+r[1] ,  o[2]])
    points.SetPoint(3, [o[0]        , o[1]         ,  o[2]+r[2]])
    
    polyLine0 = vtk.vtkPolyLine()
    polyLine0.GetPointIds().SetNumberOfIds(2)
    polyLine0.GetPointIds().SetId(0,0)
    polyLine0.GetPointIds().SetId(1,1)
        
    polyLine1 = vtk.vtkPolyLine()
    polyLine1.GetPointIds().SetNumberOfIds(2)
    polyLine1.GetPointIds().SetId(0,0)
    polyLine1.GetPointIds().SetId(1,2)
    
    polyLine2 = vtk.vtkPolyLine()
    polyLine2.GetPointIds().SetNumberOfIds(2)
    polyLine2.GetPointIds().SetId(0,0)
    polyLine2.GetPointIds().SetId(1,3)
    
    cells0 = vtk.vtkCellArray()
    cells0.InsertNextCell(polyLine0)
    cells0.InsertNextCell(polyLine1)
    cells0.InsertNextCell(polyLine2)
    
    polyData = vtk.vtkPolyData()
    polyData.SetPoints(points)
    polyData.SetLines(cells0)
    
    ac=[]
        
    
    ac.append(CreateSphere(points.GetPoint(0), 0.05, [1, 1, 1]))
    ac.append(CreateSphere(points.GetPoint(1), 0.1, [1, 0, 0]))
    ac.append(CreateSphere(points.GetPoint(2), 0.1, [0, 1, 0]))
    ac.append(CreateSphere(points.GetPoint(3), 0.1, [0, 0, 1]))
        
                
    ac.append(CreateActor(polyData))
    return ac

def RequestInformation():
    import vtk

    ############# Get I/O #############

    # Get the two inputs, and the output
    polyDataA = self.GetInputDataObject(0, 0)
    polyDataB = self.GetInputDataObject(0, 1)
    pdo = self.GetPolyDataOutput()

    # If only one input is given, raise an exception
    if polyDataA is None or polyDataB is None:
        raise Exception("\nThis filter takes 2 inputs:\n"
                        "Point Cloud Data files: pc_HHMMSSDD_NNN.vtk\n"
                        "Pose Data file: pc_HHMMSSDD_poses.vtk\n"
                        "Note that ParaView groups all the Point Cloud Data files in one\n")

    # Initialize vtkPolyData for point cloud data (PC) and pose data (P)
    polyData_PC = vtk.vtkPolyData()
    polyData_P = vtk.vtkPolyData()

    if polyDataA.GetFieldData().GetArray("timestamp") is not None and \
            polyDataB.GetPointData().GetArray("timestamp") is not None:
        pointCloudPortIndex = 0
    else:
        if polyDataB.GetFieldData().GetArray("timestamp") is not None and \
                polyDataA.GetPointData().GetArray("timestamp") is not None:
            pointCloudPortIndex = 1
        else:   # If none of the configuration above is met, raise an exception
            raise Exception("\nOne or both of the inputs don't have a \"timestamp\" Point/Field Data\n"
                            "Is this data coming from the \"Paraview Tango Recorder\" app ?\n"
                            "The input that ends with \'_poses.vtk\" must have a \"timestamp\" PointData\n"
                            "The input that ends with \'*.vtk\" must have a \"timestamp\" FieldData\n")

    def setOutputTimesteps ( algorithm , timesteps ):
        "helper routine to set timestep information"
        executive = algorithm . GetExecutive ()
        outInfo = executive . GetOutputInformation (0)
        outInfo.Remove ( executive.TIME_STEPS () )
        for timestep in timesteps :
            outInfo . Append ( executive . TIME_STEPS () , timestep )

            outInfo . Remove ( executive . TIME_RANGE () )
            outInfo . Append ( executive . TIME_RANGE () , timesteps [0])
            outInfo . Append ( executive . TIME_RANGE () , timesteps [ -1])

    def getInputTimesteps( algorithm, portindex):
        "helper routine to set timestep information"
        executive = algorithm . GetExecutive ()
        inInfo = executive . GetInputInformation (0, portindex)
        return inInfo.Get(executive.TIME_STEPS())

    myrange = getInputTimesteps(self, pointCloudPortIndex)
    setOutputTimesteps(self, myrange)

    def __init__(self, parent = None):
        super(VTKFrame, self).__init__(parent)

        self.vtkWidget = QVTKRenderWindowInteractor(self)
        vl = QtGui.QVBoxLayout(self)
        vl.addWidget(self.vtkWidget)
        vl.setContentsMargins(0, 0, 0, 0)
 
        self.ren = vtk.vtkRenderer()
        self.ren.SetBackground(0.1, 0.2, 0.4)
        self.vtkWidget.GetRenderWindow().AddRenderer(self.ren)
        self.iren = self.vtkWidget.GetRenderWindow().GetInteractor()

        points = vtk.vtkPoints()
        points.InsertNextPoint(0.0, 0.0, 0.0)
        points.InsertNextPoint(1.0, 0.0, 0.0)
        points.InsertNextPoint(0.0, 1.0, 0.0)

        pointsPolydata = vtk.vtkPolyData()
        pointsPolydata.SetPoints(points)

        vertexFilter = vtk.vtkVertexGlyphFilter()
        vertexFilter.SetInputConnection(pointsPolydata.GetProducerPort())
        vertexFilter.Update()

        polydata = vtk.vtkPolyData()
        polydata.ShallowCopy(vertexFilter.GetOutput())

        # Setup colors
        colors = vtk.vtkUnsignedCharArray()
        colors.SetNumberOfComponents(3)
        colors.SetName("Colors")
        colors.InsertNextTupleValue((255, 0, 0))
        colors.InsertNextTupleValue((0, 255, 0))
        colors.InsertNextTupleValue((0, 0, 255))

        polydata.GetPointData().SetScalars(colors)
 
        # Create a mapper
        mapper = vtk.vtkPolyDataMapper()
        mapper.SetInputConnection(polydata.GetProducerPort())
 
        # Create an actor
        actor = vtk.vtkActor()
        actor.SetMapper(mapper)
        actor.GetProperty().SetPointSize(5)
 
        self.ren.AddActor(actor)
        self.ren.ResetCamera()

        self._initialized = False

    def ReadPDB(file_name):
        pdb = vtk.vtkPDBReader()
        pdb.SetFileName(file_name)
        pdb.SetHBScale(1.0)
        pdb.SetBScale(1.0)
        pdb.Update()

        sphere = vtk.vtkSphereSource()
        sphere.SetCenter(0, 0, 0)
        sphere.SetRadius(1)

        glyph = vtk.vtkGlyph3D()
        glyph.SetInputConnection(pdb.GetOutputPort())
        glyph.SetSourceConnection(sphere.GetOutputPort())
        glyph.SetOrient(1)
        glyph.SetColorMode(1)
        glyph.SetScaleMode(2)
        glyph.SetScaleFactor(.25)
        glyph.Update()

        tube = vtk.vtkTubeFilter()
        tube.SetInputConnection(pdb.GetOutputPort())
        tube.SetNumberOfSides(6)
        tube.CappingOff()
        tube.SetRadius(0.2)
        tube.SetVaryRadius(0)
        tube.SetRadiusFactor(10)
        tube.Update()

        tubeMesh = vtk.vtkPolyData()
        tubeMesh.ShallowCopy(tube.GetOutput())
        N = tubeMesh.GetNumberOfPoints()

        rgb_colors = tubeMesh.GetPointData().GetArray("rgb_colors")
        if rgb_colors is not None:
            if rgb_colors.GetNumberOfComponents() == 3:
                for i in range(N):
                    rgb_colors.SetTupleValue(i, (127, 127, 127))

        appendFilter = vtk.vtkAppendPolyData()
        appendFilter.AddInputConnection(glyph.GetOutputPort())
        try:
            appendFilter.AddInputData(tubeMesh)
        except:
            appendFilter.AddInput(tubeMesh)
        appendFilter.Update()

        polyData = vtk.vtkPolyData()
        polyData.ShallowCopy(appendFilter.GetOutput())
        return polyData

def CreatePlanarCrossSectionPolyDataFromFile(file):

  with open(file, 'r') as f:
    read_data = f.read()

  tokens = string.split(read_data)

  offset = 2
  
  planeAppender = vtk.vtkAppendPolyData()
  outlineAppender = vtk.vtkAppendPolyData()

  # Iterate over separate pieces in the file
  while True:
    if (offset >= len(tokens)):
      break
    pointsInPiece = int(tokens[offset])
    
    newPoints = vtk.vtkPoints()
    newPoints.SetNumberOfPoints(pointsInPiece)
    
    for ptId in xrange(pointsInPiece):
      x = float(tokens[ptId*3 + 0 + offset + 1])
      y = float(tokens[ptId*3 + 1 + offset + 1])
      z = float(tokens[ptId*3 + 2 + offset + 1])
      newPoints.SetPoint(ptId, x, y, z)
    
    offset = offset + 3*pointsInPiece + 1
    
    polygon = vtk.vtkPolyData()
    polygon.SetPoints(newPoints)
    polygon.Allocate(pointsInPiece)
    polygon.InsertNextCell(vtk.VTK_POLYGON, pointsInPiece, range(pointsInPiece))

    triFilter = vtk.vtkTriangleFilter()
    triFilter.SetInputData(polygon)
   
    planeAppender.AddInputConnection(triFilter.GetOutputPort())
    
    outline = vtk.vtkPolyData()
    outline.SetPoints(newPoints)
    outline.Allocate(pointsInPiece)
    outline.InsertNextCell(vtk.VTK_POLY_LINE, pointsInPiece, range(pointsInPiece))
    outlineAppender.AddInputData(outline)
    
  planeAppender.Update()
  outlineAppender.Update()

  return (planeAppender.GetOutput(), outlineAppender.GetOutput())

def make_pData_periodic (np, x, a, lattice):
    pos = vtk.vtkPoints()
    diameter = vtk.vtkDoubleArray()
    diameter.SetNumberOfComponents(1)
    # primary cell
    for i in range(np):
        pos.InsertNextPoint(x[i*3], x[i*3+1], x[i*3+2])
        if a != []:
            diameter.InsertNextTuple1(2.0*a[i])
        else:
            diameter.InsertNextTuple1(2.0)
    # image cells
    for ix in range(-1,2):
        for iy in range(-1,2):
            for iz in range(-1,2):
                if ix == 0 and iy == 0 and iz == 0: continue
                for i in range(np):
                    pos.InsertNextPoint(x[i*3  ]+float(ix)*lattice[0],
                                        x[i*3+1]+float(iy)*lattice[1],
                                        x[i*3+2]+float(iz)*lattice[2])
                    if a != []:
                        diameter.InsertNextTuple1(2.0*a[i])
                    else:
                        diameter.InsertNextTuple1(2.0)

    # first make pData containing particle coordinates
    pData = vtk.vtkPolyData()
    pData.SetPoints(pos)
    pData.GetPointData().SetScalars(diameter)

    return pData

def make_sphereActor (x, r, rgb, opacity):
    points = vtk.vtkPoints()
    points.InsertNextPoint(x[0], x[1], x[2])

    diameter = vtk.vtkDoubleArray()
    diameter.SetNumberOfComponents(1)
    diameter.InsertNextTuple1(2.0*r)

    pData = vtk.vtkPolyData()
    pData.SetPoints(points)
    pData.GetPointData().SetScalars(diameter)

    pSource = vtk.vtkSphereSource()
    pSource.SetPhiResolution(16)
    pSource.SetThetaResolution(16)

    pGlyph = vtk.vtkGlyph3D()
    pGlyph.SetSource(pSource.GetOutput())
    pGlyph.SetInput(pData)
    pGlyph.ScalingOn()
    pGlyph.SetScaleModeToScaleByScalar()

    pMapper = vtk.vtkPolyDataMapper()
    pMapper.ScalarVisibilityOff()
    pMapper.SetInput(pGlyph.GetOutput())

    pActor = vtk.vtkActor()
    pActor.SetMapper(pMapper)
    pActor.GetProperty().SetColor(rgb[0], rgb[1], rgb[2])
    pActor.GetProperty().SetOpacity(opacity)

    return pActor

def createPolyData(faces, vtList, verts, tcoords):

    points = vtk.vtkPoints()
    points.SetDataTypeToDouble()
    points.SetNumberOfPoints(len(vtList))

    tcoordArray = vtk.vtkDoubleArray()
    tcoordArray.SetName('tcoords')
    tcoordArray.SetNumberOfComponents(2)
    tcoordArray.SetNumberOfTuples(len(vtList))

    for i, vt in enumerate(vtList):
        vi, ti = vt
        xyz = verts[vi]
        uv = tcoords[ti]

        points.SetPoint(i, xyz)
        tcoordArray.SetTuple2(i, uv[0], uv[1])

    cells = vtk.vtkCellArray()

    for i, face in enumerate(faces):
        tri = vtk.vtkTriangle()
        tri.GetPointIds().SetId(0, face[0])
        tri.GetPointIds().SetId(1, face[1])
        tri.GetPointIds().SetId(2, face[2])
        cells.InsertNextCell(tri)

    polyData = vtk.vtkPolyData()
    polyData.SetPoints(points)
    polyData.SetPolys(cells)
    polyData.GetPointData().SetTCoords(tcoordArray)
    return polyData

def testSingleTriangle():

    "Single triangle"

    h = 0.1
    # create set of points
    points = vtk.vtkPoints()
    points.SetNumberOfPoints(3)
    points.SetPoint(0, [1., -1.*h/3., -1.*h/3.])
    points.SetPoint(1, [1., 2.*h/3., -1.*h/3.])
    points.SetPoint(2, [1., -1.*h/3., 2.*h/3.])

    # create vtkPolyData object
    pdata = vtk.vtkPolyData()
    pdata.SetPoints(points)
    ptIds = vtk.vtkIdList()
    ptIds.SetNumberOfIds(3)
    ptIds.SetId(0, 0)
    ptIds.SetId(1, 1)
    ptIds.SetId(2, 2)
    pdata.Allocate(1, 1)
    pdata.InsertNextCell(vtk.VTK_POLYGON, ptIds)

    for order in range(1, 6):
        lslm = PoissonSolver(pdata, max_edge_length=1000.)
        print('order = ', order)
        print('g matrix: ', lslm.getGreenMatrix())

  def copyFirstNLines(self, sourcePolyData, lineCount):
    """make a polydata with only the first N polylines"""

    polyData = vtk.vtkPolyData()
    points = vtk.vtkPoints()
    polyData.SetPoints(points)

    lines = vtk.vtkCellArray()
    polyData.SetLines(lines)

    sourcePoints = sourcePolyData.GetPoints()
    sourceLines = sourcePolyData.GetLines()
    sourceIdList = vtk.vtkIdList()
    sourceLines.InitTraversal()
    while sourceLines.GetNextCell(sourceIdList):
        pointCount = sourceIdList.GetNumberOfIds()
        idList = vtk.vtkIdList()
        for idIndex in range(pointCount):
            sourceId = sourceIdList.GetId(idIndex)
            point = sourcePoints.GetPoint(sourceId)
            id = points.InsertNextPoint(point)
            idList.InsertNextId(id)
        lines.InsertNextCell(idList)
        if lines.GetNumberOfCells() > lineCount:
            break

    return polyData

  def ProbeData(self, coordinates, name):
    """Interpolate field values at these coordinates."""

    # Initialise locator
    locator = vtk.vtkPointLocator()
    locator.SetDataSet(self.ugrid)
    locator.SetTolerance(10.0)
    locator.Update()

    # Initialise probe
    points = vtk.vtkPoints()
    points.SetDataTypeToDouble()
    ilen, jlen = coordinates.shape
    for i in range(ilen):
      points.InsertNextPoint(coordinates[i][0], coordinates[i][1], coordinates[i][2])
    polydata = vtk.vtkPolyData()
    polydata.SetPoints(points)
    probe = vtk.vtkProbeFilter()
    probe.SetInput(polydata)
    probe.SetSource(self.ugrid)
    probe.Update()

    # Generate a list invalidNodes, containing a map from invalid nodes in the
    # result to their closest nodes in the input
    valid_ids = probe.GetValidPoints()
    valid_loc = 0
    invalidNodes = []
    for i in range(ilen):
      if valid_ids.GetTuple1(valid_loc) == i:
        valid_loc += 1
      else:
        nearest = locator.FindClosestPoint([coordinates[i][0], coordinates[i][1], coordinates[i][2]])
        invalidNodes.append((i, nearest))

    # Get final updated values
    pointdata=probe.GetOutput().GetPointData()
    vtkdata=pointdata.GetArray(name)
    nc=vtkdata.GetNumberOfComponents()
    nt=vtkdata.GetNumberOfTuples()
    array = arr([vtkdata.GetValue(i) for i in range(nt * nc)])
    
    # Fix the point data at invalid nodes
    if len(invalidNodes) > 0:
      try:
        oldField = self.ugrid.GetPointData().GetArray(name)
        components = oldField.GetNumberOfComponents()
      except:
        try:
          oldField = self.ugrid.GetCellData().GetArray(name)
          components = oldField.GetNumberOfComponents()
        except:
          raise Exception("ERROR: couldn't find point or cell field data with name "+name+" in file "+self.filename+".")
      for invalidNode, nearest in invalidNodes:
        for comp in range(nc):
          array[invalidNode * nc + comp] = oldField.GetValue(nearest * nc + comp)
          
    valShape = self.GetField(name)[0].shape
    array.shape = tuple([nt] + list(valShape))
          
    return array

    def GetLineFromWidget(self, obj, event):

        if self.Type == "freehand":
            path = vtk.vtkPolyData()
            obj.GetPath(path)
        elif self.Type == "contour":
            path = self.ImageTracerWidget.GetRepresentation().GetContourRepresentationAsPolyData()

        spacing = self.Image.GetSpacing()

        translation = [0.0, 0.0, 0.0]
        translation[self.Axis] = self.SliceVOI[self.Axis * 2] * spacing[self.Axis]

        transform = vtk.vtkTransform()
        transform.Translate(translation)

        pathTransform = vtk.vtkTransformPolyDataFilter()
        pathTransform.SetInput(path)
        pathTransform.SetTransform(transform)
        pathTransform.Update()

        self.Line = pathTransform.GetOutput()

        if self.Line.GetSource():
            self.Line.GetSource().UnRegisterAllOutputs()

    def __init__(self, pointlist=[]):
        points = vtk.vtkPoints()
        cellArr = vtk.vtkCellArray()
        Colors = vtk.vtkUnsignedCharArray()
        Colors.SetNumberOfComponents(3)
        Colors.SetName("Colors")
        
        n=0
        for p in pointlist:
            vert = vtk.vtkVertex()
            points.InsertNextPoint(p.x, p.y, p.z)
            vert.GetPointIds().SetId(0,n)
            cellArr.InsertNextCell( vert )
            col = clColor(p.cc())
            Colors.InsertNextTuple3( float(255)*col[0], float(255)*col[1], float(255)*col[2] )
            n=n+1
            
        polydata= vtk.vtkPolyData()
        polydata.SetPoints(points)
        polydata.SetVerts( cellArr )
        polydata.GetPointData().SetScalars(Colors)

        polydata.Modified()
        polydata.Update()
        self.src=polydata
        self.mapper = vtk.vtkPolyDataMapper()
        self.mapper.SetInput(self.src)
        self.SetMapper(self.mapper)

 def setEdgesPolydata(self, vd):
     self.edges = []
     self.edges = vd.getEdgesGenerators()
     self.epts = vtk.vtkPoints()
     nid = 0
     lines=vtk.vtkCellArray()
     for e in self.edges:
         p1 = self.scale*e[0]
         p2 = self.scale*e[1] 
         self.epts.InsertNextPoint( p1.x, p1.y, p1.z)
         self.epts.InsertNextPoint( p2.x, p2.y, p2.z)
         line = vtk.vtkLine()
         line.GetPointIds().SetId(0,nid)
         line.GetPointIds().SetId(1,nid+1)
         nid = nid+2
         lines.InsertNextCell(line)
     
     linePolyData = vtk.vtkPolyData()
     linePolyData.SetPoints(self.epts)
     linePolyData.SetLines(lines)
     
     mapper = vtk.vtkPolyDataMapper()
     mapper.SetInput(linePolyData)
     
     self.edge_actor = vtk.vtkActor()
     self.edge_actor.SetMapper(mapper)
     self.edge_actor.GetProperty().SetColor( camvtk.cyan )
     myscreen.addActor( self.edge_actor )
     myscreen.render() 

 def createCurtain( self, **args ):
     trajectory_points = self.getTrajectoryPoints( **args )
     extent =  self.input().GetExtent() 
     spacing =  self.input().GetSpacing() 
     nStrips = extent[5] - extent[4] 
     zmax = spacing[2] * nStrips
     z_inc = zmax / nStrips
     polydata = vtk.vtkPolyData()
     stripArray = vtk.vtkCellArray()
     stripData = [ vtk.vtkIdList() for ix in range(nStrips) ] 
     points = vtk.vtkPoints() 
     for iPt in range( trajectory_points.GetNumberOfPoints() ):
         pt_coords = trajectory_points.GetPoint( iPt )
         z = 0.0
         for iLevel in range( nStrips ):
             vtkId = points.InsertNextPoint( pt_coords[0], pt_coords[1], z )
             sd = stripData[ iLevel ]
             sd.InsertNextId( vtkId )               
             sd.InsertNextId( vtkId+1 )
             z = z + z_inc 
         points.InsertNextPoint( pt_coords[0], pt_coords[1], z )
                    
     for strip in stripData:
         stripArray.InsertNextCell(strip)
         
     polydata.SetPoints( points )
     polydata.SetStrips( stripArray )
     return polydata

 def ConvertDataSetToSurface(algorithmOutputPort):
     dataSetSurfaceFilter = vtk.vtkDataSetSurfaceFilter()
     dataSetSurfaceFilter.SetInputConnection(algorithmOutputPort)
     dataSetSurfaceFilter.Update()
     polyData = vtk.vtkPolyData()
     polyData.ShallowCopy(dataSetSurfaceFilter.GetOutput())
     return polyData

def draw_lines(nodes, color):
    points = vtk.vtkPoints()
    lines = vtk.vtkCellArray()
    nodecnt = 0
    colors = vtk.vtkUnsignedCharArray()
    colors.SetNumberOfComponents(3)
    colors.SetName("Colors")
    
    for node in nodes: 
        edges = node.getedges()
       
        for edge in edges:
            x0,y0,z0 = edge[0]
            x1,y1,z1 = edge[1]

            points.InsertNextPoint(edge[0])
            points.InsertNextPoint(edge[1])

            line = vtk.vtkLine()
            line.GetPointIds().SetId(0,nodecnt)
            line.GetPointIds().SetId(1,nodecnt+1)
            lines.InsertNextCell(line)
            nodecnt += 2
            colors.InsertNextTupleValue(color)
            
    # Create a polydata to store everything in
    linesPolyData = vtk.vtkPolyData()
    # Add the points to the dataset
    linesPolyData.SetPoints(points)
    # Add the lines to the dataset
    linesPolyData.SetLines(lines)
    linesPolyData.GetCellData().SetScalars(colors)
    return linesPolyData

 def readPolyData(file_name, readerType):
     reader = readerType()
     reader.SetFileName(file_name)
     reader.Update()
     polyData = vtk.vtkPolyData()
     polyData.ShallowCopy(reader.GetOutput())
     return polyData