- 积分
- 2040
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2017-3-30
- 最后登录
- 1970-1-1
|
楼主 |
发表于 2018-4-12 21:34:23
|
显示全部楼层
%%
%改
clear;clc;clf;
%试试读取数据并绘制KC流量
path='C:\Essay\Data\Essay\ofes_vvel';
file=dir(fullfile(path,'*.mat'));
L=0.1*2*3.1415926*6377.83*cosd(18)/360;% 北纬18度0.1个经度所代表的的横向距离,单位:千米
Q=zeros(1,780);
for d=1:length(file) %读取流速数据 以便计算流量
data=load(file(d).name);
depth=data.depth;
v=data.v;
vertical=v(181,61:151,1:35);
vv=permute(vertical,[3 2 1]);
%表层流量
S(1)=vv(1,1)*L/2*(depth(2)-depth(1))/2;
for h=2:90;
S(h)=vv(1,h)*L*(depth(2)-depth(1))/2;
s=0+S(h);
end
S(91)=vv(1,91)*L/2*(depth(2)-depth(1))/2;
sss=0;
for h=13:91;
sss=sss+(S(h));
end
%中间层流量
for j=2:34; %定义为层数
M(j,1)=vv(j,1)*L/2*(depth(j+1)-depth(j-1))/2;
for k=2:90
M(j,k)=vv(j,k)*L*(depth(j+1)-depth(j-1))/2;
end
M(j,91)=vv(j,91)*L/2*(depth(j+1)-depth(j-1))/2;
end
mm=0;
for j=2:34;
mm=mm+M(j,91);
end
mmm=0;
for j=2:20;
for k=13:90;
mmm=mmm+M(j,k);
end
end
mmmm=0;
for j=21:34;
for k=14:90;
mmmm=mmmm+M(j,k);
end
end
%底层流量
B(1)=vv(35,1)*L/2*(depth(35)-depth(34))/2;
for l=2:90
B(l)=vv(35,l)*L*(depth(35)-depth(34))/2;
b=0+B(l);
end
B(91)=vv(35,91)*L/2*(depth(35)-depth(34))/2;
bbb=0;
for l=14:91;
bbb=bbb+(B(l));
end
Q1=(sss+mm+mmm+mmmm+bbb)/1000
Q(d) = Q1;
end
解决了,是循环的参数用错了。
以上是对应于需求的程序。
|
|