import os, sys sys.path.append('..') from mesh import * from coordinate_transforms.geo_reference import Geo_reference #------------------------------------------------------------- if __name__ == "__main__": bellambi_northing = 6195000.0 north_beach_northing = 6183000.0 north_beach_easting = 308500 border = 12000.0 west = 301000.0 # From the map inner_east = 310000.0 # From the map inner_west_south = 303000.0 # From the map inner_west_north = 307000.0 # From the map east = 328421.0 sea_boundary = (east+inner_east)/2 inner_center_easting = ( inner_east +inner_west_north)/2 inner_center_northing = (bellambi_northing + north_beach_northing)/2 # from 100m dem geo = Geo_reference(xllcorner = inner_center_easting, yllcorner = inner_center_northing, zone = 56) print "***********************" print "geo ref", geo print "***********************" m = Mesh(geo_reference=geo) dict = {} dict['points'] = [[west, north_beach_northing - border], #sw [west, bellambi_northing + border], #nw [sea_boundary, bellambi_northing + border], #n_sea_boundary [east, bellambi_northing + border], #ne [east, north_beach_northing - border], #se [sea_boundary, north_beach_northing - border] #s_sea_boundary ] dict['segments'] = [[0,1],[1,2],[2,3], [3,4],[4,5],[5,0], # the outer boarder [2,5] # the sea_boundary line ] m.addVertsSegs(dict) dict = {} dict['points'] = [ [inner_west_south,north_beach_northing], #sw [inner_west_north,bellambi_northing], #nw [inner_east,bellambi_northing], #ne [inner_east,north_beach_northing] #se ] dict['segments'] = [[0,1],[1,2],[2,3],[3,0]] # the inner boarder m.addVertsSegs(dict) #Even more inner points flagstaff_x = 308000 flagstaff_y = 6189000 west_x = 305000 west_y = 6185000 north_beach_northing = 6183000.0 dict = {} dict['points'] = [ [flagstaff_x,flagstaff_y], #n [west_x,west_y], #w [north_beach_easting,north_beach_northing], #s ] dict['segments'] = [[0,1],[1,2],[2,0]] # the inner boarder m.addVertsSegs(dict) factor = 1000 #low res 10000, high res 1000 low = m.addRegionEN(east-0.5, bellambi_northing + border-0.5) low.setMaxArea(100*factor) medium = m.addRegionEN(sea_boundary-0.5, bellambi_northing + border-0.5) medium.setMaxArea(10*factor) high = m.addRegionEN(inner_east-0.5,bellambi_northing-0.5) high.setMaxArea(5*factor) #doing 0.5, fac = 1000 gives degenerate error center = m.addRegionEN(flagstaff_x,flagstaff_y-0.5) center.setMaxArea(0.05*factor) #doing 0.5, fac = 1000 gives degenerate error #m.generateMesh() m.generateMesh("pzq28.0za1000000a") m.export_mesh_file("wollongong_4regions.msh") #m.export_mesh_file("wollongong_highres_a28.tsh")