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