""" Make a line of gauges from outside boundary to near centre polygon for Onslow study to look at MOST wave and compare to ANUGA output """ import project from pylab import plot, xlabel, ylabel, title, ion, axis, savefig, close, text from utilities.polygon import poly_xy from Numeric import arange from pmesh.create_mesh import convert_points_from_latlon_to_utm x_bound, y_bound = poly_xy(project.polyAll) MOST_south = project.south MOST_north = project.north MOST_east = project.east MOST_west = project.west p0 = [MOST_south, MOST_west] p1 = [MOST_south, MOST_east] p2 = [MOST_north, MOST_east] p3 = [MOST_north, MOST_west] MOSTpoly, zone = convert_points_from_latlon_to_utm([p0, p1, p2, p3]) MOSTymax = 0 MOSTymin = 7800000.0 MOSTxmax = 0 MOSTxmin = 7800000.0 for i, point in enumerate(MOSTpoly): x = point[0] y = point[1] if y > MOSTymax: MOSTymax = y if y < MOSTymin: MOSTymin = y if x > MOSTxmax: MOSTxmax = x if x < MOSTxmin: MOSTxmin = x MOSTymax = MOSTymax-3000.0 print 'MOST clipping: ', MOSTxmin, MOSTxmax, MOSTymin, MOSTymax # nominated point in centre region for Onslow x2 = 305000.0 y2 = 7605000.0 testy = 7625000.0 d0 = project.d0 d1 = project.d1 d2 = project.d2 d3 = project.d3 d = [d0, d1, d2, d3] ion() fid = open(project.gauge_comparison,'w') s = 'Easting,Northing,Name,Elevation \n' fid.write(s) count = 0 miny = 7580000.0 maxy = 7800000.0 xplot = [] yplot = [] names = [] for i, point in enumerate(d): x = point[0] y = point[1] m = (y-y2)/(x-x2) c = y2-m*x2 if m > 0: plotx = arange(MOSTxmin,MOSTxmax,1000.0) for j, xpoint in enumerate(plotx): ypoint = m*xpoint+c if ypoint > testy and ypoint <= d0[1]: #if ypoint >= testy and ypoint <= MOSTymax: if xpoint > MOSTxmin and xpoint < MOSTxmax: xplot.append(xpoint) yplot.append(ypoint) thisname = 'Point%d' %count names.append(thisname) count += 1 s = '%.2f,%.2f,%s,%.2f \n' %(xpoint,ypoint,thisname,0.0) fid.write(s) else: plotx = arange(MOSTxmin,MOSTxmax,5000.0) for j, xpoint in enumerate(plotx): ypoint = m*xpoint+c if ypoint > testy and ypoint < MOSTymax: #testy if xpoint > MOSTxmin and xpoint < MOSTxmax: xplot.append(xpoint) yplot.append(ypoint) thisname = 'Point%d' %count names.append(thisname) count += 1 s = '%.2f,%.2f,%s,%.2f \n' %(xpoint,ypoint,thisname,0.0) fid.write(s) plot(x_bound, y_bound, xplot, yplot,'b+') for i, string in enumerate(names): if xplot[i] > 240000 and xplot[i] < 340000: text(xplot[i],yplot[i],string,fontsize=7) axis([240000,340000, miny, maxy]) savefig('boundarycomparison') close('all')