diff --git a/code/A_REVUB_initialise_minimum_example.py b/code/A_REVUB_initialise_minimum_example.py
index 983922f..fcb157b 100644
--- a/code/A_REVUB_initialise_minimum_example.py
+++ b/code/A_REVUB_initialise_minimum_example.py
@@ -208,6 +208,9 @@
 # [set by user] number of turbines per hydropower plant (used for turbine use statistics in script C)
 no_turbines = parameters_hydropower_values[np.where(parameters_hydropower_list == 'no_turbines', True, False)][0].astype(int)
 
+# [set by user] percentile value used to calculate exceedance probability of delivered power
+p_exceedance = parameters_hydropower_values[np.where(parameters_hydropower_list == 'p_exceedance', True, False)][0].astype(int)
+
 
 # %% pre.4) Time series
 
diff --git a/code/B_REVUB_main_code.py b/code/B_REVUB_main_code.py
index 16bb5a7..feb0888 100644
--- a/code/B_REVUB_main_code.py
+++ b/code/B_REVUB_main_code.py
@@ -352,6 +352,10 @@
 L_unmet_BAL_frac_bymonth = np.zeros(shape = (months_yr, len(simulation_years), HPP_number))
 L_unmet_STOR_frac_bymonth = np.zeros(shape = (months_yr, len(simulation_years), HPP_number))
 
+
+##### MISCELLANEOUS PARAMETERS #####
+
+
 # [precallocate] Failure rate (fraction of time) of operation in case of prolonged dry periods
 fraction_outage_CONV = np.zeros(shape = HPP_number)
 fraction_outage_BAL = np.zeros(shape = HPP_number)
@@ -362,6 +366,11 @@
 fraction_overflow_BAL = np.zeros(shape = HPP_number)
 fraction_overflow_STOR = np.zeros(shape = HPP_number)
 
+# [preallocate] Guaranteed power (MW) per hydropower plant
+P_CONV_total_guaranteed = np.full(HPP_number,np.nan)
+P_BAL_total_guaranteed = np.full(HPP_number,np.nan)
+P_STOR_total_guaranteed = np.full(HPP_number,np.nan)
+
 
 # %% REVUB.3) Classify HPPs
 
@@ -1930,7 +1939,8 @@
 
 # %% REVUB.5) Post-processing
 
-# [initialise] use STOR equal to BAL for reservoirs where STOR not modelled
+# [initialise] use STOR equal to BAL for reservoirs where STOR not modelled (except RoR component). 
+# This step has no physical meaning and is purely meant to avoid plotting errors in script C_multiple if users want to plot scenarios where certain dams have STOR scenarios and others only BAL.
 for HPP in range(HPP_number):
     if STOR_break[HPP] == 1:
         P_STOR_hydro_stable_hourly[:,:,HPP] = P_BAL_hydro_stable_hourly[:,:,HPP]
@@ -2044,7 +2054,14 @@
         
         # [calculate] ELCC by year in MWh/year (eq. S23)
         ELCC_STOR_yearly[y,HPP] = np.sum(L_followed_STOR_hourly[hrs_year,y,HPP])
-        
+    
+    
+    # [calculate] statistics of power generation under user-defined p_exceedance criterion
+    P_CONV_total_guaranteed[HPP] = np.nanpercentile(P_CONV_hydro_stable_hourly[:,:,HPP] + P_CONV_hydro_RoR_hourly[:,:,HPP], 100 - p_exceedance[HPP])
+    P_BAL_total_guaranteed[HPP] = np.nanpercentile(P_BAL_hydro_stable_hourly[:,:,HPP] + P_BAL_hydro_flexible_hourly[:,:,HPP] + P_BAL_hydro_RoR_hourly[:,:,HPP] + P_BAL_solar_hourly[:,:,HPP] + P_BAL_wind_hourly[:,:,HPP], 100 - p_exceedance[HPP])
+    
+    if STOR_break[HPP] == 0 and option_storage == 1:
+        P_STOR_total_guaranteed[HPP] = np.nanpercentile(P_STOR_hydro_stable_hourly[:,:,HPP] + P_STOR_hydro_flexible_hourly[:,:,HPP] + P_STOR_solar_hourly[:,:,HPP] + P_STOR_wind_hourly[:,:,HPP], 100 - p_exceedance[HPP])
     
 
 # [display] signal simulation end
diff --git a/data/parameters_simulation.xlsx b/data/parameters_simulation.xlsx
index b507522..31fad66 100644
Binary files a/data/parameters_simulation.xlsx and b/data/parameters_simulation.xlsx differ
diff --git a/manual/REVUB_manual.pdf b/manual/REVUB_manual.pdf
index 6db0f2d..0711bb6 100644
Binary files a/manual/REVUB_manual.pdf and b/manual/REVUB_manual.pdf differ