source: inundation/pmesh/test_mesh.py @ 3475

Last change on this file since 3475 was 3440, checked in by duncan, 18 years ago

comments

File size: 82.0 KB
Line 
1#!/usr/bin/env python
2#
3import tempfile
4import unittest
5
6#try:
7from pmesh.mesh import *
8#except ImportError: 
9#    from mesh import *
10   
11from load_mesh.loadASCII import *
12from coordinate_transforms.geo_reference import Geo_reference
13from geospatial_data.geospatial_data import Geospatial_data
14from utilities.polygon import  is_inside_polygon ### inside_polygon
15
16class meshTestCase(unittest.TestCase):
17    def setUp(self):
18        pass
19   
20    def tearDown(self):
21        pass
22
23    def testPointDistance(self):
24        a = Point(0.0, 0.0)
25        b = Point(0.0, 10.0)
26       
27        self.failUnless( a.DistanceToPoint(b) == 10.0,
28                        'Point DistanceToPoint is wrong!')
29   
30    def testVertexDistance(self):
31        a = Vertex (0.0, 0.0)
32        b = Vertex (0.0, 10.0)
33       
34        self.failUnless( a.DistanceToPoint(b) == 10.0,
35                        'Point DistanceToPoint is wrong!')
36       
37    def testTriangle(self):
38        a = Vertex (0.0, 0.0)
39        b = Vertex (0.0, 2.0)
40        c = Vertex (2.0,0.0)
41        d = Vertex (0.0, 4.0)
42        e = Vertex (2.0, 2.0)
43        f = Vertex (4.0,0.0)
44       
45        t1 = Triangle(b,a,c)       
46        t2 = Triangle(b,c,e)     
47        t3 = Triangle(e,c,f)     
48        t4 = Triangle(d,b,e)
49        t2.setNeighbors(t3,t4,t1)
50       
51        self.failUnless( t2.neighbors[2].vertices[0] == b, 'Triangle initialisation is wrong!')
52   
53       
54    def testSegment(self):
55        a = Vertex (0.0, 0.0)
56        b = Vertex (0.0, 10.0)
57        s = Segment(a,b, tag = 20)     
58       
59        self.failUnless( s.vertices[0].DistanceToPoint(s.vertices[1]) == 10.0,
60                        'vertices in a segment are wrong')
61       
62        self.failUnless( s.tag == 20.0,
63                        'tag in a segment are wrong')
64
65    def testdeleteUserVertex(self):
66
67       
68        mesh = Mesh()
69        a = mesh.addUserVertex(0.0, 0.0)
70        b = mesh.addUserVertex (0.0, 2.0)
71        c = mesh.addUserVertex (2.0,0.0)
72       
73        s1 = mesh.addUserSegment(a,b)
74        s2 = mesh.addUserSegment(a,c)
75        s3 = mesh.addUserSegment(c,b)
76
77        mesh.deleteMeshObject (a) 
78        self.failUnless(mesh.userSegments[0] == s3,
79                        'Bad segment. ')       
80        self.failUnless(len(mesh.userSegments) ==1,
81                        'Segments not deleted.')
82        self.failUnless(len(mesh.userVertices) == 2,
83                        'Vertex not deleted.')
84       
85 
86    # FIXME add test for minAngle   
87    def testgenerateMesh(self):
88        a = Vertex (0.0, 0.0)
89        d = Vertex (0.0, 4.0)
90        f = Vertex (4.0,0.0)
91
92        s1 = Segment(a,d)
93        s2 = Segment(d,f)
94        s3 = Segment(a,f)
95
96        r1 = Region(0.3, 0.3,tag = 1.3,maxArea = .6)
97        #print r1
98        m = Mesh(userVertices=[a,d,f], userSegments=[s1,s2,s3], regions=[r1] )
99       
100        m.generateMesh("Q", maxArea = 2.1 )         
101
102        #print m
103
104        #m.plotMeshTriangle()
105
106        result = 1.414214
107        delta  = 0.00001
108       
109        self.failUnless((m.meshTriangles[1].vertices[0].x < result + delta) or
110                        (m.meshTriangles[1].vertices[0].x > result - delta),
111                        'generated mesh is wrong!')
112
113    def test_regionalMaxArea(self):
114        v0 = Vertex (0.0, 0.0)
115        v1 = Vertex (6.0, 0.0)
116        v2 = Vertex (6.0,6.0)
117        v3 = Vertex (0.0,6.0)
118
119        s1 = Segment(v0,v1)
120        s2 = Segment(v1,v2)
121        s3 = Segment(v3,v2)
122        s4 = Segment(v3,v0)
123        s5 = Segment(v2,v0)
124
125        r1 = Region(3, 1,tag = 1.3)
126        #print r1
127        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5], regions=[r1] )
128       
129        m.generateMesh("Q", maxArea = 36 )         
130
131        #m.plotMeshTriangle()
132        #print "len(m.meshTriangles)",len(m.meshTriangles)
133
134        self.failUnless(len(m.meshTriangles) == 2, 
135                        'test_regionalMaxArea 1:generated mesh is wrong!')
136       
137        ## Another test case
138        r1 = Region(3, 1,tag = 1.3)
139        r2 = Region(1, 3,tag = 1.3, maxArea = 8)
140        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
141                 regions=[r1,r2] )
142        m.generateMesh("Q", maxArea = 36 )
143       
144        self.failUnless(len(m.meshTriangles) >= 6,
145                        'testregion_with_maxarea 2: # of tris is wrong!')   
146       
147               
148        ## Another test case
149        r1 = Region(3, 1, tag = 1.3, maxArea = 8)
150        r2 = Region(1, 3, tag = 1.3, maxArea = 8)
151        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
152                 regions=[r1,r2] )
153        m.generateMesh("Q", maxArea = 36 ) 
154        #print "len(m.meshTriangles)",len(m.meshTriangles)
155       
156        self.failUnless(len(m.meshTriangles) >= 8,
157                        'testregion_with_maxarea 3: # of tris is wrong!')
158               
159               
160        ## Another test case
161        r1 = Region(3, 1, tag = 1.3 )
162        r2 = Region(1, 3, tag = 1.3, maxArea = 8)
163        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
164                 regions=[r1,r2] )
165        m.generateMesh("Q", maxArea = 8 ) 
166        self.failUnless(len(m.meshTriangles) >= 8,
167                        'testregion_with_maxarea 4: # of tris is wrong!')   
168
169       
170        ## Another test case
171        r1 = Region(3, 1,tag = 1.3, maxArea = 8)
172        r2 = Region(1, 3,tag = 1.3, maxArea = 8)
173        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
174                 regions=[r1,r2] )
175        m.generateMesh("Q", maxArea = 36,isRegionalMaxAreas = False )     
176        self.failUnless(len(m.meshTriangles) == 2, 
177                        'test_regionalMaxArea 5:generated mesh is wrong!')
178       
179        ## Another test case
180        r1 = Region(3, 1,tag = 1.3, maxArea = 8)
181        r2 = Region(1, 3,tag = 1.3, maxArea = 8)
182        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
183                 regions=[r1,r2] )
184        m.generateMesh("Q",isRegionalMaxAreas = False )
185        self.failUnless(len(m.meshTriangles) == 2, 
186                        'test_regionalMaxArea 5:generated mesh is wrong!')
187       
188    def test_generate_mesh(self):
189        v0 = Vertex (0.0, 0.0)
190        v1 = Vertex (6.0, 0.0)
191        v2 = Vertex (6.0,6.0)
192        v3 = Vertex (0.0,6.0)
193
194        s1 = Segment(v0,v1)
195        s2 = Segment(v1,v2)
196        s3 = Segment(v3,v2)
197        s4 = Segment(v3,v0)
198        s5 = Segment(v2,v0)
199
200        r1 = Region(3, 1,tag = 1.3)
201        #print r1
202        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
203                 regions=[r1] )
204       
205        m.generate_mesh(maximum_triangle_area=36,verbose=False)         
206
207        self.failUnless(len(m.meshTriangles) == 2, 
208                        'test_regionalMaxArea 1:generated mesh is wrong!')
209       
210        ## Another test case
211        r1 = Region(3, 1,tag = 1.3)
212        r2 = Region(1, 3,tag = 1.3, maxArea = 8)
213        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
214                 regions=[r1,r2] )
215        m.generate_mesh(maximum_triangle_area=36,verbose=False) 
216       
217        self.failUnless(len(m.meshTriangles) >= 6,
218                        'testregion_with_maxarea 2: # of tris is wrong!')   
219               
220        ## Another test case
221        r1 = Region(3, 1, tag = 1.3, maxArea = 8)
222        r2 = Region(1, 3, tag = 1.3, maxArea = 8)
223        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
224                 regions=[r1,r2] )
225        m.generate_mesh(maximum_triangle_area=36,verbose=False)         
226        #print "len(m.meshTriangles)",len(m.meshTriangles)
227       
228        self.failUnless(len(m.meshTriangles) >= 8,
229                        'testregion_with_maxarea 3: # of tris is wrong!')
230                         
231        ## Another test case
232        r1 = Region(3, 1, tag = 1.3 )
233        r2 = Region(1, 3, tag = 1.3, maxArea = 8)
234        m = Mesh(userVertices=[v0,v1,v2,v3], userSegments=[s1,s2,s3,s4,s5],
235                 regions=[r1,r2] )
236        m.generate_mesh(maximum_triangle_area=8,verbose=False)   
237        self.failUnless(len(m.meshTriangles) >= 8,
238                        'testregion_with_maxarea 4: # of tris is wrong!')   
239
240        ## Another test case r1 = Region(3, 1,tag = 1.3, maxArea = 8)
241        r2 = Region(1, 3,tag = 1.3, maxArea = 8)
242        m = Mesh(userVertices=[v0,v1,v2,v3],
243        userSegments=[s1,s2,s3,s4,s5], regions=[r1,r2] )
244        m.generate_mesh(verbose=False)
245        #print "en(m.meshTriangles)", len(m.meshTriangles)
246        self.failUnless(len(m.meshTriangles) >= 8,
247        'You have issues!')
248       
249    def testdeleteUserVertex(self):
250        mesh = Mesh()
251        a = mesh.addUserVertex(0.0, 0.0)
252        b = mesh.addUserVertex (0.0, 2.0)
253        c = mesh.addUserVertex (2.0,0.0)
254       
255        s1 = mesh.addUserSegment(a,b)
256        s2 = mesh.addUserSegment(a,c)
257        s3 = mesh.addUserSegment(c,b)
258
259        mesh.deleteMeshObject (s2)
260       
261        #print ",s2 in mesh.userSegments" ,s2 in mesh.userSegments
262        self.failUnless(not(s2 in mesh.userSegments),
263                        'Bad segment. ')       
264        self.failUnless(len(mesh.userSegments) ==2,
265                        'Segments not deleted.')
266        self.failUnless(len(mesh.userVertices) == 3,
267                        'Vertex deleted, instead of segment.')
268
269    def testTriangleArea(self):
270        a = Vertex (10.0, 10.0)
271        b = Vertex (10.0, 20.0)
272        c = Vertex (20.0,10.0)
273       
274        d = Vertex (-20.0, 0.0)
275        e = Vertex (-20.0, -20.0)
276        f = Vertex (20.0,-20.0)
277       
278        t1 = Triangle(b,a,c)     
279        t2 = Triangle(e,d,f)
280       
281#         print "t1", t1
282#         print "t1 area ", t1.calcArea()
283#         print "t2", t2
284#         print "t2 area ", t2.calcArea()
285        self.failUnless( t1.calcArea() == 50 and t2.calcArea() == 400, 'Triangle area is wrong!')
286    def testisUserSegmentNew (self):
287        mesh = Mesh()
288        a = mesh.addUserVertex(0.0, 0.0)
289        b = mesh.addUserVertex (0.0, 2.0)
290        c = mesh.addUserVertex (2.0,0.0)
291        d = mesh.addUserVertex (2.0,3.0)
292       
293        s1 = mesh.addUserSegment(a,b)
294        s2 = mesh.addUserSegment(a,c)
295        s3 = mesh.addUserSegment(c,b)
296
297        self.failUnless(mesh.isUserSegmentNew(a,d) ,
298                        'Segment should be new. ')
299        self.failUnless(not(mesh.isUserSegmentNew(a,b)) ,
300                        'Segment should not be new. ')
301
302
303    def testisUserSegmentNewII (self):
304        mesh = Mesh()
305        a = mesh.addUserVertex(0.0, 0.0)
306        b = mesh.addUserVertex (0.0, 2.0)
307        c = mesh.addUserVertex (2.0,0.0)
308        d = mesh.addUserVertex (2.0,3.0)
309       
310        s1 = mesh.addUserSegment(a,b)
311        s2 = mesh.addUserSegment(a,c)
312        s3 = mesh.addUserSegment(c,b)
313
314        self.failUnless(mesh.representedUserSegment(a,d) == None,
315                        'Segment should be new. ')
316        self.failUnless(mesh.representedUserSegment(a,b) == s1 ,
317                        'Segment should not be new. ')
318       
319    def testauto_segment(self):
320        p0 = Vertex (0.0, 0.0)
321        p1 = Vertex (0.0, 4.0)
322        p2 = Vertex (4.0,4.0)
323        p3 = Vertex (4.0,0.0)
324
325        s1 = Segment(p0,p1)
326       
327        m = Mesh(userVertices=[p0, p1, p2, p3], userSegments=[s1] ) 
328        m.auto_segment()
329       
330        #print 'Len', len(m.userSegments)
331        self.failUnless(len(m.getUserSegments()) == 4 ,
332                        'userSegments is wrong!')
333     
334        m.auto_segment()
335        self.failUnless(len(m.getUserSegments()) == 4 ,
336                        'userSegments is wrong!')
337     
338    def testauto_segmentII(self):
339        p1 = Vertex (3.0, 4.0)
340        p2 = Vertex (3.0,2.0)
341        p3 = Vertex (3.0,0.0)
342        p4 = Vertex (6.0, 4.0)
343        p5 = Vertex (6.0,2.0)
344        p0 = Vertex (6.0,0.0)
345
346
347        s1 = Segment(p2,p3)
348        s2 = Segment(p4,p5)
349       
350        m = Mesh(userVertices=[p0, p1, p2, p3, p4, p5],
351                 userSegments=[s1, s2])     
352
353        m.auto_segment()
354       
355        s3 = m.representedAlphaUserSegment(p3,p0)
356        self.failUnless(not (s3 == None) ,
357                        'userSegments is wrong!')
358
359       
360        s6 = m.representedAlphaUserSegment(p1,p4)       
361        self.failUnless(not (s6 == None) ,
362                        'userSegments is wrong!')
363       
364        # remove a segment, add a point, auto_segment
365        m.alphaUserSegments.remove(s3)
366        p6 = Vertex (1.0, 2.0)
367        m.userVertices.append(p6)
368       
369        m.auto_segment()
370       
371        s1_now = m.representedUserSegment(p3,p2)
372        self.failUnless(s1_now == s1 ,
373                        'userSegments is wrong!')
374       
375        s2_now = m.representedUserSegment(p5,p4)       
376        self.failUnless(s2_now == s2 ,
377                        'userSegments is wrong!')
378       
379        s3 = m.representedAlphaUserSegment(p3,p6)       
380        self.failUnless(not (s3 == None) ,
381                        'userSegments is wrong!')
382       
383        s4 = m.representedAlphaUserSegment(p3,p6)       
384        self.failUnless(not (s4 == None) ,
385                        'userSegments is wrong!')
386       
387        s5 = m.representedAlphaUserSegment(p4,p6)       
388        self.failUnless(s5 == None ,
389                        'userSegments is wrong!')
390        #print m
391       
392    def testRegions(self):
393        a = Vertex (0.0, 0.0)
394        b = Vertex (0.0, 2.0)
395        c = Vertex (2.0,0.0)
396        d = Vertex (0.0, 4.0)
397        e = Vertex (2.0, 2.0)
398        f = Vertex (4.0,0.0)
399        g = Vertex (0.0,-2.0)
400       
401        Segment.set_default_tag("")
402        s1 = Segment(a,b)
403        s2 = Segment(b,e)
404        s3 = Segment(e,c)
405        s4 = Segment(c,a)
406        s5 = Segment(b,d)
407        s6 = Segment(e,d)
408        s7 = Segment(e,f)
409        s8 = Segment(c,f)
410        s9 = Segment(g,c)
411        s10 = Segment(g,a)
412
413        r1 = Region(0.1,0.1,tag="22")
414        r2 = Region(0.1,2.1,tag="11")
415        r3 = Region(2.1,0.1)
416       
417        m = Mesh(userVertices=[a,b,c,d,e,f,g], userSegments=[s1,s2,s3,s4,s5,s6,s7,s8,s9,s10], regions=[r1,r2,r3] )
418        m.generateMesh("Q", maxArea = 2.1 )
419        #print m
420        Triangulation =  m.getTriangulation()
421        #print Triangulation[0].attribute
422        #print Triangulation[1].attribute
423        #print Triangulation[2].attribute
424        #print Triangulation[3].attribute
425        #print Triangulation[4].attribute
426       
427        self.failUnless(Triangulation[0].attribute == "" and
428                        Triangulation[1].attribute == "22" and
429                        Triangulation[2].attribute == "" and
430                        Triangulation[3].attribute == "11" and
431                        Triangulation[4].attribute == "22" ,
432                        'region attributes are wrong!')   
433
434    def test_vertexAttribs(self):
435        a = Vertex (0.0, 0.0, attributes = [12.0,2.0])
436        d = Vertex (0.0, 4.0, attributes = [9.0,7.0])
437        f = Vertex (4.0,0.0, attributes = [14.0,3.0])
438   
439        Segment.set_default_tag("")
440        s1 = Segment(a,d)
441        s2 = Segment(d,f)
442        s3 = Segment(a,f)
443     
444        r1 = Region(0.3, 0.3, tag = 88.9)
445     
446        m = Mesh(userVertices=[a,d,f], userSegments=[s1,s2,s3], regions=[r1])
447
448        m.generateMesh("Q", maxArea = 2.1)
449
450        vert = m.getMeshVertices()
451       
452        self.failUnless(vert[0].attributes == [12.0, 2.0] and
453                        vert[1].attributes == [9.0, 7.0] and
454                        vert[2].attributes == [14.0,3.0] and
455                        vert[3].attributes == [12.232233047033631, 4.4142135623730949] and
456                        vert[4].attributes == [13.0, 2.5] ,
457                        'vertex attributes are wrong!')
458
459       
460    def test_vertexAttribs2(self):
461   
462        a = Vertex (0.0, 0.0)
463        d = Vertex (0.0, 4.0)
464        f = Vertex (4.0,0.0)
465   
466        Segment.set_default_tag("")
467        s1 = Segment(a,d)
468        s2 = Segment(d,f)
469        s3 = Segment(a,f)
470     
471        r1 = Region(0.3, 0.3, tag = 88.9)
472     
473        m = Mesh(userVertices=[a,d,f], userSegments=[s1,s2,s3], regions=[r1])
474
475        m.generateMesh("Q", maxArea = 2.1 )
476
477        vert = m.getMeshVertices()
478        self.failUnless(vert[0].attributes == [] and
479                        vert[1].attributes == [] and
480                        vert[2].attributes == [] and
481                        vert[3].attributes == [] and
482                        vert[4].attributes == [],
483                        'vertex attributes are wrong!')
484
485    def test_segtag(self):
486   
487        a = Vertex (0.0, 0.0)
488        d = Vertex (0.0, 4.0)
489        f = Vertex (4.0,0.0)
490   
491        s1 = Segment(a,d,tag = 5)
492        s2 = Segment(d,f,tag = 7)
493        s3 = Segment(a,f,tag = 9)
494     
495        m = Mesh(userVertices=[a,d,f], userSegments=[s1,s2,s3])
496
497        m.generateMesh("Q", maxArea = 2.1 )
498
499        #m.export_mesh_file("from_test_mesh.tsh")
500        seg = m.getMeshSegments()
501        #print "seg",seg
502        #print "seg[0].tag"
503        #print seg[0].tag
504        #print "seg[0].tag"
505       
506        self.failUnless(seg[0].tag == 5 and
507                        seg[1].tag == 7 and
508                        seg[2].tag == 9 and
509                        seg[3].tag == 7 and
510                        seg[4].tag == 9,
511                        'seg tags are wrong')
512           
513
514    def test_segtag2(self):
515   
516        a = Vertex (0.0, 0.0)
517        d = Vertex (0.0, 4.0)
518        f = Vertex (4.0,0.0)
519        e = Vertex (1.0,1.0)
520   
521        s1 = Segment(a,d)
522        s2 = Segment(d,f)
523        s3 = Segment(a,f)
524        s4 = Segment(a,e)
525     
526        m = Mesh(userVertices=[a,d,f,e], userSegments=[s1,s2,s3,s4])
527
528        m.generateMesh("Q", maxArea = 2.1)
529
530        seg = m.getMeshSegments()
531        self.failUnless(seg[0].tag == "exterior" and
532                        seg[1].tag == "exterior" and
533                        seg[2].tag == "exterior" and
534                        seg[3].tag == "" and
535                        seg[4].tag == "exterior",
536                        '2nd seg tags are wrong')
537
538    def test_asciiFile(self):
539   
540        a = Vertex (0.0, 0.0)  #, attributes = [1.1])
541        d = Vertex (0.0, 4.0)  #, attributes = [1.2])
542        f = Vertex (4.0,0.0)  #, attributes = [1.3])
543        e = Vertex (1.0,1.0)  #, attributes = [1.4])
544   
545        s1 = Segment(a,d)
546        s2 = Segment(d,f)
547        s3 = Segment(a,f)
548        s4 = Segment(a,e)
549     
550        m = Mesh(userVertices=[a,d,f,e], userSegments=[s1,s2,s3,s4])
551
552        m.generateMesh("Q", maxArea = 2.1 )
553        seg = m.getMeshSegments()
554       
555        fileName = tempfile.mktemp(".tsh")
556        m.export_mesh_file(fileName)
557        file = open(fileName)
558        lFile = file.read().split('\n')
559        file.close()
560        os.remove(fileName)
561       
562        #print "@^@^"
563        #for l in lFile:
564        #    print l,"<"
565        #print "@^@^"
566
567        # no need to check the title again
568        #self.failUnless(lFile[0] == "5 0 # <vertex #> <x> <y> [attributes] ...Triangulation Vertices..."
569          #              ,'Ascii file is wrong, vertex title')
570        self.failUnless(lFile[1] == "0 0.0 0.0 " and #1.1 " and
571                        lFile[2] == "1 0.0 4.0 " and #1.2 " and
572                        lFile[3] == "2 4.0 0.0 " and #1.3 " and
573                        lFile[4] == "3 1.0 1.0 " and #1.4 " and
574                        lFile[5] == "4 2.0 2.0 "  #1.25 "
575                        ,
576                        'Ascii file is wrong, vertex')
577       
578        #self.failUnless(lFile[6] == "# attribute column titles ...Triangulation Vertex Titles..."
579          #              ,'Ascii file is wrong, attribute column title')
580        self.failUnless(lFile[8] == "0 3 2 4 -1 2 3  " and
581                        lFile[9] == "1 1 0 3 3 2 -1  " and
582                        lFile[10] == "2 3 4 1 -1 1 0  " and
583                        lFile[11] == "3 2 3 0 1 -1 0  "
584                        ,
585                        'Ascii file is wrong, triangle') 
586
587        self.failUnless( lFile[13] == "0 0 1 exterior" and
588                        lFile[14] == "1 1 4 exterior" and
589                        lFile[15] == "2 2 0 exterior" and
590                        lFile[16] == "3 0 3 " and
591                        lFile[17] == "4 4 2 exterior" ,
592                        'Ascii file is wrong, segment')
593       
594       # self.failUnless(lFile[18] == '4 0 # <vertex #> <x> <y> [attributes] ...Mesh Vertices...',
595        #                'Ascii file is wrong, Mesh Vertices Title')
596       
597        self.failUnless(lFile[19] == '0 0.0 0.0 ' and #1.1 ' and
598                        lFile[20] == '1 0.0 4.0 ' and #1.2 ' and
599                        lFile[21] == '2 4.0 0.0 ' and #1.3 ' and
600                        lFile[22] == '3 1.0 1.0 ' #1.4 '
601                        ,
602                        'Ascii file is wrong, Mesh Vertices II')
603       
604        self.failUnless(lFile[24] == '0 0 1 ' and
605                        lFile[25] == '1 1 2 ' and
606                        lFile[26] == '2 0 2 ' and
607                        lFile[27] == '3 0 3 '
608                        ,'Ascii file is wrong, Mesh Segments')       
609
610 
611    def test_ascii_file(self):
612   
613        a = Vertex (0.0, 0.0) #, attributes = [1.1])
614        d = Vertex (0.0, 4.0) #, attributes = [1.2])
615        f = Vertex (4.0,0.0) #, attributes = [1.3])
616        e = Vertex (1.0,1.0) #, attributes = [1.4])
617   
618        s1 = Segment(a,d)
619        s2 = Segment(d,f)
620        s3 = Segment(a,f)
621        s4 = Segment(a,e)
622     
623        m = Mesh(userVertices=[a,d,f,e], userSegments=[s1,s2,s3,s4])
624
625        m.generateMesh("Q", maxArea = 2.1 )
626
627        seg = m.getMeshSegments()
628       
629        fileName = tempfile.mktemp(".tsh")
630        m.export_mesh_file(fileName)
631        file = open(fileName)
632        lFile = file.read().split('\n')
633        file.close()
634        os.remove(fileName)
635       
636        #print "@^@^"
637        #for l in lFile:
638        #    print l,"<"
639        #print "@^@^"
640        self.failUnless(lFile[0] == "5 0 # <# of verts> <# of vert attributes>, next lines <vertex #> <x> <y> [attributes] ...Triangulation Vertices..."
641                        ,
642                        'Ascii file is wrong, vertex title')
643        self.failUnless(lFile[1] == "0 0.0 0.0 " and #1.1 " and
644                        lFile[2] == "1 0.0 4.0 " and #1.2 " and
645                        lFile[3] == "2 4.0 0.0 " and #1.3 " and
646                        lFile[4] == "3 1.0 1.0 " and #1.4 " and
647                        lFile[5] == "4 2.0 2.0 "  #1.25 "
648                        ,
649                        'Ascii file is wrong, vertex')
650       
651        self.failUnless(lFile[6] == "# attribute column titles ...Triangulation Vertex Titles..."
652                        ,
653                        'Ascii file is wrong, attribute column title')
654        self.failUnless(lFile[7] == "4 # <# of triangles>, next lines <triangle #> [<vertex #>] [<neigbouring triangle #>] [attribute of region] ...Triangulation Triangles..." and
655                        lFile[8] == "0 3 2 4 -1 2 3  " and
656                        lFile[9] == "1 1 0 3 3 2 -1  " and
657                        lFile[10] == "2 3 4 1 -1 1 0  " and
658                        lFile[11] == "3 2 3 0 1 -1 0  "
659                        ,
660                        'Ascii file is wrong, triangle') 
661
662        self.failUnless(lFile[12] == "5 # <# of segments>, next lines <segment #> <vertex #>  <vertex #> [boundary tag] ...Triangulation Segments..." and
663                        lFile[13] == "0 0 1 exterior" and
664                        lFile[14] == "1 1 4 exterior" and
665                        lFile[15] == "2 2 0 exterior" and
666                        lFile[16] == "3 0 3 " and
667                        lFile[17] == "4 4 2 exterior" ,
668                        'Ascii file is wrong, segment')
669       
670        self.failUnless(lFile[18] == '4 0 # <# of verts> <# of vert attributes>, next lines <vertex #> <x> <y> [attributes] ...Mesh Vertices...',
671                        'Ascii file is wrong, Mesh Vertices Title')
672       
673        self.failUnless(lFile[19] == '0 0.0 0.0 ' and #1.1 ' and
674                        lFile[20] == '1 0.0 4.0 ' and #1.2 ' and
675                        lFile[21] == '2 4.0 0.0 ' and #1.3 ' and
676                        lFile[22] == '3 1.0 1.0 ' #1.4 '
677                        ,
678                        'Ascii file is wrong, Mesh Vertices II')
679       
680        self.failUnless(lFile[23] == '4 # <# of segments>, next lines <segment #> <vertex #>  <vertex #> [boundary tag] ...Mesh Segments...' and
681                        lFile[24] == '0 0 1 ' and
682                        lFile[25] == '1 1 2 ' and
683                        lFile[26] == '2 0 2 ' and
684                        lFile[27] == '3 0 3 ' and
685                        lFile[28] == '0 # <# of holes>, next lines <Hole #> <x> <y> ...Mesh Holes...' and
686                        lFile[29] == '0 # <# of regions>, next lines <Region #> <x> <y> <tag>...Mesh Regions...'
687                        ,
688                        'Ascii file is wrong, Mesh Segments')       
689 
690
691    def test_thinoutVertices(self):
692
693        v1 = Vertex(-20,-20)
694        v2 = Vertex(-11,-11)
695        v3 = Vertex(-10,-10)
696        v4 = Vertex(-9,-1)
697        v5 = Vertex(-8,2)
698        v6 = Vertex(6,3)
699        v7 = Vertex(12,9)
700        v8 = Vertex(15,3)
701        v9 = Vertex(24,3)
702        m = Mesh(userVertices = [v1,v2,v3,v4,v5,v6,v7,v8,v9])
703        m.thinoutVertices(10)
704         
705        self.failUnless(v1 in m.userVertices,
706                        'test_thinoutVertices, test 1 failed')
707        self.failUnless(v3 in m.userVertices,
708                        'test_thinoutVertices, test 2 failed')
709        self.failUnless(v4 in m.userVertices,
710                        'test_thinoutVertices, test 3 failed')
711        self.failUnless(v6 in m.userVertices,
712                        'test_thinoutVertices, test 4 failed')
713        self.failUnless(v7 in m.userVertices,
714                        'test_thinoutVertices, test 5 failed')
715        self.failUnless(v9 in m.userVertices,
716                        'test_thinoutVertices, test 6 failed')
717        self.failUnless(v5 not in m.userVertices,
718                        'test_thinoutVertices, test 7 failed')
719        self.failUnless(v2 not in m.userVertices,
720                        'test_thinoutVertices, test 8 failed')
721        self.failUnless(v8 not in m.userVertices,
722                        'test_thinoutVertices, test 9 failed')
723
724    def test_same_x_y(self):
725        v = Point(7,8)
726        f = Point(7,8)
727        f.same_x_y(v)
728
729        self.failUnless(f.same_x_y(v),
730                        'same_x_y True failed')
731        e = Point(7,9)
732        self.failUnless(not f.same_x_y(e),
733                        'same_x_y False failed')
734
735    def test_import_tsh(self):
736       
737        a_att = [5.0,2.0]
738        d_att =[4.0,2.0]
739        f_att = [3.0,2.0]
740        e_att = [2.0,2.0]
741        a_xy = [0.0, 0.0]
742        a = Vertex ( a_xy[0],a_xy[1]) #, attributes =a_att)
743        d = Vertex (0.0, 4.0) #, attributes =d_att)
744        f = Vertex (4.0,0.0) #, attributes =f_att)
745        e = Vertex (1.0,1.0) #, attributes =e_att)
746   
747        s1 = Segment(a,d, tag = "50")
748        s2 = Segment(d,f, tag = "40")
749        s3 = Segment(a,f, tag = "30")
750        s4 = Segment(a,e, tag = "20")
751     
752        r1 = Region(0.3, 0.3,tag = "1.3")
753        geo = Geo_reference(8.9,8.9,65)
754        m = Mesh(userVertices=[a,d,f,e],
755                 userSegments=[s1,s2,s3,s4],
756                 regions=[r1],
757                 geo_reference=geo)
758
759        m.generateMesh("Q", maxArea = 2.1)
760        fileName = tempfile.mktemp(".tsh")
761        #print "dgs!!!"
762        #print "****************** fileName", fileName
763        m.export_mesh_file(fileName)
764        #print "******************"
765        #print "m", m
766        #print "******************"
767        m_returned = importMeshFromFile(fileName)
768        #print "m_returned",m_returned
769        #print "******************"
770        #print "****************** fileName", fileName
771        os.remove(fileName)
772        self.failUnless(0 == m.__cmp__(m_returned),
773                        'loading and saving of a mesh failed')
774        # do this when .msh supports geo_refs
775        #self.failUnless(m.geo_reference == m_returned.geo_reference,
776        #                'loading and saving of a mesh geo refs failed')
777
778    def test_import_mesh(self):
779       
780        a_att = [5.0,2.0]
781        d_att =[4.0,2.0]
782        f_att = [3.0,2.0]
783        e_att = [2.0,2.0]
784        a_xy = [0.0, 0.0]
785        a = Vertex ( a_xy[0],a_xy[1]) #, attributes =a_att)
786        d = Vertex (0.0, 4.0) #, attributes =d_att)
787        f = Vertex (4.0,0.0) #, attributes =f_att)
788        e = Vertex (1.0,1.0) #, attributes =e_att)
789   
790        s1 = Segment(a,d, tag = "50")
791        s2 = Segment(d,f, tag = "40")
792        s3 = Segment(a,f, tag = "30")
793        s4 = Segment(a,e, tag = "20")
794     
795        r1 = Region(0.3, 0.3,tag = "1.3")
796        geo = Geo_reference(65,8.9,8.9)
797        m = Mesh(userVertices=[a,d,f,e],
798                 userSegments=[s1,s2,s3,s4],
799                 regions=[r1],
800                 geo_reference=geo)
801
802        m.generateMesh("Q", maxArea = 2.1)
803        fileName = tempfile.mktemp(".msh")
804        #print "dgs!!!"
805        #print "****************** fileName", fileName
806        m.export_mesh_file(fileName)
807        #print "******************"
808        #print "m", m
809        #print "******************"
810        m_returned = importMeshFromFile(fileName)
811        #print "m_returned",m_returned
812        #print "******************"
813        #print "****************** fileName", fileName
814        os.remove(fileName)
815        #print "m.geo_reference",m.geo_reference
816        #print "m_returned.geo_reference,",m_returned.geo_reference
817        self.failUnless(0 == m.__cmp__(m_returned),
818                        'loading and saving of a mesh failed')
819        self.failUnless(m.geo_reference == m_returned.geo_reference,
820                        'loading and saving of a mesh geo refs failed')
821
822    def test_normaliseMesh(self):
823       
824        a_att = [5.0,2.0]
825        d_att =[4.0,2.0]
826        f_att = [3.0,2.0]
827        e_att = [2.0,2.0]
828        a_xy = [10.0, 10.0]
829        a = Vertex ( a_xy[0],a_xy[1], attributes =a_att)
830        d = Vertex (15.0, 10.0, attributes =d_att)
831        f = Vertex (10.0,20.0, attributes =f_att)
832        e = Vertex (15.0,20.0, attributes =e_att)
833   
834        s1 = Segment(a,d, tag = 50)
835        s2 = Segment(d,e, tag = 40)
836        s3 = Segment(e,f, tag = 30)
837        s4 = Segment(f,a, tag = 20)
838     
839        r1 = Region(0.3, 0.3,tag = 1.3)
840        m = Mesh(userVertices=[a,d,f,e],
841                 userSegments=[s1,s2,s3,s4],
842                 regions=[r1])
843        m.normaliseMesh(1,0,1)
844        [xmin, ymin, xmax, ymax] = m.boxsize()
845        [attmin, attmax] = m.maxMinVertAtt(0)
846        self.failUnless(attmin == 0.0 and attmax == 1.0,
847                        'normalise failed')
848        self.failUnless(xmin == 0.0 and ymin == 0.0 and xmax == 0.5 and ymax == 1.0,
849                        'normalise failed')
850        m.normaliseMesh(200,-100,5)
851        [xmin, ymin, xmax, ymax] = m.boxsize()
852        [attmin, attmax] = m.maxMinVertAtt(0)
853        self.failUnless(attmin == 0.0 and attmax == 5.0,
854                        'normalise failed')
855        self.failUnless(xmin == -100.0 and ymin == -100.0 and xmax == 0.0 and ymax == 100.0,
856                        'normalise failed')
857       
858    def test_exportASCIIsegmentoutlinefile(self):
859        a = Vertex (0,0)
860        b = Vertex (0,3)
861        c = Vertex (3,3)
862        d = Vertex (1,2)
863        e = Vertex (3,1)
864     
865        s1 = Segment(a,b, tag = "50")
866        s2 = Segment(b,c, tag = "40")
867        s3 = Segment(c,a, tag = "30")
868     
869        r1 = Region(2, 1,tag = "1.3")
870        h1 = Hole(1,4)
871        m = Mesh(userVertices=[a,b,c,d,e],
872                 userSegments=[s1,s2,s3],
873                 regions=[r1],
874                 holes = [h1])     
875
876        # vertex e is outside of the outline, so
877        # it is a loner and it is removed.
878        m.generateMesh("Q", maxArea = 2.1)
879        #print "mesh ***************dsg*", m
880
881        fileName = tempfile.mktemp(".tsh")
882        m.exportASCIIsegmentoutlinefile(fileName)
883       
884        m_returned = importMeshFromFile(fileName)
885
886        #print "m_returned ****",m_returned
887        #print "****************** fileName", fileName
888        os.remove(fileName)
889
890        #Trim mesh, so it should like like m_returned
891        m.meshVertices = []
892        m.meshTriangles = []
893        m.meshSegments = []
894        m.userVertices=[a,b,c]
895        #print "mesh ***************dsg*", m
896        #print "(m.__cmp__(m_returned)", m.__cmp__(m_returned)
897        self.failUnless(0 == m.__cmp__(m),
898                        'test_exportASCIIsegmentoutlinefile:loading and saving of a mesh failed')
899        # Having problems with this on linux.
900        #The ordering of the dictionary values wasn't the same as the windows
901        #returned value (verts.values())
902        #self.failUnless(0 == m.__cmp__(m_returned),
903        #                'test_exportASCIIsegmentoutlinefile:loading and saving of a mesh failed')
904       
905        self.failUnless(3 == len(m_returned.userVertices),
906                        'segmentoutlinefile:IO of a mesh failed')
907        self.failUnless(len(m.userSegments) == len(m_returned.userSegments),
908                        'segmentoutlinefile:IO of a mesh failed')
909        for i in range(len(m.userSegments)):
910            self.failUnless(m.userSegments[i].vertices[0].x ==
911                            m_returned.userSegments[i].vertices[0].x,
912                        'loading and saving of a mesh outline fialed')
913            self.failUnless(m.userSegments[i].vertices[0].y ==
914                            m_returned.userSegments[i].vertices[0].y,
915                        'loading and saving of a mesh outline fialed')
916            self.failUnless(m.userSegments[i].vertices[1].x ==
917                            m_returned.userSegments[i].vertices[1].x,
918                        'loading and saving of a mesh outline fialed')
919            self.failUnless(m.userSegments[i].vertices[1].y ==
920                            m_returned.userSegments[i].vertices[1].y,
921                        'loading and saving of a mesh outline fialed')
922
923 
924    def test_exportASCIIsegmentoutlinefile2(self):
925        a = Vertex (0,0)
926        b = Vertex (0,1)
927        c = Vertex (1,0)
928        d = Vertex (1,1)
929        e = Vertex (0.5,0.5)
930        f  = Vertex (0.6,0.6)
931     
932        s1 = Segment(a,e, tag = "50")
933        s2 = Segment(b,e, tag = "40")
934        s3 = Segment(c,e, tag = "30")
935        s4 = Segment(d,e, tag = "30")
936     
937        r1 = Region(2, 1,tag = "1.3")
938        h1 = Hole(1,4)
939        m = Mesh(userVertices=[a,b,c,d,e],
940                 userSegments=[s1,s2,s3,s4],
941                 regions=[r1],
942                 holes = [h1])     
943       
944        fileName = tempfile.mktemp(".tsh")
945        m.exportASCIIsegmentoutlinefile(fileName)
946       
947        m_returned = importMeshFromFile(fileName)
948        #print "****************** fileName", fileName
949        os.remove(fileName)
950
951        #Trim mesh, so it should look like m_returned
952        m.meshVertices = []
953        m.meshTriangles = []
954        m.meshSegments = []
955        m.userVertices=[a,e,d,b,c]
956        #print "mesh ***************dsg*", m
957        #print "(m.__cmp__(m_returned)", m.__cmp__(m_returned)
958        self.failUnless(0 == m.__cmp__(m),
959                        'loading and saving of a mesh failed')
960
961        self.failUnless(5 == len(m_returned.userVertices),
962                        'segmentoutlinefile:IO of a mesh failed')
963        self.failUnless(len(m.userSegments) == len(m_returned.userSegments),
964                        'segmentoutlinefile:IO of a mesh failed')
965        for i in range(len(m.userSegments)):
966            self.failUnless(m.userSegments[i].vertices[0].x ==
967                            m_returned.userSegments[i].vertices[0].x,
968                        'loading and saving of a mesh outline fialed')
969            self.failUnless(m.userSegments[i].vertices[0].y ==
970                            m_returned.userSegments[i].vertices[0].y,
971                        'loading and saving of a mesh outline fialed')
972            self.failUnless(m.userSegments[i].vertices[1].x ==
973                            m_returned.userSegments[i].vertices[1].x,
974                        'loading and saving of a mesh outline fialed')
975            self.failUnless(m.userSegments[i].vertices[1].y ==
976                            m_returned.userSegments[i].vertices[1].y,
977                        'loading and saving of a mesh outline fialed')
978
979
980    def test_loadxy(self):
981        """
982        To test the mesh side of loading xya files.
983        Not the loading of xya files
984        """
985        import os
986        import tempfile
987       
988        fileName = tempfile.mktemp(".xya")
989        file = open(fileName,"w")
990        file.write("elevation speed \n\
9911.0 0.0 10.0 0.0\n\
9920.0 1.0 0.0 10.0\n\
9931.0 0.0 10.4 40.0\n")
994        file.close()
995        #print fileName
996        m = importMeshFromFile(fileName)
997        os.remove(fileName)
998        self.failUnless(m.userVertices[0].x == 1.0,
999                        'loadxy, test 1 failed')
1000        self.failUnless(m.userVertices[0].y == 0.0,
1001                        'loadxy, test 2 failed')
1002        #self.failUnless(m.userVertices[0].attributes == [10.0,0.0],
1003        #                'loadxy, test 2.2 failed')
1004        self.failUnless(m.userVertices[1].x == 0.0,
1005                        'loadxy, test 3 failed')
1006        self.failUnless(m.userVertices[1].y == 1.0,
1007                        'loadxy, test 4 failed')
1008        #self.failUnless(m.userVertices[1].attributes == [0.0,10.0],
1009        #                'loadxy, test 5 failed')
1010       
1011    def test_exportPointsFile(self):
1012        a = Vertex (0,0)
1013        b = Vertex (0,3)
1014        c = Vertex (3,3)
1015        d = Vertex (1,2)
1016        e = Vertex (3,1)
1017        #f = Vertex (3,1)
1018     
1019        s1 = Segment(a,b, tag = 50)
1020        s2 = Segment(b,c, tag = 40)
1021        s3 = Segment(c,a, tag = 30)
1022     
1023        r1 = Region(2, 1,tag = 1.3)
1024        h1 = Hole(1,4)
1025        # Warning mesh can't produce this type of data structure its self
1026        m = Mesh(userVertices=[a,b,c,d,e],
1027                 userSegments=[s1,s2,s3],
1028                 regions=[r1],
1029                 holes = [h1])
1030       
1031        fileName = tempfile.mktemp(".xya")
1032        #fileName = 't.xya'
1033        #os.remove(fileName)
1034        m.exportPointsFile(fileName)
1035        file = open(fileName)
1036        lFile = file.read().split('\n')
1037        file.close()
1038
1039        os.remove(fileName)
1040        self.failUnless(lFile[0] == "" and
1041                        lFile[1] == "0,0" and
1042                        lFile[2] == "0,3" and
1043                        lFile[3] == "3,3" 
1044                        ,
1045                        'exported Ascii xya file is wrong')
1046        self.failUnless(lFile[4] == "1,2" and
1047                        lFile[5] == "3,1" 
1048                        ,
1049                        'exported Ascii xya file is wrong')
1050       
1051        # vertex e is outside of the outline, so
1052        # it is a loner and it is removed.
1053        m.generateMesh("Q", maxArea = 2.1)
1054        fileName = tempfile.mktemp(".xya")
1055        #fileName = 't.xya'
1056        #m.export_mesh_file('m.tsh')
1057        m.exportPointsFile(fileName)
1058        file = open(fileName)
1059        lFile = file.read().split('\n')
1060        file.close()
1061        os.remove(fileName)
1062       
1063        self.failUnless(lFile[0] == "" and
1064                        lFile[1] == "0.0,0.0" and
1065                        lFile[2] == "0.0,3.0" and
1066                        lFile[3] == "3.0,3.0" and
1067                        lFile[4] == "1.0,2.0"
1068                        ,
1069                        'exported Ascii xya file is wrong')
1070           
1071    def test_exportPointsFilefile2(self):
1072        m = Mesh()
1073       
1074        fileName = tempfile.mktemp(".xya")
1075        m.exportPointsFile(fileName)
1076        file = open(fileName)
1077        lFile = file.read().split('\n')
1078        file.close()
1079
1080        os.remove(fileName)
1081        #print "************* test_mesh exportPointsFilefile"
1082        #print "lFile",lFile
1083        #print "************* test_mesh exportPointsFilefile"
1084        self.failUnless(lFile[0] == "" 
1085                        ,
1086                        'exported Ascii xya file is wrong')
1087       
1088    def test_strings2ints(self):
1089        list = ["sea","river inlet","","sea","","moat"]
1090        preset = ["moat", "internal boundary"]
1091        [intlist, converter] = segment_strings2ints(list,preset )
1092        self.failUnless(intlist == [2,3 ,0 ,2 ,0 ,0 ]
1093                        ,
1094                        'test_strings2ints produces bad intlist')
1095        self.failUnless(converter == ['moat', 'internal boundary',
1096                                      'sea', 'river inlet']
1097                        ,
1098                        'test_strings2ints produces bad converter')
1099       
1100    def test_ints2strings(self):
1101        list = ["internal boundary","sea","river inlet",
1102            "","sea","","moat","internal boundary"]
1103        outlist = ['internal boundary', 'sea', 'river inlet', 'moat',
1104                   'sea', 'moat', 'moat', 'internal boundary']
1105        preset = ["moat", "internal boundary"]
1106        [intlist, converter] = segment_strings2ints(list,preset )
1107        newlist = segment_ints2strings(intlist, converter)
1108        self.failUnless(outlist == newlist
1109                        ,
1110                        'test_strings2ints produces bad intlist')
1111        self.failUnless(converter == ['moat', 'internal boundary',
1112                                      'sea', 'river inlet']
1113                        ,
1114                        'test_strings2ints produces bad converter')
1115       
1116    def test_ints2strings2(self):
1117        list = ["","",""]
1118        preset = ["moat", "internal boundary"]
1119        [intlist, converter] = segment_strings2ints(list,preset )
1120        newlist = segment_ints2strings(intlist, converter)
1121        outlist = ['moat', 'moat', 'moat']
1122        self.failUnless(outlist == newlist
1123                        ,
1124                        'test_strings2ints produces bad intlist')
1125        self.failUnless(converter == ['moat', 'internal boundary']
1126                        ,
1127                        'test_strings2ints produces bad converter')
1128
1129       
1130    def test_removeDuplicatedVertices(self):
1131        a = Vertex (0,0)
1132        a.index = 0
1133        b = Vertex (0,3)
1134        b.index = 1
1135        c = Vertex (3,3)
1136        c.index = 2
1137        d = Vertex (1,1)
1138        d.index = 3
1139        e = Vertex (3,1)
1140        e.index = 4
1141        f = Vertex (1,1)
1142        f.index = 5
1143        g = Vertex (1,1)
1144        g.index = 6
1145        inputVerts_noDups = [a,b,c,d,e]
1146       
1147        m = Mesh(userVertices=[a,b,c,d,e,f,g])
1148        counter = m.removeDuplicatedUserVertices()
1149        UserVerts = m.getUserVertices()
1150       
1151         
1152        self.failUnless(UserVerts == inputVerts_noDups,
1153                            'duplicate verts not removed')
1154        #for userVert, inputVert in map(None, UserVerts, inputVerts_noDups):
1155        #    self.failUnless(userVert.x == inputVert.x,
1156        #                    'x duplicate verts not removed')
1157        #    self.failUnless(userVert.y == inputVert.y,
1158        #                    'y duplicate verts not removed')
1159
1160       
1161    def test_ungenerateFileLoading(self):
1162       
1163        fileName = tempfile.mktemp(".txt")
1164        file = open(fileName,"w")
1165        file.write("         1       ??      ??\n\
1166       0.0       0.0\n\
1167       1.0       0.0\n\
1168       1.0       1.0\n\
1169       0.0       1.0\n\
1170       0.0       0.0\n\
1171END\n\
1172         2      ?? ??\n\
1173       10.0       10.0\n\
1174       10.0       20.0\n\
1175       20.0       20.0\n\
1176       10.0       10.0\n\
1177END\n\
1178END\n")
1179        file.close()
1180       
1181       
1182        a = Vertex (0.0, 0.0) #, attributes = [1.1])
1183        b = Vertex (0.0, 40.0) #, attributes = [1.2])
1184        c = Vertex (40.0,40.0) #, attributes = [1.3])
1185        d = Vertex (40.0,0.0) #, attributes = [1.4])
1186   
1187        s1 = Segment(a,b)
1188        s2 = Segment(b,c)
1189        s3 = Segment(c,d)
1190        s4 = Segment(d,a)
1191     
1192        m = Mesh(userVertices=[a,b,c,d], userSegments=[s1,s2,s3,s4])
1193        dict = importUngenerateFile(fileName)
1194        #os.remove(fileName)
1195
1196        tag = "DSG"
1197        Segment.set_default_tag(tag)
1198        m.addVertsSegs(dict)
1199
1200        # have to reset this , since it's a class attribute
1201        Segment.set_default_tag("")
1202           
1203        self.failUnless(len(m.userSegments) ==11,
1204                        'Wrong segment list length.')
1205        self.failUnless(len(m.userVertices) == 11,
1206                        'Wrong vertex list length.')
1207        self.failUnless(m.userSegments[10].vertices[0] == m.userVertices[10],
1208                        'bad vertex on segment.')
1209        self.failUnless(m.userSegments[10].vertices[1] == m.userVertices[8],
1210                        'Bad segment.')
1211        self.failUnless(m.userSegments[10].tag == tag,
1212                        'wrong tag.')
1213
1214        ## let's test the method
1215        a = Vertex (0.0, 0.0) #, attributes = [1.1])
1216        b = Vertex (0.0, 40.0) #, attributes = [1.2])
1217        c = Vertex (40.0,40.0) #, attributes = [1.3])
1218        d = Vertex (40.0,0.0) #, attributes = [1.4])
1219   
1220        s1 = Segment(a,b)
1221        s2 = Segment(b,c)
1222        s3 = Segment(c,d)
1223        s4 = Segment(d,a)
1224     
1225        m = Mesh(userVertices=[a,b,c,d], userSegments=[s1,s2,s3,s4])
1226
1227        tag = "DSG"       
1228        initial_tag = "PIG"
1229        Segment.set_default_tag(initial_tag)
1230        m.import_ungenerate_file(fileName, tag=tag)
1231
1232        os.remove(fileName)
1233
1234        self.failUnless(Segment.get_default_tag() == initial_tag,
1235                        'Wrong segment list length.')
1236       
1237
1238        # have to reset this , since it's a class attribute
1239        Segment.set_default_tag("")
1240           
1241        self.failUnless(len(m.userSegments) ==11,
1242                        'Wrong segment list length.')
1243        self.failUnless(len(m.userVertices) == 11,
1244                        'Wrong vertex list length.')
1245        self.failUnless(m.userSegments[10].vertices[0] == m.userVertices[10],
1246                        'bad vertex on segment.')
1247        self.failUnless(m.userSegments[10].vertices[1] == m.userVertices[8],
1248                        'Bad segment.')
1249        self.failUnless(m.userSegments[10].tag == tag,
1250                        'wrong tag.')
1251       
1252    def test_ungenerateFileLoadingII(self):
1253       
1254        fileName = tempfile.mktemp(".txt")
1255        file = open(fileName,"w")
1256        file.write("         1       ??      ??\n\
1257       0.0       0.0\n\
1258       1.0       0.0\n\
1259       1.0       1.0\n\
1260       0.0       1.0\n\
1261       0.0       0.0\n\
1262END\n\
1263         2      ?? ??\n\
1264       10.0       10.0\n\
1265       10.0       20.0\n\
1266       20.0       20.0\n\
1267END\n\
1268END\n")
1269        file.close()
1270       
1271       
1272        a = Vertex (0.0, 0.0) #, attributes = [1.1])
1273        b = Vertex (0.0, 40.0) #, attributes = [1.2])
1274        c = Vertex (40.0,40.0) #, attributes = [1.3])
1275        d = Vertex (40.0,0.0) #, attributes = [1.4])
1276   
1277        s1 = Segment(a,b)
1278        s2 = Segment(b,c)
1279        s3 = Segment(c,d)
1280        s4 = Segment(d,a)
1281     
1282        m = Mesh(userVertices=[a,b,c,d], userSegments=[s1,s2,s3,s4])
1283        dict = importUngenerateFile(fileName)
1284        #os.remove(fileName)
1285
1286        tag = "DSG"
1287        Segment.set_default_tag(tag)
1288        m.addVertsSegs(dict)
1289
1290        self.failUnless(len(m.userSegments) ==10,
1291                        'Wrong segment list length.')
1292        self.failUnless(len(m.userVertices) == 11,
1293                        'Wrong vertex list length.')
1294
1295        # Test the method
1296        a = Vertex (0.0, 0.0) #, attributes = [1.1])
1297        b = Vertex (0.0, 40.0) #, attributes = [1.2])
1298        c = Vertex (40.0,40.0) #, attributes = [1.3])
1299        d = Vertex (40.0,0.0) #, attributes = [1.4])
1300   
1301        s1 = Segment(a,b)
1302        s2 = Segment(b,c)
1303        s3 = Segment(c,d)
1304        s4 = Segment(d,a)
1305     
1306        m = Mesh(userVertices=[a,b,c,d], userSegments=[s1,s2,s3,s4])
1307        tag = "DSG"       
1308        initial_tag = "PIG"
1309        Segment.set_default_tag(initial_tag)
1310        m.import_ungenerate_file(fileName, tag=tag)
1311
1312        os.remove(fileName)
1313
1314        self.failUnless(Segment.get_default_tag() == initial_tag,
1315                        'Wrong segment list length.')
1316       
1317       
1318        self.failUnless(len(m.userSegments) ==10,
1319                        'Wrong segment list length.')
1320        self.failUnless(len(m.userVertices) == 11,
1321                        'Wrong vertex list length.')
1322       
1323        # have to reset this , since it's a class attribute
1324        Segment.set_default_tag("")
1325       
1326    def test_addVertsSegs(self):
1327        m = Mesh()
1328        Segment.set_default_tag("food")
1329        dict = {}
1330        dict['points'] = [[0.0, 0.0], [1.0, 0.0], [1.0, 1.0]]
1331        dict['segments'] = [[0, 1], [1, 2]]
1332        dict['segment_tags'] = ['','do-op']
1333        m.addVertsSegs(dict)
1334        # have to reset this , since it's a class attribute
1335        Segment.set_default_tag("")
1336
1337       
1338        self.failUnless(len(m.userSegments) ==2,
1339                        'Wrong segment list length.')
1340        self.failUnless(len(m.userVertices) == 3,
1341                        'Wrong vertex list length.')
1342        self.failUnless(m.userSegments[0].tag =='food',
1343                        'Wrong segment tag length.')
1344        self.failUnless(m.userSegments[1].tag =='do-op',
1345                        'Wrong segment tag.')
1346       
1347    def test_addVertsSegs2(self):
1348        geo = Geo_reference(56,5,10)
1349        m = Mesh(geo_reference=geo)
1350        dict = {}
1351        dict['points'] = [[2.0, 1.0], [3.0, 1.0], [2.0, 2.0]]
1352        dict['segments'] = [[0, 1], [1, 2], [2,0]]
1353        dict['segment_tags'] = ['','do-op','']
1354        m.addVertsSegs(dict)
1355
1356    def test_exportASCIImeshfile(self):
1357   
1358        #a_att = [5,2]
1359        #d_att =[4,2]
1360        #f_att = [3,2]
1361        #e_att = [2,2]
1362        a_xy = [0.0, 0.0]
1363        a = Vertex ( a_xy[0],a_xy[1]) #, attributes =a_att)
1364        d = Vertex (0.0, 4.0) #, attributes =d_att)
1365        f = Vertex (4.0,0.0) #, attributes =f_att)
1366        e = Vertex (1.0,1.0) #, attributes =e_att)
1367   
1368        s1 = Segment(a,d, tag = "50")
1369        s2 = Segment(d,f, tag = "40")
1370        s3 = Segment(a,f, tag = "30")
1371        s4 = Segment(a,e, tag = "20")
1372     
1373        r1 = Region(0.3, 0.3,tag = "1.3", maxArea = 36)
1374
1375
1376        h1 = Hole(0.2,0.6)
1377       
1378        m = Mesh(userVertices=[a,d,f,e],
1379                 userSegments=[s1,s2,s3,s4],
1380                 regions=[r1],
1381                 holes=[h1])
1382
1383        seg = m.getUserSegments()
1384        points = m.getUserVertices()
1385        holes = m.getHoles()
1386        regions = m.getRegions()
1387        fileName = tempfile.mktemp(".tsh")
1388        m.export_mesh_file(fileName)
1389        #print "***************************fileName", fileName
1390        new_m = importMeshFromFile(fileName)
1391        os.remove(fileName)
1392       
1393
1394        #print '**@@@@@******'
1395        #print "new_m",new_m
1396        #print '**@@@@@******'
1397        #print "m",m
1398        #print '**@@@@@******'
1399       
1400        self.failUnless( new_m == m,
1401                         'loadASCIITestCase failed. test new 1')
1402           
1403    def test_Mesh2MeshList(self):
1404
1405        a_att = [5,2]
1406        d_att =[4,2]
1407        f_att = [3,2]
1408        e_att = [2,2]
1409        a_xy = [0.0, 0.0]
1410        a = Vertex ( a_xy[0],a_xy[1]) #, attributes =a_att)
1411        d = Vertex (0.0, 4.0) #, attributes =d_att)
1412        f = Vertex (4.0,0.0) #, attributes =f_att)
1413        e = Vertex (1.0,1.0) #, attributes =e_att)
1414   
1415        s1 = Segment(a,d, tag = "50")
1416        s2 = Segment(d,f, tag = "40")
1417        s3 = Segment(a,f, tag = "30")
1418        s4 = Segment(a,e, tag = "20")
1419     
1420        r1 = Region(0.3, 0.3,tag = "1.3", maxArea = 45)
1421        m = Mesh(userVertices=[a,d,f,e],
1422                 userSegments=[s1,s2,s3,s4],
1423                 regions=[r1])
1424
1425        m.generateMesh("Qa2.1")
1426
1427        seg = m.getMeshSegments()
1428        points = m.getMeshVertices()
1429        dict = m.Mesh2MeshList()
1430        #print "dict",dict
1431        # test not finished...
1432 
1433    def test_Mesh2IOTriangulationDict(self):
1434
1435        a_att = [5,2]
1436        d_att =[4,2]
1437        f_att = [3,2]
1438        e_att = [2,2]
1439        a_xy = [0.0, 0.0]
1440        a = Vertex ( a_xy[0],a_xy[1] , attributes =a_att)
1441        d = Vertex (0.0, 4.0 , attributes =d_att)
1442        f = Vertex (4.0,0.0 , attributes =f_att)
1443        e = Vertex (1.0,1.0 , attributes =e_att)
1444   
1445        s1 = Segment(a,d, tag = "50")
1446        s2 = Segment(d,f, tag = "40")
1447        s3 = Segment(a,f, tag = "30")
1448        s4 = Segment(a,e, tag = "20")
1449     
1450        r1 = Region(0.3, 0.3,tag = "1.3", maxArea = 45)
1451        m = Mesh(userVertices=[a,d,f,e],
1452                 userSegments=[s1,s2,s3,s4],
1453                 regions=[r1])
1454        titles = ['ele','friction']
1455        m.attributeTitles = titles
1456        m.generateMesh("Qa2.1")
1457
1458        seg = m.getMeshSegments()
1459        verts = m.getMeshVertices()
1460        dict = m.Mesh2IOTriangulationDict()
1461        #print "dict",dict
1462       
1463        self.failUnless( dict['vertex_attribute_titles'] == titles,
1464                         'test_Mesh2IOTriangulationDict failed. test 1')
1465        answer = [a_xy,[0.0, 4.0],[4.0,0.0], [1.0,1.0], [2.0,2.0]]
1466        #print "answer",answer
1467        #print "dict['vertices']",dict['vertices']
1468       
1469        self.failUnless( dict['vertices'] == answer,
1470                         'test_Mesh2IOTriangulationDict failed. test 2')
1471
1472       
1473        for pimport,pactual,pimpatt in map(None,dict['vertices'],
1474                                           verts,dict['vertex_attributes']):
1475            #assert all_close( pimport, (pactual.x,pactual.y))
1476            self.failUnless( pimport == [pactual.x,pactual.y],
1477                        'test_Mesh2IOTriangulationDict failed. test 2.1')
1478            self.failUnless( pimpatt == pactual.attributes,
1479                        'test_Mesh2IOTriangulationDict failed. test 2.2')
1480        self.failUnless( dict['segments'][0] == [0,1],
1481                        'test_Mesh2IOTriangulationDict failed. test 3')
1482        for segimp,segactual in map(None,dict['segment_tags'],seg):
1483            self.failUnless( segimp == segactual.tag,
1484                        'test_Mesh2IOTriangulationDict failed. test 4')
1485        #print " dict['triangles'][0]", dict['triangles'][0]
1486        self.failUnless( dict['triangles'][0] == [3,2,4],
1487                        'test_Mesh2IOTriangulationDict failed. test 5')
1488        self.failUnless( dict['triangle_neighbors'][0] == [-1,2,3],
1489                        'test_Mesh2IOTriangulationDict failed. test 6')
1490        self.failUnless( dict['triangle_tags'][0] == "1.3",
1491                         'test_Mesh2IOTriangulationDict failed. test 7')
1492
1493 
1494    def test_Mesh2IODict(self):
1495
1496        a_att = [5,2]
1497        d_att =[4,2]
1498        f_att = [3,2]
1499        e_att = [2,2]
1500        a_xy = [0.0, 0.0]
1501        a = Vertex ( a_xy[0],a_xy[1] , attributes =a_att)
1502        d = Vertex (0.0, 4.0 , attributes =d_att)
1503        f = Vertex (4.0,0.0 , attributes =f_att)
1504        e = Vertex (1.0,1.0 , attributes =e_att)
1505   
1506        s1 = Segment(a,d, tag = "50")
1507        s2 = Segment(d,f, tag = "40")
1508        s3 = Segment(a,f, tag = "30")
1509        s4 = Segment(a,e, tag = "20")
1510     
1511        r1 = Region(0.3, 0.3,tag = "1.3", maxArea = 45)
1512        m = Mesh(userVertices=[a,d,f,e],
1513                 userSegments=[s1,s2,s3,s4],
1514                 regions=[r1])
1515        titles = ['ele','friction']
1516        m.attributeTitles = titles
1517        m.generateMesh("Qa2.1")
1518
1519        seg = m.getMeshSegments()
1520        verts = m.getMeshVertices()
1521        dict = m.Mesh2IODict()
1522        #print "dict",dict
1523       
1524        self.failUnless( dict['vertex_attribute_titles'] == titles,
1525                         'test_Mesh2IOTriangulationDict failed. test 1')
1526        answer = [a_xy,[0.0, 4.0],[4.0,0.0], [1.0,1.0], [2.0,2.0]]
1527        #print "answer",answer
1528        #print "dict['vertices']",dict['vertices']
1529       
1530        self.failUnless( dict['vertices'] == answer,
1531                         'test_Mesh2IOTriangulationDict failed. test 2')
1532
1533       
1534        for pimport,pactual,pimpatt in map(None,dict['vertices'],
1535                                           verts,dict['vertex_attributes']):
1536            #assert all_close( pimport, (pactual.x,pactual.y))
1537            self.failUnless( pimport == [pactual.x,pactual.y],
1538                        'test_Mesh2IODict failed. test 2.1')
1539            self.failUnless( pimpatt == pactual.attributes,
1540                        'test_Mesh2IODict failed. test 2.2')
1541        self.failUnless( dict['segments'][0] == [0,1],
1542                        'test_Mesh2IODict failed. test 3')
1543        for segimp,segactual in map(None,dict['segment_tags'],seg):
1544            self.failUnless( segimp == segactual.tag,
1545                        'test_Mesh2IODict failed. test 4')
1546        #print " dict['triangles'][0]", dict['triangles'][0]
1547        self.failUnless( dict['triangles'][0] == [3,2,4],
1548                        'test_Mesh2IODict failed. test 5')
1549        self.failUnless( dict['triangle_neighbors'][0] == [-1,2,3],
1550                        'test_Mesh2IODict failed. test 6')
1551        self.failUnless( dict['triangle_tags'][0] == "1.3",
1552                         'test_Mesh2IODict failed. test 7')
1553
1554        seg = m.getUserSegments()
1555        points = m.getUserVertices()
1556        holes = m.getHoles()
1557        regions = m.getRegions()
1558       
1559        for pimport,pactual,pimpatt in map(None,dict['points'],points,dict['point_attributes']):
1560            self.failUnless( pimport == [pactual.x,pactual.y],
1561                        'test_Mesh2IODict failed. test 1')
1562            self.failUnless( pimpatt == pactual.attributes,
1563                        'test_Mesh2IODict failed. test 1.1')
1564        self.failUnless( dict['outline_segments'][0] == [0,1],
1565                        'test_Mesh2IODict failed. test 3')
1566        for segimp,segactual in map(None,dict['outline_segment_tags'],seg):
1567            self.failUnless( segimp == segactual.tag,
1568                        'test_Mesh2IODict failed. test 4')
1569        for holeimp,holeactual in map(None,dict['holes'],holes):
1570            self.failUnless( holeimp == [holeactual.x,holeactual.y],
1571                        'test_Mesh2IODict failed. test 5')
1572       
1573        for regimp,regactual,regattimp, regmaxarea in map(None,dict['regions'],regions, dict['region_tags'], dict['region_max_areas']):
1574            self.failUnless( regimp == [regactual.x,regactual.y],
1575                        'loadASCIITestCase failed. test 6')
1576            self.failUnless( regattimp == regactual.getTag(),
1577                        'loadASCIITestCase failed. test 7')
1578            self.failUnless( regmaxarea == regactual.getMaxArea(),
1579                        'loadASCIITestCase failed. test 7')
1580   
1581           
1582       
1583    def test_Mesh2IOOutlineDict(self):
1584
1585        a_att = [5,2]
1586        d_att =[4,2]
1587        f_att = [3,2]
1588        e_att = [2,2]
1589        a_xy = [0.0, 0.0]
1590        a = Vertex ( a_xy[0],a_xy[1] , attributes =a_att)
1591        d = Vertex (0.0, 4.0 , attributes =d_att)
1592        f = Vertex (4.0,0.0 , attributes =f_att)
1593        e = Vertex (1.0,1.0 , attributes =e_att)
1594   
1595        s1 = Segment(a,d, tag = "50")
1596        s2 = Segment(d,f, tag = "40")
1597        s3 = Segment(a,f, tag = "30")
1598        s4 = Segment(a,e, tag = "20")
1599     
1600        r1 = Region(0.3, 0.3,tag = "1.3", maxArea = 45)
1601        m = Mesh(userVertices=[a,d,f,e],
1602                 userSegments=[s1,s2,s3,s4],
1603                 regions=[r1])
1604        titles = ['ele','friction']
1605        m.attributeTitles = titles
1606        m.generateMesh("Qa2.1")
1607
1608        seg = m.getMeshSegments()
1609        verts = m.getMeshVertices()
1610        dict = m.Mesh2IOOutlineDict()
1611       
1612        seg = m.getUserSegments()
1613        points = m.getUserVertices()
1614        holes = m.getHoles()
1615        regions = m.getRegions()
1616       
1617        for pimport,pactual,pimpatt in map(None,dict['points'],points,dict['point_attributes']):
1618            self.failUnless( pimport == [pactual.x,pactual.y],
1619                        'loadASCIITestCase failed. test 1')
1620            self.failUnless( pimpatt == pactual.attributes,
1621                        'loadASCIITestCase failed. test 1.1')
1622        self.failUnless( dict['outline_segments'][0] == [0,1],
1623                        'loadASCIITestCase failed. test 3')
1624        for segimp,segactual in map(None,dict['outline_segment_tags'],seg):
1625            self.failUnless( segimp == segactual.tag,
1626                        'loadASCIITestCase failed. test 4')
1627        for holeimp,holeactual in map(None,dict['holes'],holes):
1628            self.failUnless( holeimp == [holeactual.x,holeactual.y],
1629                        'loadASCIITestCase failed. test 5')
1630        #for regimp,regactual in map(None,dict['regions'],regions):
1631         #   self.failUnless( [regimp[0],regimp[1]]==[regactual.x,regactual.y],
1632          #              'loadASCIITestCase failed. test 6')
1633           # self.failUnless( regimp[2] == regactual.getTag(),
1634            #            'loadASCIITestCase failed. test 7')
1635            #self.failUnless( regimp[3] == regactual.getMaxArea(),
1636             #           'loadASCIITestCase failed. test 7')
1637
1638           
1639        for regimp,regactual,regattimp, regmaxarea in map(None,dict['regions'],regions, dict['region_tags'], dict['region_max_areas']):
1640            self.failUnless( regimp == [regactual.x,regactual.y],
1641                        'loadASCIITestCase failed. test 6')
1642            self.failUnless( regattimp == regactual.getTag(),
1643                        'loadASCIITestCase failed. test 7')
1644            self.failUnless( regmaxarea == regactual.getMaxArea(),
1645                        'loadASCIITestCase failed. test 7')
1646
1647
1648#___________beginning of Peters tests
1649
1650    def test_set_stuff(self):
1651        """
1652        Documentation
1653        """
1654        #making a test mesh
1655        p0=[0.,2.]
1656        p1=[1.,2.]
1657        p2=[0.,1.]
1658        p3=[1.,1.]
1659        p4=[0.,0.]
1660        p5=[2.,0.]
1661        p6=[-1.,1.]
1662        point_list = [p0,p1,p2,p3,p4,p5,p6]
1663
1664        a0=[0]
1665        a1=[0]
1666        a2=[100]
1667        a3=[0]
1668        a4=[0]
1669        a5=[0]
1670        a6=[0]
1671        attribute=[a0,a1,a2,a3,a4,a5,a6]
1672       
1673        t0=[0,3,1]
1674        t1=[0,2,3]
1675        t2=[2,4,3]
1676        t3=[4,5,3]
1677        t4=[1,3,5]
1678        t5=[2,6,4]
1679
1680        n0=[4,-1,2]
1681        n1=[2,0,-1]
1682        n2=[3,1,5]
1683        n3=[4,2,-1]
1684        n4=[3,-1,0]
1685        n5=[-1,2,-1]
1686
1687        tri_list = [t0,t1,t2,t3,t4,t5]
1688        n_list = [n0,n1,n2,n3,n4,n5]
1689        for i in range(6):
1690            for j in (0,1,2):
1691                a=attribute[tri_list[i][j]]
1692                tri_list[i][j]=point_list[tri_list[i][j]]
1693                tri_list[i][j]=Vertex(tri_list[i][j][0]\
1694                                      ,tri_list[i][j][1],a)
1695            neighbours=n_list[i]
1696            tri_list[i]=Triangle(tri_list[i][0],\
1697                                 tri_list[i][1],tri_list[i][2]\
1698                                ,neighbors=neighbours)
1699
1700        #testing selectAll
1701        mesh = Mesh()
1702        mesh.attributeTitles=['attrib']
1703
1704        mesh.meshTriangles=tri_list
1705
1706        mesh.selectAllTriangles()
1707        A=mesh.sets[mesh.setID['All']]
1708        assert list_comp(tri_list,A)
1709
1710       #testing threshold
1711        mesh = Mesh()
1712        mesh.attributeTitles=['attrib']
1713
1714        mesh.meshTriangles=tri_list
1715        mesh.selectAllTriangles()
1716        mesh.threshold('All',min=30,max=35,attribute_name = 'attrib')
1717        A = [tri_list[1],tri_list[2],tri_list[5]]
1718        B = mesh.sets[mesh.setID['All']]
1719        assert list_comp(A,B)
1720       
1721
1722        A = [tri_list[3],tri_list[2],tri_list[5]]
1723        assert not list_comp(A,B)
1724
1725        #testing
1726
1727    def test_Discretised_Tuple_Set_rounding(self):
1728        #This is the hardest bit of DST
1729
1730        tol = 0.1
1731        a=Discretised_Tuple_Set(p_rel=1,t_rel= tol)
1732        m = 0.541
1733        m_up = 0.6
1734        m_down = 0.5
1735        assert m_up == a.round_up_rel(m)
1736        assert m_down == a.round_down_rel(m)
1737
1738        tol = 0.1
1739        a=Discretised_Tuple_Set(p_rel=1,t_rel = tol)
1740        m = 0.539
1741        m_up = 0.5
1742        m_down = 0.5
1743        assert m_up == a.round_up_rel(m)
1744        assert m_down == a.round_down_rel(m)
1745
1746        tol = 0.5
1747        a=Discretised_Tuple_Set(p_rel=1,t_rel = tol)
1748
1749
1750        m = 0.6
1751        m_up = 0.7
1752        m_down = 0.5
1753        assert m_up == a.round_up_rel(m)
1754        assert m_down == a.round_down_rel(m)
1755
1756        m = 0.599
1757        m_up = 0.6
1758        m_down = 0.5
1759        assert m_up == a.round_up_rel(m)
1760        assert m_down == a.round_down_rel(m)
1761
1762    def test_Discretised_Tuple_Set_get(self):
1763       
1764        tol = 0.25
1765        a=Discretised_Tuple_Set(p_rel=1,t_rel = tol)
1766        b = (1.1,1.1)
1767        a.append(b)
1768        list = [(1.2,1.),(1.,1.),(1.,1.2),(1.2,1.2)]
1769        for key in list:
1770            assert a[key][0]==b
1771            assert len(a[key])==1
1772       
1773        c = (2.1,1.)
1774        a.append(c)
1775        assert a[(2.,1.)][0]==c
1776        assert a[(2.2,1.)][0]==c
1777
1778    def test_mapped_Discretised_Tuple_Set(self):
1779
1780        def map(sequence):
1781            return [len(sequence)]
1782
1783        tol = 0.5
1784        a=Mapped_Discretised_Tuple_Set(map,p_rel=1,t_rel = tol)
1785        b = range(20)
1786        a.append(b)
1787        assert b in a[range(17)] 
1788        assert b in a[range(22)]
1789
1790        tol = 0.01
1791        a=Mapped_Discretised_Tuple_Set(map,p_rel=1,t_rel = tol)
1792        b = range(20)
1793        a.append(b)
1794        assert b in a[range(20)] 
1795        assert b in a[range(19)] 
1796        assert not range(17) in a
1797
1798#___________end of Peters tests
1799
1800    def test_add_region_from_polygon(self):
1801        m=Mesh()
1802        region = m.add_region_from_polygon([[0,0],[1,0],[0,1]],
1803                                  max_triangle_area = 88)
1804        self.failUnless(len(m.regions)==1,
1805                        'FAILED!')
1806        self.failUnless(region.getMaxArea()==88,
1807                        'FAILED!')
1808        self.failUnless(len(m.getUserSegments())==3,
1809                        'FAILED!')
1810        self.failUnless(len(m.userVertices)==3,
1811                        'FAILED!')
1812       
1813    def test_add_region_from_polygon2(self):
1814        m=Mesh()
1815        m.add_region_from_polygon([[0,0],[1,0],[1,1],[0,1]],
1816                               {'tagin':[0,1],'bom':[2]},
1817                                  max_triangle_area=10)
1818        self.failUnless(len(m.regions)==1,
1819                        'FAILED!')
1820        segs = m.getUserSegments()
1821        self.failUnless(len(segs)==4,
1822                        'FAILED!')
1823        self.failUnless(len(m.userVertices)==4,
1824                        'FAILED!') 
1825        self.failUnless(segs[0].tag=='tagin',
1826                        'FAILED!') 
1827        self.failUnless(segs[1].tag=='tagin',
1828                        'FAILED!') 
1829         
1830        self.failUnless(segs[2].tag=='bom',
1831                        'FAILED!') 
1832        self.failUnless(segs[3].tag=='',
1833                        'FAILED!') 
1834       
1835    def test_add_region_from_polygon3(self):
1836        x=-500
1837        y=-1000
1838        m=Mesh(geo_reference=Geo_reference(56,x,y))
1839
1840        # These are the absolute values
1841        polygon_absolute = [[0,0],[1,0],[1,1],[0,1]]
1842       
1843        x_p = -10
1844        y_p = -40
1845        geo_ref_poly = Geo_reference(56, x_p, y_p)
1846        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
1847       
1848        poly_point = m.add_region_from_polygon(polygon,
1849                                               {'tagin':[0,1],'bom':[2]},
1850                                               geo_reference=geo_ref_poly,
1851                                               max_triangle_area=10)
1852        # poly_point values are relative to the mesh geo-ref
1853        # make them absolute
1854        self.failUnless(is_inside_polygon([poly_point.x+x,poly_point.y+y],
1855                                       polygon_absolute, closed = False),
1856                        'FAILED!')
1857               
1858        self.failUnless(len(m.regions)==1,
1859                        'FAILED!')
1860        segs = m.getUserSegments()
1861        self.failUnless(len(segs)==4,
1862                        'FAILED!')
1863        self.failUnless(len(m.userVertices)==4,
1864                        'FAILED!') 
1865        self.failUnless(segs[0].tag=='tagin',
1866                        'FAILED!') 
1867        self.failUnless(segs[1].tag=='tagin',
1868                        'FAILED!') 
1869         
1870        self.failUnless(segs[2].tag=='bom',
1871                        'FAILED!') 
1872        self.failUnless(segs[3].tag=='',
1873                        'FAILED!')
1874        verts = m.getUserVertices()
1875        #print "User verts",verts
1876        #print 'polygon',polygon
1877        #vert values are relative
1878        for point,new_point in map(None,polygon,verts):
1879            point_x = point[0] + geo_ref_poly.get_xllcorner()
1880            new_point_x = new_point.x + m.geo_reference.get_xllcorner()
1881            point_y = point[1] + geo_ref_poly.get_yllcorner()
1882            #print "new_point.y",new_point.y
1883            #print "m.geo_ref.get_yllcorner()",m.geo_reference.get_yllcorner()
1884            new_point_y = new_point.y + m.geo_reference.get_yllcorner()
1885            #print "point_y",point_y
1886            #print "new_point_y",new_point_y
1887           
1888            self.failUnless(point_x == new_point_x, ' failed')
1889            self.failUnless(point_y == new_point_y, ' failed')
1890           
1891         
1892    def test_add_region_from_polygon4(self):
1893        x=50000
1894        y=1000
1895        m=Mesh(geo_reference=Geo_reference(56,x,y))
1896        polygon = [[0,0],[1,0],[1,1],[0,1]]
1897       
1898        m.add_region_from_polygon(polygon,
1899                               {'tagin':[0,1],'bom':[2]},
1900                                  max_triangle_area=10)
1901        self.failUnless(len(m.regions)==1,
1902                        'FAILED!')
1903        segs = m.getUserSegments()
1904        self.failUnless(len(segs)==4,
1905                        'FAILED!')
1906        self.failUnless(len(m.userVertices)==4,
1907                        'FAILED!') 
1908        self.failUnless(segs[0].tag=='tagin',
1909                        'FAILED!') 
1910        self.failUnless(segs[1].tag=='tagin',
1911                        'FAILED!') 
1912         
1913        self.failUnless(segs[2].tag=='bom',
1914                        'FAILED!') 
1915        self.failUnless(segs[3].tag=='',
1916                        'FAILED!')
1917        verts = m.getUserVertices()
1918        #print "User verts",verts
1919        #print 'polygon',polygon
1920        #vert values are relative
1921        for point,new_point in map(None,polygon,verts):
1922            point_x = point[0] 
1923            new_point_x = new_point.x + m.geo_reference.get_xllcorner()
1924            #print "point_x",point_x
1925            #print "new_point_x",new_point_x
1926            point_y = point[1] 
1927            new_point_y = new_point.y + m.geo_reference.get_yllcorner()
1928           
1929            self.failUnless(point_x == new_point_x, ' failed')
1930            self.failUnless(point_y == new_point_y, ' failed')
1931
1932
1933    def test_add_hole_from_polygon(self):
1934        x=-500
1935        y=-1000
1936        m=Mesh(geo_reference=Geo_reference(56,x,y))
1937
1938        # These are the absolute values
1939        polygon_absolute = [[0,0],[1,0],[1,1],[0,1]]
1940       
1941        x_p = -10
1942        y_p = -40
1943        geo_ref_poly = Geo_reference(56, x_p, y_p)
1944        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
1945       
1946        poly_point = m.add_hole_from_polygon(polygon,
1947                                               {'tagin':[0,1],'bom':[2]},
1948                                               geo_reference=geo_ref_poly)
1949        # poly_point values are relative to the mesh geo-ref
1950        # make them absolute
1951        #print "poly_point.x+x",poly_point.x+x
1952        #print "poly_point.y+y",poly_point.y+y
1953        #print "polygon_absolute", polygon_absolute
1954        self.failUnless(is_inside_polygon([poly_point.x+x,poly_point.y+y],
1955                                       polygon_absolute, closed = False),
1956                        'FAILED!')
1957               
1958        self.failUnless(len(m.holes)==1,
1959                        'FAILED!')
1960        segs = m.getUserSegments()
1961        self.failUnless(len(segs)==4,
1962                        'FAILED!')
1963        self.failUnless(len(m.userVertices)==4,
1964                        'FAILED!') 
1965        self.failUnless(segs[0].tag=='tagin',
1966                        'FAILED!') 
1967        self.failUnless(segs[1].tag=='tagin',
1968                        'FAILED!') 
1969         
1970        self.failUnless(segs[2].tag=='bom',
1971                        'FAILED!') 
1972        self.failUnless(segs[3].tag=='',
1973                        'FAILED!')
1974        verts = m.getUserVertices()
1975        #print "User verts",verts
1976        #print 'polygon',polygon
1977        #vert values are relative
1978        for point,new_point in map(None,polygon,verts):
1979            point_x = point[0] + geo_ref_poly.get_xllcorner()
1980            new_point_x = new_point.x + m.geo_reference.get_xllcorner()
1981            point_y = point[1] + geo_ref_poly.get_yllcorner()
1982            #print "new_point.y",new_point.y
1983            #print "m.geo_ref.get_yllcorner()",m.geo_reference.get_yllcorner()
1984            new_point_y = new_point.y + m.geo_reference.get_yllcorner()
1985            #print "point_y",point_y
1986            #print "new_point_y",new_point_y
1987           
1988            self.failUnless(point_x == new_point_x, ' failed')
1989            self.failUnless(point_y == new_point_y, ' failed')
1990
1991    def test_add_circle(self):
1992        x=-500
1993        y=-1000
1994        m=Mesh(geo_reference=Geo_reference(56,x,y))
1995
1996        # These are the absolute values
1997        tag = 'hey'
1998        segment_count = 104
1999        radius = 30
2000        circle_center_absolute = [100,80]       
2001        x_p = -.666
2002        y_p = -.777
2003        geo_ref_poly = Geo_reference(56, x_p, y_p)
2004        circle_center = \
2005                geo_ref_poly.change_points_geo_ref(circle_center_absolute)
2006        circle_center = circle_center[0] #make a list of lists a list
2007        poly_point = m.add_circle(circle_center, radius, segment_count,
2008                                  tag=tag,
2009                                  region=True,
2010                                  center_geo_reference=geo_ref_poly)
2011        # poly_point values are relative to the mesh geo-ref
2012        # make them absolute
2013        #print "poly_point.x+x",poly_point.x+x
2014        #print "polygon_absolute", polygon_absolute
2015     
2016       
2017        #m.export_mesh_file("aaat.msh")
2018       
2019        self.failUnless(len(m.regions)==1,
2020                        'FAILED!')
2021        segs = m.getUserSegments()
2022        self.failUnless(len(segs)==segment_count,
2023                        'FAILED!')
2024        self.failUnless(len(m.userVertices)==segment_count,
2025                        'FAILED!') 
2026        self.failUnless(segs[0].tag==tag,
2027                        'FAILED!') 
2028        self.failUnless(segs[1].tag==tag,
2029                        'FAILED!') 
2030         
2031        verts = m.getUserVertices()
2032       
2033        #m.export_mesh_file("aaat.msh")
2034       
2035    def FIXMEtest_auto_set_geo_reference(self):
2036        x=50000
2037        y=1000
2038        m=Mesh(geo_reference=Geo_reference(56,x,y))
2039        polygon = [[0,0],[1,0],[1,1],[0,1]]
2040       
2041        m.add_region_from_polygon(polygon,
2042                               {'tagin':[0,1],'bom':[2]},
2043                                  max_triangle_area=10)
2044        m.auto_set_geo_reference()
2045       
2046 
2047    def test_duplicat_verts_are_removed(self):
2048   
2049     
2050        a = Vertex ( 0.0 ,0.0)
2051        b = Vertex (0.0, 4.0)
2052        c = Vertex (4.0,4.0)
2053        d = Vertex (4.0,0.0)
2054        e = Vertex (4.0,0.0) # duplicate point
2055   
2056        s1 = Segment(a,b, tag = "50")
2057        s2 = Segment(b,c, tag = "40")
2058        s3 = Segment(c,d, tag = "30")
2059        s4 = Segment(d,e, tag = "no where seg")
2060        s5 = Segment(e,a, tag = "20")
2061
2062       
2063        m = Mesh(userVertices=[a,b,c,d,e],
2064                 userSegments=[s1,s2,s3,s4,s5])
2065
2066        seg = m.getUserSegments()
2067        points = m.getUserVertices()
2068        holes = m.getHoles()
2069        regions = m.getRegions()
2070        #fileName = tempfile.mktemp(".tsh")
2071        #fileName = "badmesh.tsh"
2072        #m.export_mesh_file(fileName)
2073        #print "***************************fileName", fileName
2074        #new_m = importMeshFromFile(fileName)
2075        #os.remove(fileName)
2076       
2077        m.generateMesh("Q", maxArea = 2000.1 )
2078
2079        #m.export_mesh_file("from_test_mesh.tsh")
2080        seg = m.getMeshSegments()
2081        self.failUnless(4==len(seg),
2082                        'FAILED!') 
2083
2084        vert = m.getMeshVertices() 
2085        self.failUnless(4==len(vert),
2086                        'FAILED!')
2087 
2088    def test_duplicat_verts_are_removedII(self):
2089   
2090     
2091        a = Vertex ( 0.0 ,0.0)
2092        b = Vertex (0.0, 4.0)
2093        c = Vertex (4.0,4.0)
2094        d = Vertex (4.0,0.0)
2095        e = Vertex (4.0,0.0) # duplicate point
2096        f = Vertex (49.0,0.0) # unused point
2097   
2098        s1 = Segment(a,b, tag = "50")
2099        s2 = Segment(b,c, tag = "40")
2100        s3 = Segment(c,d, tag = "30")
2101        s4 = Segment(d,e, tag = "no where seg")
2102        s5 = Segment(e,a, tag = "20")
2103
2104       
2105        m = Mesh(userVertices=[a,b,c,d,e,f],
2106                 userSegments=[s1,s2,s3,s4,s5])
2107
2108        seg = m.getUserSegments()
2109        points = m.getUserVertices()
2110        holes = m.getHoles()
2111        regions = m.getRegions()
2112        #fileName = tempfile.mktemp(".tsh")
2113        #fileName = "badmesh.tsh"
2114        #m.export_mesh_file(fileName)
2115        #print "***************************fileName", fileName
2116        #new_m = importMeshFromFile(fileName)
2117        #os.remove(fileName)
2118       
2119        m.generateMesh("Q", maxArea = 2000.1 )
2120
2121        #m.export_mesh_file("from_test_mesh.tsh")
2122        seg = m.getMeshSegments()
2123        self.failUnless(4==len(seg),
2124                        'FAILED!') 
2125
2126        vert = m.getMeshVertices() 
2127        self.failUnless(4==len(vert),
2128                        'FAILED!')
2129   
2130    def test_add_vertices(self):
2131        points_ab = [[0.1,1],[0.4,.2],[7,5],[10,5]]
2132        geo =  Geo_reference(56,23,21)
2133        points = geo.change_points_geo_ref(points_ab)
2134        spat = Geospatial_data(points, geo_reference=geo)
2135       
2136        geo_mesh =  Geo_reference(56,100,200)
2137        m = Mesh(geo_reference=geo_mesh)
2138        m.add_vertices(spat)
2139
2140        vert = m.getUserVertices()
2141        #print "vert",vert
2142        self.failUnless(4==len(vert),
2143                        'FAILED!')
2144        vert= m.get_user_vertices(absolute=True)
2145        self.failUnless(vert==points_ab,
2146                        'FAILED!')
2147
2148   
2149    def test_add_vertices_more(self):
2150        points = [[0.1,1],[0.4,.2],[7,5],[10,5]]
2151        #spat = Geospatial_data(points)
2152       
2153        m = Mesh()
2154        m.add_vertices(points)
2155
2156        vert = m.getUserVertices()
2157        #print "vert",vert
2158        self.failUnless(4==len(vert),
2159                        'FAILED!')
2160        vert= m.get_user_vertices(absolute=True)
2161        self.failUnless(vert==points,
2162                        'FAILED!')
2163   
2164    def test_add_verticesII(self):
2165        points_lat_long = [[-33,152],[-35,152],[-35,150],[-33,150]]
2166       
2167        spat = Geospatial_data(data_points=points_lat_long,
2168                               points_are_lats_longs=True)
2169        points_ab = spat.get_data_points( absolute = True)
2170        m = Mesh()
2171        m.add_vertices(spat)
2172
2173        vert = m.getUserVertices()
2174        #print "vert",vert
2175        self.failUnless(4==len(vert),
2176                        'FAILED!')
2177        vert= m.get_user_vertices(absolute=True)
2178        self.failUnless(vert==points_ab,
2179                        'FAILED!')
2180
2181        spat = Geospatial_data(data_points=points_lat_long,
2182                               points_are_lats_longs=True)
2183        points_ab = spat.get_data_points( absolute = True)
2184        geo =  Geo_reference(56,400000,6000000)
2185        spat.set_geo_reference(geo)
2186        m = Mesh()
2187        m.add_vertices(spat)
2188
2189        vert = m.getUserVertices()
2190        #print "vert",vert
2191        self.failUnless(4==len(vert),
2192                        'FAILED!')
2193        vert= m.get_user_vertices(absolute=True)
2194        self.failUnless(vert==points_ab,
2195                        'FAILED!')
2196
2197        #geo =  Geo_reference(56,23,21)
2198        #points = geo.change_points_geo_ref(points_ab)
2199       
2200    def test_get_user_vertices(self):
2201        points_ab = [[0.1,1],[0.4,.2],[7,5],[10,5]]
2202        geo =  Geo_reference(56,23,21)
2203        points = geo.change_points_geo_ref(points_ab)
2204        spat = Geospatial_data(points, geo_reference=geo)
2205       
2206        geo_mesh =  Geo_reference(56,100,200)
2207        m = Mesh(geo_reference=geo_mesh)
2208        m.add_vertices(spat)
2209
2210        vert = m.getUserVertices()
2211        #print "vert",vert
2212        self.failUnless(4==len(vert),
2213                        'FAILED!')
2214        vert= m.get_user_vertices(absolute=True)
2215        self.failUnless(vert==points_ab,
2216                        'FAILED!')
2217        vert= m.get_user_vertices(absolute=False)
2218        points_new = m.geo_reference.get_absolute(vert)
2219       
2220        self.failUnless(points_ab==points_new,
2221                        'FAILED!')
2222       
2223       
2224def list_comp(A,B):
2225    yes = len(A)==len(B)
2226    for item in A:
2227        if not item in B:
2228            yes = False
2229    return yes
2230
2231
2232           
2233#-------------------------------------------------------------
2234if __name__ == "__main__":
2235    suite = unittest.makeSuite(meshTestCase,'test')
2236    #suite = unittest.makeSuite(meshTestCase,'test_generate_mesh')
2237    runner = unittest.TextTestRunner() #verbosity=2)
2238    runner.run(suite)
2239   
Note: See TracBrowser for help on using the repository browser.