-
Notifications
You must be signed in to change notification settings - Fork 14
/
Copy pathExtract_RAS.m
197 lines (169 loc) · 9.05 KB
/
Extract_RAS.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
%==============================================================================
% FluEgg -Fluvial Egg Drift Simulator
%==============================================================================
% Copyright (c) 2013 Tatiana Garcia
% This program is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License version 3 as published by
% the Free Software Foundation (currently at http://www.gnu.org/licenses/agpl.html)
% with a permitted obligation to maintain Appropriate Legal Notices.
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
% You should have received a copy of the GNU General Public License
% along with this program. If not, see <http://www.gnu.org/licenses/>.
%%%:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::%
%% Extract output from HEC-RAS %
%:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::%
%-------------------------------------------------------------------------%
% This function is used to extract data from a HEC-RAS project and produces %
% References: %
% Goodell, C.R. 2014. %
% Breaking the HEC-RAS Code: A User's Guide to Automating HEC-RAS. A User's
% Guide to Automating HEC-RAS. h2ls. Portland, OR. %
%-------------------------------------------------------------------------%
% %
%-------------------------------------------------------------------------%
% Created by : Tatiana Garcia %
% Date : March 29, 2016 %
% Last Modified : April 18, 2016 %
%-------------------------------------------------------------------------%
% Inputs:
% Outputs:
%
%:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::%
function [Riverinputfile]=Extract_RAS(strRASProject,handles,profile)
%% Creates a COM Server for the HEC-RAS Controller
% The command above depends on the version of HEC-RAS you have, in my case
% I am using version 5.0.3
try
RC = actxserver('RAS503.HECRASController');
catch
try %HECRAS 5.0.2
RC = actxserver('RAS502.HECRASController');
catch
try %HECRAS 5.0.1
RC = actxserver('RAS501.HECRASController');
catch
try %HECRAS 5.0.0
RC = actxserver('RAS500.HECRASController');
catch
end %HECRAS 5.0.0
end %HECRAS 5.0.1
end%HECRAS 5.0.2
end %HECRAS 5.0.3
%% Open the project
%strRASProject = 'D:\Asian Carp\Asian Carp_USGS_Project\Tributaries data\Sandusky River\SANDUSKY_Hec_RAS_mod\Sandusky_mod_II\BallvilleDam_Updated.prj';
RC.Project_Open(strRASProject); %open and show interface, no need to use RC.ShowRAS in Matlab
%% Define Variables
% get variables from GUI with user selection
lngRiverID = get(handles.popup_River,'Value'); % RiverID
lngReachID = get(handles.popup_Reach,'Value'); % ReachID
%%
if profile==-1
lngProfile = get(handles.popup_HECRAS_profile,'Value')-1; % Profile Number, 1 is all profiles
else
%The first profile in HEC-RAS corresponds to MAX WS, then the profile
% corresponding to initial conditions is:
lngProfile=profile+1;
end
lngUpDn = 0; % Up/Down index for nodes with multiple sections (only used for bridges)
% Output ID of Variables ( see page 247 in reference book for more details)
% lngWS_ID = 2; % The Water Surface Elevation ID is 2.
lngVelChnl_ID = 25; % The ID of the average velocity of flow for the main channel is 23.
lngHydrDepthC_ID = 128; % Hydraulic depth in channel.
lngQChannel_ID = 7; % Flow in the main channel.
lngHydrRadiusC_ID = 210; % ID for hydraulic radius in channel.
lngMannWtdChnl_ID = 45; % ID for Conveyance weighted Manning's n for the main channel.
lngLengthChnl_ID = 42; % ID for Downstream reach length of the main channel to next XS
%(unless BR in d/s, then this is the distance to the deck/roadway).
lngNum_XS = RC.Schematic_XSCount(); %Number of XS - HEC-RAS Controller will populate.
[~,~,lngNum_RS]=RC.Geometry_GetNodes(lngRiverID,lngReachID,0,0,0);%Number of nodes
% Preallocate memory
strRS = cell(lngNum_RS,1); %Array of names of the nodes-->River station name. See page 36 in book
strNodeType = strRS; %Pre-allocate array for node type
% Preallocate memory for output vectors
%sngWS = nan(lngNum_RS,1);
sngVelChnl = nan(lngNum_XS,1);
sngHydrDepthC = nan(lngNum_XS,1);
sngQChannel = nan(lngNum_XS,1);
sngHydrRadiusC = nan(lngNum_XS,1);
sngMannWtdChnl = nan(lngNum_XS,1);
sngLengthChnl = nan(lngNum_XS,1);
%% Run the current plan
%RC.Compute_CurrentPlan(0,0); %from book: = RC.Compute_CurrentPlan(lngMessages,strMessages(), True)
% RC.Compute_HideComputationWindow; %To hide Computation Window
% Uncoment the lines below for info about current plan and project
% Do not delete, we might need this in the future
%==========================================================================
%Current_Plan= RC.CurrentPlanFile()
%Current_Geometry_File= RC.CurrentGeomFile()
%[lngNumProf,strProfileName]=RC.Output_GetProfiles(0,0) %Profiles names and todatl number
%% Output Results
XC_counter=0;
% Extracts variable info from XS to Xs
for i=1:lngNum_RS
strNodeType{i}=RC.Geometry.NodeCType(lngRiverID,lngReachID,i);
if strcmp(strNodeType{i},'')% An empty strings (i.e. ' ') denotes a cross section;
XC_counter=XC_counter+1;
% Here we assing the river station name to each node
strRS{XC_counter} = RC.Geometry.NodeRS(lngRiverID,lngReachID,i);
% Extracts average velocity of flow for the main channel
sngVelChnl(XC_counter) = abs(RC.Output_NodeOutput(lngRiverID,lngReachID,i,lngUpDn,...
lngProfile,lngVelChnl_ID));
% Extracts hydraulic depth in channel
sngHydrDepthC(XC_counter) = RC.Output_NodeOutput(lngRiverID,lngReachID,i,lngUpDn,...
lngProfile,lngHydrDepthC_ID);
% Extracts flow in the main channel.
sngQChannel(XC_counter) = RC.Output_NodeOutput(lngRiverID,lngReachID,i,lngUpDn,...
lngProfile,lngQChannel_ID);
% Extracts hydraulic radius in channel.
sngHydrRadiusC(XC_counter) = RC.Output_NodeOutput(lngRiverID,lngReachID,i,lngUpDn,...
lngProfile,lngHydrRadiusC_ID);
% Extracts conveyance weighted Manning's n for the main channel.
sngMannWtdChnl(XC_counter) = RC.Output_NodeOutput(lngRiverID,lngReachID,i,lngUpDn,...
lngProfile,lngMannWtdChnl_ID);
% Extracts conveyance weighted Manning's n for the main channel.
sngLengthChnl(XC_counter) = RC.Output_NodeOutput(lngRiverID,lngReachID,i,lngUpDn,...
lngProfile,lngLengthChnl_ID);
% Extracts Water Surface Elevation at each XS
%sngWS(i) = RC.Output_NodeOutput(lngRiverID,lngReachID,i,lngUpDn,...
% lngProfile,lngWS_ID);
end%if is a cross section
end %for
sngLengthChnl(end)=0; %There is not data of length channel for the last cell
%% Generate Rivirinputfile dataset for the FluEgg model
Riverinputfile=Generate_Rivirinputfile();
try
%Kill Hec-ras
!taskkill /im ras.exe
RC.QuitRAS
catch
end
delete(RC);
%% :::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
function Riverinputfile=Generate_Rivirinputfile()
CumlDistance=[0;cumsum(sngLengthChnl(1:end-1))/1000];%In km
CumlDistance=[((CumlDistance(1:end-1)+CumlDistance(2:end))/2); CumlDistance(end)]; %FluEgg cells are located between 2 XS
CellNumber=(1:1:lngNum_XS)';
ks=(8.1.*sngMannWtdChnl.*sqrt(9.81)).^6;
Ustar=abs(sngVelChnl./(8.1*((sngHydrRadiusC./ks).^(1/6))));
%% Temperature choice
Temperature=str2double(get(handles.Const_Temp(2),'String'));
if isnan(Temperature) %If user didn't input temperature
ed = errordlg('Please input temperature','Error');
set(ed, 'WindowStyle', 'modal');
uiwait(ed);
return
end
Riverinputfile=[CellNumber CumlDistance sngHydrDepthC sngQChannel...
sngVelChnl zeros(lngNum_XS,1) zeros(lngNum_XS,1) Ustar ...
Temperature*ones(lngNum_XS,1)]; % default temperature 22C
end
%% :::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
end
%:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::%
%% <<<<<<<<<<<<<<<<<<<<<<<<< END OF FUNCTION >>>>>>>>>>>>>>>>>>>>>>>>>>>>%%
%:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::%
%% Notes:
% When there is not data the controller puts a very large number: 3.39999995214436e+38