source: trunk/anuga_core/source/anuga/pmesh/test_mesh.py @ 7753

Last change on this file since 7753 was 7711, checked in by hudson, 14 years ago

Refactored geometry classes to live in their own folder.

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