diff --git a/Analysis_solid_gauging.py b/Analysis_solid_gauging.py
index a4de37bb83524f0f0d9d09f8878d6f0ec9590286..87924635ce165f4144a2e8626f2b00e6d7a032aa 100644
--- a/Analysis_solid_gauging.py
+++ b/Analysis_solid_gauging.py
@@ -40,6 +40,7 @@ path_R = 'C:/Users/jessica.laible/AppData/Local/Programs/R/R-4.2.1/bin/x64/Rscri
# vertical distance above the bottom (BD)
position_diver = 0.12 # (m)
HADCP_option = True # Choose if HADCP section or not
+add_plots = False # Choose if additional plots should be made
# Add paths to system
cwd = os.getcwd()
@@ -246,11 +247,11 @@ if choice[3] == 0:
for i in range(len(analysis_data)):
meas = analysis_data.iloc[i,:]
if meas['Sampler'] == 'BD':
- sand_flux_point_i = sediment_flux_BD(outpath_figures, meas, choice, unit)
+ sand_flux_point_i = sediment_flux_BD(outpath_figures, meas, choice, unit, add_plots)
sand_flux_point.append(sand_flux_point_i)
concentration_sand_point.append(np.nan)
elif meas['Sampler'] == 'BD_C':
- sand_flux_point_i = sediment_flux_BD_C(outpath_figures, meas, choice, unit)
+ sand_flux_point_i = sediment_flux_BD_C(outpath_figures, meas, choice, unit, add_plots)
sand_flux_point.append(sand_flux_point_i)
concentration_sand_point.append(np.nan)
concentration_fines.append(np.nan)
@@ -297,12 +298,12 @@ elif choice[3] == 1:
for i in range(len(analysis_data)):
meas = analysis_data.iloc[i,:]
if meas['Sampler'] == 'BD':
- sand_flux_point_i, c_sand_point_g_l_from_flux_i = sediment_flux_BD(outpath_figures, meas, choice, unit)
- sand_flux_point.append(sand_flux_point_i)
+ sand_flux_point_i, c_sand_point_g_l_from_flux_i = sediment_flux_BD(outpath_figures, meas, choice, unit, add_plots)
+ sand_flux_point.append(sand_flux_point_i)
concentration_sand_point.append(c_sand_point_g_l_from_flux_i)
concentration_fines.append(np.nan)
elif meas['Sampler'] == 'BD_C':
- sand_flux_point_i, c_sand_point_g_l_from_flux_i = sediment_flux_BD_C(outpath_figures, meas, choice, unit)
+ sand_flux_point_i, c_sand_point_g_l_from_flux_i = sediment_flux_BD_C(outpath_figures, meas, choice, uni, add_plotst)
sand_flux_point.append(sand_flux_point_i)
concentration_sand_point.append(c_sand_point_g_l_from_flux_i)
concentration_fines.append(np.nan)
@@ -396,11 +397,11 @@ if choice[2] == 1 and choice[3] == 0 and number_gsd >= 2: # if no ADCP-data avai
#%% Uncertainty estimation
analysis_data, U_F, uncertainty, u_p, u_param, maxpost, maxpost_alpha, maxpost_lnCr = uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
- choice, v2_d_index, summary_ISO, summary_fine, map_data, cwd, path_R)
+ choice, v2_d_index, summary_ISO, summary_fine, map_data, cwd, path_R, add_plots)
#%% SDC method
analysis_data, summary_SDC, Conc_SDC_cell_export, sand_flux_SDC_kg_s_export = SDC_method(analysis_data, map_data, outpath,
- outpath_figures, sampling_date, maxpost, maxpost_alpha, maxpost_lnCr)
+ outpath_figures, sampling_date, maxpost, maxpost_alpha, maxpost_lnCr, add_plots)
#%% Compare results from different methods
@@ -419,7 +420,7 @@ if choice[3] == 1:
ax=ax, width=0.6, rot= 0 )
ax.set_xlim(-0.5,4.5)
- ax.set_ylabel('$\mathregular{\overline{C}_{sand}}$ (g/l)', fontsize = 16, weight = 'bold')
+ ax.set_ylabel('$\mathregular{\overline{C_{sand}}}$ (g/l)', fontsize = 16, weight = 'bold')
ax.set_xlabel('Method', fontsize = 16, weight = 'bold')
ax.tick_params(axis='both', which='major', labelsize = 14)
fig.tight_layout()
diff --git a/classes/SDC_fines_method.py b/classes/SDC_fines_method.py
index 632b9a7d78e91e599cb68cf18d5fe75f1b828dac..ea5e12e89628b7b53d0e748d755e041e9f74cb6c 100644
--- a/classes/SDC_fines_method.py
+++ b/classes/SDC_fines_method.py
@@ -55,11 +55,11 @@ def SDC_fines_method(analysis_data, map_data, outpath, outpath_figures, sampling
most_frequent_sampler= analysis_data['Sampler'].value_counts().idxmax()
# Fontsizes
- fontsize_axis = 14
- fontsize_legend = 14
- fontsize_legend_title = 14
- fontsize_text = 14
- fontsize_ticks = 12
+ fontsize_axis = 16
+ fontsize_legend = 16
+ fontsize_legend_title = 16
+ fontsize_text = 16
+ fontsize_ticks = 14
#%% SDC_fines profile
stats_SDC_fines_profile = []
@@ -272,16 +272,14 @@ def SDC_fines_method(analysis_data, map_data, outpath, outpath_figures, sampling
# Add color bar and axis labels
cb = fig.colorbar(c, pad=0.02, shrink=0.7, anchor=(0, 0.1))
cb.ax.set_ylabel(canvas.tr('$\mathregular{C_{fine}}$ (g/l)'), weight = 'bold')
- cb.ax.yaxis.label.set_fontsize(12)
- cb.ax.tick_params(labelsize=12)
+ cb.ax.yaxis.label.set_fontsize(16)
+ cb.ax.tick_params(labelsize=14)
fig.ax.invert_yaxis()
fig.ax.fill_between(x_fill, np.nanmax(track_depth)+2, depth_fill, color='w') # below bathy
fig.ax.plot(track, track_depth, color='k', linewidth=1.5) # bathy
- fig.ax.set_xlabel(canvas.tr('Distance (m)'), fontsize = fontsize_axis, weight = 'bold')
- fig.ax.set_ylabel(canvas.tr('Depth (m)'), fontsize = fontsize_axis, weight = 'bold')
- fig.ax.xaxis.label.set_fontsize(12)
- fig.ax.yaxis.label.set_fontsize(12)
+ fig.ax.set_xlabel(canvas.tr('Distance (m)'), fontsize = fontsize_axis+6, weight = 'bold')
+ fig.ax.set_ylabel(canvas.tr('Depth (m)'), fontsize = fontsize_axis+6, weight = 'bold')
fig.ax.tick_params(axis='both', direction='in', bottom=True, top=True, left=True, right=True)
fig.ax.set_ylim(top=0, bottom=np.nanmax(track_depth)+0.3)
lower_limit = track[0]-0.5
@@ -298,10 +296,10 @@ def SDC_fines_method(analysis_data, map_data, outpath, outpath_figures, sampling
legend = fig.ax.legend(loc = 'upper right', fontsize=fontsize_legend,
title = 'Sampler', bbox_to_anchor=(1.17, 1))
plt.setp(legend.get_title(), fontsize = fontsize_legend_title)
- fig.ax.text(1.9, 0.12, '$\overline{C}_{fine}$ = ' + str(np.round(Conc_mean_SDC_fines,3)) + ' g/l',
+ fig.ax.text(2.2, 0.1, '$\overline{C_{\mathregular{fine}}}$ = ' + str(np.round(Conc_mean_SDC_fines,3)) + ' g/l',
transform = ax.transAxes, fontsize = fontsize_text)
- fig.ax.text(2.7, 0.12, '$\Phi_{total}$ = ' + str(np.round(total_fine_flux_SDC_kg_s,0)) + ' kg/s', transform = ax.transAxes,
- fontsize = fontsize_text)
+ # fig.ax.text(2.7, 0.12, '$\Phi_{total}$ = ' + str(np.round(total_fine_flux_SDC_kg_s,0)) + ' kg/s', transform = ax.transAxes,
+ # fontsize = fontsize_text)
canvas.draw()
show_figure(fig)
@@ -663,7 +661,6 @@ def SDC_fines_method(analysis_data, map_data, outpath, outpath_figures, sampling
cb.ax.set_yticklabels(['0.01','0.1','0.2','0.4','0.6','0.8','1.2','2.0','3.0', '4.0', '5.0','6.0','7.0', '10.0'],ha='right')
cb.ax.yaxis.set_tick_params(pad=35)
cb.ax.set_ylabel(canvas.tr('$\mathregular{C_{fine}/C_{sand}}$'), fontsize = fontsize_axis, weight = 'bold')
- cb.ax.yaxis.label.set_fontsize(12)
cb.ax.tick_params(labelsize=12)
fig.ax.invert_yaxis()
fig.ax.fill_between(x_fill, np.nanmax(track_depth)+2, depth_fill, color='w') # below bathy
@@ -671,8 +668,6 @@ def SDC_fines_method(analysis_data, map_data, outpath, outpath_figures, sampling
fig.ax.set_xlabel(canvas.tr('Distance (m)'), fontsize = fontsize_axis, weight = 'bold')
fig.ax.set_ylabel(canvas.tr('Depth (m)'), fontsize = fontsize_axis, weight = 'bold')
- fig.ax.xaxis.label.set_fontsize(12)
- fig.ax.yaxis.label.set_fontsize(12)
fig.ax.tick_params(axis='both', direction='in', bottom=True, top=True, left=True, right=True)
fig.ax.set_ylim(top=0, bottom=np.nanmax(track_depth)+0.3)
lower_limit = track[0]-0.5
@@ -689,7 +684,7 @@ def SDC_fines_method(analysis_data, map_data, outpath, outpath_figures, sampling
legend = fig.ax.legend(loc = 'upper right', fontsize=fontsize_legend,
title = 'Sampler', bbox_to_anchor=(1.17, 1))
plt.setp(legend.get_title(), fontsize = fontsize_legend_title)
- fig.ax.text(2, 0.12, '$\overline{C_{fine}/C_{sand}}$ = ' + str(np.round(mean_ratio,3)),
+ fig.ax.text(2, 0.1, '$\overline{C_{fine}/C_{sand}}$ = ' + str(np.round(mean_ratio,3)),
transform = ax.transAxes, fontsize = fontsize_text)
canvas.draw()
diff --git a/classes/SDC_method.py b/classes/SDC_method.py
index 272dcbe7964661719c74cba8ceb48e94b80d12cc..5ebecbb1cf7d690887e8c216b711215262a718fc 100644
--- a/classes/SDC_method.py
+++ b/classes/SDC_method.py
@@ -32,7 +32,7 @@ from Common_functions import show_figure
import matplotlib.colors as colors
-def SDC_method(analysis_data, map_data, outpath, outpath_figures, sampling_date, maxpost, maxpost_alpha, maxpost_lnCr):
+def SDC_method(analysis_data, map_data, outpath, outpath_figures, sampling_date, maxpost, maxpost_alpha, maxpost_lnCr, add_plots):
#%%
abscissa_values = (analysis_data['Abscissa']).unique().tolist()
@@ -167,130 +167,36 @@ def SDC_method(analysis_data, map_data, outpath, outpath_figures, sampling_date,
figname = '_Conc_profiles_log_SDC_unc_article' + most_frequent_sampler
fig.savefig(outpath_figures + sampling_date + figname + '.png', dpi = 400, bbox_inches='tight')
-
-
- #%% Exponential profile
- # stats_vertical_profile = []
- # x_ranges_vertical = []
-
- # for i in range(len(abscissa_values)):
- # x = np.array(analysis_data['z_h'][analysis_data['Abscissa'] == abscissa_values[i]])
- # Y = np.log(np.array(analysis_data['Concentration_sand_g_l'][analysis_data['Abscissa'] == abscissa_values[i]]))
- # array_x_vertical = np.linspace(1, 0, num=100)
- # x_ranges_vertical.append(array_x_vertical)
-
- # res = stats.linregress(x, Y)
- # # correct slope (increasing C with depth)
- # res = list(res)
- # if res[0] >= 0:
- # res[0] = -0.2
- # res[1] = np.max(Y)
- # res[2] = np.nan
- # print('Slope and intercept of the SDC sand profile corrected for vertical i = ' + str(i))
- # stats_vertical_profile.append(res)
-
- # x_ranges_vertical = pd.DataFrame(x_ranges_vertical)
- # stats_vertical_profile = pd.DataFrame(stats_vertical_profile)
- # stats_vertical_profile.columns = ['slope', 'intercept','rvalue', 'pvalue', 'stderr']
- # stats_vertical_profile['Cr_vertical'] = np.exp(stats_vertical_profile['intercept'])
- # stats_vertical_profile['Crh_vertical'] = stats_vertical_profile['Cr_vertical']/depth_verticals
- # stats_vertical_profile['alphah_vertical'] = -(stats_vertical_profile['slope'])
- # stats_vertical_profile['alpha_vertical'] = stats_vertical_profile['alphah_vertical']/depth_verticals
- # stats_vertical_profile['R2'] = stats_vertical_profile['rvalue']**2
-
- # Conc_vertical_profile = pd.DataFrame([np.exp(stats_vertical_profile['intercept'][i])*np.exp(stats_vertical_profile['slope'][i]*x_ranges_vertical.iloc[i,:])
- # for i in range(len(stats_vertical_profile))])
-
-
-#%% Plot all vertical profiles in one graph
- # fig, ax = plt.subplots(figsize=(10, 6), dpi = 100)
-
- # for i in range(len(abscissa_values)):
- # ax.plot(Conc_vertical_profile.iloc[i,:], x_ranges_vertical.iloc[i,:],
- # linestyle = '-', linewidth = 2, color = colorss[i],
- # label = str(abscissa_values[i]))
-
- # ax.plot(analysis_data['Concentration_sand_g_l'][analysis_data['Abscissa']== abscissa_values[i]],
- # analysis_data['z_h'][analysis_data['Abscissa']== abscissa_values[i]],
- # color = colorss[i], linestyle = ' ',
- # markersize = 9, marker = markerss[i],
- # markeredgewidth = 0.2, markeredgecolor = 'black')
- # # ax.text(0.02, 1, 'a)', fontsize = 14)
- # legend = ax.legend(fontsize = fontsize_legend, title = 'Abscissa',
- # loc = 'upper right', facecolor = 'white', framealpha = 1, bbox_to_anchor = (1.2,1))
- # plt.setp(legend.get_title(), fontsize=fontsize_legend_title)
-
- # ax.set_ylim(-0.05,1.05)
- # ax.set_xlim(0)
- # ax.set_ylabel('z/h (-)', fontsize = fontsize_axis, weight = 'bold')
- # ax.set_xlabel('Concentration (g/l)', fontsize = fontsize_axis, weight = 'bold')
- # ax.grid(linewidth = 0.2)
- # ax.tick_params(axis='both', which='major', labelsize = fontsize_ticks)
- # fig.tight_layout()
- # figname = '_Conc_profiles_SDC' + most_frequent_sampler
- # fig.savefig(outpath_figures + sampling_date + figname + '.png', dpi = 200, bbox_inches='tight')
-
- # # log-scale
- # cmap = plt.cm.get_cmap('nipy_spectral')
- # colorss = cmap(np.linspace(0,1,len(abscissa_values)))
-
- # fig, ax = plt.subplots(figsize=(10, 6), dpi = 100)
-
- # for i in range(len(abscissa_values)):
- # ax.plot(Conc_vertical_profile.iloc[i,:], x_ranges_vertical.iloc[i,:],
- # linestyle = '-', linewidth = 2, color = colorss[i],
- # label = str(abscissa_values[i]))
-
- # ax.plot(analysis_data['Concentration_sand_g_l'][analysis_data['Abscissa']== abscissa_values[i]],
- # analysis_data['z_h'][analysis_data['Abscissa']== abscissa_values[i]],
- # color = colorss[i], linestyle = ' ',
- # markersize = 9, marker = markerss[i],
- # markeredgewidth = 0.2, markeredgecolor = 'black')
-
- # legend = ax.legend(fontsize = fontsize_legend, title = 'Abscissa',
- # loc = 'upper right', facecolor = 'white', framealpha = 1, bbox_to_anchor = (1.2,1))
- # plt.setp(legend.get_title(), fontsize=fontsize_legend_title)
-
- # ax.set_xscale('log')
- # ax.set_ylim(-0.05,1.05)
- # #ax.set_xlim(0)
- # ax.set_ylabel('z/h (-)', fontsize = fontsize_axis, weight = 'bold')
- # ax.set_xlabel('Concentration (g/l)', fontsize = fontsize_axis, weight = 'bold')
- # ax.grid(linewidth = 0.2)
- # ax.tick_params(axis='both', which='major', labelsize = fontsize_ticks)
- # fig.tight_layout()
- # figname = '_Conc_profiles_log_SDC_' + most_frequent_sampler
- # fig.savefig(outpath_figures + sampling_date + figname + '.png', dpi = 200, bbox_inches='tight')
-
#%% Plot vertical profiles for each abscissa
- for i in range(len(abscissa_values)):
- fig, ax = plt.subplots(figsize=(10, 6), dpi = 100)
- ax.invert_yaxis()
- ax.plot(maxpost_g_l.iloc[:,i], x_ranges_vertical[:],
- linestyle = '-', linewidth = 2, color = colorss[i], label = 'Abscissa ' + str(i +1))
-
- ax.plot(analysis_data['Concentration_sand_g_l'][analysis_data['Abscissa']== abscissa_values[i]],
- analysis_data['z_h'][analysis_data['Abscissa']== abscissa_values[i]],
- color = colorss[i], linestyle = ' ',
- markersize = 9, marker = markerss[i],
- markeredgewidth = 0.2, markeredgecolor = 'black')
-
- a = np.round(maxpost_alpha[i], 2)
- b = np.round(maxpost_lnCr[i], 2)
- # ax.text(0.75, 0.92, '$y=%3.7s*e^{%3.7sx}$'%(a, b),
- # transform=ax.transAxes, fontsize = fontsize_text)
- # ax.text(0.75, 0.85, '${R²}$ = ' + str(np.round(stats_vertical_profile['R2'][i], 2)),
- # transform=ax.transAxes, fontsize = fontsize_text)
-
- ax.set_ylim(-0.05,1.05)
- #ax.set_xlim(0)
- ax.set_ylabel('z/h (-)', fontsize = 16, weight = 'bold')
- ax.set_xlabel('Concentration (g/l)', fontsize = 16, weight = 'bold')
- ax.grid(linewidth = 0.2)
- ax.tick_params(axis='both', which='major', labelsize = 14)
- fig.tight_layout()
- figname = sampling_date + '_Conc_profiles_SDC_'
- fig.savefig(outpath_figures + figname + str(abscissa_values[i]) + '_' + most_frequent_sampler + '.png', dpi = 200, bbox_inches='tight')
+ if add_plots == True:
+ for i in range(len(abscissa_values)):
+ fig, ax = plt.subplots(figsize=(10, 6), dpi = 100)
+ ax.invert_yaxis()
+ ax.plot(maxpost_g_l.iloc[:,i], x_ranges_vertical[:],
+ linestyle = '-', linewidth = 2, color = colorss[i], label = 'Abscissa ' + str(i +1))
+
+ ax.plot(analysis_data['Concentration_sand_g_l'][analysis_data['Abscissa']== abscissa_values[i]],
+ analysis_data['z_h'][analysis_data['Abscissa']== abscissa_values[i]],
+ color = colorss[i], linestyle = ' ',
+ markersize = 9, marker = markerss[i],
+ markeredgewidth = 0.2, markeredgecolor = 'black')
+
+ a = np.round(maxpost_alpha[i], 2)
+ b = np.round(maxpost_lnCr[i], 2)
+ # ax.text(0.75, 0.92, '$y=%3.7s*e^{%3.7sx}$'%(a, b),
+ # transform=ax.transAxes, fontsize = fontsize_text)
+ # ax.text(0.75, 0.85, '${R²}$ = ' + str(np.round(stats_vertical_profile['R2'][i], 2)),
+ # transform=ax.transAxes, fontsize = fontsize_text)
+
+ ax.set_ylim(-0.05,1.05)
+ #ax.set_xlim(0)
+ ax.set_ylabel('z/h (-)', fontsize = 16, weight = 'bold')
+ ax.set_xlabel('Concentration (g/l)', fontsize = 16, weight = 'bold')
+ ax.grid(linewidth = 0.2)
+ ax.tick_params(axis='both', which='major', labelsize = 14)
+ fig.tight_layout()
+ figname = sampling_date + '_Conc_profiles_SDC_'
+ fig.savefig(outpath_figures + figname + str(abscissa_values[i]) + '_' + most_frequent_sampler + '.png', dpi = 200, bbox_inches='tight')
#%% # Define MAP grid
cell_height = map_data.depth_cells_border[1][0]-map_data.depth_cells_border[0][0]
diff --git a/classes/Sed_flux_big_bend.py b/classes/Sed_flux_big_bend.py
index bfd946d5c9c4acc89806cfd0cc5e2b8714f55c08..54a1cac95424ca6c6fe9e542a2ea4e7d6510347e 100644
--- a/classes/Sed_flux_big_bend.py
+++ b/classes/Sed_flux_big_bend.py
@@ -49,12 +49,12 @@ def sed_flux_big_bend(outpath_figures) :
# Parametres
BD_grain_size_steps = [75, 90, 100, 110, 130, 150, 210]
name = '_big_bend'
- title = name.replace('_', ' ')
- title = title[1:].capitalize()
+ # title = name.replace('_', ' ')
+ # title = title[1:].capitalize()
colors = ['tab:blue','tab:orange', 'tab:green', 'tab:red', 'tab:purple', 'tab:brown', 'tab:pink']
# Define and plot alpha
- # plt.figure(figsize=(9,6))
+ plt.figure(figsize=(9,6))
k = 0
for k in range(len(BD_grain_size_steps)) :
x_range = np.arange(globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][0],
@@ -63,40 +63,56 @@ def sed_flux_big_bend(outpath_figures) :
globals()['pars' + str(BD_grain_size_steps[k]) + str(name)], globals()['cov' + str(BD_grain_size_steps[k]) + str(name)] = curve_fit(
f=exponential, xdata=globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)],
ydata=globals()['y_alpha' + str(BD_grain_size_steps[k]) + str(name)], maxfev=5000)
- # plt.plot(x_range_extrap, exponential(x_range_extrap,
- # globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='--',
- # color=colors[k])
- # plt.plot(x_range, exponential(x_range,
- # globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='-',
- # color=colors[k], label=f' {BD_grain_size_steps[k]}')
- globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = np.round(exponential(x_range,
+ plt.plot(x_range_extrap, exponential(x_range_extrap,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='--',
+ color=colors[k])
+ plt.plot(x_range, exponential(x_range,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='-',
+ color=colors[k], label=f' {BD_grain_size_steps[k]}')
+
+ extrapolate_values = np.round(exponential(x_range,
globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), 3)
- # plt.plot([], [], linestyle='--', color='black', label='extrapolation')
+
+ globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = extrapolate_values
+
+ plt.plot([], [], linestyle='--', color='black', label='extrapolation')
- # plt.legend(title='$D_{50}$ (\u03BCm)')
+ plt.legend(title='$D_{50}$ (\u03BCm)', fontsize = 16, bbox_to_anchor = (1.05, 1))
# plt.grid(linewidth=0.2)
- # plt.xlim(0.4, 1.8)
- # plt.ylim(0.6, 2.5)
- # plt.xlabel('Velocity m/s')
- # plt.ylabel(r'$ \alpha $ (-)')
+ plt.xlim(0.4, 1.8)
+ plt.ylim(0.6, 2.5)
+ plt.xlabel('Velocity (m/s)', fontsize = 16, weight = 'bold')
+ plt.ylabel(r'$\mathregular{\alpha}$ (-)', fontsize = 16, weight = 'bold')
# plt.title(str(title) + ' nozzle')
- # figname = 'Big_bend'
- # plt.savefig(outpath_figures + figname + '.png', dpi=400, bbox_inches='tight')
+ figname = 'Big_bend'
+ plt.savefig(outpath_figures + figname + '.png', dpi=400, bbox_inches='tight')
# Define values out of lower range (velocity < 0.5 m/s)
# Define values out of upper range (velocity > 2.4 m/s)
# Small straight nozzle
globals()['alpha' + str(name)] = []
+ alpha_temp = []
k = 0
- for k in range(len(BD_grain_size_steps)) :
- globals()['alpha' + str(BD_grain_size_steps[k]) + '_lower_outrange' + str(name)] = np.array(
- [np.min(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])] * 120)
- globals()['alpha' + str(BD_grain_size_steps[k]) + '_upper_outrange' + str(name)] = np.array(
- [np.max(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])] * 120)
- globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = np.concatenate(
- (globals()['alpha' + str(BD_grain_size_steps[k]) + '_lower_outrange' + str(name)], globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)]
- , globals()['alpha' + str(BD_grain_size_steps[k]) + '_upper_outrange' + str(name)]))
- globals()['alpha' + str(name)].append(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])
-
- sed_flux_big_bend.alpha_big_bend = globals()['alpha' + str(name)]
\ No newline at end of file
+ for k in range(len(BD_grain_size_steps)) :
+ x_range = np.arange(globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][0],
+ globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][-1], 0.01)
+ extrapolate_values = np.round(exponential(x_range,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), 3)
+
+ vel_possiblerange = np.round(np.arange(0.01, 4.01, 0.01), 3)
+ idx_start = np.where(vel_possiblerange == np.round(x_range[0], 2))[0]
+ idx_end = np.where(vel_possiblerange == np.round(x_range[-1], 2))[0]
+
+ alpha_start = [extrapolate_values[0]] * (idx_start[0] - 1)
+ alpha_end = [extrapolate_values[-1]] * (len(vel_possiblerange) - idx_end[0])
+
+ alpha_k = np.concatenate([alpha_start, extrapolate_values, alpha_end])
+ alpha_temp.append(alpha_k)
+
+
+ sed_flux_big_bend.alpha_big_bend = alpha_temp
+ # sed_flux_big_bend.alpha_big_bend = globals()['alpha' + str(name)]
+
+ # test = pd.DataFrame([x_range, extrapolate_values])
+
diff --git a/classes/Sed_flux_big_straight.py b/classes/Sed_flux_big_straight.py
index 9e9675491cb4a67742c3c646eda84b61220b75da..80cedad1228e91669ed85910b7c0cf6abfc9321d 100644
--- a/classes/Sed_flux_big_straight.py
+++ b/classes/Sed_flux_big_straight.py
@@ -19,6 +19,7 @@ along with this program. If not, see <https://www.gnu.org/licenses/>.
import numpy as np
from scipy.optimize import curve_fit
+import matplotlib.pyplot as plt
#############################################
@@ -47,12 +48,12 @@ def sed_flux_big_straight(outpath_figures) :
# Parametres
BD_grain_size_steps = [75, 90, 100, 110, 130, 150, 210]
name = '_big_straight'
- title = name.replace('_', ' ')
- title = title[1:].capitalize()
-
+ # title = name.replace('_', ' ')
+ # title = title[1:].capitalize()
+ colors = ['tab:blue','tab:orange', 'tab:green', 'tab:red', 'tab:purple', 'tab:brown', 'tab:pink']
+
# Define and plot alpha
-
- # plt.figure(figsize=(9,6))
+ plt.figure(figsize=(9,6))
k = 0
for k in range(len(BD_grain_size_steps)) :
x_range = np.arange(globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][0],
@@ -61,26 +62,26 @@ def sed_flux_big_straight(outpath_figures) :
globals()['pars' + str(BD_grain_size_steps[k]) + str(name)], globals()['cov' + str(BD_grain_size_steps[k]) + str(name)] = curve_fit(
f=exponential, xdata=globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)],
ydata=globals()['y_alpha' + str(BD_grain_size_steps[k]) + str(name)], maxfev=5000)
- # plt.plot(x_range_extrap, exponential(x_range_extrap,
- # globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='--',
- # color=colors[k])
- # plt.plot(x_range, exponential(x_range,
- # globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='-',
- # color=colors[k], label=f' {BD_grain_size_steps[k]}')
+ plt.plot(x_range_extrap, exponential(x_range_extrap,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='--',
+ color=colors[k])
+ plt.plot(x_range, exponential(x_range,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='-',
+ color=colors[k], label=f' {BD_grain_size_steps[k]}')
globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = np.round(exponential(x_range,
globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), 3)
- # plt.plot([], [], linestyle='--', color='black', label='extrapolation')
+ plt.plot([], [], linestyle='--', color='black', label='extrapolation')
- # plt.legend(title='$D_{50}$ (\u03BCm)')
+ plt.legend(title='$D_{50}$ (\u03BCm)', fontsize = 16, bbox_to_anchor = (1.05, 1))
# plt.grid(linewidth=0.2)
- # plt.xlim(0.4, 1.8)
- # plt.ylim(0.6, 2.5)
- # plt.xlabel('Velocity m/s')
- # plt.ylabel(r'$ \alpha $ (-)')
- # plt. title(str(title) + ' nozzle')
+ plt.xlim(0.4, 1.8)
+ plt.ylim(0.6, 2.5)
+ plt.xlabel('Velocity (m/s)', fontsize = 16, weight = 'bold')
+ plt.ylabel(r'$\mathregular{\alpha}$ (-)', fontsize = 16, weight = 'bold')
+ # plt.title(str(title) + ' nozzle')
- # figname = 'Big_straight'
- # # plt.savefig(outpath_figures + figname + '.png', dpi=400, bbox_inches='tight')
+ figname = 'Big_straight'
+ plt.savefig(outpath_figures + figname + '.png', dpi=400, bbox_inches='tight')
# Define values out of lower range (velocity < 0.5 m/s)
# x_lower_outrange_big = np.arange(0, 0.5, 0.01)
@@ -88,15 +89,23 @@ def sed_flux_big_straight(outpath_figures) :
# Define values out of upper range (velocity > 2.4 m/s)
# Small straight nozzle
globals()['alpha' + str(name)] = []
+ alpha_temp = []
k = 0
- for k in range(len(BD_grain_size_steps)) :
- globals()['alpha' + str(BD_grain_size_steps[k]) + '_lower_outrange' + str(name)] = np.array(
- [np.min(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])] * 120)
- globals()['alpha' + str(BD_grain_size_steps[k]) + '_upper_outrange' + str(name)] = np.array(
- [np.max(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])] * 120)
- globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = np.concatenate(
- (globals()['alpha' + str(BD_grain_size_steps[k]) + '_lower_outrange' + str(name)], globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)]
- , globals()['alpha' + str(BD_grain_size_steps[k]) + '_upper_outrange' + str(name)]))
- globals()['alpha' + str(name)].append(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])
-
- sed_flux_big_straight.alpha_big_straight = globals()['alpha' + str(name)]
+ for k in range(len(BD_grain_size_steps)) :
+ x_range = np.arange(globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][0],
+ globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][-1], 0.01)
+ extrapolate_values = np.round(exponential(x_range,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), 3)
+
+ vel_possiblerange = np.round(np.arange(0.01, 4.01, 0.01), 3)
+ idx_start = np.where(vel_possiblerange == np.round(x_range[0], 2))[0]
+ idx_end = np.where(vel_possiblerange == np.round(x_range[-1], 2))[0]
+
+ alpha_start = [extrapolate_values[0]] * (idx_start[0] - 1)
+ alpha_end = [extrapolate_values[-1]] * (len(vel_possiblerange) - idx_end[0])
+
+ alpha_k = np.concatenate([alpha_start, extrapolate_values, alpha_end])
+ alpha_temp.append(alpha_k)
+
+
+ sed_flux_big_straight.alpha_big_straight = alpha_temp
\ No newline at end of file
diff --git a/classes/Sed_flux_small_bend.py b/classes/Sed_flux_small_bend.py
index 327b00d7e7a630f3e4a02617ef7f065dbd4256ea..c6a98c6737169a90f02a79b53c8fdfceb45f6cc2 100644
--- a/classes/Sed_flux_small_bend.py
+++ b/classes/Sed_flux_small_bend.py
@@ -50,13 +50,12 @@ def sed_flux_small_bend(outpath_figures) :
# Parametres
BD_grain_size_steps = [75, 90, 100, 110, 130, 150, 210]
name = '_small_bend'
- title = name.replace('_', ' ')
- title = title[1:].capitalize()
+ # title = name.replace('_', ' ')
+ # title = title[1:].capitalize()
colors = ['tab:blue','tab:orange', 'tab:green', 'tab:red', 'tab:purple', 'tab:brown', 'tab:pink']
# Define and plot alpha
-
- # plt.figure(figsize=(9,6))
+ plt.figure(figsize=(9,6))
k = 0
for k in range(len(BD_grain_size_steps)) :
x_range = np.arange(globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][0],
@@ -65,26 +64,24 @@ def sed_flux_small_bend(outpath_figures) :
globals()['pars' + str(BD_grain_size_steps[k]) + str(name)], globals()['cov' + str(BD_grain_size_steps[k]) + str(name)] = curve_fit(
f=exponential, xdata=globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)],
ydata=globals()['y_alpha' + str(BD_grain_size_steps[k]) + str(name)], maxfev=5000)
- # plt.plot(x_range_extrap, exponential(x_range_extrap,
- # globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='--',
- # color=colors[k])
- # plt.plot(x_range, exponential(x_range,
- # globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='-',
- # color=colors[k], label=f' {BD_grain_size_steps[k]}')
+ plt.plot(x_range_extrap, exponential(x_range_extrap,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='--',
+ color=colors[k])
+ plt.plot(x_range, exponential(x_range,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='-',
+ color=colors[k], label=f' {BD_grain_size_steps[k]}')
globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = np.round(exponential(x_range,
globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), 3)
- # plt.plot([], [], linestyle='--', color='black', label='extrapolation')
+ plt.plot([], [], linestyle='--', color='black', label='extrapolation')
- # plt.legend(title='$D_{50}$ (\u03BCm)')
- # plt.grid(linewidth=0.2)
- # plt.xlim(0.4, 2.8)
- # plt.ylim(0.6, 2.5)
- # plt.xlabel('Velocity m/s')
- # plt.ylabel(r'$ \alpha $ (-)')
- # plt.title(str(title) + ' nozzle')
+ # plt.legend(title='$D_{50}$ (\u03BCm)', fontsize = 16, bbox_to_anchor = (1.219, 0.3))
+ plt.xlim(0.4, 1.8)
+ plt.ylim(0.6, 2.5)
+ plt.xlabel('Velocity (m/s)', fontsize = 16, weight = 'bold')
+ plt.ylabel(r'$\mathregular{\alpha}$ (-)', fontsize = 16, weight = 'bold')
- # figname = 'Small_bend'
- # # plt.savefig(outpath_figures + figname + '.png', dpi=400, bbox_inches='tight')
+ figname = 'Small_bend'
+ plt.savefig(outpath_figures + figname + '.png', dpi=400, bbox_inches='tight')
# Define values out of lower range (velocity < 0.6 m/s)
x_lower_outrange_small = np.arange(0, 0.6, 0.01)
@@ -94,12 +91,22 @@ def sed_flux_small_bend(outpath_figures) :
# Small straight nozzle
globals()['alpha' + str(name)] = []
k = 0
+ alpha_temp= []
for k in range(len(BD_grain_size_steps)) :
- globals()['alpha' + str(BD_grain_size_steps[k]) + '_upper_outrange' + str(name)] = np.array(
- [np.max(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])] * 120)
- globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = np.concatenate(
- (alpha_lower_outrange_small, globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)]
- , globals()['alpha' + str(BD_grain_size_steps[k]) + '_upper_outrange' + str(name)]))
- globals()['alpha' + str(name)].append(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])
-
- sed_flux_small_bend.alpha_small_bend = globals()['alpha' + str(name)]
+ x_range = np.arange(globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][0],
+ globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][-1], 0.01)
+ extrapolate_values = np.round(exponential(x_range,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), 3)
+
+ vel_possiblerange = np.round(np.arange(0.01, 4.01, 0.01), 3)
+ idx_start = np.where(vel_possiblerange == np.round(x_range[0], 2))[0]
+ idx_end = np.where(vel_possiblerange == np.round(x_range[-1], 2))[0]
+
+ alpha_start = [extrapolate_values[0]] * (idx_start[0] - 1)
+ alpha_end = [extrapolate_values[-1]] * (len(vel_possiblerange) - idx_end[0])
+
+ alpha_k = np.concatenate([alpha_start, extrapolate_values, alpha_end])
+ alpha_temp.append(alpha_k)
+
+
+ sed_flux_small_bend.alpha_small_bend = alpha_temp
\ No newline at end of file
diff --git a/classes/Sed_flux_small_straight.py b/classes/Sed_flux_small_straight.py
index c81f40e0407fd80d820e1faa0859e86bf7ac6d89..8319479c182ef9f8b5e356ac8e3312a1dd28de4c 100644
--- a/classes/Sed_flux_small_straight.py
+++ b/classes/Sed_flux_small_straight.py
@@ -19,6 +19,7 @@ along with this program. If not, see <https://www.gnu.org/licenses/>.
import numpy as np
from scipy.optimize import curve_fit
+import matplotlib.pyplot as plt
#############################################
@@ -44,17 +45,15 @@ def exponential(x, a, b, c):
return a*x**3+b*x**2+c
def sed_flux_small_straight(outpath_figures) :
-
# Parametres
BD_grain_size_steps = [75, 90, 100, 110, 130, 150, 210]
name = '_small_straight'
- title = name.replace('_', ' ')
- title = title[1:].capitalize()
+ # title = name.replace('_', ' ')
+ # title = title[1:].capitalize()
colors = ['tab:blue','tab:orange', 'tab:green', 'tab:red', 'tab:purple', 'tab:brown', 'tab:pink']
# Define and plot alpha
-
- # plt.figure(figsize=(9,6))
+ plt.figure(figsize=(9,6))
k = 0
for k in range(len(BD_grain_size_steps)) :
x_range = np.arange(globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][0],
@@ -63,26 +62,24 @@ def sed_flux_small_straight(outpath_figures) :
globals()['pars' + str(BD_grain_size_steps[k]) + str(name)], globals()['cov' + str(BD_grain_size_steps[k]) + str(name)] = curve_fit(
f=exponential, xdata=globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)],
ydata=globals()['y_alpha' + str(BD_grain_size_steps[k]) + str(name)], maxfev=5000)
- # plt.plot(x_range_extrap, exponential(x_range_extrap,
- # globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='--',
- # color=colors[k])
- # plt.plot(x_range, exponential(x_range,
- # globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='-',
- # color=colors[k], label=f' {BD_grain_size_steps[k]}')
+ plt.plot(x_range_extrap, exponential(x_range_extrap,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='--',
+ color=colors[k])
+ plt.plot(x_range, exponential(x_range,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), linestyle='-',
+ color=colors[k], label=f' {BD_grain_size_steps[k]}')
globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = np.round(exponential(x_range,
globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), 3)
- # plt.plot([], [], linestyle='--', color='black', label='extrapolation')
+ plt.plot([], [], linestyle='--', color='black', label='extrapolation')
- # plt.legend(title='$D_{50}$ (\u03BCm)')
- # plt.grid(linewidth=0.2)
- # plt.xlim(0.4, 2.8)
- # plt.ylim(0.6, 2.5)
- # plt.xlabel('Velocity m/s')
- # plt.ylabel(r'$ \alpha $ (-)')
- # plt.title(str(title) + ' nozzle')
+ # plt.legend(title='$D_{50}$ (\u03BCm)', fontsize = 16, bbox_to_anchor = (1.05, 1))
+ plt.xlim(0.4, 1.8)
+ plt.ylim(0.6, 2.5)
+ plt.xlabel('Velocity (m/s)', fontsize = 16, weight = 'bold')
+ plt.ylabel(r'$\mathregular{\alpha}$ (-)', fontsize = 16, weight = 'bold')
- # figname = 'Small_straight'
- # # plt.savefig(outpath_figures + figname + '.png', dpi=400, bbox_inches='tight')
+ figname = 'Small_straight'
+ plt.savefig(outpath_figures + figname + '.png', dpi=400, bbox_inches='tight')
# Define values out of lower range (velocity < 0.6 m/s)
@@ -93,12 +90,22 @@ def sed_flux_small_straight(outpath_figures) :
# Small straight nozzle
globals()['alpha' + str(name)] = []
k = 0
+ alpha_temp= []
for k in range(len(BD_grain_size_steps)) :
- globals()['alpha' + str(BD_grain_size_steps[k]) + '_upper_outrange' + str(name)] = np.array(
- [np.max(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])] * 120)
- globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)] = np.concatenate(
- (alpha_lower_outrange_small, globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)]
- , globals()['alpha' + str(BD_grain_size_steps[k]) + '_upper_outrange' + str(name)]))
- globals()['alpha' + str(name)].append(globals()['alpha' + str(BD_grain_size_steps[k]) + str(name)])
-
- sed_flux_small_straight.alpha_small_straight = globals()['alpha' + str(name)]
+ x_range = np.arange(globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][0],
+ globals()['x_vel' + str(BD_grain_size_steps[k]) + str(name)][-1], 0.01)
+ extrapolate_values = np.round(exponential(x_range,
+ globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][0],globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][1], globals()['pars' + str(BD_grain_size_steps[k]) + str(name)][2]), 3)
+
+ vel_possiblerange = np.round(np.arange(0.01, 4.01, 0.01), 3)
+ idx_start = np.where(vel_possiblerange == np.round(x_range[0], 2))[0]
+ idx_end = np.where(vel_possiblerange == np.round(x_range[-1], 2))[0]
+
+ alpha_start = [extrapolate_values[0]] * (idx_start[0] - 1)
+ alpha_end = [extrapolate_values[-1]] * (len(vel_possiblerange) - idx_end[0])
+
+ alpha_k = np.concatenate([alpha_start, extrapolate_values, alpha_end])
+ alpha_temp.append(alpha_k)
+
+
+ sed_flux_small_straight.alpha_small_straight = alpha_temp
\ No newline at end of file
diff --git a/classes/Sediment_flux_BD_single.py b/classes/Sediment_flux_BD_single.py
index c9de9f1a8fff543f6b21a37d81dcb54102938ca7..2e820647d35ff52ebcf7732bfa43a1a8964ee6f0 100644
--- a/classes/Sediment_flux_BD_single.py
+++ b/classes/Sediment_flux_BD_single.py
@@ -28,7 +28,7 @@ from Sed_flux_big_straight import sed_flux_big_straight
from Sed_flux_big_bend import sed_flux_big_bend
from Functions import find_nearest, get_sec
-def sediment_flux_BD(outpath_figures, meas, choice, unit) :
+def sediment_flux_BD(outpath_figures, meas, choice, unit, add_plots) :
# %%################### Calculate sediment flux ###########################
unit_mass = unit.unit_masses
@@ -45,18 +45,18 @@ def sediment_flux_BD(outpath_figures, meas, choice, unit) :
BD_grain_size_steps = [75, 90, 100, 110, 130, 150, 210]
- # Create Alpha graphs
+ # Create Alpha relations
sed_flux_small_straight(outpath_figures)
+ sed_flux_small_bend(outpath_figures)
+ sed_flux_big_straight(outpath_figures)
+ sed_flux_big_bend(outpath_figures)
+
alpha_small_straight = sed_flux_small_straight.alpha_small_straight
-
- sed_flux_small_bend(outpath_figures)
alpha_small_bend = sed_flux_small_bend.alpha_small_bend
-
- sed_flux_big_straight(outpath_figures)
alpha_big_straight = sed_flux_big_straight.alpha_big_straight
-
- sed_flux_big_bend(outpath_figures)
alpha_big_bend = sed_flux_big_bend.alpha_big_bend
+
+
# %%################### Determine alpha and nozzle surface for samples #########
@@ -90,7 +90,7 @@ def sediment_flux_BD(outpath_figures, meas, choice, unit) :
if pd.isna(meas['Sampling_duration_s']):
meas['Sampling_duration_s'] = meas['Filling_time_field']
-
+ print(alpha)
# %%################### Calculate point solid flux #######################
if unit_mass == 'mg':
sand_flux_point_i = alpha * meas['Dry_mass_sand'] * \
diff --git a/classes/Uncertainty_estimation.py b/classes/Uncertainty_estimation.py
index 01bb71adf3904f4f2f4e06a231d84495c7ddcfeb..9f8d5f4954fb7f3156074b86cd69f6ad9ae94a77 100644
--- a/classes/Uncertainty_estimation.py
+++ b/classes/Uncertainty_estimation.py
@@ -24,7 +24,7 @@ import os
from scipy import stats
def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
- choice, v2_d_index, summary_ISO, summary_fine, map_data, cwd, path_R):
+ choice, v2_d_index, summary_ISO, summary_fine, map_data, cwd, path_R, add_plots):
#%%
@@ -159,48 +159,51 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
for i in range(len(abscissa_values)):
vv = analysis_data[analysis_data['Abscissa'] == abscissa_values[i]]
# Measured parameters
- if choice[2] == 1 and choice[3] == 1:
- if np.count_nonzero(~np.isnan(np.array(vv['D50']))) >= 2:
- d50_vertical = np.mean(vv['D50'])*1e-6
- d90_bed = 0.0004
- velocity_vert = vv['Mean_velocity_vertical'].iloc[0]
- h = vv['max_depth_vertical'].iloc[0]
- # roughness height
- kst = 0.05
- ks = 2* d90_bed
- # kinematic viscosity of water
- nu = 0.000001
- # acceleration of the gravity
- g = 9.81
- # relative density
- rho = 2.65
- # sedimentological diameter
- ds = d50_vertical*(g*(rho-1)/nu**2)**(1/3) # use mean d50 per vertical
- # settling velocity
- ws = (nu/d50_vertical)*((10.36**2 + 1.049 * ds**3)**0.5 - 10.36)
- # Critical bed shear stress
- theta_cr = 0.3/(1 + 1.2*ds) + 0.055*(1-np.exp(-0.02*ds))
- # van Karman constant
- kappa = 0.41
- # skin friction coefficient
- fc = 2*(kappa/(1+np.log(ks/30*h)))**2
- # Shields parameter
- theta = 0.5*fc*velocity_vert**2/((rho-1)*d50_vertical*g)
- #turbulant Schmidt number
- sigma = 1
- # skin friction coefficient for total shear velocity
- fc2 = 2*(kappa/(1+np.log(kst/30/h)))**2
- # Total shear velocity
- shear_vel = velocity_vert*np.sqrt(fc2/2)
+ if choice[2] == 1 and choice[3] == 1:
+ alpha_est = stats_vertical_profile['slope'][i]
+ lnCr_est = stats_vertical_profile['intercept'][i]
+ # if np.count_nonzero(~np.isnan(np.array(vv['D50']))) >= 3:
+ # d50_vertical = np.mean(vv['D50'])*1e-6
+ # d90_bed = 0.0006
+ # velocity_vert = vv['Mean_velocity_vertical'].iloc[0]
+ # h = vv['max_depth_vertical'].iloc[0]
+ # # roughness height
+ # kst = 0.05
+ # ks = 2* d90_bed
+ # # kinematic viscosity of water
+ # nu = 0.000001
+ # # acceleration of the gravity
+ # g = 9.81
+ # # relative density
+ # rho = 2.65
+ # # sedimentological diameter
+ # ds = d50_vertical*(g*(rho-1)/nu**2)**(1/3) # use mean d50 per vertical
+ # # settling velocity
+ # ws = (nu/d50_vertical)*((10.36**2 + 1.049 * ds**3)**0.5 - 10.36)
+ # # Critical bed shear stress
+ # theta_cr = 0.3/(1 + 1.2*ds) + 0.055*(1-np.exp(-0.02*ds))
+ # # van Karman constant
+ # kappa = 0.41
+ # # skin friction coefficient
+ # fc = 2*(kappa/(1+np.log(ks/30*h)))**2
+ # # Shields parameter
+ # theta = 0.5*fc*velocity_vert**2/((rho-1)*d50_vertical*g)
+ # #turbulant Schmidt number
+ # sigma = 1
+ # # skin friction coefficient for total shear velocity
+ # fc2 = 2*(kappa/(1+np.log(kst/30/h)))**2
+ # # Total shear velocity
+ # shear_vel = velocity_vert*np.sqrt(fc2/2)
- # Estimated a priori parameters
- alpha_est = -6*ws/(sigma*kappa*shear_vel*h)
- Cr = 1.5*1000*theta*np.exp(-0.2*ds)*np.exp(-4.5*theta_cr/theta) # Reference concentration Cr
- lnCr_est = np.log(Cr)
+ # # Estimated a priori parameters
+ # alpha_est = -6*ws/(sigma*kappa*shear_vel*h)
+ # Cr = 1.5*1000*theta*np.exp(-0.2*ds)*np.exp(-4.5*theta_cr/theta) # Reference concentration Cr
+ # lnCr_est = np.log(Cr)
+ # print(i, np.round(d50_vertical,5), np.round(ds,5), np.round(theta_cr,5),np.round(fc,5), np.round(theta,5))
- else:
- alpha_est = stats_vertical_profile['slope'][i]
- lnCr_est = stats_vertical_profile['intercept'][i]
+ # else:
+ # alpha_est = stats_vertical_profile['slope'][i]
+ # lnCr_est = stats_vertical_profile['intercept'][i]
else:
alpha_est = stats_vertical_profile['slope'][i]
lnCr_est = stats_vertical_profile['intercept'][i]
@@ -233,111 +236,7 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
path_info = pd.DataFrame([outpath, sampling_date, most_frequent_sampler])
path_info.to_csv(cwd + '\\Info_BaM.csv', sep = ';', index = False)
-
- #%% 2. Uncertainty due to the vertical integration using BaM! (Renard et al. 2006; Mansanarez et al. 2019)
- # # Estimate a priori parameters for the model
- # a_priori_est = []
- # for i in range(len(abscissa_values)):
- # vv = analysis_data[analysis_data['Abscissa'] == abscissa_values[i]]
- # # Measured parameters
- # if choice[2] == 1:
- # if np.count_nonzero(~np.isnan(np.array(vv['D50']))) >= 2:
- # d50_vertical = np.mean(vv['D50'])*1e-6
- # else:
- # d50_vertical = 0.00013
- # else:
- # d50_vertical = 0.00013
- # d90_bed = 0.0004
- # #u_d50_vertical = 0.2
- # # velocity per vertical
- # if choice[3] == 1:
- # velocity_vert = vv['Mean_velocity_vertical'].iloc[0]
- # else:
- # velocity_vert = 1
- # #u_velocity_vert = 0.05
- # # water depth
- # if choice[3] == 1:
- # h = vv['max_depth_vertical'].iloc[0]
- # else:
- # h = 1
- # u_h = 0.05
- # # roughness height
- # kst = 0.05
- # #u_kst = 0.5*kst
- # ks = 2* d90_bed
- # #u_ks = 0.3 *ks
-
-
- # # Fixed parameters (errors u are ±)
- # # kinematic viscosity of water
- # nu = 0.000001
- # #u_nu = 0.05
- # # acceleration of the gravity
- # g = 9.81
- # #u_g = 0.01
- # # relative density
- # rho = 2.65
- # #u_rho = 0.02
- # # sedimentological diameter
- # ds = d50_vertical*(g*(rho-1)/nu**2)**(1/3) # use mean d50 per vertical
- # #u_ds = ds*(1/9*((u_g/g)**2 + (u_rho/rho)**2+4*(u_nu/nu)**2)+ (u_d50_vertical/d50_vertical)**2)
- # # settling velocity
- # ws = (nu/d50_vertical)*((10.36**2 + 1.049 * ds**3)**0.5 - 10.36)
- # u_ws = 0.05
- # # Critical bed shear stress
- # theta_cr = 0.3/(1 + 1.2*ds) + 0.055*(1-np.exp(-0.02*ds))
- # #u_theta_cr = 0.1
- # # van Karman constant
- # kappa = 0.41
- # u_kappa = 0
- # # skin friction coefficient
- # fc = 2*(kappa/(1+np.log(ks/30*h)))**2
- # #u_fc = 0.05
- # # Shields parameter
- # theta = 0.5*fc*velocity_vert**2/((rho-1)*d50_vertical*g)
- # #u_theta = 0.05
- # #turbulant Schmidt number
- # sigma = 1
- # u_sigma = 0.2
- # # skin friction coefficient for total shear velocity
- # fc2 = 2*(kappa/(1+np.log(kst/30/h)))**2
- # #u_fc2 = 0.05 *fc2
- # # Reference concentration Cr
- # Cr = 1.5*1000*theta*np.exp(-0.2*ds)*np.exp(-4.5*theta_cr/theta)
- # # u_Cr = 0.3 * 1000 * Cr
-
- # # Total shear velocity
- # shear_vel = velocity_vert*np.sqrt(fc2/2)
- # u_shear_vel = 0.05
-
- # # Estimated a priori parameters
- # alpha_est = -6*ws/(sigma*kappa*shear_vel*h)
- # #u_alpha_est = np.sqrt((u_ws**2/ws**2 + u_sigma**2/sigma**2 + u_kappa**2/kappa**2 + u_shear_vel**2/shear_vel**2 + u_h**2/h**2)*alpha_est**2)
- # u_alpha_est = np.sqrt(u_ws**2 + u_sigma**2 + u_kappa**2 + u_shear_vel**2 + u_h**2)
-
- # lnCr_est = np.log(Cr)
- # u_lnCr_est = 0.5
-
- # # Prepare export data
- # a_pri_est = pd.DataFrame([vv['Abscissa_Name'].iloc[0], alpha_est, u_alpha_est, lnCr_est, u_lnCr_est]).transpose()
- # a_priori_est.append(a_pri_est)
-
- # # Export a priori parameters for RBaM simulation
- # a_priori_est = pd.concat(a_priori_est)
- # a_priori_est = a_priori_est.reset_index(drop = True)
- # a_priori_inf = pd.Series([0]*len(a_priori_est))
- # a_priori_inf.loc[0] = sampling_date
- # a_priori_inf.loc[1] = most_frequent_sampler
- # # a_priori_inf.loc[2] = 2 # std for lnCr
- # # a_priori_inf.loc[3] = 2 # std for alpha
- # a_priori_est = pd.concat([(a_priori_inf), a_priori_est], axis = 1)
- # a_priori_est.columns = ['Info','Abscissa_Name', 'alpha_est', 'u_alpha_est', 'lnCr_est', 'u_lnCr_est']
- # a_priori_est.to_csv(outpath + sampling_date + '_Estimated_a_priori_params_RBaM_' + most_frequent_sampler +'.csv', sep = ';', index=False)
-
- # path_info = pd.DataFrame([outpath, sampling_date, most_frequent_sampler])
- # path_info.to_csv(cwd + '\\Info_BaM.csv', sep = ';', index = False)
-
#%% Perform BaM simulations in R
#### Run code in Rstudio
@@ -423,17 +322,7 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
# s = np.random.normal(mu, sigma, len(sampling_heights_interp))
s = np.random.normal(mu, sigma, len(results_cooking))
unc_sample.append(s)
- unc_sample = pd.DataFrame(unc_sample).transpose()
-
-
- # ##
- # unc_sample = []
- # for i in range(len(results_cooking)):
- # mu, sigma = 0, vv['ulnCs'][2]
- # #mu, sigma = 0, results_cooking.iloc[i,2]
- # s = np.random.normal(mu, sigma, len(sampling_heights_interp))
- # unc_sample.append(s)
- # unc_sample = pd.DataFrame(unc_sample)
+ unc_sample = pd.DataFrame(unc_sample).transpose()
# Calculate point concentration with error
unc_sample = conc_sampl * unc_sample ##
@@ -484,6 +373,8 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
# u_rem - remnant uncertainty due to the vertical integration
u_rem = [np.sqrt(u_p[i]**2 - u_param[i]**2)
for i in range(len(u_param))]
+ import math
+ u_rem = [0 if math.isnan(x) else x for x in u_rem]
u2_rem_final = np.sum([Q_j[j]**2/Q_total**2 *u_rem[j]**2
for j in range(len(abscissa_values))])
U2_rem = u2_rem_final*4
@@ -497,6 +388,8 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
maxpost = pd.DataFrame(maxpost).transpose()
maxpost.columns = abscissa_values
+ # maxpost = maxpost.iloc[::-1]
+ maxpost.reset_index(drop = True, inplace = True)
#%% Plot profiles for param
# plots only the last vertical
@@ -548,17 +441,20 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
elinewidth = 1, capsize = 2,
marker = 'D', color = 'black', markersize = 7, lw = 0,
label = 'Samples', zorder = 10)
- # add maxpost
- idx_max_post = np.argmax(results_cooking['logPost'])
- max_post = [np.exp(results_cooking['lnCr'][idx_max_post]-results_cooking['alpha'][idx_max_post]*h[i])
- for i in range(len(h))]
- ax.plot(max_post,h,
+ # # add maxpost
+ # idx_max_post = np.argmax(results_cooking['logPost'])
+ # max_post = [np.exp(results_cooking['lnCr'][idx_max_post]-results_cooking['alpha'][idx_max_post]*h[i])
+ # for i in range(len(h))]
+ # ax.plot(max_post,h,
+ # color= 'darkorange', lw = 2, label = r'$\mathregular{C_{n_0}(z)}$')
+
+ ax.plot(maxpost.iloc[:,l],h,
color= 'darkorange', lw = 2, label = r'$\mathregular{C_{n_0}(z)}$')
ax.legend(loc = 'upper right', fontsize = 18)
ax.set_xlabel(r'C (mg l$^{\rm{-1}}$)', fontsize = 20, weight = 'bold')
ax.set_ylabel('z (m)', fontsize = 20, weight = 'bold')
- ax.set_xlim(0,200)
+ ax.set_xlim(0,2000)
ax.set_ylim(0,h[-1])
ax.tick_params(axis='both', which='major', labelsize = 18)
@@ -567,50 +463,49 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
-
- # Plot overline(C(z))
- fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
-
- ax.hist(mean_c, bins = 50, color= 'teal')
- ax.vlines(mean_conc, 0,50,
- color = 'darkred', lw = 2, label = r'$\mathregular{\overline{C}}$')
-
- ax.legend(loc = 'upper right', fontsize = 18)
- ax.set_xlabel(r'$\mathregular{\overline{C}}$ (mg/l)', fontsize = 20, weight = 'bold')
- ax.set_ylabel('Frequency', fontsize = 20, weight = 'bold')
- ax.tick_params(axis='both', which='major', labelsize = 18)
- ax.set_xlim(0,200)
-
- fig.tight_layout()
- figname = '_hist_overline_C_vertical_' + str(l+1) + '_'
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
-
-
- # Plot overline(lnC(z)) # c)
- fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
-
- ax.hist(mean_c_log, bins = 50, alpha = 0.8, color= 'teal')
- ax.vlines(mean_conc_log, 0,150,
- color = 'darkred', lw = 2, label = r'$\mathregular{\overline{ln(\overline{C})}}$')
- ax.vlines(mean_conc_log+ std_param, 0,150,
- color = 'indigo', lw = 3,ls = '--', label = r'$\mathregular{\overline{ln(\overline{C})} \pm u_{p,param}}$')
- ax.vlines(mean_conc_log-std_param, 0,150,
- color = 'indigo', lw = 3,ls = '--')
-
- ax.legend(loc = 'upper right', fontsize = 18, framealpha = 1)
- ax.set_xlabel(r'$\mathregular{ln(C)}$ (mg l$^{\rm{-1}}$)', fontsize = 20, weight = 'bold')
- ax.set_ylabel('Frequency', fontsize = 20, weight = 'bold')
- ax.tick_params(axis='both', which='major', labelsize = 18)
- ax.set_xlim(1,5)
- ax.set_ylim(0,38)
+ if add_plots == True:
+ # Plot overline(C(z))
+ fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
- fig.tight_layout()
- figname = '_hist_ln_overline_C_vertical_' + str(l+1) + '_'
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
-
-
+ ax.hist(mean_c, bins = 50, color= 'teal')
+ ax.vlines(mean_conc, 0,50,
+ color = 'darkred', lw = 2, label = r'$\mathregular{\overline{C}}$')
+
+ ax.legend(loc = 'upper right', fontsize = 18)
+ ax.set_xlabel(r'$\mathregular{\overline{C}}$ (mg/l)', fontsize = 20, weight = 'bold')
+ ax.set_ylabel('Frequency', fontsize = 20, weight = 'bold')
+ ax.tick_params(axis='both', which='major', labelsize = 18)
+ ax.set_xlim(0,200)
+
+ fig.tight_layout()
+ figname = '_hist_overline_C_vertical_' + str(l+1) + '_'
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
+
+
+ # Plot overline(lnC(z)) # c)
+ fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
+
+ ax.hist(mean_c_log, bins = 50, alpha = 0.8, color= 'teal')
+ ax.vlines(mean_conc_log, 0,150,
+ color = 'darkred', lw = 2, label = r'$\mathregular{\overline{ln(\overline{C})}}$')
+ ax.vlines(mean_conc_log+ std_param, 0,150,
+ color = 'indigo', lw = 3,ls = '--', label = r'$\mathregular{\overline{ln(\overline{C})} \pm u_{p,param}}$')
+ ax.vlines(mean_conc_log-std_param, 0,150,
+ color = 'indigo', lw = 3,ls = '--')
+
+ ax.legend(loc = 'upper right', fontsize = 18, framealpha = 1)
+ ax.set_xlabel(r'$\mathregular{ln(C)}$ (mg l$^{\rm{-1}}$)', fontsize = 20, weight = 'bold')
+ ax.set_ylabel('Frequency', fontsize = 20, weight = 'bold')
+ ax.tick_params(axis='both', which='major', labelsize = 18)
+ ax.set_xlim(1,5)
+ ax.set_ylim(0,38)
+
+ fig.tight_layout()
+ figname = '_hist_ln_overline_C_vertical_' + str(l+1) + '_'
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
+
#%% Plot profiles for struc
# plots only the last vertical
@@ -653,7 +548,7 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
- #%% Plot mean C natural (depth-integrated)
+ # Plot mean C natural (depth-integrated)
fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
ax.vlines(mean_c_rem[0], 0,h[-1], color= 'darkgoldenrod', lw =0.5, alpha = 0.3,
@@ -710,79 +605,77 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
+ if add_plots == True:
+ # Plot mean lnC (depth-integrated)
+ fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
+
+ ax.vlines(mean_c_rem_log[0], 0,h[-1], color= 'teal',lw = 1, alpha = 0.3,
+ label = r'$\mathregular{ln(\overline{C_{mod,n}})}$')
+ for i in range(len(mean_c_rem_log)):
+ ax.vlines(mean_c_rem_log[i], 0,h[-1], color= 'teal', lw = 1, alpha = 0.3)
+ ax.vlines(mean_conc_rem_log, 0,h[-1],
+ color = 'midnightblue', lw = 2, label = r'$\mathregular{ln(\overline{C_{mod}})}$')
+ ax.vlines(mean_conc_rem_log-std_conc_rem, 0,h[-1],
+ color = 'mediumorchid', lw = 2, ls = '--', label = r'$\mathregular{ln(\overline{C_{mod}}) \pm u_p}$')
+ ax.vlines(mean_conc_rem_log+std_conc_rem, 0,h[-1],
+ color = 'mediumorchid', lw = 2, ls = '--')
+
+ ax.legend(loc = 'upper right', fontsize = 18, framealpha = 1)
+ ax.set_xlabel(r'$\mathregular{ln(C)}$ (mg/l)', fontsize = 20, weight = 'bold')
+ ax.set_ylabel('z (m)', fontsize = 20, weight = 'bold')
+ ax.tick_params(axis='both', which='major', labelsize = 18)
+ ax.set_ylim(0,h[-1])
+ ax.set_xlim(1,5)
+
+ fig.tight_layout()
+ figname = '_ln_overline_C_mod_vertical_' + str(l+1) + '_'
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
- # Plot mean lnC (depth-integrated)
- fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
-
- ax.vlines(mean_c_rem_log[0], 0,h[-1], color= 'teal',lw = 1, alpha = 0.3,
- label = r'$\mathregular{ln(\overline{C_{mod,n}})}$')
- for i in range(len(mean_c_rem_log)):
- ax.vlines(mean_c_rem_log[i], 0,h[-1], color= 'teal', lw = 1, alpha = 0.3)
- ax.vlines(mean_conc_rem_log, 0,h[-1],
- color = 'midnightblue', lw = 2, label = r'$\mathregular{ln(\overline{C_{mod}})}$')
- ax.vlines(mean_conc_rem_log-std_conc_rem, 0,h[-1],
- color = 'mediumorchid', lw = 2, ls = '--', label = r'$\mathregular{ln(\overline{C_{mod}}) \pm u_p}$')
- ax.vlines(mean_conc_rem_log+std_conc_rem, 0,h[-1],
- color = 'mediumorchid', lw = 2, ls = '--')
-
- ax.legend(loc = 'upper right', fontsize = 18, framealpha = 1)
- ax.set_xlabel(r'$\mathregular{ln(C)}$ (mg/l)', fontsize = 20, weight = 'bold')
- ax.set_ylabel('z (m)', fontsize = 20, weight = 'bold')
- ax.tick_params(axis='both', which='major', labelsize = 18)
- ax.set_ylim(0,h[-1])
- ax.set_xlim(1,5)
-
- fig.tight_layout()
- figname = '_ln_overline_C_mod_vertical_' + str(l+1) + '_'
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
-
-
- # Plot overline(C(z))
- fig, ax = plt.subplots(figsize=(6, 6), dpi = 100)
-
- ax.hist(mean_c, bins = 50, color= 'teal')
- ax.vlines(mean_conc, 0,170,
- color = 'darkred', lw = 2, label = r'$\mathregular{\overline{C_{mod}}}$')
-
- ax.legend(loc = 'upper right', fontsize = 18)
- ax.set_xlabel(r'C (mg l$^{\rm{-1}}$)', fontsize = 20, weight = 'bold')
- ax.set_ylabel('Frequency', fontsize = 20, weight = 'bold')
- ax.tick_params(axis='both', which='major', labelsize = 18)
- ax.set_xlim(0,200)
- ax.set_ylim(0,50)
-
- fig.tight_layout()
- figname = '_hist_overline_C_mod_vertical_' + str(l+1) + '_'
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
-
-
- # Plot hist ln(overline(C(z))) # f)
- fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
-
- ax.hist(mean_c_log, bins = 50, color= 'peru', alpha = 0.8)
- ax.vlines(mean_conc_log, 0,150,
- color = 'firebrick', lw = 2, label = r'$\mathregular{\overline{ln(\overline{C_{mod}})}}$')
- ax.vlines(mean_conc_log+ std_conc_rem, 0,150,
- color = 'mediumorchid', lw = 3,ls = '--', label = r'$\mathregular{\overline{ln(\overline{C_{mod}})} \pm u_p}$')
- ax.vlines(mean_conc_log-std_conc_rem, 0,150,
- color = 'mediumorchid', lw = 3,ls = '--')
-
- ax.legend(loc = 'upper right', fontsize = 18, framealpha = 1)
- ax.set_xlabel(r'$\mathregular{ln(C)}$ (mg l$^{\rm{-1}}$)', fontsize = 20, weight = 'bold')
- ax.set_ylabel('Frequency', fontsize = 20, weight = 'bold')
- ax.tick_params(axis='both', which='major', labelsize = 18)
- ax.set_xlim(1,5)
- ax.set_ylim(0,37)
-
- fig.tight_layout()
- figname = '_hist_ln_overline_C_mod_vertical_' + str(l+1) + '_'
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
- fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
-
-
+
+ # Plot overline(C(z))
+ fig, ax = plt.subplots(figsize=(6, 6), dpi = 100)
+ ax.hist(mean_c, bins = 50, color= 'teal')
+ ax.vlines(mean_conc, 0,170,
+ color = 'darkred', lw = 2, label = r'$\mathregular{\overline{C_{mod}}}$')
+
+ ax.legend(loc = 'upper right', fontsize = 18)
+ ax.set_xlabel(r'C (mg l$^{\rm{-1}}$)', fontsize = 20, weight = 'bold')
+ ax.set_ylabel('Frequency', fontsize = 20, weight = 'bold')
+ ax.tick_params(axis='both', which='major', labelsize = 18)
+ ax.set_xlim(0,200)
+ ax.set_ylim(0,50)
+
+ fig.tight_layout()
+ figname = '_hist_overline_C_mod_vertical_' + str(l+1) + '_'
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
+
+
+ # Plot hist ln(overline(C(z))) # f)
+ fig, ax = plt.subplots(figsize=(6, 8), dpi = 100)
+
+ ax.hist(mean_c_log, bins = 50, color= 'peru', alpha = 0.8)
+ ax.vlines(mean_conc_log, 0,150,
+ color = 'firebrick', lw = 2, label = r'$\mathregular{\overline{ln(\overline{C_{mod}})}}$')
+ ax.vlines(mean_conc_log+ std_conc_rem, 0,150,
+ color = 'mediumorchid', lw = 3,ls = '--', label = r'$\mathregular{\overline{ln(\overline{C_{mod}})} \pm u_p}$')
+ ax.vlines(mean_conc_log-std_conc_rem, 0,150,
+ color = 'mediumorchid', lw = 3,ls = '--')
+
+ ax.legend(loc = 'upper right', fontsize = 18, framealpha = 1)
+ ax.set_xlabel(r'$\mathregular{ln(C)}$ (mg l$^{\rm{-1}}$)', fontsize = 20, weight = 'bold')
+ ax.set_ylabel('Frequency', fontsize = 20, weight = 'bold')
+ ax.tick_params(axis='both', which='major', labelsize = 18)
+ ax.set_xlim(1,5)
+ ax.set_ylim(0,37)
+
+ fig.tight_layout()
+ figname = '_hist_ln_overline_C_mod_vertical_' + str(l+1) + '_'
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 300, bbox_inches='tight')
+ fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.eps', dpi = 300, bbox_inches='tight')
+
#%% 3. Uncertainty u_m due to lateral integration
# Based on nomograph (Guy & Norman 1970, based on Hubbell 1960)
@@ -895,9 +788,7 @@ def uncertainty_analysis(analysis_data, outpath, outpath_figures, sampling_date,
figname = '_Contribution_to_U_F_article'
fig.savefig(outpath_figures + '\\' + str(sampling_date) + figname + most_frequent_sampler + '.png', dpi = 400, bbox_inches='tight')
- #%% Export data
- maxpost = maxpost.iloc[::-1]
- maxpost.reset_index(drop = True, inplace = True)
+ #%% Export data
uncertainty = pd.DataFrame([u_sys, u_m, np.sqrt(u2_p_final),np.sqrt(U2_param), np.sqrt(U2_rem), U_C_perc, U_Q, U_F]).transpose()
uncertainty.columns = ['u\'_sys', 'u\'_m', 'u\'_p', 'u\'_p,param','u\'_p,struc', 'U\'_C', 'U\'_Q', 'U\'F']
diff --git a/classes/__pycache__/SDC_fines_method.cpython-38.pyc b/classes/__pycache__/SDC_fines_method.cpython-38.pyc
index 5bac8b637606f61d3e9cd5bd8bacff5e10710f44..80a58e1e2853650abf28a69cad2bb56eb3e53327 100644
Binary files a/classes/__pycache__/SDC_fines_method.cpython-38.pyc and b/classes/__pycache__/SDC_fines_method.cpython-38.pyc differ
diff --git a/classes/__pycache__/SDC_method.cpython-38.pyc b/classes/__pycache__/SDC_method.cpython-38.pyc
index 88350a5e69c7e26f5c1a2845ba65b6656994b6f7..e21284ed46ad6d8d6d5e38799777e8bb6debdccd 100644
Binary files a/classes/__pycache__/SDC_method.cpython-38.pyc and b/classes/__pycache__/SDC_method.cpython-38.pyc differ
diff --git a/classes/__pycache__/Sed_flux_big_bend.cpython-38.pyc b/classes/__pycache__/Sed_flux_big_bend.cpython-38.pyc
index 5d43f5b1c76eb827b2e3056999c6ad1f2a63fb38..349b62f06455ff9a5dd8d3ddacd254475dea8b71 100644
Binary files a/classes/__pycache__/Sed_flux_big_bend.cpython-38.pyc and b/classes/__pycache__/Sed_flux_big_bend.cpython-38.pyc differ
diff --git a/classes/__pycache__/Sed_flux_big_straight.cpython-38.pyc b/classes/__pycache__/Sed_flux_big_straight.cpython-38.pyc
index d0996c9e2e14e7d1e851ffbcc439097d75c6f1ad..a732de2f7bbee1ac7ba6fd66e7caac1a5eb158d1 100644
Binary files a/classes/__pycache__/Sed_flux_big_straight.cpython-38.pyc and b/classes/__pycache__/Sed_flux_big_straight.cpython-38.pyc differ
diff --git a/classes/__pycache__/Sed_flux_small_bend.cpython-38.pyc b/classes/__pycache__/Sed_flux_small_bend.cpython-38.pyc
index e21ab33c99001233a000905bd44d9d36852de4f3..7083782b3d1303d6d3318af406eb9e43e1cbfca8 100644
Binary files a/classes/__pycache__/Sed_flux_small_bend.cpython-38.pyc and b/classes/__pycache__/Sed_flux_small_bend.cpython-38.pyc differ
diff --git a/classes/__pycache__/Sed_flux_small_straight.cpython-38.pyc b/classes/__pycache__/Sed_flux_small_straight.cpython-38.pyc
index a7cdf704a6c62245f3c3043161970dc4425af8f6..fd75831099edbb3c8bbcf23a4f984cbd877fdf50 100644
Binary files a/classes/__pycache__/Sed_flux_small_straight.cpython-38.pyc and b/classes/__pycache__/Sed_flux_small_straight.cpython-38.pyc differ
diff --git a/classes/__pycache__/Sediment_flux_BD_single.cpython-38.pyc b/classes/__pycache__/Sediment_flux_BD_single.cpython-38.pyc
index 2bf8666b0d6313ba8499e4f30e6f8fec227cfa3d..4dd3580ea4f24ed20c83f5db5defc6376c03d905 100644
Binary files a/classes/__pycache__/Sediment_flux_BD_single.cpython-38.pyc and b/classes/__pycache__/Sediment_flux_BD_single.cpython-38.pyc differ
diff --git a/classes/__pycache__/Uncertainty_estimation.cpython-38.pyc b/classes/__pycache__/Uncertainty_estimation.cpython-38.pyc
index d28a92c22b93675b2f981ff17ee60a4e38eebbe5..38ba114b8e76cf1fd8106e8c7ca91e5954b8f5e0 100644
Binary files a/classes/__pycache__/Uncertainty_estimation.cpython-38.pyc and b/classes/__pycache__/Uncertainty_estimation.cpython-38.pyc differ
diff --git a/path_files.txt b/path_files.txt
index 533de3e668237ac8b88a08c88fd465b082316dbf..5eeb52337555f22bcc1743eef500adbc1085c449 100644
--- a/path_files.txt
+++ b/path_files.txt
@@ -1 +1 @@
-C:\Users\jessica.laible\Documents\Mesures\Manips&Sites\Grenoble_Campus\20230510_Campus
\ No newline at end of file
+C:\Users\jessica.laible\Documents\Mesures\Manips&Sites\Grenoble_Campus\20221129_Campus
\ No newline at end of file