Package PyFoam :: Package Applications :: Module DisplayBlockMesh
[hide private]
[frames] | no frames]

Source Code for Module PyFoam.Applications.DisplayBlockMesh

  1  #!/usr/bin/env python 
  2   
  3  import sys 
  4   
  5  from PyFoam.RunDictionary.ParsedBlockMeshDict import ParsedBlockMeshDict 
  6  from PyFoam.Applications.PyFoamApplication import PyFoamApplication 
  7  from PyFoam.Error import error,warning 
  8   
9 -def doImports():
10 try: 11 global Tkinter 12 import Tkinter 13 global vtk 14 try: 15 from paraview import vtk 16 print "Using VTK implementation from Paraview" 17 except ImportError: 18 print "Trying system-VTK" 19 import vtk 20 global vtkTkRenderWindowInteractor 21 from vtk.tk.vtkTkRenderWindowInteractor import vtkTkRenderWindowInteractor 22 except ImportError,e: 23 error("Error while importing modules:",e)
24
25 -class DisplayBlockMesh(PyFoamApplication):
26 - def __init__(self):
27 description=""" 28 Reads the contents of a blockMeshDict-file and displays the 29 vertices as spheres (with numbers). The blocks are sketched by 30 lines. One block can be seceted with a slider. It will be 31 displayed as a green cube with the local directions x1,x2 and 32 x3. Also a patch that is selected by a slider will be sketched 33 by blue squares 34 """ 35 PyFoamApplication.__init__(self,description=description,usage="%prog [options] <blockMeshDict>",interspersed=True,nr=1)
36
37 - def run(self):
38 doImports() 39 40 self.renWin = vtk.vtkRenderWindow() 41 self.ren = vtk.vtkRenderer() 42 self.renWin.AddRenderer(self.ren) 43 self.renWin.SetSize(600, 600) 44 self.ren.SetBackground(0.7, 0.7, 0.7) 45 self.ren.ResetCamera() 46 self.cam = self.ren.GetActiveCamera() 47 48 self.axes = vtk.vtkCubeAxesActor2D() 49 self.axes.SetCamera(self.ren.GetActiveCamera()) 50 51 self.undefinedActor=vtk.vtkTextActor() 52 self.undefinedActor.GetPositionCoordinate().SetCoordinateSystemToNormalizedDisplay() 53 self.undefinedActor.GetPositionCoordinate().SetValue(0.05,0.2) 54 self.undefinedActor.GetTextProperty().SetColor(1.,0.,0.) 55 self.undefinedActor.SetInput("") 56 57 try: 58 self.readFile() 59 except Exception,e: 60 warning("While reading",self.parser.getArgs()[0],"this happened:",e) 61 raise e 62 63 self.ren.ResetCamera() 64 65 self.root = Tkinter.Tk() 66 self.root.withdraw() 67 68 # Create the toplevel window 69 self.top = Tkinter.Toplevel(self.root) 70 self.top.title("blockMesh-Viewer") 71 self.top.protocol("WM_DELETE_WINDOW", self.quit) 72 73 # Create some frames 74 self.f1 = Tkinter.Frame(self.top) 75 self.f2 = Tkinter.Frame(self.top) 76 77 self.f1.pack(side="top", anchor="n", expand=1, fill="both") 78 self.f2.pack(side="bottom", anchor="s", expand="f", fill="x") 79 80 # Create the Tk render widget, and bind the events 81 self.rw = vtkTkRenderWindowInteractor(self.f1, width=400, height=400, rw=self.renWin) 82 self.rw.pack(expand="t", fill="both") 83 84 self.blockHigh=Tkinter.IntVar() 85 self.blockHigh.set(-1) 86 87 self.oldBlock=-1 88 self.blockActor=None 89 self.blockTextActor=None 90 91 self.patchHigh=Tkinter.IntVar() 92 self.patchHigh.set(-1) 93 94 self.oldPatch=-1 95 self.patchActor=None 96 self.patchTextActor=vtk.vtkTextActor() 97 self.patchTextActor.GetPositionCoordinate().SetCoordinateSystemToNormalizedDisplay() 98 self.patchTextActor.GetPositionCoordinate().SetValue(0.05,0.1) 99 self.patchTextActor.GetTextProperty().SetColor(0.,0.,0.) 100 self.patchTextActor.SetInput("Patch: <none>") 101 102 self.scroll=Tkinter.Scale(self.f2,orient='horizontal', 103 from_=-1,to=len(self.blocks)-1,resolution=1,tickinterval=1, 104 label="Block (-1 is none)", 105 variable=self.blockHigh,command=self.colorBlock) 106 107 self.scroll.pack(side="top", expand="t", fill="x") 108 109 self.scroll2=Tkinter.Scale(self.f2,orient='horizontal', 110 from_=-1,to=len(self.patches.keys())-1,resolution=1,tickinterval=1, 111 label="Patch (-1 is none)", 112 variable=self.patchHigh,command=self.colorPatch) 113 114 self.scroll2.pack(side="top", expand="t", fill="x") 115 116 self.f3 = Tkinter.Frame(self.f2) 117 self.f3.pack(side="bottom", anchor="s", expand="f", fill="x") 118 119 self.b1 = Tkinter.Button(self.f3, text="Quit", command=self.quit) 120 self.b1.pack(side="left", expand="t", fill="x") 121 self.b2 = Tkinter.Button(self.f3, text="Reread blockMeshDict", command=self.reread) 122 self.b2.pack(side="left", expand="t", fill="x") 123 124 self.root.update() 125 126 self.iren = self.renWin.GetInteractor() 127 self.istyle = vtk.vtkInteractorStyleSwitch() 128 129 self.iren.SetInteractorStyle(self.istyle) 130 self.istyle.SetCurrentStyleToTrackballCamera() 131 132 self.addProps() 133 134 self.iren.Initialize() 135 self.renWin.Render() 136 self.iren.Start() 137 138 self.root.mainloop()
139
140 - def readFile(self):
141 self.blockMesh=ParsedBlockMeshDict(self.parser.getArgs()[0]) 142 143 self.vertices=self.blockMesh.vertices() 144 self.vActors=[None]*len(self.vertices) 145 146 self.blocks=self.blockMesh.blocks() 147 self.patches=self.blockMesh.patches() 148 149 self.vRadius=self.blockMesh.typicalLength()/50 150 151 for i in range(len(self.vertices)): 152 self.addVertex(i) 153 154 self.setAxes() 155 156 self.undefined=[] 157 158 for i in range(len(self.blocks)): 159 self.addBlock(i) 160 161 for a in self.blockMesh.arcs(): 162 self.makeArc(a) 163 164 if len(self.undefined)>0: 165 self.undefinedActor.SetInput("Undefined vertices: "+str(self.undefined)) 166 else: 167 self.undefinedActor.SetInput("")
168
169 - def addUndefined(self,i):
170 if not i in self.undefined: 171 self.undefined.append(i)
172
173 - def addProps(self):
174 self.ren.AddViewProp(self.axes) 175 self.ren.AddActor2D(self.patchTextActor) 176 self.ren.AddActor2D(self.undefinedActor)
177
178 - def addPoint(self,coord,factor=1):
179 sphere=vtk.vtkSphereSource() 180 sphere.SetRadius(self.vRadius*factor) 181 sphere.SetCenter(coord) 182 mapper=vtk.vtkPolyDataMapper() 183 mapper.SetInputConnection(sphere.GetOutputPort()) 184 actor = vtk.vtkActor() 185 actor.SetMapper(mapper) 186 self.ren.AddActor(actor) 187 188 return actor
189
190 - def addVertex(self,index):
191 coord=self.vertices[index] 192 self.vActors[index]=self.addPoint(coord) 193 text=vtk.vtkVectorText() 194 text.SetText(str(index)) 195 tMapper=vtk.vtkPolyDataMapper() 196 tMapper.SetInput(text.GetOutput()) 197 tActor = vtk.vtkFollower() 198 tActor.SetMapper(tMapper) 199 tActor.SetScale(2*self.vRadius,2*self.vRadius,2*self.vRadius) 200 tActor.AddPosition(coord[0]+self.vRadius,coord[1]+self.vRadius,coord[2]+self.vRadius) 201 tActor.SetCamera(self.cam) 202 tActor.GetProperty().SetColor(1.0,0.,0.) 203 self.ren.AddActor(tActor)
204
205 - def addLine(self,index1,index2):
206 try: 207 c1=self.vertices[index1] 208 c2=self.vertices[index2] 209 except: 210 if index1>=len(self.vertices): 211 self.addUndefined(index1) 212 if index2>=len(self.vertices): 213 self.addUndefined(index2) 214 return None 215 line=vtk.vtkLineSource() 216 line.SetPoint1(c1) 217 line.SetPoint2(c2) 218 mapper=vtk.vtkPolyDataMapper() 219 mapper.SetInputConnection(line.GetOutputPort()) 220 actor = vtk.vtkActor() 221 actor.SetMapper(mapper) 222 self.ren.AddActor(actor) 223 return actor
224
225 - def makeDirection(self,index1,index2,label):
226 try: 227 c1=self.vertices[index1] 228 c2=self.vertices[index2] 229 except: 230 return None,None 231 line=vtk.vtkLineSource() 232 line.SetPoint1(c1) 233 line.SetPoint2(c2) 234 tube=vtk.vtkTubeFilter() 235 tube.SetRadius(self.vRadius*0.5) 236 tube.SetNumberOfSides(10) 237 tube.SetInput(line.GetOutput()) 238 text=vtk.vtkVectorText() 239 text.SetText(label) 240 tMapper=vtk.vtkPolyDataMapper() 241 tMapper.SetInput(text.GetOutput()) 242 tActor = vtk.vtkFollower() 243 tActor.SetMapper(tMapper) 244 tActor.SetScale(self.vRadius,self.vRadius,self.vRadius) 245 tActor.AddPosition((c1[0]+c2[0])/2+self.vRadius,(c1[1]+c2[1])/2+self.vRadius,(c1[2]+c2[2])/2+self.vRadius) 246 tActor.SetCamera(self.cam) 247 tActor.GetProperty().SetColor(0.0,0.,0.) 248 return tube.GetOutput(),tActor
249
250 - def makeSpline(self,lst):
251 points = vtk.vtkPoints() 252 for i in range(len(lst)): 253 v=lst[i] 254 points.InsertPoint(i,v[0],v[1],v[2]) 255 spline=vtk.vtkParametricSpline() 256 spline.SetPoints(points) 257 spline.ClosedOff() 258 splineSource=vtk.vtkParametricFunctionSource() 259 splineSource.SetParametricFunction(spline) 260 mapper=vtk.vtkPolyDataMapper() 261 mapper.SetInputConnection(splineSource.GetOutputPort()) 262 actor = vtk.vtkActor() 263 actor.SetMapper(mapper) 264 self.ren.AddActor(actor)
265
266 - def makeArc(self,data):
267 try: 268 self.makeSpline([self.vertices[data[0]],data[1],self.vertices[data[2]]]) 269 except: 270 if data[0]>=len(self.vertices): 271 self.addUndefined(data[0]) 272 if data[2]>=len(self.vertices): 273 self.addUndefined(data[2]) 274 275 self.addPoint(data[1],factor=0.5)
276
277 - def makeFace(self,lst):
278 points = vtk.vtkPoints() 279 side = vtk.vtkCellArray() 280 side.InsertNextCell(len(lst)) 281 for i in range(len(lst)): 282 try: 283 v=self.vertices[lst[i]] 284 except: 285 self.addUndefined(lst[i]) 286 return None 287 points.InsertPoint(i,v[0],v[1],v[2]) 288 side.InsertCellPoint(i) 289 result=vtk.vtkPolyData() 290 result.SetPoints(points) 291 result.SetPolys(side) 292 293 return result
294
295 - def addBlock(self,index):
296 b=self.blocks[index] 297 298 self.addLine(b[ 0],b[ 1]) 299 self.addLine(b[ 3],b[ 2]) 300 self.addLine(b[ 7],b[ 6]) 301 self.addLine(b[ 4],b[ 5]) 302 self.addLine(b[ 0],b[ 3]) 303 self.addLine(b[ 1],b[ 2]) 304 self.addLine(b[ 5],b[ 6]) 305 self.addLine(b[ 4],b[ 7]) 306 self.addLine(b[ 0],b[ 4]) 307 self.addLine(b[ 1],b[ 5]) 308 self.addLine(b[ 2],b[ 6]) 309 self.addLine(b[ 3],b[ 7])
310
311 - def setAxes(self):
312 append=vtk.vtkAppendPolyData() 313 for a in self.vActors: 314 if a!=None: 315 append.AddInput(a.GetMapper().GetInput()) 316 self.axes.SetInput(append.GetOutput())
317 318 319 # Define a quit method that exits cleanly.
320 - def quit(self):
321 self.root.quit()
322
323 - def reread(self):
324 self.ren.RemoveAllViewProps() 325 self.patchActor=None 326 self.blockActor=None 327 self.blockTextActor=None 328 self.addProps() 329 self.readFile() 330 331 tmpBlock=int(self.blockHigh.get()) 332 if not tmpBlock<len(self.blocks): 333 tmpBlock=len(self.blocks)-1 334 self.scroll.config(to=len(self.blocks)-1) 335 self.blockHigh.set(tmpBlock) 336 self.colorBlock(tmpBlock) 337 338 tmpPatch=int(self.patchHigh.get()) 339 if not tmpPatch<len(self.patches.keys()): 340 tmpPatch=len(self.patches.keys())-1 341 self.scroll2.config(to=len(self.patches.keys())-1) 342 self.patchHigh.set(tmpPatch) 343 self.colorPatch(tmpPatch) 344 345 self.renWin.Render()
346
347 - def colorBlock(self,value):
348 newBlock=int(value) 349 if self.oldBlock>=0 and self.blockActor!=None: 350 self.ren.RemoveActor(self.blockActor) 351 for ta in self.blockTextActor: 352 self.ren.RemoveActor(ta) 353 self.blockActor=None 354 self.blockTextActor=None 355 if newBlock>=0: 356 append=vtk.vtkAppendPolyData() 357 append2=vtk.vtkAppendPolyData() 358 b=self.blocks[newBlock] 359 append.AddInput(self.makeFace([b[0],b[1],b[2],b[3]])) 360 append.AddInput(self.makeFace([b[4],b[5],b[6],b[7]])) 361 append.AddInput(self.makeFace([b[0],b[1],b[5],b[4]])) 362 append.AddInput(self.makeFace([b[3],b[2],b[6],b[7]])) 363 append.AddInput(self.makeFace([b[0],b[3],b[7],b[4]])) 364 append.AddInput(self.makeFace([b[1],b[2],b[6],b[5]])) 365 d,t1=self.makeDirection(b[0],b[1],"x1") 366 append.AddInput(d) 367 self.ren.AddActor(t1) 368 d,t2=self.makeDirection(b[0],b[3],"x2") 369 append.AddInput(d) 370 self.ren.AddActor(t2) 371 d,t3=self.makeDirection(b[0],b[4],"x3") 372 append.AddInput(d) 373 self.ren.AddActor(t3) 374 self.blockTextActor=(t1,t2,t3) 375 mapper=vtk.vtkPolyDataMapper() 376 mapper.SetInputConnection(append.GetOutputPort()) 377 self.blockActor = vtk.vtkActor() 378 self.blockActor.SetMapper(mapper) 379 self.blockActor.GetProperty().SetColor(0.,1.,0.) 380 self.blockActor.GetProperty().SetOpacity(0.3) 381 self.ren.AddActor(self.blockActor) 382 383 self.oldBlock=newBlock 384 self.renWin.Render()
385
386 - def colorPatch(self,value):
387 newPatch=int(value) 388 if self.oldPatch>=0 and self.patchActor!=None: 389 self.ren.RemoveActor(self.patchActor) 390 self.patchActor=None 391 self.patchTextActor.SetInput("Patch: <none>") 392 if newPatch>=0: 393 name=self.patches.keys()[newPatch] 394 subs=self.patches[name] 395 append=vtk.vtkAppendPolyData() 396 for s in subs: 397 append.AddInput(self.makeFace(s)) 398 mapper=vtk.vtkPolyDataMapper() 399 mapper.SetInputConnection(append.GetOutputPort()) 400 self.patchActor = vtk.vtkActor() 401 self.patchActor.SetMapper(mapper) 402 self.patchActor.GetProperty().SetColor(0.,0.,1.) 403 self.patchActor.GetProperty().SetOpacity(0.3) 404 self.ren.AddActor(self.patchActor) 405 self.patchTextActor.SetInput("Patch: "+name) 406 407 self.oldPatch=newPatch 408 self.renWin.Render()
409