{{{#!python
# -*- coding: utf-8 -*-
"""
Created on Sun Mar 23 11:17:34 2014

@author: fuji
"""

import vtk
import numpy as np

# Strctured Grid
sGrid = vtk.vtkStructuredGrid()

with open('officeFlow.txt', 'r') as CFDfile:
    nx, ny, nz = CFDfile.readline().split()
    dims = [int(nx),int(ny),int(nz)]

    sGrid.SetDimensions(dims[0],dims[1],dims[2])
    nodes = vtk.vtkPoints()
    nodes.SetNumberOfPoints(dims[0] * dims[1] * dims[2])

    scalArray = vtk.vtkDoubleArray()
    scalArray.SetNumberOfComponents(1)
    scalArray.SetName("Pressure")
    scalArray.SetNumberOfTuples(dims[0]*dims[1]*dims[2])

    velArray = vtk.vtkDoubleArray()
    velArray.SetNumberOfComponents(3)
    velArray.SetName("Velocity")
    velArray.SetNumberOfTuples(dims[0]*dims[1]*dims[2])

    for k in range(dims[2]):
        kOffset = k * dims[0] * dims[1]
        for j in range(dims[1]):
            jOffset = j * dims[0]
            for i in range(dims[0]):
                offset = i + jOffset + kOffset
                x,y,z,u,v,w,p = CFDfile.readline().split()
                nodes.SetPoint(offset,[float(x), float(y), float(z)])
                vel = np.array([float(u), float(v), float(w)])
                scalArray.SetValue(offset,float(p))
                velArray.SetTupleValue(offset,vel)

    sGrid.SetPoints(nodes)
    sGrid.GetPointData().SetVectors(velArray)
    sGrid.GetPointData().SetScalars(scalArray)
#
#
writer = vtk.vtkXMLStructuredGridWriter()
writer.SetInput(sGrid)

writer.SetFileName("officeFlow.vts")
writer.SetDataModeToAscii()
writer.Write()
}}}