source: anuga_work/development/steve/get_triangle_data.py @ 3929

Last change on this file since 3929 was 3514, checked in by duncan, 19 years ago

Hi all,
I'm doing a change in the anuga structure, moving the code to

\anuga_core\source\anuga

After you have done an svn update, the PYTHONPATH has to be changed to;
PYTHONPATH = anuga_core/source/

This is part of changes required to make installation of anuga quicker and reducing the size of our sandpits.

If any imports are broken, try fixing them. With adding anuga. to them for example. If this seems to have really broken things, email/phone me.

Cheers
Duncan

File size: 3.1 KB
Line 
1import sys
2from os import sep
3sys.path.append('..'+sep+'pyvolution')
4
5
6import least_squares
7
8alpha = 1
9mesh_file = 'merimbula_10785.tsh'
10point_file = 'meri0.xya'
11mesh_output_file = 'merimbula_10785_%g.tsh'%alpha
12
13print mesh_output_file
14
15from shallow_water import Domain
16from pmesh2domain import pmesh_to_domain_instance
17from anuga.pyvolution.util import file_function, Polygon_function, read_polygon
18from Numeric import zeros, Float, maximum, minimum
19
20#-------
21# Domain
22filename = mesh_output_file
23print 'Creating domain from', filename
24domain = pmesh_to_domain_instance(filename, Domain)
25print "Number of triangles = ", len(domain)
26
27
28
29
30
31
32triangles = domain.triangles
33
34vertices = domain.coordinates
35
36tri_vertices = domain.vertex_coordinates
37
38Bed = domain.quantities['elevation']
39
40bed = Bed.vertex_values
41
42from visual import *
43
44f = frame()
45model = faces(frame=f)
46
47
48max_x = max(vertices[:,0])
49min_x = min(vertices[:,0])
50max_y = max(vertices[:,1])
51min_y = min(vertices[:,1])
52range_x = max_x - min_x
53range_y = max_y - min_y
54
55range_xy = max(range_x, range_y)
56
57print max_x, min_x, max_y, min_y , range_x, range_y, range_xy
58
59print 'Calculating triangles'
60
61pos     = zeros( (6*len(domain),3), Float)
62colour  = zeros( (6*len(domain),3), Float)
63normals = zeros( (6*len(domain),3), Float)
64
65for i in range( len(domain) ):
66#    v0 = vector(2.0*(vertices[triangles[i,0],0]-min_x)/range_xy-1.0, \
67#            2.0*(vertices[triangles[i,0],1]-min_y)/range_xy-1.0,bed[i,0]/50.0)
68#    v1 = vector(2.0*(vertices[triangles[i,1],0]-min_x)/range_xy-1.0, \
69#            2.0*(vertices[triangles[i,1],1]-min_y)/range_xy-1.0,bed[i,1]/50.0)
70#    v2 = vector(2.0*(vertices[triangles[i,2],0]-min_x)/range_xy-1.0, \
71#            2.0*(vertices[triangles[i,2],1]-min_y)/range_xy-1.0,bed[i,2]/50.0)
72
73
74    v0 = vector(2.0*(tri_vertices[i,0]-min_x)/range_xy-1.0, \
75            2.0*(tri_vertices[i,1]-min_y)/range_xy-1.0,bed[i,0]/50.0)
76    v1 = vector(2.0*(tri_vertices[i,2]-min_x)/range_xy-1.0, \
77            2.0*(tri_vertices[i,3]-min_y)/range_xy-1.0,bed[i,1]/50.0)
78    v2 = vector(2.0*(tri_vertices[i,4]-min_x)/range_xy-1.0, \
79            2.0*(tri_vertices[i,5]-min_y)/range_xy-1.0,bed[i,2]/50.0)
80
81
82    try:
83        normal = norm( cross(v1-v0, v2-v0) )
84    except:
85        normal = vector(0,0,1)
86
87    pos[6*,:] = v0
88    pos[6*i+1,:] = v1
89    pos[6*i+2,:] = v2
90    pos[6*i+3,:] = v0
91    pos[6*i+4,:] = v2
92    pos[6*i+5,:] = v1
93
94    colour[6*,:] = color.blue
95    colour[6*i+1,:] = color.blue
96    colour[6*i+2,:] = color.blue
97    colour[6*i+3,:] = color.blue
98    colour[6*i+4,:] = color.blue
99    colour[6*i+5,:] = color.blue
100
101    normals[6*,:] = normal
102    normals[6*i+1,:] = normal
103    normals[6*i+2,:] = normal
104    normals[6*i+3,:] = -normal
105    normals[6*i+4,:] = -normal
106    normals[6*i+5,:] = -normal
107
108    #FacetedTriangle( v0, v1, v2 )
109    #print v0
110    #print v1
111    #print v2
112
113print 'Drawing plot'
114model.pos    = pos
115model.color  = colour
116model.normal = normals
117
118#v0 = vector(0.0,0.0,0.0)
119#v1 = vector(1.0,0.0,0.1)
120#v2 = vector(1.0,1.0,0.1)
121
122
Note: See TracBrowser for help on using the repository browser.