source: anuga_core/source/anuga/visualiser/realtime.py @ 7452

Last change on this file since 7452 was 7452, checked in by steve, 15 years ago

Committing latest parallel code

File size: 5.9 KB
Line 
1#from  import Float, zeros, shape
2import numpy as num
3from Tkinter import Button, E, Tk, W
4from threading import Event
5from visualiser import Visualiser
6from vtk import vtkCellArray, vtkPoints, vtkPolyData
7
8class RealtimeVisualiser(Visualiser):
9    """A VTK-powered realtime visualiser which runs in its own thread.
10    In addition to the functions provided by the standard visualiser,
11    the following additional functions are provided:
12
13    update() - Sync the visualiser to the current state of the model.
14    Should be called inside the evolve loop.
15
16    evolveFinished() - Clean up synchronisation constructs that tie the
17    visualiser to the evolve loop. Call this after the evolve loop finishes
18    to ensure a clean shutdown.
19    """
20    def __init__(self, source):
21        """The source parameter is assumed to be a Domain.
22        """
23        Visualiser.__init__(self, source)
24
25        self.running = True
26
27        self.xmin = None
28        self.xmax = None
29        self.ymin = None
30        self.ymax = None
31        self.zmin = None
32        self.zmax = None
33
34        # Synchronisation Constructs
35        self.sync_idle = Event()
36        self.sync_idle.clear()
37        self.sync_unpaused = Event()
38        self.sync_unpaused.set()
39        self.sync_redrawReady = Event()
40        self.sync_redrawReady.clear()
41       
42    def setup_grid(self):
43        self.vtk_cells = vtkCellArray()
44        triangles = self.source.get_triangles()
45        N_tri = len(self.source)
46        verticies = self.source.get_vertex_coordinates()
47        N_vert = len(verticies)
48
49        # Also build vert_index - a list of the x & y values of each vertex
50        self.vert_index = num.zeros((N_vert,2), num.float)
51        for n in range(N_tri):
52            self.vtk_cells.InsertNextCell(3)
53            for v in range(3):
54                self.vert_index[n * 3 + v] = verticies[n * 3 + v]
55                self.vtk_cells.InsertCellPoint(n * 3 + v)
56
57    def update_height_quantity(self, quantityName, dynamic=True):
58        N_vert = len(self.source.get_vertex_coordinates())
59        qty_index = num.zeros(N_vert, num.float)
60        triangles = self.source.get_triangles()
61        vertex_values, _ = self.source.get_quantity(quantityName).get_vertex_values(xy=False, smooth=False)
62
63        for n in range(N_vert):
64            qty_index[n] = vertex_values[n]
65
66        points = vtkPoints()
67        for v in range(N_vert):
68            points.InsertNextPoint(self.vert_index[v][0],
69                                   self.vert_index[v][1],
70                                   qty_index[v] * self.height_zScales[quantityName]
71                                   + self.height_offset[quantityName])
72            if self.xmin == None or self.xmin > self.vert_index[v][0]:
73                self.xmin = self.vert_index[v][0]
74            if self.xmax == None or self.xmax < self.vert_index[v][0]:
75                self.xmax = self.vert_index[v][0]
76            if self.ymin == None or self.ymin > self.vert_index[v][1]:
77                self.ymin = self.vert_index[v][1]
78            if self.ymax == None or self.ymax < self.vert_index[v][1]:
79                self.ymax = self.vert_index[v][1]
80            if self.zmin == None or self.zmin > qty_index[v] * self.height_zScales[quantityName] + self.height_offset[quantityName]:
81                self.zmin = qty_index[v] * self.height_zScales[quantityName] + self.height_offset[quantityName]
82            if self.zmax == None or self.zmax < qty_index[v] * self.height_zScales[quantityName] + self.height_offset[quantityName]:
83                self.zmax = qty_index[v] * self.height_zScales[quantityName] + self.height_offset[quantityName]
84
85        polydata = self.vtk_polyData[quantityName] = vtkPolyData()
86        polydata.SetPoints(points)
87        polydata.SetPolys(self.vtk_cells)
88
89    def get_3d_bounds(self):
90        return [self.xmin, self.xmax, self.ymin, self.ymax, self.zmin, self.zmax]
91       
92    def build_quantity_dict(self):
93        triangles = self.source.get_triangles()
94        quantities = {}
95        for q in self.source.get_quantity_names():
96            quantities[q], _ = self.source.get_quantity(q).get_vertex_values(xy=False)
97        return quantities
98
99    def setup_gui(self):
100        Visualiser.setup_gui(self)
101        self.tk_pauseResume = Button(self.tk_controlFrame, text="Pause", command=self.pauseResume)
102        self.tk_pauseResume.grid(row=1, column=0, sticky=E+W)
103
104    def pauseResume(self):
105        if self.sync_unpaused.isSet():
106            self.sync_unpaused.clear()
107            self.tk_pauseResume.config(text="Resume")
108        else:
109            self.sync_unpaused.set()
110            self.tk_pauseResume.config(text="Pause")
111
112    def shutdown(self):
113        Visualiser.shutdown(self)
114        self.running = False
115        self.sync_idle.set()
116        self.sync_unpaused.set()
117
118    def redraw(self):
119        if self.running and self.sync_unpaused.isSet():
120            self.sync_redrawReady.wait()
121            self.sync_redrawReady.clear()
122            self.redraw_quantities()
123            self.sync_idle.set()
124        Visualiser.redraw(self)
125
126    def update(self,pause=False):
127        """Sync the visualiser to the domain. Call this in the evolve loop."""
128           
129        if self.running:
130            self.sync_redrawReady.set()
131            self.sync_idle.wait()
132            self.sync_idle.clear()
133            self.sync_unpaused.wait()
134
135        if pause and self.running:
136            if self.sync_unpaused.isSet():
137                self.sync_unpaused.clear()
138                self.tk_pauseResume.config(text="Resume")
139               
140                self.sync_redrawReady.set()
141                self.sync_idle.wait()
142                self.sync_idle.clear()
143                self.sync_unpaused.wait()
144           
145        return self.running
146
147    def evolveFinished(self):
148        """Stop the visualiser from waiting on signals from the evolve loop.
149        Call this just after the evolve loop to ensure a clean shutdown."""
150        self.running = False
151        self.sync_redrawReady.set()
Note: See TracBrowser for help on using the repository browser.