Updated processing script to track method of road drowning. Added new…

… script to track elevation change in B3D domains

data/Ocracoke_init_data/Ocracoke-CASCADE-parameters.yaml

@@ -46,15 +46,15 @@ Shrub_ON: 0
 SprayDist: 170
 StormSeries: []
 StormStart: 2
-TMAX: 23
+TMAX: 126
 TideAmp: 1.2
 TimeFruit: 5.0
 UprootLimit: -0.2
 beta: 0.04
 disp_mu: -0.721891
 disp_sigma: 1.5
 dune_file: C:\Users\frank\PycharmProjects\CASCADE\data\Ocracoke_init_data\dunes\Sample_1_dune.npy
-elevation_file: C:\Users\frank\PycharmProjects\CASCADE\data\Ocracoke_init_data\elevations\Sample_1_topography.npy
+elevation_file: C:\Users\frank\PycharmProjects\CASCADE\data\Ocracoke_init_data\elevations\Topography_2019\Domain_49_topography_2019.npy
 enable_sandbags: false
 growth_param_file: growthparam_1000dam.npy
 k_sf: 19324.753176079226
@@ -65,4 +65,4 @@ rmin: 0.55
 s_sf_eq: 0.01698590034494604
 sandbag_elevation: 1.8
 sandbag_need: false
-storm_file: C:\Users\frank\PycharmProjects\CASCADE\data\Ocracoke_init_data\storms\Hindcast_Storms\OCR_1974_2022_Final_Hindcast_Storms.npy
+storm_file: C:\Users\frank\PycharmProjects\CASCADE\data\Ocracoke_init_data\storms\Synthetic_Storms\Ten_Percent_Storms\OCR_Future_StormList_99_10.npy

scripts/ocracoke_ms/Future_Runs_Test.py

@@ -5,18 +5,18 @@
 from scipy import stats as st
 import os
 
-Change_Rates = np.loadtxt('C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Revised_Offshore_Datum\\All_Shoreline_Change_Rates.csv',skiprows=1,delimiter=',')
-Subset_Change_Rates = np.loadtxt('C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Revised_Offshore_Datum\\All_Annual_Change_Rates.csv',skiprows=1,delimiter=',')
-
 #os.chdir('C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Model Runs\\Future Runs')
 os.chdir('E:\\Model_Runs')
 
 
 Save_Path = 'C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Model Runs\\Summary_Values\\'
 
-Management_Name = ['Nourishment']#['Status_Quo','Natural', 'Nourishment']
-RSLR_Rate = ['IL','I','IH']
+Management_Name = ['Status_Quo','Natural', 'Nourishment']
+#RSLR_Rate = ['IL','I','IH']
+RSLR_Rate = ['IL','IH']
 Sink_Name = ['Erosional_Sink','Accretional_Sink']
+#Sink_Name = ['Accretional_Sink']
+
 Storm_Level = 'Baseline'
 
 Base_Name_List = []
@@ -49,6 +49,7 @@ def Process_Batch(Base_Name,
     buffer_length = 15
     All_EP_Change = []
     All_Roadway_Abandonment = []
+    All_Abandonment_Reason = []
     Island_Drowning = {}
     Years_Modeled_List = []
     Relocation_TS = []
@@ -76,11 +77,12 @@ def Process_Batch(Base_Name,
         All_EP_Change.append(copy.deepcopy(shoreline_change))
         All_Shoreline_Positions.append(copy.deepcopy(all_shoreline_positions))
         # Calculate Roadway Abandonment metrics
-        roadway_abandonment = Calculate_Roadway_Abandonmet(cascade=Cascade_List,
+        roadway_abandonment,abandonmet_reason = Calculate_Roadway_Abandonmet(cascade=Cascade_List,
                                                            years_modeled=Model_Run_Year,
                                                            buffer_length=buffer_length,
                                                            number_barrier3d_models = number_barrier3d_models)
         All_Roadway_Abandonment.append(copy.deepcopy(roadway_abandonment))
+        All_Abandonment_Reason.append(copy.deepcopy(abandonmet_reason))
         # Calculate roadway relocation
         roadway_relocation, relocation_frequency = Calculate_Roadway_Relocation(cascade=Cascade_List,
                                                                                 years_modeled=Model_Run_Year,
@@ -148,7 +150,8 @@ def Process_Batch(Base_Name,
                        'All_Shoreline_Positions':All_Shoreline_Positions,
                        'Total_Nourishment_Volume':Total_Volume_TS,
                        'All_Nourishment_TS':All_Nourishment_TS,
-                       'All_Overwash_volume':Total_Overwash_TS
+                       'All_Overwash_volume':Total_Overwash_TS,
+                       'Roadway_Abandonment_Reason':All_Abandonment_Reason,
                        }
 
     All_Values_Data_Frame = pd.DataFrame(All_Values_Dict)
@@ -291,6 +294,7 @@ def Calculate_Nourishment_Volume(cascade,years_modeled,buffer_length,number_barr
     Nourishment_Data = cascade.nourishments
     All_Nourishment_TS = []
     Total_Nourishment_Per_Grid = []
+    Nourishment_Event_Volume = cascade.nourishment_volume
     for l in range(buffer_length, (number_barrier3d_models - buffer_length)):
         All_Nourishment_TS.append(copy.deepcopy(Nourishment_Data[l].nourishment_volume_TS))
         Total_Nourishment_Per_Grid.append(np.sum(Nourishment_Data[l].nourishment_volume_TS))
@@ -342,15 +346,29 @@ def Calculate_Island_Interior_Width_Change(cascade, years_modeled, buffer_length
 def Calculate_Roadway_Abandonmet(cascade, years_modeled, buffer_length, number_barrier3d_models):
     # Find times the roadway broke and save the year that it did
     Road_Drowning_Years = []
+    Road_Drowning_Reason = []
     for m in range(buffer_length, (number_barrier3d_models - buffer_length - 1)):
         Road_Data = cascade.roadways[m]
+
+        # Record reason for road drowning
+        if Road_Data.relocation_break == 1:
+            Road_Drowning_Reason.append(copy.deepcopy('No relocation room'))
+        elif Road_Data.drown_break == True:
+            Road_Drowning_Reason.append(copy.deepcopy('Too much water'))
+        elif cascade.b3d_break == 1:
+            Road_Drowning_Reason.append(copy.deepcopy('Island Drowning'))
+        else:
+            Road_Drowning_Reason.append(copy.deepcopy('NA'))
+
+        # Record the year of roadway abandonment
         if Road_Data.relocation_break == 1 or Road_Data.drown_break == True:
             Road_Drowning_Years.append(copy.deepcopy(Road_Data.time_index))
         elif Road_Data.drown_break == int(0) and Road_Data.time_index == 1:
             Road_Drowning_Years.append(copy.deepcopy(1))
         else:
             Road_Drowning_Years.append(copy.deepcopy(years_modeled))
-    return(Road_Drowning_Years)
+        x = 10
+    return(Road_Drowning_Years,Road_Drowning_Reason)
 
 def Calculate_Roadway_Relocation(cascade, years_modeled,buffer_length, number_barrier3d_models):
     Relocations = []

scripts/ocracoke_ms/Hindcast_Test_Script.py

@@ -219,9 +219,9 @@
 plt.legend(loc='lower left')
 plt.title('Shoreline Change: 1974-1997')
 plt.ylabel('Shoreline Change Rate (m/yr)')
-plt.xlabel('Domain')
+plt.xlabel('Model Domain')
 plt.tight_layout()
-plt.savefig(fname=('C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Figures\\Hindcast_Shoreline_1974_1997.png'))
+plt.savefig(fname=('C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Figures\\Hindcast_Shoreline_1974_1997.eps'),format='eps')
 plt.show()
 
 # Plot all the different OCR storm intensities
@@ -238,11 +238,28 @@
 plt.legend(loc='lower left')
 plt.title('Shoreline Change: 1997-2020')
 plt.ylabel('Shoreline Change Rate (m/yr)')
-plt.xlabel('Domain')
+plt.xlabel('Model Domain')
 plt.tight_layout()
-plt.savefig(fname=('C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Figures\\Hindcast_Shoreline_1997_2020.png'))
+plt.savefig(fname=('C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Figures\\Hindcast_Shoreline_1997_2020.eps'),format='eps')
 plt.show()
 
+
+ax = plt.gca()
+ax.set_ylim([-10, 10])
+plt.plot(domain_nums, LRR_74_97, label = 'Calibration',color='#1f77b4')
+plt.plot(domain_nums, LRR_97_20, label = 'Test',color='#ff7f0e')
+plt.axhline(y = 0, color = 'k', linestyle = '--')
+plt.legend(loc='lower left')
+plt.title('Historic Shoreline Change: 1974-2020')
+plt.ylabel('Shoreline Change Rate (m/yr)')
+plt.xlabel('Model Domain')
+plt.tight_layout()
+plt.savefig(fname=('C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Figures\\Hindcast_Shoreline_1974_2020.eps'),format='eps')
+plt.show()
+
+
+
+
 print('Hello Purr')
 '''
 #plt.plot(domain_nums, LRR_97_20, label = 'Historic Change',color='grey')

scripts/ocracoke_ms/Test_Changes.py

@@ -19,7 +19,7 @@
 Save_Path = 'C:\\Users\\frank\\OneDrive - University of North Carolina at Chapel Hill\\Chapter 3\\Model Runs\\Summary_Values\\'
 
 run_name_batch = []
-run_name_batch.append('OCR_IH_Nourishment_S1_Erosional_Sink_Test_N_1_Cells')
+run_name_batch.append('OCR_IL_Nourishment_S1_10_Accretional_Sink')
 run_name_batch.append('OCR_IH_Nourishment_S1_Erosional_Sink_Test_N_2_Cells')
 run_name_batch.append('OCR_IH_Nourishment_S1_Erosional_Sink_Test_N_60P')
 run_name_batch.append('OCR_IH_Natural_S1_Erosional_Sink')

scripts/ocracoke_ms/Test_Elevation_Change.py

@@ -0,0 +1,53 @@
+#Test changes and graphs for 1 Cascade instance
+
+import copy
+
+import numpy as np
+import pandas as pd
+from scipy import stats as st
+import os
+
+os.chdir('E:\\Model_Runs')
+
+run_name ='OCR_IH_Status_Quo_S84_Accretional_Sink'
+
+output = np.load(run_name + ".npz", allow_pickle=True)
+cascade = output["cascade"]
+cascade = cascade[0]
+b3d = cascade.barrier3d
+ny = np.size(b3d)
+
+Buffer_Domains = 15
+
+Domains_of_Interest = range(Buffer_Domains,(len(b3d)-Buffer_Domains))
+
+
+for k in Domains_of_Interest:
+    Temp_B3D = b3d[k]
+    Shoreline_Change_TS = Temp_B3D._ShorelineChangeTS
+    #for k in range(0,Temp_B3D.time_index):
+    initial_elev = Temp_B3D.DomainTS[0]
+    final_elev = Temp_B3D.DomainTS[Temp_B3D.time_index-1]
+    distance_traveled = int(abs(np.sum(Shoreline_Change_TS[0:int(Temp_B3D.time_index)])))
+
+    blank_cells = np.full((distance_traveled,int(Temp_B3D.BarrierLength)),-0.3)
+
+    updated_final_array = np.concatenate((blank_cells,final_elev),axis=0)
+
+    if len(initial_elev) > len(updated_final_array):
+        print('Add Extra Cells to final array')
+        len_difference = len(initial_elev) - len(updated_final_array)
+        extra_cells = np.full((len_difference,int(Temp_B3D.BarrierLength)),-0.3)
+        updated_final_array = np.concatenate((updated_final_array,extra_cells),axis=0)
+    elif len(initial_elev) < len(updated_final_array):
+        print('Add extra cells to initial array')
+        len_difference = len(updated_final_array) - len(initial_elev)
+        extra_cells = np.full((len_difference,int(Temp_B3D.BarrierLength)),-0.3)
+        initial_elev = np.concatenate((initial_elev,extra_cells),axis=0)
+    else:
+        print('They have equal lengths')
+
+    elev_difference = np.subtract(updated_final_array,initial_elev)
+
+    c = 20
+x = 10

scripts/ocracoke_ms/plot_script.py

@@ -374,10 +374,10 @@ def plot_ElevAnimation_CASCADE(
 
 
 os.chdir('C:\\Users\\frank\\PycharmProjects\\CASCADE\\Run_output')
-run_name = "OCR_IH_Status_Quo_S29_Erosional_Sink"
+run_name = "OCR_IL_Natural_S42_10_Accretional_Sink"
 
 name_prefix = run_name
-nt_run = 67
+nt_run = 120
 number_barrier3d_models = 5
 Run_Marsh_Dynamics = False
 
@@ -408,7 +408,7 @@ def plot_ElevAnimation_CASCADE(
     beach_management_ny=beach_management_ny,
     roadway_management_ny=roadway_management_ny,
     y_lim=None,
-    z_lim=None,
+    z_lim=6,
     fig_size=None,
     fig_eps=False,
     Model_Grids_Of_Interest = Model_Grids_Of_Interest,

scripts/ocracoke_ms/run_forward_batch_forward_simulations.py

@@ -17,26 +17,26 @@
 start_year = 2024
 
 # Set scenario type
-status_quo = True
+status_quo = False
 nourishment = True
 
 if nourishment == True:
     Management_name = '_Nourishment_'
     nourish_beach = True
     nourishment_time = None
-    nourishment_volume = 120
+    nourishment_volume = 100
 elif status_quo == True:
     Management_name = "_Status_Quo_"
     nourishment_time = None
     nourish_beach = False
-    nourishment_volume = 120
+    nourishment_volume = 100
 else:
     Management_name = '_Natural_'
     nourishment_time = None
     nourish_beach = False
-    nourishment_volume = 120
+    nourishment_volume = 100
 
-Storms = 'Baseline'
+Storms = '10'
 
 # RSLR Data
 RSLR_Type = 'IH'
@@ -64,8 +64,11 @@
 
 Sink_Options = ['Accretional_Sink','Erosional_Sink']
 
+start_num = 0
+end_num = 100
+
 run_name = []
-for snames in range(0,100):
+for snames in range(start_num,end_num):
     if Storms == 'Baseline':
         name_base = 'OCR_'+str(RSLR_Type)+str(Management_name)+'S'+str(snames)
     else:
@@ -78,10 +81,10 @@
 
 s_file = []
 if Storms == 'Baseline':
-    for storm_num in range(0,100):
+    for storm_num in range(start_num,end_num):
         s_file.append(copy.deepcopy('C:\\Users\\frank\\PycharmProjects\\CASCADE\\data\\Ocracoke_init_data\\storms\\Synthetic_Storms\\OCR_Future_StormList_'+str(storm_num)+'_baseline.npy'))
 else:
-    for storm_num in range(0, 100):
+    for storm_num in range(start_num, end_num):
         s_file.append(copy.deepcopy(
             'C:\\Users\\frank\\PycharmProjects\\CASCADE\\data\\Ocracoke_init_data\\storms\\Synthetic_Storms\\Ten_Percent_Storms\\OCR_Future_StormList_' + str(
                 storm_num) + '_10.npy'))