Added two new visualisations that use the vtk toolbox

wind_tools/VTKbased/cutoffVisualizationTool.py

@@ -0,0 +1,152 @@
+#!/usr/bin/python3
+# -*- coding: utf-8 -*-
+
+import vtk
+from PyQt5.QtWidgets import QWidget, QSlider, QLabel, QVBoxLayout
+from PyQt5.QtCore import Qt
+import numpy as np
+from vtk.qt.QVTKRenderWindowInteractor import QVTKRenderWindowInteractor
+
+
+class cutoffInterface(QWidget):
+    def __init__(self, fileLoc):
+        super().__init__()
+        self.OpacB = 3
+        self.OpacV = .1
+
+        # Start the VTK viewer and the user interface
+        self.ren, self.oTF, self.scalar_range, self.scalar_bar = instantiateVTKviewer(fileLoc)
+
+        self.initUI()
+        self.rW = self.vtkWidget.GetRenderWindow()
+        self.iren = self.rW.GetInteractor()
+        self.rW.AddRenderer(self.ren)
+
+        # Change the VTK element interaction style and show the GUI
+        self.iren.SetInteractorStyle(vtk.vtkInteractorStyleTrackballCamera())
+        self.iren.Initialize()
+        self.show()
+
+    def initUI(self):
+        layout = QVBoxLayout()
+        self.setLayout(layout)
+        self.l1 = QLabel('Cutoff value = 3 m/s')
+        self.l1.setAlignment(Qt.AlignLeft)
+        layout.addWidget(self.l1)
+
+        sldCuttoff = QSlider(Qt.Horizontal)
+        sldCuttoff.setSingleStep(.1)
+        sldCuttoff.setGeometry(30, 40, 100, 30)
+        sldCuttoff.valueChanged[int].connect(self.changeOpacityBoundary)
+        layout.addWidget(sldCuttoff)
+
+        self.l2 = QLabel('Transparancy value = 10%')
+        self.l2.setAlignment(Qt.AlignLeft)
+        layout.addWidget(self.l2)
+
+        sldTransparancy = QSlider(Qt.Horizontal)
+        sldTransparancy.setSingleStep(.1)
+        sldTransparancy.setGeometry(30, 40, 100, 30)
+        sldTransparancy.valueChanged[int].connect(self.changeOpacity)
+        layout.addWidget(sldTransparancy)
+
+        self.vtkWidget = QVTKRenderWindowInteractor()
+        layout.addWidget(self.vtkWidget)
+
+        self.setGeometry(300, 100, 1000, 800)
+        self.setWindowTitle('Slicer Interface')
+
+    def changeOpacityBoundary(self, value):
+        self.OpacB = (self.scalar_range[0] +
+                      (self.scalar_range[1]-self.scalar_range[0])*value/100)
+        self.l1.setText('Cutoff value = %.3f m/s' % self.OpacB)
+        self.updateLut()
+
+    def changeOpacity(self, value):
+        self.OpacV = value/1000
+        self.l2.setText('Transparancy value = %.1f%%' % (value/10))
+        self.updateLut()
+
+    def updateLut(self):
+        self.oTF.RemoveAllPoints()
+        self.oTF.AddPoint(self.scalar_range[0], self.OpacV)
+        self.oTF.AddPoint(self.OpacB-.01, self.OpacV)
+        self.oTF.AddPoint(self.OpacB+.01, 0)
+        self.oTF.AddPoint(self.scalar_range[1], 0)
+        self.rW.Render()
+
+
+def instantiateVTKviewer(fileLoc):
+    # Create the reader for the data
+    reader = vtk.vtkXMLImageDataReader()
+    reader.SetFileName(fileLoc)
+    reader.Update()
+    flowField = reader.GetOutput()
+    scalar_range = flowField.GetScalarRange()
+
+    # Create a custom lut, it's used both as a mapper and at the scalar_bar
+    lut = vtk.vtkLookupTable()
+    lut.SetTableRange(scalar_range)
+    lut.Build()
+
+    # create the scalar_bar
+    scalar_bar = vtk.vtkScalarBarActor()
+    scalar_bar.SetLookupTable(lut)
+    scalar_bar.SetNumberOfLabels(6)
+    scalar_bar.GetLabelTextProperty().SetFontFamilyToCourier()
+    scalar_bar.GetLabelTextProperty().SetJustificationToRight()
+    scalar_bar.GetLabelTextProperty().SetVerticalJustificationToCentered()
+    scalar_bar.GetLabelTextProperty().BoldOff()
+    scalar_bar.GetLabelTextProperty().ItalicOff()
+    scalar_bar.GetLabelTextProperty().ShadowOff()
+    scalar_bar.GetLabelTextProperty().SetColor(0, 0, 0)
+
+    # Setup rendering
+    # Create transfer mapping scalar value to opacity
+    opacityTransferFunction = vtk.vtkPiecewiseFunction()
+    opacityTransferFunction.AddPoint(0, .1)
+    opacityTransferFunction.AddPoint(2.95, .1)
+    opacityTransferFunction.AddPoint(3.05, 0)
+    opacityTransferFunction.AddPoint(7, 0)
+
+    # Create transfer mapping scalar value to color
+    colorTransferFunction = vtk.vtkColorTransferFunction()
+    for s in np.linspace(scalar_range[0], scalar_range[1], 200):
+        col = [0, 0, 0]
+        lut.GetColor(s, col)
+        colorTransferFunction.AddRGBPoint(s, col[0], col[1], col[2])
+
+    # The property describes how the data will look
+    volumeProperty = vtk.vtkVolumeProperty()
+    volumeProperty.SetColor(colorTransferFunction)
+    volumeProperty.SetScalarOpacity(opacityTransferFunction)
+    volumeProperty.SetInterpolationTypeToLinear()
+
+    volumeProperty.SetComponentWeight(0, 1.0)
+    volumeProperty.SetComponentWeight(1, 0.0)
+    volumeProperty.SetComponentWeight(2, 0.0)
+
+    # The mapper / ray cast function know how to render the data
+    volumeMapper = vtk.vtkGPUVolumeRayCastMapper()
+    volumeMapper.SetInputConnection(reader.GetOutputPort())
+
+    # The volume holds the mapper and the property and
+    # can be used to position/orient the volume
+    volume = vtk.vtkVolume()
+    volume.SetMapper(volumeMapper)
+    volume.SetProperty(volumeProperty)
+
+    # Setup camera
+    camera = vtk.vtkCamera()
+    camera.SetPosition(-800, -400, 300)
+    camera.SetFocalPoint(0, 0, 0)
+    camera.SetViewUp(0, 0, 1)
+
+    renderer = vtk.vtkRenderer()
+    renderer.AddVolume(volume)
+    renderer.AddActor(scalar_bar)
+    renderer.SetBackground(240/255, 240/255, 240/255)
+    renderer.SetActiveCamera(camera)
+    renderer.ResetCamera()
+
+    return renderer, opacityTransferFunction, scalar_range, scalar_bar

wind_tools/VTKbased/flowfieldWriter.py

@@ -0,0 +1,221 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Mon Feb 12 16:07:47 2018
+
+@author: roald
+"""
+
+
+import vtk
+import numpy as np
+from pandas import DataFrame
+
+
+def create_vti_from_df(fileName, df, spacing, dimensions, origin):
+    """Create a .vti fileformat for use with the vtk toolbox
+
+    input:
+        fileName: The fileName of the .vti file to be created
+        df: a pandas table with the columns, x,y,z,u,v,w of all relevant flow info
+        spacing: The spacing in the x, y and z direction between points
+        dimensions: a tuple or list of integers with the number of points in xyz
+        origin: the origin of the flow field, for reconstructing turbine coords
+
+    output: -
+
+    Roald Storm, 2018 """
+    imageData = vtk.vtkImageData()
+    imageData.SetDimensions(list(int(x) for x in dimensions))
+    imageData.AllocateScalars(vtk.VTK_DOUBLE, 3)
+    imageData.SetSpacing(spacing)
+    imageData.SetOrigin(origin)
+
+    # Fill every entry of the image data with "2.0"
+    i = 0
+    for z in range(int(dimensions[2])):
+        for y in range(int(dimensions[1])):
+            for x in range(int(dimensions[0])):
+                imageData.SetScalarComponentFromDouble(
+                        x, y, z, 0, df.u[i])
+                imageData.SetScalarComponentFromDouble(
+                        x, y, z, 1, df.v[i])
+                imageData.SetScalarComponentFromDouble(
+                        x, y, z, 2, df.w[i])
+                i += 1
+
+    writer = vtk.vtkXMLImageDataWriter()
+    writer.SetFileName(fileName)
+    writer.SetInputData(imageData)
+    writer.Write()
+
+
+def create_df_from_vti(fileName):
+    """Read flow array output from SOWFA
+
+
+    input:
+        filename: name of .vti file to open
+
+    output:
+        df: a pandas table with the columns, x,y,z,u,v,w of all relevant flow info
+        spacing: The spacing in the x, y and z direction between points
+        dimensions: a tuple or list of integers with the number of points in xyz
+        origin: the origin of the flow field, for reconstructing turbine coords
+
+    Roald Storm, 2018 """
+
+    reader = vtk.vtkXMLImageDataReader()
+    reader.SetFileName(fileName)
+    reader.Update()
+
+    readerOutput = reader.GetOutput()
+    dimensions = readerOutput.GetDimensions()
+    Nx, Ny, Nz = dimensions
+    spacing = readerOutput.GetSpacing()
+    origin = readerOutput.GetOrigin()
+
+    xRange = np.arange(0, Nx*spacing[0], spacing[0])
+    yRange = np.arange(0, Ny*spacing[1], spacing[1])
+    zRange = np.arange(0, Nz*spacing[2], spacing[2])
+
+    pts = np.array([(x, y, z) for z in zRange for y in yRange for x in xRange])
+    vals = np.zeros(pts.shape)
+
+    for z in range(0, Nz):
+        for y in range(0, Ny):
+            for x in range(0, Nx):
+                vals[z*Ny*Nx + y*Nx + x] = [
+                        readerOutput.GetScalarComponentAsDouble(x, y, z, 0),
+                        readerOutput.GetScalarComponentAsDouble(x, y, z, 1),
+                        readerOutput.GetScalarComponentAsDouble(x, y, z, 2)]
+
+    df = DataFrame(columns=('x', 'y', 'z', 'u', 'v', 'w'))
+
+    df['u'] = vals[:, 0]
+    df['v'] = vals[:, 1]
+    df['w'] = vals[:, 2]
+    df['x'] = pts[:, 0]
+    df['y'] = pts[:, 1]
+    df['z'] = pts[:, 2]
+
+    return df, spacing, dimensions, origin
+
+
+def create_df_from_vtk_poly(fileName):
+    """Read flow array output from SOWFA
+
+
+    input:
+        filename: name of .vtk file to open with a 2D flowfield slice specified
+        as PolyData
+
+    output:
+        df: a pandas table with the columns, x,y,z,u,v,w of all relevant flow info
+        spacing: The spacing in the x, y and z direction between points
+        dimensions: a tuple or list of integers with the number of points in xyz
+        origin: the origin of the flow field, for reconstructing turbine coords
+
+    Roald Storm, 2018 """
+
+    reader = vtk.vtkGenericDataObjectReader()
+    reader.SetFileName(fileName)
+    reader.Update()
+    readerOutput = reader.GetOutput()
+
+    if not reader.IsFilePolyData():
+        raise ValueError('This .vtk file does not have PolyData')
+
+    if readerOutput.GetCell(0).GetPoints().GetData().GetNumberOfTuples() != 4:
+        raise ValueError('This .vtk file is not a 2d flowfield slice')
+
+    ncells = readerOutput.GetNumberOfCells()
+    cellData = readerOutput.GetCellData()
+
+    pts = np.zeros([ncells, 3])
+    vals = np.zeros([ncells, 3])
+    data = cellData.GetArray(0)
+
+    for i in range(data.GetNumberOfTuples()):
+        vals[i] = data.GetTuple(i)
+        pts[i] = [min(i) for i in zip(*[readerOutput.GetCell(i).GetPoints().GetData().GetTuple(x) for x in range(4)])]
+
+    df = DataFrame(columns=('x', 'y', 'z', 'u', 'v', 'w'))
+
+    df['u'] = vals[:, 0]
+    df['v'] = vals[:, 1]
+    df['w'] = vals[:, 2]
+    df['x'] = pts[:, 0]
+    df['y'] = pts[:, 1]
+    df['z'] = pts[:, 2]
+
+    df = df.sort_values(by=['z', 'y', 'x']).reset_index(drop=True)
+
+    dimensions = (len(df[(df['y'] == 0)].index),
+                  len(df[(df['x'] == 0)].index),
+                  1)
+    spacing = (df[df['x'] != 0].iloc[0].loc['x'],
+               df[df['y'] != 0].iloc[0].loc['y'],
+               0)
+    origin = (0.0, 0.0, 0.0)
+
+    return df, spacing, dimensions, origin
+
+
+def create_df_from_vtk_unstructuredgrid(fileName):
+    """Read flow array output from SOWFA
+
+
+    input:
+        filename: name of .vtk file to open with 3D flowfield specified as an
+        unstructured grid
+
+    output:
+        df: a pandas table with the columns, x,y,z,u,v,w of all relevant flow info
+        spacing: The spacing in the x, y and z direction between points
+        dimensions: a tuple or list of integers with the number of points in xyz
+        origin: the origin of the flow field, for reconstructing turbine coords
+
+    Roald Storm, 2018 """
+
+    reader = vtk.vtkGenericDataObjectReader()
+    reader.SetFileName(fileName)
+    reader.Update()
+    readerOutput = reader.GetOutput()
+
+    if not reader.IsFileUnstructuredGrid():
+        raise ValueError('This .vtk file does not have Unstructured Grid Data')
+
+    if readerOutput.GetCell(0).GetPoints().GetData().GetNumberOfTuples() != 8:
+        raise ValueError('This .vtk file is not a 3d flowfield slice')
+
+    ncells = readerOutput.GetNumberOfCells()
+    cellData = readerOutput.GetCellData()
+
+    pts = np.zeros([ncells, 3])
+    vals = np.zeros([ncells, 3])
+    data = cellData.GetArray(1)
+
+    for i in range(data.GetNumberOfTuples()):
+        vals[i] = data.GetTuple(i)
+        pts[i] = [min(i) for i in zip(*[readerOutput.GetCell(i).GetPoints().GetData().GetTuple(x) for x in range(8)])]
+
+    df = DataFrame(columns=('x', 'y', 'z', 'u', 'v', 'w'))
+
+    df['u'] = vals[:, 0]
+    df['v'] = vals[:, 1]
+    df['w'] = vals[:, 2]
+    df['x'] = pts[:, 0]
+    df['y'] = pts[:, 1]
+    df['z'] = pts[:, 2]
+
+    df = df.sort_values(by=['z', 'y', 'x']).reset_index(drop=True)
+
+    dimensions = (len(df[(df['y'] == 0) & (df['z'] == 0)].index),
+                  len(df[(df['x'] == 0) & (df['z'] == 0)].index),
+                  len(df[(df['x'] == 0) & (df['y'] == 0)].index))
+    spacing = (df[df['x'] != 0].iloc[0].loc['x'],
+               df[df['y'] != 0].iloc[0].loc['y'],
+               df[df['z'] != 0].iloc[0].loc['z'])
+    origin = (0.0, 0.0, 0.0)
+
+    return df, spacing, dimensions, origin