From 30c002afc30f30073a926e3beaf9161770679086 Mon Sep 17 00:00:00 2001 From: Roald Storm Date: Mon, 12 Feb 2018 16:56:59 +0100 Subject: [PATCH 1/4] Added two tools to visualize flowData and one script to export a pandas dataframe specifying flowdata to a .vti file. --- .../VTKbased/cutoffVisualizationTool.py | 176 +++++++++++++++++ wind_tools/VTKbased/flowfieldWriter.py | 45 +++++ wind_tools/VTKbased/slicerTool.py | 178 ++++++++++++++++++ 3 files changed, 399 insertions(+) create mode 100644 wind_tools/VTKbased/cutoffVisualizationTool.py create mode 100644 wind_tools/VTKbased/flowfieldWriter.py create mode 100644 wind_tools/VTKbased/slicerTool.py diff --git a/wind_tools/VTKbased/cutoffVisualizationTool.py b/wind_tools/VTKbased/cutoffVisualizationTool.py new file mode 100644 index 0000000..1dd62a3 --- /dev/null +++ b/wind_tools/VTKbased/cutoffVisualizationTool.py @@ -0,0 +1,176 @@ +#!/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 + + +class cutoffInterface(QWidget): + def __init__(self, fileLoc): + super().__init__() + + self.OpacB = 3 + self.OpacV = .1 + + # Start the VTK viewer and the user interface + self.renderWindow, self.oTF, self.scalar_range, self.scalar_bar = instantiateVTKviewer(fileLoc) + self.initUI() + + # Make the vtk application interactive as well + iren = vtk.vtkRenderWindowInteractor() + iren.SetInteractorStyle(vtk.vtkInteractorStyleTrackballCamera()) + iren.SetRenderWindow(self.renderWindow) + + # create the scalar_bar_widget + scalar_bar_widget = vtk.vtkScalarBarWidget() + scalar_bar_widget.SetInteractor(iren) + scalar_bar_widget.SetScalarBarActor(self.scalar_bar) + scalar_bar_widget.On() + + iren.Initialize() + iren.Start() + + 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.setGeometry(300, 300, 250, 200) + self.setWindowTitle('Slicer Interface') + self.show() + + 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.renderWindow.Render() + + +def instantiateVTKviewer(fileLoc): + # Create the standard renderer, render window and interactor + ren = vtk.vtkRenderer() + renWin = vtk.vtkRenderWindow() + renWin.AddRenderer(ren) + iren = vtk.vtkRenderWindowInteractor() + iren.SetRenderWindow(renWin) + + # Create the reader for the data + reader = vtk.vtkXMLImageDataReader() + reader.SetFileName(fileLoc) + reader.Update() + flowField = reader.GetOutput() + scalar_range = flowField.GetScalarRange() +# scalar_range = (2.0, scalar_range[1]) + + # 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.SetOrientationToHorizontal() + 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.SetBackground(1, 1, 1) + renderer.SetActiveCamera(camera) + renderer.ResetCamera() + + # Set renderingwindow and render for the first time + renderWindow = vtk.vtkRenderWindow() + renderWindow.AddRenderer(renderer) + renderWindow.Render() + renderWindow.AddObserver('AbortCheckEvent', CheckAbort) + + return renderWindow, opacityTransferFunction, scalar_range, scalar_bar + + +def CheckAbort(obj, event): + if obj.GetEventPending() != 0: + obj.SetAbortRender(1) + +#cutoffInterface('C:/Users/roald/Dropbox/cases_for_roald/sowfaCases/CT_test_U8_pitch_00_tilt_+10/array.mean/20600.4/array.mean0D_UAvg.vtk') diff --git a/wind_tools/VTKbased/flowfieldWriter.py b/wind_tools/VTKbased/flowfieldWriter.py new file mode 100644 index 0000000..9e1f6d9 --- /dev/null +++ b/wind_tools/VTKbased/flowfieldWriter.py @@ -0,0 +1,45 @@ +# -*- coding: utf-8 -*- +""" +Created on Mon Feb 12 16:07:47 2018 + +@author: roald +""" + + +import vtk + + +def create_vti_from_df(df_flow, spacing, dimensions, fileName): + """Create a .vti fileformat for use with the vtk toolbox + + input: + df_flow: a flow dataframe from flow_field, can be SOWFA or FLORIS + spacing: The spacing in the x, y and z direction between points + dimensions: a tuple or ist of integers with the number of points in xyz + fileName: The fileName of the .vti file to be created + + 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) + + # 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_flow.u[i]) + imageData.SetScalarComponentFromDouble( + x, y, z, 1, df_flow.v[i]) + imageData.SetScalarComponentFromDouble( + x, y, z, 2, df_flow.w[i]) + i += 1 + + writer = vtk.vtkXMLImageDataWriter() + writer.SetFileName(fileName) + writer.SetInputData(imageData) + writer.Write() diff --git a/wind_tools/VTKbased/slicerTool.py b/wind_tools/VTKbased/slicerTool.py new file mode 100644 index 0000000..fa9d00a --- /dev/null +++ b/wind_tools/VTKbased/slicerTool.py @@ -0,0 +1,178 @@ +#!/usr/bin/python3 +# -*- coding: utf-8 -*- + +import vtk +from PyQt5.QtWidgets import QWidget, QSlider, QRadioButton, QGridLayout +from PyQt5.QtCore import Qt + + +class slicerInterface(QWidget): + def __init__(self, fileLoc): + super().__init__() + + # Start the VTK viewer and the user interface + (self.rW, self.renderer, self.scalar_bar, self.dims, self.cellSizes, + self.reader, self.lut, self.scalar_range) = instantiateVTKviewer2(fileLoc) + self.initUI() + + # Instantiate the dynamical slice that will be changed on user input + self.dynActor = makePlaneAt((1, 0, 0), (0, 0, 0), + self.reader, self.lut, self.scalar_range) + self.renderer.AddActor(self.dynActor) + self.rButtonMap = {'X': 0, 'Y': 1, 'Z': 2} + self.sliceVal = 0 + self.ax = 0 + + # Make the vtk application interactive as well + iren = vtk.vtkRenderWindowInteractor() + iren.SetInteractorStyle(vtk.vtkInteractorStyleTrackballCamera()) + iren.SetRenderWindow(self.rW) + + # create the scalar_bar_widget + scalar_bar_widget = vtk.vtkScalarBarWidget() + scalar_bar_widget.SetInteractor(iren) + scalar_bar_widget.SetScalarBarActor(self.scalar_bar) + scalar_bar_widget.On() + + iren.Initialize() + iren.Start() + + def initUI(self): + Grid = QGridLayout() + self.setLayout(Grid) + + rButtonX = QRadioButton('X') + rButtonX.setChecked(True) + rButtonX.toggled.connect(self.on_radio_button_toggled) + Grid.addWidget(rButtonX, 0, 0) + + rButtonY = QRadioButton('Y') + rButtonY.toggled.connect(self.on_radio_button_toggled) + Grid.addWidget(rButtonY, 0, 1) + + rButtonZ = QRadioButton('Z') + rButtonZ.toggled.connect(self.on_radio_button_toggled) + Grid.addWidget(rButtonZ, 0, 2) + + sld = QSlider(Qt.Horizontal) + sld.setSingleStep(.1) + sld.setGeometry(30, 40, 100, 30) + sld.valueChanged[int].connect(self.changeSliderValue) + Grid.addWidget(sld, 1, 0, 1, 3) + + self.setGeometry(300, 300, 250, 150) + self.setWindowTitle('Slicer Interface') + self.show() + + def changeSliderValue(self, value): + self.sliceVal = value/100 + self.changeSlice() + + def on_radio_button_toggled(self): + self.ax = self.rButtonMap[self.sender().text()] + self.changeSlice() + + def changeSlice(self): + self.renderer.RemoveActor(self.dynActor) + d = [0, 0, 0] + d[self.ax] = 1 + l = [0, 0, 0] + l[self.ax] = self.dims[self.ax]*self.cellSizes[self.ax]*self.sliceVal + self.dynActor = makePlaneAt(d, l, self.reader, + self.lut, self.scalar_range) + self.renderer.AddActor(self.dynActor) + self.rW.Render() + + +def instantiateVTKviewer2(fileLoc): + # Read the source file. + reader = vtk.vtkXMLImageDataReader() + reader.SetFileName(fileLoc) + reader.Update() + flowField = reader.GetOutput() + scalar_range = flowField.GetScalarRange() + + # Create a custom lut. The lut is used both at the mapper and at the + # scalar_bar + lut = vtk.vtkLookupTable() + lut.SetTableRange(scalar_range) + lut.Build() + + dims = flowField.GetDimensions() + cellSizes = flowField.GetSpacing() + + # create the scalar_bar + scalar_bar = vtk.vtkScalarBarActor() + scalar_bar.SetOrientationToHorizontal() + scalar_bar.SetLookupTable(lut) + scalar_bar.SetNumberOfLabels(8) + 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) + + # Create the three planes at the back sides of the volume + cutActor1 = makePlaneAt((1, 0, 0), (dims[0]*cellSizes[0]*.99, 0, 0), + reader, lut, scalar_range) + cutActor2 = makePlaneAt((0, 1, 0), (0, dims[1]*cellSizes[1]*.99, 0), + reader, lut, scalar_range) + cutActor3 = makePlaneAt((0, 0, 1), (0, 0, dims[2]*cellSizes[2]*.01), + reader, lut, scalar_range) + + # Setup camera + camera = vtk.vtkCamera() + camera.SetPosition(-800, -400, 300) + camera.SetFocalPoint(0, 0, 0) + camera.SetViewUp(0, 0, 1) + + # Setup lights, One for each plane to prevent any kind of shadow artefacts + light1 = vtk.vtkLight() + light1.SetPosition(0, 0, 0) + light1.SetFocalPoint(1, 0, 0) + light2 = vtk.vtkLight() + light2.SetPosition(0, 0, 0) + light2.SetFocalPoint(0, 1, 0) + light3 = vtk.vtkLight() + light3.SetPosition(0, 0, 1000) + light3.SetFocalPoint(0, 0, 1) + + # Setup rendering + renderer = vtk.vtkRenderer() + renderer.AddActor(cutActor1) + renderer.AddActor(cutActor2) + renderer.AddActor(cutActor3) + renderer.SetBackground(1, 1, 1) + renderer.SetActiveCamera(camera) + renderer.ResetCamera() + renderer.AddLight(light1) + renderer.AddLight(light2) + renderer.AddLight(light3) + + # Set renderingwindow and render for the first time + renderWindow = vtk.vtkRenderWindow() + renderWindow.AddRenderer(renderer) + renderWindow.Render() + return (renderWindow, renderer, scalar_bar, dims, cellSizes, + reader, lut, scalar_range) + + +def makePlaneAt(d, l, reader, lut, scalar_range): + plane = vtk.vtkPlane() + plane.SetOrigin(l) + plane.SetNormal(d) + + planeCut = vtk.vtkCutter() + planeCut.SetInputConnection(reader.GetOutputPort()) + planeCut.SetCutFunction(plane) + + cutMapper = vtk.vtkPolyDataMapper() + cutMapper.SetInputConnection(planeCut.GetOutputPort()) + cutMapper.SetScalarRange(scalar_range) + cutMapper.SetLookupTable(lut) + + Actor = vtk.vtkActor() + Actor.SetMapper(cutMapper) + return Actor From d7680512d79c6acf2b338d507baf7b1e266c76d3 Mon Sep 17 00:00:00 2001 From: Roald Storm Date: Tue, 13 Feb 2018 12:28:51 +0100 Subject: [PATCH 2/4] Made an additional function that gets a pandas dataframe from the .vti format, also stored the oriring in the .vti files --- .../VTKbased/cutoffVisualizationTool.py | 2 - wind_tools/VTKbased/flowfieldWriter.py | 70 +++++++++++++++++-- wind_tools/flow/flow_field.py | 16 +++-- 3 files changed, 73 insertions(+), 15 deletions(-) diff --git a/wind_tools/VTKbased/cutoffVisualizationTool.py b/wind_tools/VTKbased/cutoffVisualizationTool.py index 1dd62a3..d83ac80 100644 --- a/wind_tools/VTKbased/cutoffVisualizationTool.py +++ b/wind_tools/VTKbased/cutoffVisualizationTool.py @@ -172,5 +172,3 @@ def instantiateVTKviewer(fileLoc): def CheckAbort(obj, event): if obj.GetEventPending() != 0: obj.SetAbortRender(1) - -#cutoffInterface('C:/Users/roald/Dropbox/cases_for_roald/sowfaCases/CT_test_U8_pitch_00_tilt_+10/array.mean/20600.4/array.mean0D_UAvg.vtk') diff --git a/wind_tools/VTKbased/flowfieldWriter.py b/wind_tools/VTKbased/flowfieldWriter.py index 9e1f6d9..4592e74 100644 --- a/wind_tools/VTKbased/flowfieldWriter.py +++ b/wind_tools/VTKbased/flowfieldWriter.py @@ -7,16 +7,19 @@ import vtk +import numpy as np +from pandas import DataFrame -def create_vti_from_df(df_flow, spacing, dimensions, fileName): +def create_vti_from_df(fileName, df, spacing, dimensions, origin): """Create a .vti fileformat for use with the vtk toolbox input: - df_flow: a flow dataframe from flow_field, can be SOWFA or FLORIS - spacing: The spacing in the x, y and z direction between points - dimensions: a tuple or ist of integers with the number of points in xyz 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: - @@ -25,6 +28,7 @@ def create_vti_from_df(df_flow, spacing, dimensions, fileName): 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 @@ -32,14 +36,66 @@ def create_vti_from_df(df_flow, spacing, dimensions, fileName): for y in range(int(dimensions[1])): for x in range(int(dimensions[0])): imageData.SetScalarComponentFromDouble( - x, y, z, 0, df_flow.u[i]) + x, y, z, 0, df.u[i]) imageData.SetScalarComponentFromDouble( - x, y, z, 1, df_flow.v[i]) + x, y, z, 1, df.v[i]) imageData.SetScalarComponentFromDouble( - x, y, z, 2, df_flow.w[i]) + 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 diff --git a/wind_tools/flow/flow_field.py b/wind_tools/flow/flow_field.py index 02e5fbe..f49ef4d 100644 --- a/wind_tools/flow/flow_field.py +++ b/wind_tools/flow/flow_field.py @@ -20,10 +20,11 @@ def get_flow_file(case_folder): """Given a case folder, find the flow file - input: case_folder: name of case folder + input + case_folder: name of case folder output: - flow_file: full path name of flow file""" + flow_file: full path name of flow file""" array_folder = os.path.join(case_folder,'array.mean') time_folder = os.path.join(array_folder,os.listdir(array_folder)[0]) @@ -40,7 +41,9 @@ def read_flow_frame_SOWFA(filename): input: filename: name of flow array to open output: - df: a pandas table with the columns, x,y,z,u,v,w of all relavent flow info + 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 Paul Fleming, 2018 """ @@ -57,6 +60,7 @@ def read_flow_frame_SOWFA(filename): if 'ORIGIN' in read_data: origin = tuple([float(d) for d in read_data.rstrip().split(' ')[1:]]) + dimensions = tuple(int(x) for x in dimensions) # Set up x, y, z as lists xRange = np.arange(0,dimensions[0]*spacing[0],spacing[0]) yRange = np.arange(0,dimensions[1]*spacing[1],spacing[1]) @@ -75,11 +79,11 @@ def get_flow_frame_FLORIS(floris): """Read flow array output from SOWFA - input: floris: a floris model which has already been run to extract flow from + input: + floris: a floris model which has already been run to extract flow from output: - df: a pandas table with the columns, x,y,z,u,v,w of all relavent flow info - origin: the origin of the flow field, for reconstructing turbine coords + df: a pandas table with the columns, x,y,z,u,v,w of all relevant flow info Paul Fleming, 2018 """ From 1e3aa54fe15e4d955567b0d11cc2559660789cc9 Mon Sep 17 00:00:00 2001 From: Roald Storm Date: Wed, 14 Feb 2018 17:11:36 +0100 Subject: [PATCH 3/4] succesfully integrated vtk with qt5 which is callable from jupyter --- .../VTKbased/cutoffVisualizationTool.py | 58 +++++--------- wind_tools/VTKbased/slicerTool.py | 78 +++++++++---------- wind_tools/flow/flow_field.py | 10 +-- 3 files changed, 60 insertions(+), 86 deletions(-) diff --git a/wind_tools/VTKbased/cutoffVisualizationTool.py b/wind_tools/VTKbased/cutoffVisualizationTool.py index d83ac80..378efd8 100644 --- a/wind_tools/VTKbased/cutoffVisualizationTool.py +++ b/wind_tools/VTKbased/cutoffVisualizationTool.py @@ -5,32 +5,27 @@ 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.renderWindow, self.oTF, self.scalar_range, self.scalar_bar = instantiateVTKviewer(fileLoc) - self.initUI() - - # Make the vtk application interactive as well - iren = vtk.vtkRenderWindowInteractor() - iren.SetInteractorStyle(vtk.vtkInteractorStyleTrackballCamera()) - iren.SetRenderWindow(self.renderWindow) + self.ren, self.oTF, self.scalar_range, self.scalar_bar = instantiateVTKviewer(fileLoc) - # create the scalar_bar_widget - scalar_bar_widget = vtk.vtkScalarBarWidget() - scalar_bar_widget.SetInteractor(iren) - scalar_bar_widget.SetScalarBarActor(self.scalar_bar) - scalar_bar_widget.On() + self.initUI() + self.rW = self.vtkWidget.GetRenderWindow() + self.iren = self.rW.GetInteractor() + self.rW.AddRenderer(self.ren) - iren.Initialize() - iren.Start() + # 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() @@ -55,9 +50,11 @@ def initUI(self): sldTransparancy.valueChanged[int].connect(self.changeOpacity) layout.addWidget(sldTransparancy) - self.setGeometry(300, 300, 250, 200) + self.vtkWidget = QVTKRenderWindowInteractor() + layout.addWidget(self.vtkWidget) + + self.setGeometry(300, 100, 1000, 800) self.setWindowTitle('Slicer Interface') - self.show() def changeOpacityBoundary(self, value): self.OpacB = (self.scalar_range[0] + @@ -76,24 +73,16 @@ def updateLut(self): self.oTF.AddPoint(self.OpacB-.01, self.OpacV) self.oTF.AddPoint(self.OpacB+.01, 0) self.oTF.AddPoint(self.scalar_range[1], 0) - self.renderWindow.Render() + self.rW.Render() def instantiateVTKviewer(fileLoc): - # Create the standard renderer, render window and interactor - ren = vtk.vtkRenderer() - renWin = vtk.vtkRenderWindow() - renWin.AddRenderer(ren) - iren = vtk.vtkRenderWindowInteractor() - iren.SetRenderWindow(renWin) - # Create the reader for the data reader = vtk.vtkXMLImageDataReader() reader.SetFileName(fileLoc) reader.Update() flowField = reader.GetOutput() scalar_range = flowField.GetScalarRange() -# scalar_range = (2.0, scalar_range[1]) # Create a custom lut, it's used both as a mapper and at the scalar_bar lut = vtk.vtkLookupTable() @@ -102,7 +91,6 @@ def instantiateVTKviewer(fileLoc): # create the scalar_bar scalar_bar = vtk.vtkScalarBarActor() - scalar_bar.SetOrientationToHorizontal() scalar_bar.SetLookupTable(lut) scalar_bar.SetNumberOfLabels(6) scalar_bar.GetLabelTextProperty().SetFontFamilyToCourier() @@ -156,19 +144,9 @@ def instantiateVTKviewer(fileLoc): renderer = vtk.vtkRenderer() renderer.AddVolume(volume) - renderer.SetBackground(1, 1, 1) + renderer.AddActor(scalar_bar) + renderer.SetBackground(240/255, 240/255, 240/255) renderer.SetActiveCamera(camera) renderer.ResetCamera() - # Set renderingwindow and render for the first time - renderWindow = vtk.vtkRenderWindow() - renderWindow.AddRenderer(renderer) - renderWindow.Render() - renderWindow.AddObserver('AbortCheckEvent', CheckAbort) - - return renderWindow, opacityTransferFunction, scalar_range, scalar_bar - - -def CheckAbort(obj, event): - if obj.GetEventPending() != 0: - obj.SetAbortRender(1) + return renderer, opacityTransferFunction, scalar_range, scalar_bar diff --git a/wind_tools/VTKbased/slicerTool.py b/wind_tools/VTKbased/slicerTool.py index fa9d00a..2010319 100644 --- a/wind_tools/VTKbased/slicerTool.py +++ b/wind_tools/VTKbased/slicerTool.py @@ -4,6 +4,7 @@ import vtk from PyQt5.QtWidgets import QWidget, QSlider, QRadioButton, QGridLayout from PyQt5.QtCore import Qt +from vtk.qt.QVTKRenderWindowInteractor import QVTKRenderWindowInteractor class slicerInterface(QWidget): @@ -11,31 +12,26 @@ def __init__(self, fileLoc): super().__init__() # Start the VTK viewer and the user interface - (self.rW, self.renderer, self.scalar_bar, self.dims, self.cellSizes, - self.reader, self.lut, self.scalar_range) = instantiateVTKviewer2(fileLoc) + (self.ren, self.scalar_bar, self.dims, self.cellSizes, self.reader, + self.lut, self.scalar_range) = instantiateVTKviewer2(fileLoc) + self.initUI() + self.rW = self.vtkWidget.GetRenderWindow() + self.iren = self.rW.GetInteractor() + self.rW.AddRenderer(self.ren) # Instantiate the dynamical slice that will be changed on user input self.dynActor = makePlaneAt((1, 0, 0), (0, 0, 0), self.reader, self.lut, self.scalar_range) - self.renderer.AddActor(self.dynActor) + self.ren.AddActor(self.dynActor) self.rButtonMap = {'X': 0, 'Y': 1, 'Z': 2} self.sliceVal = 0 self.ax = 0 - # Make the vtk application interactive as well - iren = vtk.vtkRenderWindowInteractor() - iren.SetInteractorStyle(vtk.vtkInteractorStyleTrackballCamera()) - iren.SetRenderWindow(self.rW) - - # create the scalar_bar_widget - scalar_bar_widget = vtk.vtkScalarBarWidget() - scalar_bar_widget.SetInteractor(iren) - scalar_bar_widget.SetScalarBarActor(self.scalar_bar) - scalar_bar_widget.On() - - iren.Initialize() - iren.Start() + # Change the VTK element interaction style and show the GUI + self.iren.SetInteractorStyle(vtk.vtkInteractorStyleTrackballCamera()) + self.iren.Initialize() + self.show() def initUI(self): Grid = QGridLayout() @@ -60,9 +56,11 @@ def initUI(self): sld.valueChanged[int].connect(self.changeSliderValue) Grid.addWidget(sld, 1, 0, 1, 3) - self.setGeometry(300, 300, 250, 150) + self.vtkWidget = QVTKRenderWindowInteractor() + Grid.addWidget(self.vtkWidget, 2, 0, 1, 3) + + self.setGeometry(300, 100, 1000, 800) self.setWindowTitle('Slicer Interface') - self.show() def changeSliderValue(self, value): self.sliceVal = value/100 @@ -73,14 +71,14 @@ def on_radio_button_toggled(self): self.changeSlice() def changeSlice(self): - self.renderer.RemoveActor(self.dynActor) - d = [0, 0, 0] - d[self.ax] = 1 - l = [0, 0, 0] - l[self.ax] = self.dims[self.ax]*self.cellSizes[self.ax]*self.sliceVal - self.dynActor = makePlaneAt(d, l, self.reader, + self.ren.RemoveActor(self.dynActor) + Nor = [0, 0, 0] + Nor[self.ax] = 1 + Ori = [0, 0, 0] + Ori[self.ax] = self.dims[self.ax]*self.cellSizes[self.ax]*self.sliceVal + self.dynActor = makePlaneAt(Nor, Ori, self.reader, self.lut, self.scalar_range) - self.renderer.AddActor(self.dynActor) + self.ren.AddActor(self.dynActor) self.rW.Render() @@ -92,8 +90,8 @@ def instantiateVTKviewer2(fileLoc): flowField = reader.GetOutput() scalar_range = flowField.GetScalarRange() - # Create a custom lut. The lut is used both at the mapper and at the - # scalar_bar + # Create a custom LookUpTable. + # The lut is used both at the mapper and at the scalar_bar lut = vtk.vtkLookupTable() lut.SetTableRange(scalar_range) lut.Build() @@ -103,7 +101,6 @@ def instantiateVTKviewer2(fileLoc): # create the scalar_bar scalar_bar = vtk.vtkScalarBarActor() - scalar_bar.SetOrientationToHorizontal() scalar_bar.SetLookupTable(lut) scalar_bar.SetNumberOfLabels(8) scalar_bar.GetLabelTextProperty().SetFontFamilyToCourier() @@ -115,11 +112,14 @@ def instantiateVTKviewer2(fileLoc): scalar_bar.GetLabelTextProperty().SetColor(0, 0, 0) # Create the three planes at the back sides of the volume - cutActor1 = makePlaneAt((1, 0, 0), (dims[0]*cellSizes[0]*.99, 0, 0), + origin = flowField.GetOrigin() + box = [o + d*s for o, d, s in zip(origin, dims, cellSizes)] + + cutActor1 = makePlaneAt((1, 0, 0), (box[0]*.99, 0, 0), reader, lut, scalar_range) - cutActor2 = makePlaneAt((0, 1, 0), (0, dims[1]*cellSizes[1]*.99, 0), + cutActor2 = makePlaneAt((0, 1, 0), (0, box[1]*.99, 0), reader, lut, scalar_range) - cutActor3 = makePlaneAt((0, 0, 1), (0, 0, dims[2]*cellSizes[2]*.01), + cutActor3 = makePlaneAt((0, 0, 1), (0, 0, box[2]*.01), reader, lut, scalar_range) # Setup camera @@ -144,25 +144,21 @@ def instantiateVTKviewer2(fileLoc): renderer.AddActor(cutActor1) renderer.AddActor(cutActor2) renderer.AddActor(cutActor3) - renderer.SetBackground(1, 1, 1) + renderer.AddActor(scalar_bar) + renderer.SetBackground(240/255, 240/255, 240/255) renderer.SetActiveCamera(camera) renderer.ResetCamera() renderer.AddLight(light1) renderer.AddLight(light2) renderer.AddLight(light3) - # Set renderingwindow and render for the first time - renderWindow = vtk.vtkRenderWindow() - renderWindow.AddRenderer(renderer) - renderWindow.Render() - return (renderWindow, renderer, scalar_bar, dims, cellSizes, - reader, lut, scalar_range) + return (renderer, scalar_bar, dims, cellSizes, reader, lut, scalar_range) -def makePlaneAt(d, l, reader, lut, scalar_range): +def makePlaneAt(Nor, Ori, reader, lut, scalar_range): plane = vtk.vtkPlane() - plane.SetOrigin(l) - plane.SetNormal(d) + plane.SetOrigin(Ori) + plane.SetNormal(Nor) planeCut = vtk.vtkCutter() planeCut.SetInputConnection(reader.GetOutputPort()) diff --git a/wind_tools/flow/flow_field.py b/wind_tools/flow/flow_field.py index f49ef4d..032b759 100644 --- a/wind_tools/flow/flow_field.py +++ b/wind_tools/flow/flow_field.py @@ -16,7 +16,7 @@ import pandas as pd -def get_flow_file(case_folder): +def get_flow_file(case_folder, fileType = ".vtk"): """Given a case folder, find the flow file @@ -28,8 +28,8 @@ def get_flow_file(case_folder): array_folder = os.path.join(case_folder,'array.mean') time_folder = os.path.join(array_folder,os.listdir(array_folder)[0]) - flow_file = os.path.join(time_folder,os.listdir(time_folder)[0]) - flow_file = os.path.join(time_folder,os.listdir(time_folder)[0]) + fileNames = [fi for fi in os.listdir(time_folder) if fi.endswith(fileType)] + flow_file = os.path.join(time_folder,fileNames[0]) return flow_file @@ -56,11 +56,11 @@ def read_flow_frame_SOWFA(filename): if 'SPACING' in read_data: spacing = tuple([float(d) for d in read_data.rstrip().split(' ')[1:]]) if 'DIMENSIONS' in read_data: - dimensions = tuple([float(d) for d in read_data.rstrip().split(' ')[1:]]) + dims = tuple([float(d) for d in read_data.rstrip().split(' ')[1:]]) + dimensions = tuple(int(x) for x in dims) if 'ORIGIN' in read_data: origin = tuple([float(d) for d in read_data.rstrip().split(' ')[1:]]) - dimensions = tuple(int(x) for x in dimensions) # Set up x, y, z as lists xRange = np.arange(0,dimensions[0]*spacing[0],spacing[0]) yRange = np.arange(0,dimensions[1]*spacing[1],spacing[1]) From b605e1f804399d878d794ef404085b59f17f6157 Mon Sep 17 00:00:00 2001 From: Roald Storm Date: Mon, 26 Feb 2018 16:11:22 +0100 Subject: [PATCH 4/4] Added the option to read .vtk files that use polydata or unstructured grid data --- wind_tools/VTKbased/flowfieldWriter.py | 120 +++++++++++++++++++++++++ 1 file changed, 120 insertions(+) diff --git a/wind_tools/VTKbased/flowfieldWriter.py b/wind_tools/VTKbased/flowfieldWriter.py index 4592e74..e459121 100644 --- a/wind_tools/VTKbased/flowfieldWriter.py +++ b/wind_tools/VTKbased/flowfieldWriter.py @@ -99,3 +99,123 @@ def create_df_from_vti(fileName): 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