/[escript]/trunk/pyvisi/pyvisi/renderers/vtk/isosurface_plot.py
ViewVC logotype

Contents of /trunk/pyvisi/pyvisi/renderers/vtk/isosurface_plot.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 422 - (show annotations)
Fri Jan 6 03:10:54 2006 UTC (13 years, 5 months ago) by cochrane
File MIME type: text/x-python
File size: 16522 byte(s)
Implemented plotting of escript Data objects.  This only works with the vtk
renderer module at this stage.


1 # Copyright (C) 2004-2005 Paul Cochrane
2 #
3 # This program is free software; you can redistribute it and/or
4 # modify it under the terms of the GNU General Public License
5 # as published by the Free Software Foundation; either version 2
6 # of the License, or (at your option) any later version.
7 #
8 # This program is distributed in the hope that it will be useful,
9 # but WITHOUT ANY WARRANTY; without even the implied warranty of
10 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 # GNU General Public License for more details.
12 #
13 # You should have received a copy of the GNU General Public License
14 # along with this program; if not, write to the Free Software
15 # Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
16
17 # $Id: isosurface_plot.py,v 1.2 2006/01/03 08:22:11 paultcochrane Exp $
18
19 """
20 Class and functions associated with a pyvisi IsosurfacePlot objects
21 """
22
23 # generic imports
24 from pyvisi.renderers.vtk.common import debugMsg
25 import Numeric
26 import os
27 import copy
28
29 # module specific imports
30 from pyvisi.renderers.vtk.plot import Plot
31
32 __revision__ = '$Revision: 1.2 $'
33
34 class IsosurfacePlot(Plot):
35 """
36 Isosurface plot
37 """
38 def __init__(self, scene):
39 """
40 Initialisation of the IsosurfacePlot class
41
42 @param scene: The Scene to render the plot in
43 @type scene: Scene object
44 """
45 debugMsg("Called IsosurfacePlot.__init__()")
46 Plot.__init__(self, scene)
47
48 self.renderer = scene.renderer
49 self.renderer.addToInitStack("# IsosurfacePlot.__init__()")
50
51 # labels and stuff
52 self.title = None
53 self.xlabel = None
54 self.ylabel = None
55 self.zlabel = None
56
57 # how many contours?
58 self.numContours = 5
59
60 # contour range
61 self.contMin = None
62 self.contMax = None
63
64 # default values for fname, format and scalars
65 self.fname = None
66 self.format = None
67 self.scalars = None
68
69 self.escriptData = False
70 self.otherData = False
71
72 # add the plot to the scene
73 scene.add(self)
74
75 def setData(self, *dataList, **options):
76 """
77 Set data to the plot
78
79 @param dataList: List of data to set to the plot
80 @type dataList: tuple
81
82 @param options: Dictionary of keyword options to the method
83 @type options: dict
84
85 @param fname: the name of the input vtk file
86 @type fname: string
87
88 @param format: the format of the input vtk file ('vtk' or 'vtk-xml')
89 @type format: string
90
91 @param scalars: the name of the scalar data in the vtk file to use
92 @type scalars: string
93 """
94 debugMsg("Called setData() in IsosurfacePlot()")
95
96 self.renderer.runString("# IsosurfacePlot.setData()")
97
98 # process the options, if any
99 ## fname
100 if options.has_key('fname'):
101 fname = options['fname']
102 else:
103 fname = None
104 ## format
105 if options.has_key('format'):
106 format = options['format']
107 else:
108 format = None
109 ## scalars
110 if options.has_key('scalars'):
111 scalars = options['scalars']
112 else:
113 scalars = None
114
115 # we want to pass this info around
116 self.fname = fname
117 self.format = format
118 self.scalars = scalars
119
120 # reset the default values for shared info
121 self.escriptData = False
122 self.otherData = False
123
124 # do a sanity check on the inputs
125 if len(dataList) == 0 and fname is None:
126 raise ValueError, \
127 "You must specify a data list or an input filename"
128
129 if len(dataList) != 0 and fname is not None:
130 raise ValueError, \
131 "You cannot specify a data list as well as an input file"
132
133 if fname is not None and scalars is None:
134 debugMsg("No scalars specified; using default in vtk")
135
136 if fname is not None and format is None:
137 raise ValueError, "You must specify an input file format"
138
139 if fname is None and format is not None:
140 raise ValueError, "Format specified but no input filename"
141
142 # if have just a data list, check the objects passed in to see if
143 # they are escript data objects or not
144 if len(dataList) != 0:
145 for obj in dataList:
146 try:
147 obj.convertToNumArray()
148 # ok, we've got escript data, set the flag
149 self.escriptData = True
150 except AttributeError:
151 self.otherData = True
152
153 # if we have both escript data and other data, barf as can't handle
154 # that yet
155 if self.escriptData and self.otherData:
156 raise TypeError, \
157 "Sorry can't handle both escript and other data yet"
158
159 # now generate the code for the case when we have just escript data
160 # passed into setData()
161 if self.escriptData:
162 # get the relevant bits of data
163 if len(dataList) == 1:
164 # only one data variable, will need to get the domain from it
165 escriptZ = dataList[0]
166 escriptX = escriptZ.getDomain().getX()
167 else:
168 errorString = \
169 "Expecting only 1 element in data list. I got %d" \
170 % len(dataList)
171 raise ValueError, errorString
172
173 # convert the data to numarrays
174 fieldData = escriptZ.convertToNumArray()
175 domainData = escriptX.convertToNumArray()
176
177 xData = domainData[:,0]
178 yData = domainData[:,1]
179 zData = domainData[:,2]
180
181 # now check the dimensionality of the data and the grid; make
182 # sure everything looks right before doing anything more
183 if len(xData.shape) != 1:
184 raise ValueError, "xData not 1D. I got %d dims" % \
185 len(xData.shape)
186
187 if len(yData.shape) != 1:
188 raise ValueError, "yData not 1D. I got %d dims" % \
189 len(yData.shape)
190
191 if len(zData.shape) != 1:
192 raise ValueError, "zData not 1D. I got %d dims" % \
193 len(zData.shape)
194
195 if len(fieldData.shape) != 1:
196 raise ValueError, "fieldData not 1D. I got %d dims" % \
197 len(fieldData.shape)
198
199 # now check that the length of the various vectors is correct
200 dataLen = fieldData.shape[0]
201 if xData.shape[0] != dataLen:
202 raise ValueError, \
203 "xData length doesn't agree with fieldData length"
204
205 if yData.shape[0] != dataLen:
206 raise ValueError, \
207 "yData length doesn't agree with fieldData length"
208
209 if zData.shape[0] != dataLen:
210 raise ValueError, \
211 "zData length doesn't agree with fieldData length"
212
213 # it looks like everything is ok... share the data around
214 ### the x data
215 self.renderer.renderDict['_x'] = copy.deepcopy(xData)
216
217 ### the y data
218 self.renderer.renderDict['_y'] = copy.deepcopy(yData)
219
220 ### the z data
221 self.renderer.renderDict['_z'] = copy.deepcopy(zData)
222
223 ### the field data
224 self.renderer.renderDict['_field'] = copy.deepcopy(fieldData)
225
226 ### construct the grid, data and points arrays
227
228 # first the points array
229 evalString = "_points = vtk.vtkPoints()\n"
230 evalString += "_points.SetNumberOfPoints(%d)\n" % dataLen
231 evalString += "for _i in range(%d):\n" % dataLen
232 evalString += " _points.SetPoint(_i, _x[_i], _y[_i], _z[_i])\n"
233 self.renderer.runString(evalString)
234
235 # now the data array
236 evalString = "_data = vtk.vtkFloatArray()\n"
237 evalString += "_data.SetNumberOfValues(%d)\n" % dataLen
238 evalString += "_data.SetNumberOfTuples(1)\n"
239 evalString += "for _i in range(%d):\n" % dataLen
240 evalString += " _data.SetTuple1(_i, _field[_i])\n"
241 self.renderer.runString(evalString)
242
243 # now make the grid
244 evalString = "_grid = vtk.vtkUnstructuredGrid()\n"
245 evalString += "_grid.SetPoints(_points)\n"
246 evalString += "_grid.GetPointData().SetScalars(_data)\n"
247 self.renderer.runString(evalString)
248
249 # do the delaunay 3D to get proper isosurfaces
250 evalString = "_del3D = vtk.vtkDelaunay3D()\n"
251 evalString += "_del3D.SetInput(_grid)\n"
252 evalString += "_del3D.SetOffset(2.5)\n"
253 evalString += "_del3D.SetTolerance(0.001)\n"
254 evalString += "_del3D.SetAlpha(0.0)\n"
255 self.renderer.runString(evalString)
256
257 elif self.otherData:
258
259 raise ImplementationError, "Can't process plain array data yet"
260
261 # run the stuff for when we're reading from file
262 if fname is not None:
263
264 # check to see if the file exists
265 if not os.path.exists(fname):
266 raise SystemError, "File %s not found" % fname
267
268 if format == 'vtk':
269 # read old-style vtk files
270 evalString = "_reader = vtk.vtkUnstructuredGridReader()\n"
271 elif format == 'vtk-xml':
272 # read vtk xml files
273 evalString = "_reader = vtk.vtkXMLUnstructuredGridReader()\n"
274 else:
275 # barf
276 raise ValueError, "Unknown format. I got %s" % format
277
278 evalString += "_reader.SetFileName(\"%s\")\n" % fname
279 evalString += "_reader.Update()\n"
280 evalString += "_grid = _reader.GetOutput()\n"
281 self.renderer.runString(evalString)
282
283 # need to do a delaunay 3D here to get decent looking isosurfaces
284 evalString = "_del3D = vtk.vtkDelaunay3D()\n"
285 evalString += "_del3D.SetInput(_grid)\n"
286 evalString += "_del3D.SetOffset(2.5)\n"
287 evalString += "_del3D.SetTolerance(0.001)\n"
288 evalString += "_del3D.SetAlpha(0.0)"
289 self.renderer.runString(evalString)
290
291 # get the model centre and bounds
292 evalString = "_centre = _grid.GetCenter()\n"
293 evalString += "_bounds = _grid.GetBounds()\n"
294 self.renderer.runString(evalString)
295
296 # set up a contour filter
297 evalString = "_cont = vtk.vtkContourGrid()\n"
298 evalString += "_cont.SetInput(_del3D.GetOutput())\n"
299
300 # if contMin and contMax are or aren't set then handle the different
301 # situations
302 if self.contMin is not None and self.contMax is not None:
303 evalString += "_cont.GenerateValues(%d, %f, %f)\n" % \
304 (self.numContours, self.contMin, self.contMax)
305 elif self.contMin is not None and self.contMax is None:
306 evalString += "(_contMin, _contMax) = _grid."
307 evalString += "GetPointData().GetScalars().GetRange()\n"
308 evalString += "_cont.GenerateValues(%d, %f, _contMax)\n" % \
309 (self.numContours, self.contMin)
310 elif self.contMin is None and self.contMax is not None:
311 evalString += "(_contMin, _contMax) = _grid."
312 evalString += "GetPointData().GetScalars().GetRange()\n"
313 evalString += "_cont.GenerateValues(%d, _contMin, %f)\n" % \
314 (self.numContours, self.contMax)
315 elif self.contMin is None and self.contMax is None:
316 evalString += "(_contMin, _contMax) = _grid."
317 evalString += "GetPointData().GetScalars().GetRange()\n"
318 evalString += "_cont.GenerateValues(%d, _contMin, _contMax)\n" % \
319 (self.numContours)
320 else:
321 # barf, really shouldn't have got here
322 raise ValueError, \
323 "Major problems in IsosurfacePlot: contMin and contMax"
324
325 evalString += "_cont.GenerateValues(5, 0.25, 0.75)\n"
326 evalString += "_cont.ComputeScalarsOn()"
327 self.renderer.runString(evalString)
328
329 def render(self):
330 """
331 Does IsosurfacePlot object specific (pre)rendering stuff
332 """
333 debugMsg("Called IsosurfacePlot.render()")
334
335 self.renderer.runString("# IsosurfacePlot.render()")
336
337 # set up the mapper
338 evalString = "_mapper = vtk.vtkDataSetMapper()\n"
339 evalString += "_mapper.SetInput(_cont.GetOutput())\n"
340 evalString += "_mapper.ScalarVisibilityOn()"
341 self.renderer.runString(evalString)
342
343 # set up the actor
344 evalString = "_actor = vtk.vtkActor()\n"
345 evalString += "_actor.SetMapper(_mapper)"
346 self.renderer.runString(evalString)
347
348 # add to the renderer
349 evalString = "_renderer.AddActor(_actor)"
350 self.renderer.runString(evalString)
351
352 # set up the text properties for nice text
353 evalString = "_font_size = 18\n"
354 evalString += "_textProp = vtk.vtkTextProperty()\n"
355 evalString += "_textProp.SetFontSize(_font_size)\n"
356 evalString += "_textProp.SetFontFamilyToArial()\n"
357 evalString += "_textProp.BoldOff()\n"
358 evalString += "_textProp.ItalicOff()\n"
359 evalString += "_textProp.ShadowOff()\n"
360 evalString += "_textProp.SetColor(0.0, 0.0, 0.0)"
361 self.renderer.runString(evalString)
362
363 # if a title is set, put it in here
364 if self.title is not None:
365 # make a title
366 evalString = "_title = vtk.vtkTextMapper()\n"
367 evalString += "_title.SetInput(\"%s\")\n" % self.title
368
369 # make the title text use the text properties
370 evalString += "_titleProp = _title.GetTextProperty()\n"
371 evalString += "_titleProp.ShallowCopy(_textProp)\n"
372 evalString += "_titleProp.SetJustificationToCentered()\n"
373 evalString += "_titleProp.SetVerticalJustificationToTop()\n"
374
375 # make the actor for the title
376 evalString += "_titleActor = vtk.vtkTextActor()\n"
377 evalString += "_titleActor.SetMapper(_title)\n"
378 evalString += "_titleActor.GetPositionCoordinate()."
379 evalString += "SetCoordinateSystemToNormalizedDisplay()\n"
380 evalString += "_titleActor.GetPositionCoordinate()."
381 evalString += "SetValue(0.5, 0.95)"
382 self.renderer.runString(evalString)
383
384 # add to the renderer
385 evalString = "_renderer.AddActor(_titleActor)"
386 self.renderer.runString(evalString)
387
388 # put an outline around the data
389 evalString = "_outline = vtk.vtkOutlineSource()\n"
390 evalString += "_outline.SetBounds(_bounds)\n"
391
392 # make its mapper
393 evalString += "_outlineMapper = vtk.vtkPolyDataMapper()\n"
394 evalString += "_outlineMapper.SetInput(_outline.GetOutput())\n"
395
396 # make its actor
397 evalString += "_outlineActor = vtk.vtkActor()\n"
398 evalString += "_outlineActor.SetMapper(_outlineMapper)\n"
399 evalString += "_outlineActor.GetProperty().SetColor(0,0,0)"
400 self.renderer.runString(evalString)
401
402 # add to the renderer
403 evalString = "_renderer.AddActor(_outlineActor)"
404 self.renderer.runString(evalString)
405
406 # make a lookup table for the colour map and invert it (colours look
407 # better when it's inverted)
408 evalString = "_lut = vtk.vtkLookupTable()\n"
409 evalString += "_refLut = vtk.vtkLookupTable()\n"
410 evalString += "_lut.Build()\n"
411 evalString += "_refLut.Build()\n"
412 evalString += "for _j in range(256):\n"
413 evalString += " _lut.SetTableValue(_j, "
414 evalString += "_refLut.GetTableValue(255-_j))"
415 self.renderer.runString(evalString)
416
417 # add some axes
418 evalString = "_axes = vtk.vtkCubeAxesActor2D()\n"
419 evalString += "_axes.SetInput(_grid)\n"
420 evalString += "_axes.SetCamera(_renderer.GetActiveCamera())\n"
421 evalString += "_axes.SetLabelFormat(\"%6.4g\")\n"
422 evalString += "_axes.SetFlyModeToOuterEdges()\n"
423 evalString += "_axes.SetFontFactor(0.8)\n"
424 evalString += "_axes.SetAxisTitleTextProperty(_textProp)\n"
425 evalString += "_axes.SetAxisLabelTextProperty(_textProp)\n"
426 ### xlabel
427 if self.xlabel is not None:
428 evalString += "_axes.SetXLabel(\"%s\")\n" % self.xlabel
429 else:
430 evalString += "_axes.SetXLabel(\"\")\n"
431 ### ylabel
432 if self.ylabel is not None:
433 evalString += "_axes.SetYLabel(\"%s\")\n" % self.ylabel
434 else:
435 evalString += "_axes.SetYLabel(\"\")\n"
436 ### zlabel
437 if self.zlabel is not None:
438 evalString += "_axes.SetZLabel(\"%s\")\n" % self.zlabel
439 else:
440 evalString += "_axes.SetZLabel(\"\")\n"
441 evalString += "_axes.SetNumberOfLabels(5)\n"
442 evalString += "_axes.GetProperty().SetColor(0,0,0)"
443 self.renderer.runString(evalString)
444
445 # add to the renderer
446 evalString = "_renderer.AddActor(_axes)"
447 self.renderer.runString(evalString)
448
449 # play around with lighting
450 evalString = "_light1 = vtk.vtkLight()\n"
451 evalString += "_light1.SetFocalPoint(_centre)\n"
452 evalString += "_light1.SetPosition(_centre[0]-_bounds[1], "
453 evalString += "_centre[1]-_bounds[3], _centre[2]+_bounds[5])\n"
454 evalString += "_renderer.AddLight(_light1)\n"
455 evalString += "_light2 = vtk.vtkLight()\n"
456 evalString += "_light2.SetFocalPoint(_centre)\n"
457 evalString += "_light2.SetPosition(_centre[0]+_bounds[1], "
458 evalString += "_centre[1]+_bounds[3], _centre[2]-_bounds[5])\n"
459 evalString += "_renderer.AddLight(_light2)"
460 self.renderer.runString(evalString)
461
462 # this shouldn't be here!!!!
463 self.renderer.runString("_renderer.SetBackground(1,1,1)")
464
465 return
466
467
468 # vim: expandtab shiftwidth=4:

  ViewVC Help
Powered by ViewVC 1.1.26