Index: trunk/LMDZ.COMMON/libf/evolution/deftank/README =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/deftank/README (revision 3458) +++ trunk/LMDZ.COMMON/libf/evolution/deftank/README (revision 3460) @@ -49,6 +49,12 @@ Bash script file to clean the folder after a PEM simulation. -# output_layering.py: +# multiple_exec.sh: + Bash script to execute multiple scripts in subdirectories, useful to launch multiple simulations at once. + +# visu_layering.py: Python script file to output the stratification data from the "startpem.nc" files. + +# visu_evol_layering.py: + Python script file to output the stratification data over time from the "startpem.nc" files. Note: Index: trunk/LMDZ.COMMON/libf/evolution/deftank/visu_evol_layering.py =================================================================== --- trunk/LMDZ.COMMON/libf/evolution/deftank/visu_evol_layering.py (revision 3458) +++ trunk/LMDZ.COMMON/libf/evolution/deftank/visu_evol_layering.py (revision 3460) @@ -13,5 +13,5 @@ folder_path = "starts" base_name = "restartpem" -dz = 1 # Discrization step of the reference grid for the elevation [m] +dz = 0.1 # Discrization step of the reference grid for the elevation [m] infofile = 'info_PEM.txt' ######################### @@ -39,4 +39,6 @@ ngrid = Dataset(os.path.join(folder_path,base_name + "1.nc"),'r').dimensions['physical_points'].size # ngrid is the same for all files nslope = Dataset(os.path.join(folder_path,base_name + "1.nc"),'r').dimensions['nslope'].size # nslope is the same for all files + longitude = Dataset(os.path.join(folder_path,base_name + "1.nc"),'r').variables['longitude'][:] + latitude = Dataset(os.path.join(folder_path,base_name + "1.nc"),'r').variables['latitude'][:] for i in range(1,nfile + 1): file_path = os.path.join(folder_path,base_name + str(i) + ".nc") @@ -53,17 +55,17 @@ max_top_elevation = max(max_top_elevation,np.max(ds.variables[slope_var_name][:])) - print(f" > number of files = {nfile}") - print(f" > ngrid = {ngrid}") - print(f" > nslope = {nslope}") - print(f" > max(nb_str_max) = {max_nb_str}") - print(f" > max(top_elevation) = {max_top_elevation}") + print(f"> number of files = {nfile}") + print(f"> ngrid = {ngrid}") + print(f"> nslope = {nslope}") + print(f"> max(nb_str_max) = {max_nb_str}") + print(f"> max(top_elevation) = {max_top_elevation}") # Concatenate stratif variables with dimension 'nb_str_max' along the "Time" dimension stratif_data = [] stratif_heights = np.zeros((ngrid,nfile,nslope,max_nb_str)) - stratif_co2ice = -1.*np.ones((ngrid,nfile,nslope,max_nb_str)) - stratif_h2oice = -1.*np.ones((ngrid,nfile,nslope,max_nb_str)) - stratif_dust = -1.*np.ones((ngrid,nfile,nslope,max_nb_str)) - stratif_air = -1.*np.ones((ngrid,nfile,nslope,max_nb_str)) + stratif_co2ice = np.zeros((ngrid,nfile,nslope,max_nb_str)) + stratif_h2oice = np.zeros((ngrid,nfile,nslope,max_nb_str)) + stratif_dust = np.zeros((ngrid,nfile,nslope,max_nb_str)) + stratif_air = np.zeros((ngrid,nfile,nslope,max_nb_str)) for var_name in datasets[0].variables: if 'top_elevation' in var_name: @@ -108,5 +110,5 @@ stratif_data = [stratif_heights,stratif_co2ice,stratif_h2oice,stratif_dust,stratif_air] - return stratif_data, max_top_elevation + return stratif_data, max_top_elevation, longitude, latitude ### Function to interpolate the stratification data on a reference grid @@ -114,22 +116,17 @@ # Define the reference ref_grid ref_grid = np.arange(0,max_top_elevation,dz) - print(f" > Number of ref_grid points = {len(ref_grid)}") - - # Extrapolate the top_elevation of strata - for i in range(np.shape(stratif_data)[1]): - for j in range(np.shape(stratif_data)[2]): - for k in range(np.shape(stratif_data)[3]): - if len(np.where(stratif_data[0][i,j,k,:] == 0)[0]) > 1: # If there are trailing zeros in top_elevation - i0 = np.where(stratif_data[0][i,j,k,:] == 0)[0][1] - stratif_data[0][i,j,k,i0:] = np.linspace(stratif_data[0][i,j,k,i0 - 1],max_top_elevation,np.shape(stratif_data)[4] - i0) + print(f"> Number of ref_grid points = {len(ref_grid)}") # Interpolate the strata properties on the ref_grid - gridded_stratif_data = np.zeros((np.shape(stratif_data)[0] - 1,np.shape(stratif_data)[1],np.shape(stratif_data)[2],np.shape(stratif_data)[3],len(ref_grid))) + gridded_stratif_data = -1.*np.ones((np.shape(stratif_data)[0] - 1,np.shape(stratif_data)[1],np.shape(stratif_data)[2],np.shape(stratif_data)[3],len(ref_grid))) for iprop in range(1,np.shape(stratif_data)[0]): for i in range(np.shape(stratif_data)[1]): for j in range(np.shape(stratif_data)[2]): for k in range(np.shape(stratif_data)[3]): - f = interpolate.interp1d(stratif_data[0][i,j,k,:],stratif_data[iprop][i,j,k,:],kind = 'next',fill_value = "extrapolate") - gridded_stratif_data[iprop - 1,i,j,k,:] = f(ref_grid) + i_hmax = np.max(np.nonzero(stratif_data[0][i,j,k,:])) + hmax = stratif_data[0][i,j,k,i_hmax] + i_zmax = np.searchsorted(ref_grid,hmax,'left') + f = interpolate.interp1d(stratif_data[0][i,j,k,:i_hmax + 1],stratif_data[iprop][i,j,k,:i_hmax + 1],kind = 'next')#,fill_value = "extrapolate") + gridded_stratif_data[iprop - 1,i,j,k,:i_zmax] = f(ref_grid[:i_zmax]) return ref_grid, gridded_stratif_data @@ -153,5 +150,5 @@ ### Processing -stratif_data, max_top_elevation = process_files(folder_path,base_name) +stratif_data, max_top_elevation, longitude, latitude = process_files(folder_path,base_name) ref_grid, gridded_stratif_data = interpolate_data(stratif_data,max_top_elevation,dz) date_time = read_infofile(infofile)