Python vtk.vtkXMLUnstructuredGridWriter() Examples

def createDEM_v2():

    ds = xr.open_dataset('../data/output/Peru_20160601-20180530_comp4.nc')
    
    points = vtk.vtkPoints()
    quad = vtk.vtkQuad()
    cells = vtk.vtkCellArray()
    
    numPoints = ds.south_north.size*ds.west_east.size
   
    print('Write points \n')
    for i,j in product(ds.south_north.values,ds.west_east.values):
            points.InsertNextPoint(ds.lat.isel(south_north=i,west_east=j), ds.lon.isel(south_north=i,west_east=j), ds.HGT.sel(south_north=i,west_east=j).values/6370000.0)
    
    print('Write cells \n')
    for idx in range(points.GetNumberOfPoints()-ds.west_east.size):
        if (idx%ds.west_east.size != 0):
            quad.GetPointIds().SetId(0,idx)
            quad.GetPointIds().SetId(1,idx+1)
            quad.GetPointIds().SetId(2,idx+ds.west_east.size+1)
            quad.GetPointIds().SetId(3,idx+ds.west_east.size)
            cells.InsertNextCell(quad)

    print('Create unstructured grid \n') 
    grid = vtk.vtkUnstructuredGrid()
    grid.SetPoints(points)
    grid.SetCells(vtk.VTK_QUAD, cells)
    
    writer = vtk.vtkXMLUnstructuredGridWriter()
    writer.SetFileName('cosipy.vtu')
    writer.SetInputData(grid)
    writer.Write() 

def write_updated_mesh(mesh, new_vtu_file_path):
    """ Function to write the new VTU file

    Function 'write_updated_mesh' crete new VTU file with utdated value given
    by 'mesh' and save at 'new_vtu_file_path'

    Args:
        mesh (vtkhelpers object instance): Python instance of SU2 result file
                                           with added force and normal vectors
        new_vtu_file_path (str): New VTU file path

    """

    # To write .vtk file
    #writer = vtk.vtkUnstructuredGridWriter()
    #writer.SetFileType(0)

    # To write .vtu file
    writer = vtk.vtkXMLUnstructuredGridWriter()

    try:
        source = mesh.GetOutput()
    except AttributeError:
        source = mesh
    writer.SetInputData(source)
    writer.SetFileName(new_vtu_file_path)
    writer.Update()


# TODO: maybe create some exteral function to cope with SU2Mesh, get coord, get marker ... 

def _save_unstructured_grid(filename, vtkUnstructGrid, directory=''):
        """Saves a VTK unstructured grid file (vtu) for an already generated
        :class:`pyvista.UnstructuredGrid` object.

        Parameters
        ----------

        filename : str
            path to the output vtk file or just its name if directory is specified

        directory : str
            directory where the UBC GIF file lives

        """
        if not isinstance(vtkUnstructGrid, _vtk.vtkUnstructuredGrid):
            raise RuntimeError('`_save_unstructured_grid` can only handle `vtkUnstructuredGrid` objects. `%s` is not supported.' % vtkUnstructGrid.__class__)
        # Check the extension of the filename
        fname = os.path.join(directory, filename)
        ext = os.path.splitext(fname)[1]
        if ext is '':
            fname = fname + '.vtu'
        elif ext not in '.vtu':
            raise IOError('{:s} is an incorrect extension, has to be .vtu'.format(ext))
        # Make the writer
        vtuWriteFilter = _vtkUnstWriter()
        if float(_vtk_version.split('.')[0]) >= 6:
            vtuWriteFilter.SetInputDataObject(vtkUnstructGrid)
        else:
            vtuWriteFilter.SetInput(vtkUnstructGrid)
        vtuWriteFilter.SetFileName(fname)
        # Write the file
        vtuWriteFilter.Update() 

def _get_writer(filetype, filename):
    import vtk

    if filetype in "vtk-ascii":
        logging.warning("VTK ASCII files are only meant for debugging.")
        writer = vtk.vtkUnstructuredGridWriter()
        writer.SetFileTypeToASCII()
    elif filetype == "vtk-binary":
        writer = vtk.vtkUnstructuredGridWriter()
        writer.SetFileTypeToBinary()
    elif filetype == "vtu-ascii":
        logging.warning("VTU ASCII files are only meant for debugging.")
        writer = vtk.vtkXMLUnstructuredGridWriter()
        writer.SetDataModeToAscii()
    elif filetype == "vtu-binary":
        writer = vtk.vtkXMLUnstructuredGridWriter()
        writer.SetDataModeToBinary()
    elif filetype == "xdmf2":
        writer = vtk.vtkXdmfWriter()
    elif filetype == "xdmf3":
        writer = vtk.vtkXdmf3Writer()
    else:
        assert filetype == "exodus", "Unknown file type '{}'.".format(filename)
        writer = vtk.vtkExodusIIWriter()
        # if the mesh contains vtkmodeldata information, make use of it
        # and write out all time steps.
        writer.WriteAllTimeStepsOn()

    return writer 

def createDEM_v1():

    ds = xr.open_dataset('../data/output/Peru_20160601-20180530_comp4.nc')
    
    points = vtk.vtkPoints()
    
    numPoints = ds.south_north.size*ds.west_east.size
   
    print('Write points \n')
    for i,j in product(ds.south_north.values,ds.west_east.values):
            points.InsertNextPoint(ds.lat.isel(south_north=i,west_east=j), ds.lon.isel(south_north=i,west_east=j), ds.HGT.isel(south_north=i,west_east=j).values/6370000.0)
    
    print('Create unstructured grid \n') 
    polydata = vtk.vtkPolyData()
    polydata.SetPoints(points)

    delaunay = vtk.vtkDelaunay2D()
    delaunay.SetInputData(polydata)
    delaunay.Update()

#    subdivision = vtk.vtkButterflySubdivisionFilter()
#    subdivision.SetInputConnection(delaunay.GetOutputPort())
#    subdivision.Update()

    #smoother = vtk.vtkWindowedSincPolyDataFilter()
    #smoother.SetInputConnection(delaunay.GetOutputPort())
    #smoother.SetNumberOfIterations(5)
    #smoother.BoundarySmoothingOff()
    #smoother.FeatureEdgeSmoothingOff()
    #smoother.SetFeatureAngle(120.0)
    #smoother.SetPassBand(.001)
    #smoother.NonManifoldSmoothingOff()
    #smoother.NormalizeCoordinatesOff()
    #smoother.Update()

    appendFilter = vtk.vtkAppendFilter()
    appendFilter.AddInputData(delaunay.GetOutput())
    appendFilter.Update()

    unstructuredGrid = vtk.vtkUnstructuredGrid()
    unstructuredGrid.ShallowCopy(appendFilter.GetOutput())

    writer = vtk.vtkXMLUnstructuredGridWriter()
    writer.SetFileName('cosipy.vtu')
    writer.SetInputData(unstructuredGrid)
    writer.Write() 

def add_scalar(var, timestamp):

    vtkFile = vtk.vtkXMLUnstructuredGridReader()
    vtkFile.SetFileName('cosipy.vtu')
    vtkFile.Update()
    
    # Find cellId by coordinates
    pointLocator =  vtk.vtkPointLocator()
    pointLocator.SetDataSet(vtkFile.GetOutput())
    pointLocator.BuildLocator()
    
    ds = xr.open_dataset('../data/output/Peru_20160601-20180530_comp4.nc')
    ds = ds.sel(time=timestamp)
    
    ds_sub = ds[var].stack(x=['south_north','west_east']) 
    ds_sub = ds_sub.dropna(dim='x')
    lats = ds_sub.x.lat.values
    lons = ds_sub.x.lon.values
    data = ds_sub.values
    print(lats)

    numPoints = vtkFile.GetOutput().GetNumberOfPoints()
    scalar = np.empty(numPoints)
    scalar[:] = np.nan

    interpField = numpy_support.numpy_to_vtk(scalar)
    interpField.SetName(var)
    vtkFile.GetOutput().GetPointData().AddArray(interpField)
    vtkFile.Update()

    print('Write points \n')
    for i in np.arange(len(data)):
        # Get height
        alt = ds.HGT.sel(lat=lats[i],lon=lons[i]).values/6370000.0
        
        pointId = vtk.mutable(0) 
        Id = vtk.vtkIdList()
        pointId = pointLocator.FindClosestPoint([lons[i],lats[i],alt])
        vtkFile.GetOutput().GetPointData().GetArray(var).InsertTuple1(pointId,data[i])

    writer = vtk.vtkXMLUnstructuredGridWriter()
    writer.SetFileName('cosipy.vtu')
    writer.SetInputData(vtkFile.GetOutput())
    writer.Write()

    #plotSurface(vtkFile)