setplot.py.html | |
Source file: setplot.py | |
Directory: /Users/rjl/git/clawpack/geoclaw/examples/multi-layer/plane_wave | |
Converted: Wed Dec 28 2016 at 23:12:17 using clawcode2html | |
This documentation file will not reflect any later changes in the source file. |
""" Set up the plot figures, axes, and items to be done for each frame. This module is imported by the plotting routines and then the function setplot is called to set the plot parameters. """ from __future__ import absolute_import def setplot(plotdata=None, bathy_location=0.15, bathy_angle=0.0, bathy_left=-1.0, bathy_right=-0.2): """Setup the plotting data objects. Input: plotdata, an instance of pyclaw.plotters.data.ClawPlotData. Output: a modified version of plotdata. returns plotdata object """ import os import numpy as np import matplotlib.pyplot as plt from clawpack.visclaw import geoplot, gaugetools import clawpack.clawutil.data as clawutil import clawpack.amrclaw.data as amrclaw import clawpack.geoclaw.data import clawpack.geoclaw.multilayer.plot as ml_plot if plotdata is None: from clawpack.visclaw.data import ClawPlotData plotdata = ClawPlotData() # Load data from output clawdata = clawutil.ClawInputData(2) clawdata.read(os.path.join(plotdata.outdir,'claw.data')) amrdata = amrclaw.AmrclawInputData(clawdata) amrdata.read(os.path.join(plotdata.outdir,'amr.data')) geodata = clawpack.geoclaw.data.GeoClawData() geodata.read(os.path.join(plotdata.outdir,'geoclaw.data')) multilayer_data = clawpack.geoclaw.data.MultilayerData() multilayer_data.read(os.path.join(plotdata.outdir,'multilayer.data')) def transform_c2p(x,y,x0,y0,theta): return ((x+x0)*np.cos(theta) - (y+y0)*np.sin(theta), (x+x0)*np.sin(theta) + (y+y0)*np.cos(theta)) def transform_p2c(x,y,x0,y0,theta): return ( x*np.cos(theta) + y*np.sin(theta) - x0, -x*np.sin(theta) + y*np.cos(theta) - y0) # Setup bathymetry reference lines with open(os.path.join(plotdata.outdir,"bathy_geometry.data"), 'r') as bathy_geometry_file: bathy_location = float(bathy_geometry_file.readline()) bathy_angle = float(bathy_geometry_file.readline()) x = [0.0,0.0] y = [0.0,1.0] x1,y1 = transform_c2p(x[0],y[0],bathy_location, 0.0,bathy_angle) x2,y2 = transform_c2p(x[1],y[1],bathy_location, 0.0,bathy_angle) if abs(x1 - x2) < 10**-3: x = [x1, x1] y = [clawdata.lower[1], clawdata.upper[1]] else: m = (y1 - y2) / (x1 - x2) x[0] = (clawdata.lower[1] - y1) / m + x1 y[0] = clawdata.lower[1] x[1] = (clawdata.upper[1] - y1) / m + x1 y[1] = clawdata.upper[1] ref_lines = [((x[0],y[0]),(x[1],y[1]))] plotdata.clearfigures() plotdata.save_frames = False # ======================================================================== # Generic helper functions # ======================================================================== def pcolor_afteraxes(current_data): bathy_ref_lines(current_data) gauge_locations(current_data) def contour_afteraxes(current_data): # gauge_locations(current_data) # m_to_km_labels() plt.hold(True) pos = -80.0 * (23e3 / 180) + 500e3 - 5e3 plt.plot([pos,pos],[-300e3,300e3],'b',[pos-5e3,pos-5e3],[-300e3,300e3],'y') plt.hold(False) wind_contours(current_data) bathy_ref_lines(current_data) def profile_afteraxes(current_data): pass def bathy_ref_lines(current_data): plt.hold(True) for ref_line in ref_lines: x1 = ref_line[0][0] y1 = ref_line[0][1] x2 = ref_line[1][0] y2 = ref_line[1][1] plt.plot([x1,x2],[y1,y2],'y--',linewidth=1) plt.hold(False) def gauge_locations(current_data,gaugenos='all'): plt.hold(True) gaugetools.plot_gauge_locations(current_data.plotdata, \ gaugenos=gaugenos, format_string='kx', add_labels=True) plt.hold(False) # ======================================================================== # Axis limits #xlimits = [amrdata.xlower,amrdata.xupper] xlimits = [-0.5,0.5] #ylimits = [amrdata.ylower,amrdata.yupper] ylimits = [-0.5,0.5] eta = [multilayer_data.eta[0],multilayer_data.eta[1]] top_surface_limits = [eta[0]-0.03,eta[0]+0.03] internal_surface_limits = [eta[1]-0.015,eta[1]+0.015] # top_surface_limits = [eta[0]-0.3,eta[0]+0.3] # internal_surface_limits = [eta[1]-0.15,eta[1]+0.15] top_speed_limits = [0.0,0.1] internal_speed_limits = [0.0,0.03] # Single layer test limits # top_surface_limits = [eta[0]-2.5,eta[0]+2.5] # top_speed_limits = [0.0,6.0] # ======================================================================== # Surface Elevations # ======================================================================== plotfigure = plotdata.new_plotfigure(name='Surface', figno=0) plotfigure.show = True plotfigure.kwargs = {'figsize':(14,4)} # Top surface plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Top Surface' plotaxes.axescmd = 'subplot(1,2,1)' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = pcolor_afteraxes ml_plot.add_surface_elevation(plotaxes,1,bounds=top_surface_limits) # ml_plot.add_surface_elevation(plotaxes,1,bounds=[-0.06,0.06]) # ml_plot.add_surface_elevation(plotaxes,1) ml_plot.add_land(plotaxes) # Bottom surface plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Internal Surface' plotaxes.axescmd = 'subplot(1,2,2)' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = pcolor_afteraxes # ml_plot.add_surface_elevation(plotaxes,2,bounds=[-300-0.5,-300+0.5]) ml_plot.add_surface_elevation(plotaxes,2,bounds=internal_surface_limits) # ml_plot.add_surface_elevation(plotaxes,2) ml_plot.add_land(plotaxes) # ======================================================================== # Depths # ======================================================================== plotfigure = plotdata.new_plotfigure(name='Depths', figno=42) plotfigure.show = False plotfigure.kwargs = {'figsize':(14,4)} # Top surface plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Top Layer Depth' plotaxes.axescmd = 'subplot(1,2,1)' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = pcolor_afteraxes ml_plot.add_layer_depth(plotaxes,1,bounds=[-0.1,1.1]) ml_plot.add_land(plotaxes) # Bottom surface plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Bottom Layer Depth' plotaxes.axescmd = 'subplot(1,2,2)' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = pcolor_afteraxes ml_plot.add_layer_depth(plotaxes,2,bounds=[-0.1,0.7]) ml_plot.add_land(plotaxes) # ======================================================================== # Water Speed # ======================================================================== plotfigure = plotdata.new_plotfigure(name='speed', figno=1) plotfigure.show = True plotfigure.kwargs = {'figsize':(14,4)} # Top layer speed plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Currents - Top Layer' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.axescmd = 'subplot(1,2,1)' plotaxes.afteraxes = pcolor_afteraxes # add_speed(plotaxes,1,bounds=[0.00,0.2]) ml_plot.add_speed(plotaxes,1,bounds=top_speed_limits) # add_speed(plotaxes,1) ml_plot.add_land(plotaxes) # Bottom layer speed plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Currents - Bottom Layer' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.axescmd = 'subplot(1,2,2)' plotaxes.afteraxes = pcolor_afteraxes # add_speed(plotaxes,2,bounds=[0.0,1e-10]) ml_plot.add_speed(plotaxes,2,bounds=internal_speed_limits) # add_speed(plotaxes,2) ml_plot.add_land(plotaxes) # Individual components plotfigure = plotdata.new_plotfigure(name='speed_components',figno=401) plotfigure.show = False plotfigure.kwargs = {'figsize':(14,14)} # Top layer plotaxes = plotfigure.new_plotaxes() plotaxes.title = "X-Velocity - Top Layer" plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.axescmd = 'subplot(2,2,1)' plotaxes.afteraxes = pcolor_afteraxes # add_x_velocity(plotaxes,1,bounds=[-1e-10,1e-10]) ml_plot.add_x_velocity(plotaxes,1) ml_plot.add_land(plotaxes) plotaxes = plotfigure.new_plotaxes() plotaxes.title = "Y-Velocity - Top Layer" plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.axescmd = 'subplot(2,2,2)' plotaxes.afteraxes = pcolor_afteraxes # add_y_velocity(plotaxes,1,bounds=[-0.000125,0.000125]) ml_plot.add_y_velocity(plotaxes,1) ml_plot.add_land(plotaxes) # Bottom layer plotaxes = plotfigure.new_plotaxes() plotaxes.title = "X-Velocity - Bottom Layer" plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.axescmd = 'subplot(2,2,3)' plotaxes.afteraxes = pcolor_afteraxes # add_x_velocity(plotaxes,2,bounds=[-1e-10,1e-10]) ml_plot.add_x_velocity(plotaxes,2) ml_plot.add_land(plotaxes) plotaxes = plotfigure.new_plotaxes() plotaxes.title = "Y-Velocity - Bottom Layer" plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.axescmd = 'subplot(2,2,4)' plotaxes.afteraxes = pcolor_afteraxes # add_y_velocity(plotaxes,2,bounds=[-0.8e-6,.8e-6]) ml_plot.add_y_velocity(plotaxes,2) ml_plot.add_land(plotaxes) # ======================================================================== # Profile Plots # ======================================================================== plotfigure = plotdata.new_plotfigure(name='profile', figno=4) plotfigure.show = False # Top surface plotaxes = plotfigure.new_plotaxes() plotaxes.xlimits = xlimits plotaxes.ylimits = [-2000,20] plotaxes.title = "Profile of depth" plotaxes.afteraxes = profile_afteraxes slice_index = 30 # Internal surface def bathy_profile(current_data): return current_data.x[:,slice_index], b(current_data)[:,slice_index] def lower_surface(current_data): if multilayer_data.init_type == 2: return current_data.x[:,slice_index], eta2(current_data)[:,slice_index] elif multilayer_data.init_type == 6: return current_data.y[slice_index,:], eta2(current_data)[slice_index,:] def upper_surface(current_data): if multilayer_data.init_type == 2: return current_data.x[:,slice_index], eta1(current_data)[:,slice_index] elif multilayer_data.init_type == 6: return current_data.y[slice_index,:], eta1(current_data)[slice_index,:] def top_speed(current_data): if multilayer_data.init_type == 2: return current_data.x[:,slice_index], water_u1(current_data)[:,slice_index] elif multilayer_data.init_type == 6: return current_data.y[slice_index,:], water_u1(current_data)[slice_index,:] def bottom_speed(current_data): if multilayer_data.init_type == 2: return current_data.x[:,slice_index], water_u2(current_data)[:,slice_index] elif multilayer_data.init_type == 6: return current_data.y[slice_index,:], water_u2(current_data)[slice_index,:] # Bathy plotitem = plotaxes.new_plotitem(plot_type='1d_from_2d_data') plotitem.map_2d_to_1d = bathy_profile plotitem.plot_var = 0 plotitem.amr_plotstyle = ['-','+','x'] plotitem.color = 'k' plotitem.show = True # Internal Interface plotitem = plotaxes.new_plotitem(plot_type='1d_from_2d_data') plotitem.map_2d_to_1d = lower_surface plotitem.plot_var = 7 plotitem.amr_plotstyle = ['-','+','x'] plotitem.color = 'b' plotitem.show = True # Upper Interface plotitem = plotaxes.new_plotitem(plot_type='1d_from_2d_data') plotitem.map_2d_to_1d = upper_surface plotitem.plot_var = 6 plotitem.amr_plotstyle = ['-','+','x'] plotitem.color = (0.2,0.8,1.0) plotitem.show = True # ======================================================================== # Combined Profile Plot # ======================================================================== # add_combined_profile_plot(plotdata,0.25,direction='y',figno=120) # add_combined_profile_plot(plotdata,0.8,direction='y',figno=121) # # add_velocities_profile_plot(plotdata,0.25,direction='y',figno=130) # add_velocities_profile_plot(plotdata,0.8,direction='y',figno=131) # ======================================================================== # Bathy Profile # ======================================================================== plotfigure = plotdata.new_plotfigure(name='bathy_profile',figno=20) plotfigure.show = False plotaxes = plotfigure.new_plotaxes() plotaxes.xlimits = xlimits plotaxes.title = "Bathymetry Profile" plotaxes.scaled = 'equal' plotitem = plotaxes.new_plotitem(plot_type='2d_imshow') plotitem.plot_var = ml_plot.b plotitem.imshow_cmin = -4000 plotitem.imshow_cmax = 10 plotitem.add_colorbar = True plotitem.amr_celledges_show = [0,0,0] plotitem.amr_patchedges_show = [1,1,1] plotitem.show = True # ======================================================================== # Figures for momentum # ======================================================================== plotfigure = plotdata.new_plotfigure(name='x-momentum', figno=13) plotfigure.show = False # Set up for axes in this figure: plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'X-Velocity' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = pcolor_afteraxes # Water # plotitem = plotaxes.new_plotitem(plot_type='2d_pcolor') # # plotitem.plot_var = geoplot.surface # plotitem.plot_var = water_u # plotitem.pcolor_cmap = colormaps.make_colormap({1.0:'r',0.5:'w',0.0:'b'}) # # plotitem.pcolor_cmin = -1.e-10 # # plotitem.pcolor_cmax = 1.e-10 # # plotitem.pcolor_cmin = -2.5 # -3.0 # # plotitem.pcolor_cmax = 2.5 # 3.0 # plotitem.add_colorbar = True # plotitem.amr_celledges_show = [0,0,0] # plotitem.amr_patchedges_show = [1,1,1] # Land plotitem = plotaxes.new_plotitem(plot_type='2d_pcolor') plotitem.show = True plotitem.plot_var = geoplot.land plotitem.pcolor_cmap = geoplot.land_colors plotitem.pcolor_cmin = 0.0 plotitem.pcolor_cmax = 80.0 plotitem.add_colorbar = False plotitem.amr_celledges_show = [0,0,0] plotitem.amr_patchedges_show = [1,1,1] plotfigure = plotdata.new_plotfigure(name='y-momentum', figno=14) plotfigure.show = False # Set up for axes in this figure: plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Y-Velocity' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = pcolor_afteraxes # Water # plotitem = plotaxes.new_plotitem(plot_type='2d_pcolor') # # plotitem.plot_var = geoplot.surface # plotitem.plot_var = water_v # plotitem.pcolor_cmap = colormaps.make_colormap({1.0:'r',0.5:'w',0.0:'b'}) # # plotitem.pcolor_cmin = -1.e-10 # # plotitem.pcolor_cmax = 1.e-10 # # plotitem.pcolor_cmin = -2.5 # -3.0 # # plotitem.pcolor_cmax = 2.5 # 3.0 # plotitem.add_colorbar = True # plotitem.amr_celledges_show = [0,0,0] # plotitem.amr_patchedges_show = [1,1,1] # Land ml_plot.add_land(plotaxes) # ======================================================================== # Contour plot for surface # ======================================================================== plotfigure = plotdata.new_plotfigure(name='contour_surface',figno=15) plotfigure.show = False plotfigure.kwargs = {'figsize':(14,4)} # Set up for axes in this figure: # Top Surface plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Top Surface' plotaxes.axescmd = 'subplot(1,2,1)' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = contour_afteraxes ml_plot.add_surface_elevation(plotaxes,plot_type='contour',surface=1,bounds=[-2.5,-1.5,-0.5,0.5,1.5,2.5]) ml_plot.add_land(plotaxes,plot_type='contour') # Internal Surface plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Internal Surface' plotaxes.axescmd = 'subplot(1,2,2)' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = contour_afteraxes ml_plot.add_surface_elevation(plotaxes,plot_type='contour',surface=2,bounds=[-2.5,-1.5,-0.5,0.5,1.5,2.5]) ml_plot.add_land(plotaxes,plot_type='contour') # ======================================================================== # Contour plot for speed # ======================================================================== plotfigure = plotdata.new_plotfigure(name='contour_speed',figno=16) plotfigure.show = False plotfigure.kwargs = {'figsize':(14,4)} # Set up for axes in this figure: plotaxes = plotfigure.new_plotaxes() plotaxes.title = 'Current' plotaxes.scaled = True plotaxes.xlimits = xlimits plotaxes.ylimits = ylimits plotaxes.afteraxes = contour_afteraxes # Surface plotitem = plotaxes.new_plotitem(plot_type='2d_contour') plotitem.plot_var = ml_plot.water_speed_depth_ave plotitem.kwargs = {'linewidths':1} # plotitem.contour_levels = [1.0,2.0,3.0,4.0,5.0,6.0] plotitem.contour_levels = [0.5,1.5,3,4.5,6.0] plotitem.amr_contour_show = [1,1,1] plotitem.amr_celledges_show = [0,0,0] plotitem.amr_patchedges_show = [1,1,1] plotitem.amr_contour_colors = 'k' # plotitem.amr_contour_colors = ['r','k','b'] # color on each level # plotitem.amr_grid_bgcolor = ['#ffeeee', '#eeeeff', '#eeffee'] plotitem.show = True # Land plotitem = plotaxes.new_plotitem(plot_type='2d_contour') plotitem.plot_var = geoplot.land plotitem.contour_nlevels = 40 plotitem.contour_min = 0.0 plotitem.contour_max = 100.0 plotitem.amr_contour_colors = ['g'] # color on each level plotitem.amr_patch_bgcolor = ['#ffeeee', '#eeeeff', '#eeffee'] plotitem.amr_celledges_show = 0 plotitem.amr_patchedges_show = 0 plotitem.show = True # ======================================================================== # Vorticity Plot # ======================================================================== # plotfigure = plotdata.new_plotfigure(name='vorticity',figno=17) # plotfigure.show = False # plotaxes = plotfigure.new_plotaxes() # plotaxes.title = "Vorticity" # plotaxes.scaled = True # plotaxes.xlimits = xlimits # plotaxes.ylimits = ylimits # plotaxes.afteraxes = pcolor_afteraxes # # # Vorticity # plotitem = plotaxes.new_plotitem(plot_type='2d_imshow') # plotitem.plot_var = 9 # plotitem.imshow_cmap = plt.get_cmap('PRGn') # # plotitem.pcolor_cmap = plt.get_cmap('PuBu') # # plotitem.pcolor_cmin = 0.0 # # plotitem.pcolor_cmax = 6.0 # plotitem.imshow_cmin = -1.e-2 # plotitem.imshow_cmax = 1.e-2 # plotitem.add_colorbar = True # plotitem.amr_celledges_show = [0,0,0] # plotitem.amr_patchedges_show = [1] # # # Land # plotitem = plotaxes.new_plotitem(plot_type='2d_pcolor') # plotitem.plot_var = geoplot.land # plotitem.pcolor_cmap = geoplot.land_colors # plotitem.pcolor_cmin = 0.0 # plotitem.pcolor_cmax = 80.0 # plotitem.add_colorbar = False # plotitem.amr_celledges_show = [0,0,0] # ======================================================================== # Figures for gauges # ======================================================================== # Top plotfigure = plotdata.new_plotfigure(name='Surface & topo', figno=300, \ type='each_gauge') plotfigure.show = True plotfigure.clf_each_gauge = True # Set up for axes in this figure: plotaxes = plotfigure.new_plotaxes() plotaxes.xlimits = [0.0, 1.0] plotaxes.ylimits = top_surface_limits plotaxes.title = 'Top Surface' # Plot surface as blue curve: plotitem = plotaxes.new_plotitem(plot_type='1d_plot') plotitem.plot_var = 6 plotitem.plotstyle = 'b-' # Bottom plotfigure = plotdata.new_plotfigure(name='Bottom Surface Gauge', type='each_gauge') plotfigure.show = True plotfigure.clf_each_gauge = True # Set up for axes in this figure: plotaxes = plotfigure.new_plotaxes() plotaxes.xlimits = [0.0,1.0] # plotaxes.ylimits = internal_surface_limits plotaxes.title = 'Bottom Surface' # Plot surface as blue curve: plotitem = plotaxes.new_plotitem(plot_type='1d_plot') plotitem.plot_var = 7 plotitem.plotstyle = 'b-' # Plot topo as green curve: # plotitem = plotaxes.new_plotitem(plot_type='1d_plot') # plotitem.plot_var = geoplot.gaugetopo # plotitem.plotstyle = 'g+' #----------------------------------------- # Parameters used only when creating html and/or latex hardcopy # e.g., via pyclaw.plotters.frametools.printframes: plotdata.printfigs = True # print figures plotdata.print_format = 'png' # file format plotdata.print_framenos = 'all' # list of frames to print plotdata.print_fignos = 'all' # list of figures to print plotdata.html = True # create html files of plots? plotdata.latex = False # create latex file of plots? plotdata.latex_figsperline = 2 # layout of plots plotdata.latex_framesperline = 1 # layout of plots plotdata.latex_makepdf = False # also run pdflatex? plotdata.parallel = True # make multiple frame png's at once return plotdata