source: anuga_core/source/anuga/pmesh/test_mesh.py @ 6172

Last change on this file since 6172 was 6172, checked in by rwilson, 15 years ago

Removed code that *might* have stomped on the Numeric 'num' module.

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