Post processing for bank position timeseries

Adding post processing to output csv file with time series of bank position for each scenario

Functions/BankPos.m

@@ -0,0 +1,27 @@
+function BankY = BankPos(Dist, Bed, Radius, WL)
+%BANKPOS   compute outside bank waters edge position
+if Radius > 0
+    % Left hand bend, right bank is outside bank
+    BankCell = find(Bed<WL, 1, 'last');
+    if BankCell < length(Dist)
+        BankY = Dist(BankCell) + ...
+                      (WL - Bed(BankCell)) * ...
+                      (Dist(BankCell+1) - Dist(BankCell)) / ...
+                      (Bed(BankCell+1) - Bed(BankCell));
+    else
+        BankY = Dist(end);
+    end
+else
+    % Right hand bend, left bank is outside bank
+    BankCell = find(Bed<WL, 1, 'first');
+    if BankCell > 1
+        BankY = Dist(BankCell) - ...
+                      (WL - Bed(BankCell)) * ...
+                      (Dist(BankCell) - Dist(BankCell-1)) / ...
+                      (Bed(BankCell-1) - Bed(BankCell));
+    else
+        BankY = Dist(1);
+    end
+end    
+
+end

Functions/BankPosError.m

@@ -45,32 +45,4 @@
     ErrorSign = -sign(BankError);
 end
 AbsBankError = abs(BankError);
-end
-
-function BankY = BankPos(Dist, Bed, Radius, WL)
-% Private function to compute outside bank waters edge position
-if Radius > 0
-    % Left hand bend, right bank is outside bank
-    BankCell = find(Bed<WL, 1, 'last');
-    if BankCell < length(Dist)
-        BankY = Dist(BankCell) + ...
-                      (WL - Bed(BankCell)) * ...
-                      (Dist(BankCell+1) - Dist(BankCell)) / ...
-                      (Bed(BankCell+1) - Bed(BankCell));
-    else
-        BankY = Dist(end);
-    end
-else
-    % Right hand bend, left bank is outside bank
-    BankCell = find(Bed<WL, 1, 'first');
-    if BankCell > 1
-        BankY = Dist(BankCell) - ...
-                      (WL - Bed(BankCell)) * ...
-                      (Dist(BankCell) - Dist(BankCell-1)) / ...
-                      (Bed(BankCell-1) - Bed(BankCell));
-    else
-        BankY = Dist(1);
-    end
-end    
-
 end

Functions/BankPosTS.m

@@ -0,0 +1,24 @@
+function [Times,BankCoord] = BankPosTS(SnapshotFolder, WL, Radius)
+%BANKPOSTS   Extract timeseries of bank position from XChannel snapshots
+
+% get list of snapshot files
+Snapshots = dir(fullfile(SnapshotFolder,'*.out'));
+Snapshots = struct2table(Snapshots);
+Snapshots = Snapshots.name;
+NoOfFiles = size(Snapshots,1);
+
+% extract times from filename
+Times = regexp(Snapshots,'\d+','match','once');
+Times = str2double(Times);
+
+% Loop through files calculating bank position
+BankCoord = nan(NoOfFiles,1);
+for FileNo = 1:NoOfFiles
+    Geometry = csvread(fullfile(SnapshotFolder,Snapshots{FileNo}));
+    BankCoord(FileNo) = BankPos(Geometry(:,1), Geometry(:,2), Radius, WL);
+end
+
+[Times,SortIndex] = sort(Times);
+BankCoord = BankCoord(SortIndex);
+end
+

OutputBankPosTS.m

@@ -0,0 +1,27 @@
+% Post processing routine to loop through all results in outputs folder and
+% output time series of bank position.
+
+addpath Functions
+
+OutputsFolder = 'Outputs\';
+WL = 212.5;
+Radius = 1;
+
+Folders = struct2table(dir(OutputsFolder));
+Folders = Folders.name(Folders.isdir);
+Folders(ismember(Folders,{'.','..'})) = [];
+Scenario = str2double(regexp(Folders,'\d+','match','once'));
+[~,SortIndex] = sort(Scenario);
+Folders = Folders(SortIndex);
+
+for FolderNo = 1:size(Folders)
+    [BankPosition.ModelTime,BankPosition.(Folders{FolderNo,1})] = ...
+        BankPosTS(fullfile(OutputsFolder,Folders{FolderNo}, ...
+                           ['snapshots_',Folders{FolderNo}]), WL, Radius);
+    BankPosition.(Folders{FolderNo,1}) = ...
+        BankPosition.(Folders{FolderNo,1}) - ...
+        BankPosition.(Folders{FolderNo,1})(1);
+end
+
+BankPosition = struct2table(BankPosition);
+writetable(BankPosition,'Outputs\BankPosition.csv')