/[escript]/trunk/pyvisi/examples/ballPlotExample.py
ViewVC logotype

Contents of /trunk/pyvisi/examples/ballPlotExample.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 605 - (show annotations)
Mon Mar 20 07:16:27 2006 UTC (13 years, 7 months ago) by gross
File MIME type: text/x-python
File size: 14419 byte(s)
pyvisi with vtk runs noe on the altix
1 # Copyright (C) 2004 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: ballPlotExample.py,v 1.1 2005/05/25 05:22:32 paultcochrane Exp $
18
19 ## @file ballPlotExample.py
20
21 """
22 Example of plotting spheres with pyvisi
23 """
24
25 # what plotting method are we using?
26 method = 'pyvisi'
27
28 import os, sys
29 import random
30
31 # set up some data to plot
32 from Numeric import *
33
34 # the three axes in space
35 # this will give us 10 particles (_not_ 1000)
36 x = arange(10, typecode=Float)
37 y = arange(10, typecode=Float)
38 z = arange(10, typecode=Float)
39
40 # 3D position information
41 posArray = []
42 for i in range(len(x)):
43 for j in range(len(y)):
44 for k in range(len(z)):
45 posArray.append( (x[i], y[j], z[k]) )
46
47 # radius information
48 random.seed()
49 radiiArray = zeros(len(x)*len(y)*len(z), typecode=Float)
50 for i in range(len(x)*len(y)*len(z)):
51 radiiArray[i] = random.random()*0.8
52
53 # tag information
54 random.seed()
55 tagsArray = zeros(len(x)*len(y)*len(z), typecode=Int)
56 for i in range(len(x)*len(y)*len(z)):
57 tagsArray[i] = int(random.random()*10)
58
59 # plot it using one of the three methods
60 if method == 'pyvisi':
61
62 # example code for how a user would write a script in pyvisi
63 from esys.pyvisi import * # base level visualisation stuff
64 # import the objects to render the scene using the specific renderer
65 #from esys.pyvisi.renderers.gnuplot import * # gnuplot
66 from esys.pyvisi.renderers.vtk import * # vtk
67
68 # define the scene object
69 # a Scene is a container for all of the kinds of things you want to put
70 # into your plot for instance, images, meshes, arrow/vector/quiver plots,
71 # contour plots, spheres etc.
72 scene = Scene()
73
74 # create a BallPlot object
75 plot = BallPlot(scene)
76
77 # add some helpful info to the plot
78 plot.title = 'Example ball plot'
79
80 # assign some data to the plot
81 # one way of doing it
82 # (tags indirectly determine colour of the spheres in the plot)
83 plot.setData(points=posArray, radii=radiiArray, tags=tagsArray)
84
85 # render the scene
86 scene.render(pause=True, interactive=True)
87
88 # without specifying a tags array input
89 plot.setData(points=posArray, radii=radiiArray)
90 # render the scene
91 scene.render(pause=True, interactive=True)
92
93 # another way loading an old style-vtk file
94 plot.setData(fname="cp_test_0.vtk",
95 format="vtk",
96 radii="radius",
97 colors="particleTag")
98
99 # render the scene to screen
100 scene.render(pause=True, interactive=True)
101
102 # another way loading a vtk xml file
103 plot.setData(fname="cp_test_0.xml",
104 format="vtk-xml",
105 radii="radius",
106 colors="particleTag")
107
108 # render the scene to screen
109 scene.render(pause=True, interactive=True)
110
111 # save the scene out to file
112 ## as png
113 plot.setData(fname="cp_test_0.xml",
114 format="vtk-xml",
115 radii="radius",
116 colors="particleTag")
117 scene.save(fname="ballPlotExample.png", format="png")
118 ## as postscript
119 plot.setData(fname="cp_test_0.xml",
120 format="vtk-xml",
121 radii="radius",
122 colors="particleTag")
123 scene.save(fname="ballPlotExample.ps", format="PS")
124
125 elif method == 'vtk':
126 #### original vtk code
127
128 import vtk
129 import os
130
131 # we're loading a file; make sure it exists first
132 if os.path.exists("cp_test_0.xml"):
133 pass
134 else:
135 raise IOError, "File not found"
136
137 # create the reader of the file
138 _reader = vtk.vtkXMLUnstructuredGridReader()
139 _reader.SetFileName("cp_test_0.xml")
140 _reader.Update()
141
142 # read the output into an unstructured grid
143 _grid = _reader.GetOutput()
144
145 # grab the radius data for the radii of the balls
146 _radii = _grid.GetPointData().GetScalars("radius")
147
148 # grab the tag data and use for colouring the balls
149 _tags = _grid.GetPointData().GetScalars("particleTag")
150
151 # work out dynamically the number of different tags so that can use this
152 # information to automatically set the scalar range for colouring
153 _numPoints = _tags.GetNumberOfTuples()
154 _valueDict = {}
155 for i in range(_numPoints):
156 _tagValue = _tags.GetValue(i)
157 _valueDict[_tagValue] = 1
158
159 _numTags = len(_valueDict.keys())
160
161 _tagValues = _valueDict.keys()
162 _tagValues.sort()
163
164 # count the number of tags, and make an evenly spaced array of points
165 # between zero and one, then use these as the scalars to colour by
166 _scaledTags = vtk.vtkFloatArray()
167 _scaledTags.SetNumberOfTuples(_numPoints)
168 _scaledTags.SetNumberOfComponents(1)
169 _scaledTags.SetName("tags")
170 for i in range(_numPoints):
171 _tagValue = _tags.GetValue(i)
172 for j in range(_numTags):
173 if _tagValues[j] == _tagValue:
174 _scaledTags.InsertTuple1(i, float(j)/float(_numTags-1))
175
176 # in vtk 4.2 have to set up an array of two components to get
177 # the data through the glyph object to the mapper so create this now
178 _data = vtk.vtkFloatArray()
179 _data.SetNumberOfComponents(3)
180 _data.SetNumberOfTuples(_radii.GetNumberOfTuples())
181 _data.CopyComponent(0, _radii, 0)
182 _data.CopyComponent(1, _tags, 0)
183 _data.CopyComponent(2, _scaledTags, 0)
184 _data.SetName("data")
185
186 # add the data array to the grid (again, only specific to vtk 4.2 and 4.4)
187 _grid.GetPointData().AddArray(_data)
188
189 # make the data the active scalars
190 _grid.GetPointData().SetActiveScalars("data")
191
192 # to make sphere glyphs need a sphere source
193 _sphere = vtk.vtkSphereSource()
194 _sphere.SetRadius(1.0) # set to 0.5 by default in vtk, we need it to be 1.0
195 _sphere.SetThetaResolution(5) # chunky, but will be good in large models
196 _sphere.SetPhiResolution(5)
197
198 # the spheres are 3D glyphs so set that up now
199 _glyph = vtk.vtkGlyph3D()
200 _glyph.ScalingOn()
201 _glyph.SetScaleModeToScaleByScalar() # scale by the radii scalars
202 _glyph.SetColorModeToColorByScalar() # colour by the tag scalars
203 _glyph.SetScaleFactor(1.0) # just in case
204 _glyph.SetInput(_grid)
205 _glyph.SetSource(_sphere.GetOutput()) # grab the 3D glyph to use
206 _glyph.ClampingOff() # so we even see the tiny spheres
207
208 # set up a stripper (this will speed up rendering a lot)
209 _stripper = vtk.vtkStripper()
210 _stripper.SetInput(_glyph.GetOutput())
211
212 # set up the mapper
213 _mapper = vtk.vtkPolyDataMapper()
214 _mapper.SetInput(_stripper.GetOutput())
215 _mapper.ScalarVisibilityOn()
216 _mapper.ColorByArrayComponent("data", 2)
217 _mapper.SetScalarRange(0, 1)
218
219 # set up the actor
220 _actor = vtk.vtkActor()
221 _actor.SetMapper(_mapper)
222
223 # text properties
224 _font_size = 14
225 _textProp = vtk.vtkTextProperty()
226 _textProp.SetFontSize(_font_size)
227 _textProp.SetFontFamilyToArial()
228 _textProp.BoldOff()
229 _textProp.ItalicOff()
230 _textProp.ShadowOff()
231
232 # add a title
233 _titleMapper = vtk.vtkTextMapper()
234 _title = "Example ball plot"
235 _titleMapper.SetInput(_title)
236
237 _titleProp = _titleMapper.GetTextProperty()
238 _titleProp.ShallowCopy(_textProp)
239 _titleProp.SetJustificationToCentered()
240 _titleProp.SetVerticalJustificationToTop()
241
242 # set up the text actor
243 _titleActor = vtk.vtkTextActor()
244 _titleActor.SetMapper(_titleMapper)
245 _titleActor.GetPositionCoordinate().SetCoordinateSystemToNormalizedDisplay()
246 _titleActor.GetPositionCoordinate().SetValue(0.5, 0.95)
247
248 # set up the renderer and the render window
249 _ren = vtk.vtkRenderer()
250 _renWin = vtk.vtkRenderWindow()
251 _renWin.SetSize(640, 480)
252 _renWin.AddRenderer(_ren)
253
254 # add the actor
255 _ren.AddActor(_actor)
256 _ren.AddActor(_titleActor)
257
258 # render the scene
259 _iren = vtk.vtkRenderWindowInteractor()
260 _iren.SetRenderWindow(_renWin)
261 _iren.Initialize()
262 _renWin.Render()
263 _iren.Start()
264
265 # convert the render window to an image
266 _win2imgFilter = vtk.vtkWindowToImageFilter()
267 _win2imgFilter.SetInput(_renWin)
268
269 # save the image to file
270 _outWriter = vtk.vtkPNGWriter()
271 _outWriter.SetInput(_win2imgFilter.GetOutput())
272 _outWriter.SetFileName("ballPlotExample.png")
273 _outWriter.Write()
274
275 elif method == 'povray':
276 ### povray code
277
278 # load the data from the vtk file (yes, I know this is a bit dodgy)
279 import vtk
280
281 # create the reader of the file
282 _reader = vtk.vtkXMLUnstructuredGridReader()
283 _reader.SetFileName("cp_test_0.xml")
284 #_reader.SetFileName("/home/cochrane/raid2/vis4people/steffen/frame_0.xml")
285 _reader.Update()
286
287 # read the output into an unstructured grid
288 _grid = _reader.GetOutput()
289
290 _modelCentre = _grid.GetCenter()
291 _xMin, _xMax, _yMin, _yMax, _zMin, _zMax = _grid.GetBounds()
292
293 # grab the points where the data sit
294 _vtkPoints = _grid.GetPoints()
295
296 # grab the radius data for the radii of the balls
297 _vtkRadii = _grid.GetPointData().GetScalars("radius")
298
299 # grab the tag data and use for colouring the balls
300 _vtkTags = _grid.GetPointData().GetScalars("particleTag")
301
302 # work out dynamically the number of different tags so that can use this
303 # information to automatically set the scalar range for colouring
304 _numPoints = _vtkTags.GetNumberOfTuples()
305 _valueDict = {}
306 for i in range(_numPoints):
307 _tagValue = _vtkTags.GetValue(i)
308 _valueDict[_tagValue] = 1
309
310 _numTags = len(_valueDict.keys())
311
312 _tagValues = _valueDict.keys()
313 _tagValues.sort()
314
315 # count the number of tags, and make an evenly spaced array of points
316 # between zero and one, then use these as the scalars to colour by
317 _vtkScaledTags = vtk.vtkFloatArray()
318 _vtkScaledTags.SetNumberOfTuples(_numPoints)
319 _vtkScaledTags.SetNumberOfComponents(1)
320 _vtkScaledTags.SetName("tags")
321 for i in range(_numPoints):
322 _tagValue = _vtkTags.GetValue(i)
323 for j in range(_numTags):
324 if _tagValues[j] == _tagValue:
325 _vtkScaledTags.InsertTuple1(i, float(j)/float(_numTags-1))
326
327 # use vtk to generate the colour map, will have to do this myself at
328 # some point
329 _lut = vtk.vtkLookupTable()
330 _lut.Build()
331
332 _red = zeros(_numPoints, typecode=Float)
333 _green = zeros(_numPoints, typecode=Float)
334 _blue = zeros(_numPoints, typecode=Float)
335 for i in range(_numPoints):
336 _red[i], _green[i], _blue[i] = _lut.GetColor(_vtkScaledTags.GetValue(i))
337
338 # now convert the information we want (radii, colours, positions) into
339 # array objects so that I can play with them as per normal in python
340
341 # note: this is an inefficient way to do this, I can do it in one loop,
342 # but this way makes the meaning of the code a lot clearer
343
344 ### the points
345 _xData = zeros(_numPoints, typecode=Float)
346 _yData = zeros(_numPoints, typecode=Float)
347 _zData = zeros(_numPoints, typecode=Float)
348 for i in range(_numPoints):
349 _xData[i], _yData[i], _zData[i] = _vtkPoints.GetPoint(i)
350
351 ### the radii
352 _radii = zeros(_numPoints, typecode=Float)
353 for i in range(_numPoints):
354 _radii[i] = _vtkRadii.GetValue(i)
355
356 ### the tags
357 _scaledTags = zeros(_numPoints, typecode=Float)
358 _tags = zeros(_numPoints, typecode=Int)
359 for i in range(_numPoints):
360 _scaledTags[i] = _vtkScaledTags.GetValue(i)
361 _tags[i] = _vtkTags.GetValue(i)
362
363 ### generate the pov file
364
365 # open the pov file to write to
366 pov = open("ballPlotExample.pov", "w")
367
368 # the include files to add
369 pov.write("#include \"shapes.inc\"\n")
370 pov.write("#include \"colors.inc\"\n")
371
372 # the camera
373 pov.write("camera {\n")
374 pov.write(" location <%f, %f, -100>\n" %
375 (_modelCentre[0], _modelCentre[1]))
376 pov.write(" direction <0, 0, 2>\n")
377 pov.write(" up <0, 1, 0>\n")
378 pov.write(" right <4/3, 0, 0>\n")
379 pov.write(" look_at <%f, %f, %f>\n" %
380 (_modelCentre[0], _modelCentre[1], _modelCentre[2]))
381 pov.write("}\n")
382
383 # the light source
384 pov.write("light_source {\n")
385 pov.write(" <%f, %f, -300>\n" % (_modelCentre[0], _modelCentre[1]))
386 pov.write(" colour White\n")
387 pov.write("}\n")
388
389 # the spheres
390 for i in range(_numPoints):
391 pov.write("sphere { \n")
392 pov.write(" <%f, %f, %f>, %f\n" %
393 (_xData[i], _yData[i], _zData[i], _radii[i]))
394 pov.write(" pigment {\n")
395 if _tags[i] != 20:
396 pov.write(" colour red %f green %f blue %f\n" %
397 (_red[i], _green[i], _blue[i]))
398 else:
399 pov.write(" rgbt <%f, %f, %f, 0.90>\n" %
400 (_red[i], _green[i], _blue[i]))
401 pov.write(" }\n")
402 pov.write("}\n")
403
404 # put a title on it
405 pov.write("text {\n")
406 pov.write(" ttf \"timrom.ttf\" \"Example ball plot\" 0.1, 0\n")
407 pov.write(" pigment {\n")
408 pov.write(" colour White\n")
409 pov.write(" }\n")
410 pov.write(" scale <3, 3, 1>\n")
411 pov.write(" translate <%f, %f, 0>\n" %
412 (_modelCentre[0]-10, 1.2*_yMax))
413 pov.write("}\n")
414
415 # close the file
416 pov.close()
417
418 ### generate the ini file
419
420 # open the ini file to write to
421 ini = open("ballPlotExample.ini", "w")
422
423 # the output resolution
424 ini.write("Width=640\n")
425 ini.write("Height=480\n")
426
427 # anti-aliasing settings
428 ini.write("Antialias=on\n")
429
430 # generate png files
431 ini.write("Output_File_Type=N\n")
432
433 # the name of the input pov file
434 ini.write("Input_File_Name=ballPlotExample.pov\n")
435
436 # pause when done
437 ini.write("Pause_When_Done=on\n")
438
439 # close the file
440 ini.close()
441
442 # run povray on the file
443 os.system("povray ballPlotExample.ini")
444
445
446 else:
447 print "Eeek! What plotting method am I supposed to use???"
448
449 # vim: expandtab shiftwidth=4:
450

  ViewVC Help
Powered by ViewVC 1.1.26