[195] | 1 | from visual import * |
---|
| 2 | |
---|
| 3 | def normal(v1,v2): |
---|
| 4 | """Safe computation of normalised cross product |
---|
| 5 | """ |
---|
| 6 | try: |
---|
| 7 | n = norm(cross(v1, v2)) |
---|
| 8 | except ZeroDivisionError: |
---|
| 9 | print v1, v2 |
---|
| 10 | raise Exception |
---|
| 11 | #n = vector(0,0,0) |
---|
| 12 | |
---|
[1295] | 13 | return n |
---|
[195] | 14 | |
---|
| 15 | |
---|
| 16 | class Triangle: |
---|
| 17 | def __init__(self,v0,v1,v2,color=(1,1,1), frame=None): |
---|
| 18 | """Create one flat triangular panel with two sides |
---|
| 19 | """ |
---|
| 20 | |
---|
| 21 | #Outward normal |
---|
[1295] | 22 | |
---|
[195] | 23 | self.normal=normal(v1-v0, v2-v0) |
---|
| 24 | |
---|
| 25 | self.panel = faces(frame = frame) |
---|
[1295] | 26 | self.backpanel = faces(frame = frame) |
---|
[195] | 27 | for v in (v0,v1,v2): |
---|
| 28 | self.panel.append(pos=v, normal = self.normal, color=color) |
---|
| 29 | |
---|
| 30 | for v in (v2,v1,v0): |
---|
[1295] | 31 | self.backpanel.append(pos=v, normal = -self.normal, color=color) |
---|
[195] | 32 | |
---|
[1295] | 33 | |
---|
[195] | 34 | def move(self,v): |
---|
| 35 | """Move panel in direction given by vector v |
---|
| 36 | """ |
---|
| 37 | |
---|
| 38 | for vertex in self.panel.pos: |
---|
[1295] | 39 | vertex += v |
---|
[195] | 40 | |
---|
| 41 | def update_height(self,d,fix_baseline=False): |
---|
| 42 | """Change height of triangle by displacing the third vertex by |
---|
| 43 | d in the direction perpendicular to the baseline (v1-v0) and the |
---|
| 44 | outward normal vector |
---|
| 45 | """ |
---|
| 46 | |
---|
| 47 | v0 = self.panel.pos[0] |
---|
| 48 | v1 = self.panel.pos[1] |
---|
[1295] | 49 | v2 = self.panel.pos[2] |
---|
| 50 | |
---|
[195] | 51 | n = normal(v1-v0, self.normal) |
---|
| 52 | |
---|
| 53 | if fix_baseline: |
---|
| 54 | self.panel.pos[2] -= d*n |
---|
| 55 | else: |
---|
[1295] | 56 | self.panel.pos[:2] += d*n |
---|
[195] | 57 | |
---|
| 58 | def update_color(self, c): |
---|
| 59 | """Change color to c |
---|
| 60 | """ |
---|
| 61 | |
---|
| 62 | self.panel.color = c |
---|
| 63 | |
---|
[1295] | 64 | |
---|
[195] | 65 | def set_vertexheights(self, heights, floor_heights = None): |
---|
| 66 | from Numeric import zeros, Float |
---|
| 67 | |
---|
| 68 | from config import minimum_allowed_height as hmin |
---|
| 69 | if floor_heights is None: |
---|
| 70 | floor_heights = zeros(heights.shape, Float) |
---|
| 71 | |
---|
[1295] | 72 | all_vertices_below_threshold = True |
---|
[195] | 73 | for k in range(3): |
---|
| 74 | w = heights[k] |
---|
| 75 | z = floor_heights[k] |
---|
| 76 | |
---|
| 77 | if w-z >= hmin: |
---|
| 78 | all_vertices_below_threshold = False |
---|
[1295] | 79 | |
---|
[195] | 80 | vertex = self.panel.pos[k] |
---|
| 81 | vertex[2] = w |
---|
| 82 | |
---|
| 83 | vertex = self.backpanel.pos[2-k] |
---|
| 84 | vertex[2] = w |
---|
| 85 | |
---|
| 86 | |
---|
[1295] | 87 | #Update color to visualise dry areas |
---|
| 88 | if all_vertices_below_threshold: |
---|
[195] | 89 | self.panel.color = self.bottom_color |
---|
[1295] | 90 | self.backpanel.color = self.bottom_color |
---|
[195] | 91 | else: |
---|
| 92 | self.panel.color = self.top_color |
---|
| 93 | self.backpanel.color = self.top_color |
---|
| 94 | |
---|
| 95 | |
---|
[1295] | 96 | #update normal |
---|
[195] | 97 | v0 = self.panel.pos[0] |
---|
| 98 | v1 = self.panel.pos[1] |
---|
| 99 | v2 = self.panel.pos[2] |
---|
| 100 | |
---|
| 101 | |
---|
| 102 | n = normal(v1-v0, v2-v0) |
---|
| 103 | self.panel.normal=n |
---|
| 104 | self.backpanel.normal=-n |
---|
| 105 | |
---|
| 106 | |
---|
| 107 | |
---|
[1295] | 108 | |
---|
[195] | 109 | def create_surface(domain): |
---|
| 110 | """Link surface of domains to their visual representations |
---|
| 111 | """ |
---|
| 112 | |
---|
| 113 | fr = frame() #Default frame to contain all objects |
---|
| 114 | s=0 |
---|
| 115 | |
---|
[773] | 116 | Q = domain.quantities['stage'] |
---|
[195] | 117 | try: |
---|
| 118 | Z = domain.quantities['elevation'] |
---|
| 119 | bed = True |
---|
| 120 | except: |
---|
| 121 | bed = False |
---|
[1295] | 122 | |
---|
[195] | 123 | N = Q.vertex_values.shape[0] |
---|
| 124 | |
---|
| 125 | domain.visuals = [] |
---|
| 126 | for i in range(N): |
---|
| 127 | |
---|
| 128 | z0 = Q.vertex_values[i, 0] |
---|
| 129 | z1 = Q.vertex_values[i, 1] |
---|
[1295] | 130 | z2 = Q.vertex_values[i, 2] |
---|
[195] | 131 | |
---|
| 132 | x0, y0 = domain.get_vertex_coordinate(i,0) |
---|
| 133 | x1, y1 = domain.get_vertex_coordinate(i,1) |
---|
[1295] | 134 | x2, y2 = domain.get_vertex_coordinate(i,2) |
---|
| 135 | |
---|
[195] | 136 | V0 = vector(x0, y0, z0) |
---|
| 137 | V1 = vector(x1, y1, z1) |
---|
| 138 | V2 = vector(x2, y2, z2) |
---|
| 139 | |
---|
| 140 | |
---|
| 141 | #Top surface |
---|
| 142 | c0 = 0.1 |
---|
| 143 | c1 = 0.4 |
---|
| 144 | c2 = 0.99 |
---|
| 145 | if s: |
---|
| 146 | c2 = 0.7*c2 #To show triangles in slightly different shades |
---|
| 147 | |
---|
| 148 | s = 1-s |
---|
| 149 | |
---|
| 150 | col = (c0,c1,c2) |
---|
[1295] | 151 | |
---|
[195] | 152 | visual_top = Triangle(V0,V1,V2,color=col,frame=fr) |
---|
| 153 | visual_top.top_color = col |
---|
| 154 | |
---|
[1295] | 155 | |
---|
[195] | 156 | #Bottom surface |
---|
| 157 | #v0, v1, v2 = volume.vertices |
---|
| 158 | #v0 = vector(v0.x, v0.y, v0.field_values[index]) |
---|
| 159 | #v1 = vector(v1.x, v1.y, v1.field_values[index]) |
---|
| 160 | #v2 = vector(v2.x, v2.y, v2.field_values[index]) |
---|
| 161 | if bed: |
---|
| 162 | z0 = Z.vertex_values[i, 0] |
---|
| 163 | z1 = Z.vertex_values[i, 1] |
---|
| 164 | z2 = Z.vertex_values[i, 2] |
---|
| 165 | else: |
---|
| 166 | z0 = z1 = z2 = 0.0 |
---|
| 167 | |
---|
| 168 | V0 = vector(x0, y0, z0) |
---|
| 169 | V1 = vector(x1, y1, z1) |
---|
[1295] | 170 | V2 = vector(x2, y2, z2) |
---|
[195] | 171 | |
---|
| 172 | c0 = 0.3 |
---|
| 173 | c1 = 0.3 |
---|
| 174 | c2 = 0.3 |
---|
| 175 | if s: |
---|
| 176 | c2 = 0.4*c2 #To show triangles in slightly different shades |
---|
| 177 | |
---|
[1295] | 178 | col = (c0,c1,c2) |
---|
[195] | 179 | visual_bottom = Triangle(V0, V1, V2, color=col,frame=fr) |
---|
| 180 | visual_top.bottom_color=col |
---|
[1295] | 181 | |
---|
[195] | 182 | domain.visuals.append( (visual_top, visual_bottom) ) |
---|
[1295] | 183 | |
---|
[195] | 184 | update(domain) |
---|
[1295] | 185 | #print 'Scale', scene.scale |
---|
[195] | 186 | |
---|
| 187 | |
---|
[1295] | 188 | def update(domain): |
---|
[195] | 189 | """Update vertex heights. |
---|
[1295] | 190 | The argument index refers to which conserved quantity to visualise. |
---|
[195] | 191 | If domain.smooth is set True, vertex values will be averaged |
---|
| 192 | yielding a smoother surface. |
---|
| 193 | """ |
---|
| 194 | |
---|
| 195 | from Numeric import array |
---|
| 196 | |
---|
| 197 | |
---|
[773] | 198 | Q = domain.quantities['stage'] |
---|
[195] | 199 | N = Q.vertex_values.shape[0] |
---|
[1295] | 200 | |
---|
[195] | 201 | #print scene.forward |
---|
| 202 | #FIXME: Use smoother from pyvolution instead |
---|
| 203 | if domain.smooth: |
---|
| 204 | #Get all average point values |
---|
| 205 | vertex_heights = {} |
---|
| 206 | |
---|
| 207 | for k in range(N): |
---|
| 208 | for i in range(3): |
---|
[324] | 209 | vertex = domain.triangles[k, i] |
---|
[195] | 210 | if vertex_heights.has_key(vertex): |
---|
| 211 | vertex_heights[vertex].append( |
---|
| 212 | Q.vertex_values[k, i]) |
---|
| 213 | else: |
---|
| 214 | vertex_heights[vertex] = [] |
---|
| 215 | vertex_heights[vertex].append( |
---|
| 216 | Q.vertex_values[k, i]) |
---|
| 217 | |
---|
| 218 | |
---|
| 219 | |
---|
[1295] | 220 | for k in range(N): |
---|
[195] | 221 | if domain.smooth: |
---|
| 222 | #The averages |
---|
| 223 | x = zeros(3, Float) |
---|
[1295] | 224 | for i in range(3): |
---|
[324] | 225 | vertex = domain.triangles[k, i] |
---|
[195] | 226 | A = array(vertex_heights[vertex]) |
---|
| 227 | x[i] = sum(A)/len(A) |
---|
[1295] | 228 | else: |
---|
[195] | 229 | #The true values |
---|
| 230 | x = [Q.vertex_values[k, 0], |
---|
| 231 | Q.vertex_values[k, 1], |
---|
| 232 | Q.vertex_values[k, 2]] |
---|
| 233 | |
---|
[1295] | 234 | |
---|
[195] | 235 | #Do it |
---|
| 236 | floor_heights = array([pos[2] for pos in domain.visuals[k][1].panel.pos]) |
---|
[1295] | 237 | |
---|
[195] | 238 | domain.visuals[k][0].set_vertexheights(x, floor_heights) |
---|
| 239 | |
---|
| 240 | |
---|
| 241 | |
---|
| 242 | |
---|
| 243 | scene.width = 1000 |
---|
| 244 | scene.height = 800 |
---|
| 245 | |
---|
| 246 | #Original |
---|
| 247 | scene.center = (0.5,0.5,0) |
---|
[1295] | 248 | scene.forward = vector(0.0, 0.5, -0.5) |
---|
[195] | 249 | |
---|
| 250 | #Temporary (for bedslope) |
---|
[1295] | 251 | #scene.forward = vector(0.0006, 0.7, -0.03) |
---|
[195] | 252 | |
---|
| 253 | |
---|
| 254 | #Temporary for hackett - begin |
---|
| 255 | #scene.autoscale = 0 |
---|
| 256 | #scene.scale = (0.002, 0.002, 0.01) #Scale z so that countours stand out more |
---|
| 257 | #scene.center = (300.0,500.0,-10) |
---|
| 258 | #Temporary for hackett - end |
---|
| 259 | |
---|
[1295] | 260 | |
---|
[195] | 261 | scene.ambient = 0.4 |
---|
| 262 | scene.lights = [(0.6, 0.3, 0.2), (0.1, -0.5, 0.4), (-0.1, 0.1, -0.4), |
---|
| 263 | (-0.2, 0.2, 0.1)] |
---|
| 264 | |
---|
| 265 | |
---|
| 266 | |
---|
| 267 | |
---|