Changeset 2276 for inundation/pmesh


Ignore:
Timestamp:
Jan 25, 2006, 12:11:06 PM (19 years ago)
Author:
duncan
Message:

adding new pmesh.mesh interface

Location:
inundation/pmesh
Files:
2 added
3 edited

Legend:

Unmodified
Added
Removed
  • inundation/pmesh/create_gong.py

    r1086 r2276  
    3333    m = Mesh(geo_reference=geo)
    3434   
    35     dict = {}
    36     dict['points'] = [[west, north_beach_northing - border],  #sw
     35    factor = 1000 #low res 10000, high res 1000
     36   
     37    poly = [[west, north_beach_northing - border],  #sw
    3738                         [west, bellambi_northing + border],  #nw
    3839                         [sea_boundary,
     
    4344                          north_beach_northing - border]    #s_sea_boundary
    4445                         ]
    45     dict['segments'] = [[0,1],[1,2],[2,3],
    46                            [3,4],[4,5],[5,0], # the outer boarder
    47                            [2,5]   # the sea_boundary line
    48                            ]
     46    tags = {'wall':[0,1,2,3,4,5]}
     47
     48    m.add_region_from_polygon(poly, tags, max_triangle_area=factor*1000)
     49    poly =  [
     50        [inner_west_south,north_beach_northing], #sw
     51        [inner_west_north,bellambi_northing],   #nw
     52        [inner_east,bellambi_northing], #ne
     53        [inner_east,north_beach_northing] #se
     54        ]
     55    m.add_region_from_polygon(poly, max_triangle_area=factor*100)
    4956   
    50     m.addVertsSegs(dict)
    51 
    52     dict = {}
    53     dict['points'] = [
    54                          [inner_west_south,north_beach_northing], #sw
    55                          [inner_west_north,bellambi_northing],   #nw
    56                          [inner_east,bellambi_northing], #ne
    57                          [inner_east,north_beach_northing] #se
    58                          ]
    59     dict['segments'] = [[0,1],[1,2],[2,3],[3,0]] # the inner boarder
    60     m.addVertsSegs(dict)
    6157
    6258    #Even more inner points
     
    6864    north_beach_northing = 6183000.0
    6965   
    70     dict = {}
    71     dict['points'] = [
     66    poly = [
    7267                         [flagstaff_x,flagstaff_y], #n
    7368                         [west_x,west_y],   #w
    7469                         [north_beach_easting,north_beach_northing], #s
    7570                         ]
    76     dict['segments'] = [[0,1],[1,2],[2,0]] # the inner boarder
    77     m.addVertsSegs(dict)
     71    m.add_region_from_polygon(poly, max_triangle_area=factor*10)
    7872   
    79     factor = 1000 #low res 10000, high res 1000
    80     low = m.addRegionEN(east-0.5, bellambi_northing + border-0.5)
    81     low.setMaxArea(100*factor)
    82    
    83     medium = m.addRegionEN(sea_boundary-0.5, bellambi_northing + border-0.5)
    84     medium.setMaxArea(10*factor)
    85    
    86     high = m.addRegionEN(inner_east-0.5,bellambi_northing-0.5)
    87     high.setMaxArea(5*factor) #doing 0.5, fac = 1000 gives degenerate error
    88    
    89     center = m.addRegionEN(flagstaff_x,flagstaff_y-0.5)
    90     center.setMaxArea(0.05*factor) #doing 0.5, fac = 1000 gives degenerate error
    9173    #m.generateMesh()
    9274    m.generateMesh("pzq28.0za1000000a")
  • inundation/pmesh/mesh.py

    r2261 r2276  
    604604   
    605605    def add_region(self, x,y, geo_reference=None):
    606         [x,y] = self.geo_reference.change_points_geo_ref([x,y],
     606        """
     607        adds a point, which represents a region.
     608
     609        The point data can have it's own geo_refernece.
     610        If geo_reference is None the data is asumed to be absolute
     611        """
     612        [[x,y]] = self.geo_reference.change_points_geo_ref([x,y],
    607613                                                 points_geo_ref=geo_reference)
    608614        return self._addRegion(x, y)
     
    614620   
    615621    def add_region_from_polygon(self, polygon, tags=None,
    616                              max_area=None, geo_reference=None):
     622                                max_triangle_area=None, geo_reference=None):
    617623        """
    618624        Add a polygon with tags to the current mesh, as a region.
    619625        The maxArea of the region can be specified.
    620626
    621         If a geo_ref is given, this is used.
     627        If a geo_reference of the polygon points is given, this is used.
    622628        If not;
    623629        The x,y info is assumed to be Easting and Northing, absolute,
     
    632638       
    633639        """
    634         #FIXME: take into account georef on mesh side
    635         polygon = self.geo_reference.change_points_geo_ref(polygon,
    636                                                  points_geo_ref=geo_reference)
    637         #print "polygon - should be relative to mesh geo_ref",polygon
     640        #take into account georef on mesh side
     641        #polygon = self.geo_reference.change_points_geo_ref(polygon,
     642                         #                        points_geo_ref=geo_reference)
     643        #print "polygon - should be relative to poly geo_ref",polygon
     644
     645        #get absolute values
     646        if geo_reference is not None:
     647            polygon = geo_reference.get_absolute(polygon)
     648        # polygon is now absolute
     649        #print "polygon  should be absolute",polygon
    638650        #create points, segs and tags
    639651        region_dict = {}
     
    666678   
    667679        self.addVertsSegs(region_dict) #this is assuming absolute geos
    668    
    669         #get inner point
    670         inner_point = point_in_polygon(polygon)
    671         inner = self.add_region(inner_point[0], inner_point[1],
    672                                geo_reference=self.geo_reference)
    673 
    674         if max_area is not None:
    675             inner.setMaxArea(max_area)
     680
     681        inner = None
     682       
     683        if max_triangle_area is not None:
     684            #get inner point - absolute values
     685            inner_point = point_in_polygon(polygon)
     686            inner = self.add_region(inner_point[0], inner_point[1],
     687                                    geo_reference=None)
     688            inner.setMaxArea(max_triangle_area)
     689           
    676690        return inner
    677691       
     
    728742        return s
    729743       
    730     def generateMesh(self, mode = None, maxArea = None,
     744    def generate_mesh(self,
     745                      maximum_triangle_area=None,
     746                      minimum_triangle_angle=None,
     747                      verbose=False):
     748        if verbose is True:
     749            silent = ''
     750        else:
     751            silent = 'Q'
     752        self.generateMesh(mode = silent +"pzq"+str(minimum_triangle_angle)
     753                                  +"a"+str(maximum_triangle_area)
     754                                  +"a")
     755        #The last a is so areas for regions will be used
     756       
     757    def generateMesh(self, mode = None, maxArea = None, minAngle= None,
    731758                     isRegionalMaxAreas = True):
    732759        """
     
    743770            self.mode = mode
    744771       
    745         if not re.match('p',self.mode):
     772        if not self.mode.find('p'):
    746773            self.mode += 'p' #p - read a planar straight line graph.
    747774            #there must be segments to use this switch
    748775            # TODO throw an aception if there aren't seg's
    749776            # it's more comlex than this.  eg holes
    750         if not re.match('z',self.mode):
     777        if not self.mode.find('z'):
    751778            self.mode += 'z' # z - Number all items starting from zero
    752779                             # (rather than one)
    753         if not re.match('n',self.mode):
     780        if self.mode.find('n'):
    754781            self.mode += 'n' # n - output a list of neighboring triangles
    755         if not re.match('A',self.mode):
     782        if not self.mode.find('A'):
    756783            self.mode += 'A' # A - output region attribute list for triangles
    757         if not re.match('V',self.mode) and not re.match('Q',self.mode):
     784
     785        if not self.mode.find('V') and not self.mode.find('Q'):
    758786            self.mode += 'V' # V - output info about what Triangle is doing
    759787       
     788        if not self.mode.find('q') and minAngle is not None:
     789            #   print "**********8minAngle******** ",minAngle
     790            min_angle = 'q' + str(minAngle)
     791            self.mode += min_angle # z - Number all items starting from zero
     792                             # (rather than one)
    760793        if maxArea != None:
    761794            self.mode += 'a' + str(maxArea)
    762 
     795        #FIXME (DSG-DSG) This isn't explained.  I don't think it works
     796        # well with  maxArea = None as well.
    763797        if isRegionalMaxAreas:
    764798            self.mode += 'a'
     
    12551289
    12561290        Assume the values are in Eastings and Northings, with no reference
    1257         point
     1291        point. eg absolute
    12581292        """
    12591293        if not outlineDict.has_key('segment_tags'):
  • inundation/pmesh/test_mesh.py

    r2261 r2276  
    66from load_mesh.loadASCII import *
    77from coordinate_transforms.geo_reference import Geo_reference
     8from utilities.polygon import inside_polygon
    89
    910class meshTestCase(unittest.TestCase):
    1011    def setUp(self):
    1112        pass
    12        
    13 
    14 
     13   
    1514    def tearDown(self):
    1615        pass
     
    7978       
    8079 
    81        
     80    # FIXME add test for minAngle   
    8281    def testgenerateMesh(self):
    8382        a = Vertex (0.0, 0.0)
     
    16571656    def test_add_region_from_polygon(self):
    16581657        m=Mesh()
    1659         m.add_region_from_polygon([[0,0],[1,0],[0,1]])
     1658        region = m.add_region_from_polygon([[0,0],[1,0],[0,1]],
     1659                                  max_triangle_area = 88)
    16601660        self.failUnless(len(m.regions)==1,
     1661                        'FAILED!')
     1662        self.failUnless(region.getMaxArea()==88,
    16611663                        'FAILED!')
    16621664        self.failUnless(len(m.getUserSegments())==3,
     
    16681670        m=Mesh()
    16691671        m.add_region_from_polygon([[0,0],[1,0],[1,1],[0,1]],
    1670                                {'tagin':[0,1],'bom':[2]})
     1672                               {'tagin':[0,1],'bom':[2]},
     1673                                  max_triangle_area=10)
    16711674        self.failUnless(len(m.regions)==1,
    16721675                        'FAILED!')
     
    16871690       
    16881691    def test_add_region_from_polygon3(self):
    1689         x=0
    1690         y=0
     1692        x=-500
     1693        y=-1000
    16911694        m=Mesh(geo_reference=Geo_reference(56,x,y))
    1692         polygon = [[0,0],[1,0],[1,1],[0,1]]
    1693         x_p = 1000
    1694         y_p = 4000
     1695
     1696        # These are the absolute values
     1697        polygon_absolute = [[0,0],[1,0],[1,1],[0,1]]
     1698       
     1699        x_p = -10
     1700        y_p = -40
    16951701        geo_ref_poly = Geo_reference(56, x_p, y_p)
    1696         polygon = geo_ref_poly.change_points_geo_ref(polygon)
    1697         #print "polygon", polygon
    1698         m.add_region_from_polygon(polygon,
    1699                                {'tagin':[0,1],'bom':[2]},
    1700                                   geo_reference=geo_ref_poly)
     1702        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
     1703       
     1704        poly_point = m.add_region_from_polygon(polygon,
     1705                                               {'tagin':[0,1],'bom':[2]},
     1706                                               geo_reference=geo_ref_poly,
     1707                                               max_triangle_area=10)
     1708        # poly_point values are relative to the mesh geo-ref
     1709        # make them absolute
     1710        self.failUnless(inside_polygon([poly_point.x+x,poly_point.y+y],
     1711                                       polygon_absolute, closed = False),
     1712                        'FAILED!')
     1713               
    17011714        self.failUnless(len(m.regions)==1,
    17021715                        'FAILED!')
     
    17231736            new_point_x = new_point.x + m.geo_reference.get_xllcorner()
    17241737            point_y = point[1] + geo_ref_poly.get_yllcorner()
     1738            #print "new_point.y",new_point.y
     1739            #print "m.geo_ref.get_yllcorner()",m.geo_reference.get_yllcorner()
    17251740            new_point_y = new_point.y + m.geo_reference.get_yllcorner()
    17261741            #print "point_y",point_y
     
    17311746           
    17321747         
    1733         #REJIG
    17341748    def test_add_region_from_polygon4(self):
    1735         x=0
    1736         y=0
     1749        x=50000
     1750        y=1000
    17371751        m=Mesh(geo_reference=Geo_reference(56,x,y))
    17381752        polygon = [[0,0],[1,0],[1,1],[0,1]]
    17391753       
    17401754        m.add_region_from_polygon(polygon,
    1741                                {'tagin':[0,1],'bom':[2]})
     1755                               {'tagin':[0,1],'bom':[2]},
     1756                                  max_triangle_area=10)
    17421757        self.failUnless(len(m.regions)==1,
    17431758                        'FAILED!')
Note: See TracChangeset for help on using the changeset viewer.