Postprocessing tools (#216)

* Re-ordered the existing AeoLiS simulation examples

The examples by deVries2023 are merged with the other existing examples and a new README file is added, covering all available examples.

* Updated plotting functions with improved documentation

* Add files via upload

Author: bartvanwesten

aeolis/examples/README.txt

@@ -0,0 +1,77 @@
+This folder contains a variation of AeoLiS simulation setups.
+
+
+1D - (simple) One-dimensional examples:
+
+1D/case1_small_waves
+<description>
+
+1D/case2_larger_waves
+<description>
+
+1D/case3_erosion_avalanching
+<description>
+
+1D/case4_surface_moisture
+<description>
+
+
+2D - (simple) Two-dimensional examples
+
+2D/Barchan_dune
+<description>
+
+2D/Parabolic_dune
+<description>
+
+2D/Carrara_Coast
+<description>
+
+
+deVries2023 - Simulations describing ...
+
+deVries2023/Run 1 - Fetch effects and spatial variability in supply
+<description>
+
+deVries2023/Run 2 - Fetch effects and wind directionality
+<description>
+
+deVries2023/Run 3 - Fetch effects and temporal variability in supply
+<description>
+
+deVries2023/Run 4 - Predictions of decadal development
+<description>
+
+
+grainsizevariations - SHORT DESCRIPTION
+<description>
+
+
+sandengine_small_grids - Fast Sand Engine simulation
+Simulation of the subaerial Sand Engine evolution, describing supply-limiting effects like moisture and sediment sorting on a larger scale
+
+
+vanWesten2024 - Landform simulations as described in the AeoLiS landform publication (van Westen, 2024)
+
+vanWesten2024/barchan_morocco
+<description>
+
+vanWesten2024/barchan_rotate
+<description>
+
+vanWesten2024/blowout
+<description>
+
+vanWesten2024/embryo
+<description>
+
+vanWesten2024/parabolic
+<description>
+
+
+
+
+
+
+
+

tools/aeolis_output.nc

@@ -17,7 +17,7 @@
 
 def plot1d(ncfile, itransects=0, itimes=0, ifrac=0, params='zb'):
 
-    '''Generic function for plotting 2D parameters from AeoLiS netcdf results
+    '''Generic function for plotting 1D parameters from AeoLiS netcdf results
 
     Parameters
     ----------
@@ -30,17 +30,18 @@ def plot1d(ncfile, itransects=0, itimes=0, ifrac=0, params='zb'):
     ifrac : integer
         index of fraction, only used when the results has multiple dimensions ( >=4 )
     params :
-        string or list of strings describing the name of the parameter for plotting
+        string or list of strings describing the name of the parameter for plotting'
+        the parameters that can be plotted are defined in the AeoLiS input file.
 
     Returns
     -------
-    figs/ axs?
+    figs/ ax
         ...
 
     Examples
     --------
-    >>> ploted(parse_value('T'))
-        bool
+    >>> plot1d(ncfile, itransects=0, itimes=0, ifrac=0, params='zb')            # plot the bed level for the first transect and timestep 0
+    >>> plot1d(ncfile, itransects=0, itimes=-1, ifrac=0, params=['zb', 'Hs'])   # plot the bed level and wave height for the first transect and last timestep
 
     '''
 
@@ -103,17 +104,16 @@ def plot2d(ncfile, itimes=0, ifrac=0, param='zb', cmap=plt.cm.jet, clim='auto',
         Boolean for indicating plotting difference between timesteps or the actual value
     itimeDelta : integer
         index for the starting timestep in the netcdf-file when merging, only used when delta == Bool
-    quiver?
 
     Returns
     -------
-    figs/ axs?
+    figs/ ax
         ...
 
     Examples
     --------
-    >>> ploted(parse_value('T'))
-        bool
+    >>> plot2d(ncfile, itimes=0, ifrac=0, param='zb')            # plot the bed level for True first timestep
+    >>> plot2d(ncfile, itimes=-1, ifrac=0, param='zb', cmap=plt.cm.jet, clim='auto', delta=False, itimeDelta=0)   # plot the bed level difference for the first timestep
 
     '''
 
@@ -153,100 +153,3 @@ def plot2d(ncfile, itimes=0, ifrac=0, param='zb', cmap=plt.cm.jet, clim='auto',
             plt.colorbar(pc)
 
     return fig, ax
-
-
-def plot3d(ncfile, itimes, scalez=1., vangle=90., hangle=0.):
-    
-    '''Generic function for plotting 2D parameters from AeoLiS netcdf results
-
-    Parameters
-    ----------
-    ncfile : str
-        filename of the netcdf-file containing the results
-    itimes : integer or array of integers
-        index or list of indices describing the timestep-indices in the netcdf-file for plotting
-
-    Returns
-    -------
-    figs/ axs?
-        ...
-
-    Examples
-    --------
-    >>> ploted(parse_value('T'))
-        bool
-
-    '''
-    
-    rgb=255.
-    sand = [230/rgb,230/rgb,210/rgb]
-    blue = [0/rgb, 166/rgb, 214/rgb]
-    sand = [230/rgb,230/rgb,210/rgb]
-    vegetation = [120/rgb,180/rgb,50/rgb]
-    white = [255/rgb,255/rgb,255/rgb]
-    gray = [100/rgb,100/rgb,100/rgb]
-    
-        
-    if type(itimes) == int:
-        itimes = [itimes]
-        
-    for it in itimes:
-    
-        with netCDF4.Dataset(ncfile, 'r') as ds:
-            
-            # get spatial dimensions and bed levels
-            x = ds.variables['x'][:,:]
-            y = ds.variables['y'][:,:]
-            z = ds.variables['zb'][it,:,:]
-            
-            fig = plt.figure()
-            ax = plt.gca(projection='3d')
-    
-            x_scale = 1.
-            y_scale = 1. # np.max(x)/np.max(y)
-            z_scale = scalez*np.max(z)/np.max(x)
-            
-            scale=np.diag([x_scale, y_scale, z_scale, 1.0])
-            scale=scale*(1.0/scale.max())
-            scale[3,3]=1.0
-            
-            def short_proj():
-                return np.dot(Axes3D.get_proj(ax), scale)
-            
-            ax.get_proj = short_proj   
-            
-            # lay-out of axes
-            ax.set_ylabel('alongshore distance [m]')
-            ax.set_xlabel('cross-shore distance [m]')
-            ax.axes.get_xaxis().set_ticks([])
-            ax.axes.get_yaxis().set_ticks([])
-            ax.axes.get_xaxis().set_visible(False)
-            ax.axes.get_yaxis().set_visible(False)
-            
-            # perspective
-            ax.view_init(vangle, hangle)
-            ax.set_proj_type(proj_type='persp')
-            
-            # make the panes transparent           
-            for axis in [ax.xaxis, ax.yaxis, ax.zaxis]:
-                axis.set_pane_color((0/rgb,166/rgb,214/rgb,0.0))
-                axis._axinfo["grid"]['color'] =  (1,1,1,0) 
-                axis.line.set_color((1.0, 1.0, 1.0, 0.0)) 
-            
-            ax.set_axis_off()
-            pos_z = z.copy()
-            neg_z = z.copy()
-            # neg_z[(neg_z > 0.2)] = np.nan
-            pos_z[(pos_z < 0.2)] = np.nan
-            
-            # ax.plot_surface(y, x, neg_z[:,:], color=white, linewidth=0, antialiased=False, shade=False, alpha=1.0, rstride=1, cstride=1)
-            ax.plot_surface(y, x, pos_z[:,:], color=sand, linewidth=0, antialiased=False, shade=True, alpha=1.0, rstride=1, cstride=1)          
-            
-            # plt.savefig(r'c:\Users\weste_bt\AeoLiS\Examples\Parabolic figures\plottoptime' + str('{:03}'.format(time_index)) + '.png', dpi=200)
-            # plt.show()
-            # plt.close(fig)
-            
-    return fig, ax
-
-
-    

tools/visualization/plotting.py renamed to tools/postprocessing/plotting_functions.py

@@ -0,0 +1,69 @@
+from plotting_functions import plot1d, plot2d
+import matplotlib.pyplot as plt
+import os
+
+# This script is used to plot the results of the AeoLiS model
+# The functions plot1d and plot2d are defined in the plotting_functions.py file
+# Using this script as a template, you can easily start plotting some basic results generated by the AeoLiS model
+
+# Get the path to the netcdf-file containing the results of the AeoLiS model run
+# This is relative to the location of the current script
+ncfile = os.path.join(os.path.dirname(os.path.realpath(__file__)), 'aeolis_output.nc')	
+
+
+############################################################################################################
+# The function plot2d is used to plot 2D results of the AeoLiS model
+
+# - ncfile: path to the netcdf-file containing the results of the AeoLiS model run
+# - itimes: index of the time step to plot (default: 0)
+# - ifrac: grainfraction to plot (default: 0)
+# - param: parameter to plot (default: 'zb')
+# - cmap: colormap to use (default: plt.cm.jet)
+# - clim: color limits [min, max], or 'auto' will automatically determine the color limits (default: 'auto')
+# - delta: plot the difference between the current and previous time step (default: False)
+# - itimeDelta: time step difference (default: -1), only used if delta=True
+
+# Plot the bed level at the first time step
+fig, ax = plot2d(ncfile)
+ax.set_title('Bed level at the first time step')
+
+# Plot the bed level at the last time step
+fig, ax = plot2d(ncfile, param='zb', itimes=-1)
+ax.set_title('Bed level at the last time step')
+
+# Plot the bed level at the first time step with a different colormap and color limits
+fig, ax = plot2d(ncfile, param='zb', itimes=0, delta=True, itimeDelta=-1, cmap=plt.cm.RdBu, clim=[-1.5, 1.5])
+ax.set_title('Bed level change between the first and last time step with custom colormap and color limits')
+
+# Plot the wind velocity at the last timestep using the viridis colormap
+fig, ax = plot2d(ncfile, param='uw', itimes=-1, cmap=plt.cm.viridis)
+ax.set_title('Wind velocity at the last time step')
+
+# Plot the velocity threshold at the last timestep only for the 1st fraction
+fig, ax = plot2d(ncfile, param='uth', itimes=-1, ifrac=0, clim=[0, 0.4])
+ax.set_title('Velocity threshold at the last time step for the 1st fraction')
+
+
+############################################################################################################
+# The function plot1d is used to plot 1D results of the AeoLiS model of a single transect
+# Timesteps, transects and parameters can be given as a list to plot multiple results at once
+
+# - ncfile: path to the netcdf-file containing the results of the AeoLiS model run
+# - itransects: index of the transect to plot (default: 0)
+# - itimes: index of the time step to plot (default: 0)
+# - ifrac: grainfraction to plot (default: 0)
+# - params: list of parameters to plot (default: ['zb'])
+
+# Plot the bed level and water level at the first transect and last timestep
+fig, ax = plot1d(ncfile, itransects=0, itimes=-1, ifrac=0, params=['zb', 'zs'])
+ax.set_title('Bed level and water level at the first transect and last timestep')
+
+# Plot the bed level at the first, middle and last timestep for the first transect
+fig, ax = plot1d(ncfile, itransects=0, itimes=[0, 4, -1], ifrac=0, params='zb')
+ax.set_title('Bed level at the first, middle and last timestep for the first transect')
+
+# Plot the bed level for the first and last transect at the last timestep
+fig, ax = plot1d(ncfile, itransects=[0, -1], itimes=-1, ifrac=0, params='zb')
+ax.set_title('Bed level for the first and last transect at the last timestep')
+
+plt.show()