- 实习目的
熟悉垂直速度在气象中的应用,掌握垂直速度的实际编程计算。
- 实习内容
已知2020年的6h再分析资料,要素场有温度、高度、相对湿度和水平风场。利用第二种修正方案计算各日(各小组负责日期的第2日)20时区域(10ºN-60ºN, 60ºE-160ºE)在200、500和700hPa上的垂直速度,并给出相应的图。(根据实习内容和资料说明,编写计算你熟悉的语言程序(如Fortran、c和matlab等),并输出计算结果和绘图(如NCL, Matlab和grads)。)
- 实习结果分析
3.1 计算的程序
%读取.nc资料数据_%变量
u=ncread(['E:\天诊\','uwnd.2020.nc'],'uwnd');%水平风场u
v=ncread(['E:\天诊\','vwnd.2020.nc'],'vwnd');%水平风场v
% (4维)索引变量-lon,lat,level,time
lon=ncread(['E:\天诊\','air.2020.nc'],'lon');
lat=ncread(['E:\天诊\','air.2020.nc'],'lat');
level=ncread(['E:\天诊\','air.2020.nc'],'level');
time=ncread(['E:\天诊\','air.2020.nc'],'time');
%10-60,60-160——6月9日——北京20时
timez=(datenum(2020,6,9,12,0,0)-datenum('1800-01-01 00:00:0.0'))*24; %将世界12时转化为资料时间
time_z=find(time==timez);
minlon=find(lon==60);
maxlon=find(lon==160);
minlat=find(lat==10);
maxlat=find(lat==60);
%%%%%%%%%%%
latt=lat(maxlat:minlat)/180*pi;
r=6371*10^3;
a=1/(2*r);
d=2.5/180*pi;%经纬度网格距转为弧度制
level_z=find(level==200);
%转换索引下标——%200-500-700hpa——%取出所需资料(lon,lat,level,time)
for level_z=1:size(level,1)
uu=u(minlon:maxlon,maxlat:minlat,level_z,time_z);
vv=v(minlon:maxlon,maxlat:minlat,level_z,time_z);
uu=uu';
vv=vv';
%计算所有散度
for i=2:size(uu,1)-1
for j=2:size(uu,2)-1
D(i,j,level_z)=a*((uu(i,j+1)-uu(i,j-1))/(cos(latt(i))*d) -(vv(i+1,j)-vv(i-1,j))/d-2*vv(i,j)*tan(latt(i)));
end
end
end
%计算omg
omg(:,:,1)=zeros(size(uu,1)-1,size(uu,2)-1,1);
for level_z=2:size(level,1)
omg(:,:,level_z)=omg(:,:,level_z-1)+(D(:,:,level_z)+D(:,:,level_z-1))/2*(level(level_z)-level(level_z-1));
end
%%%%%%%O'Brien方案修正
omg_xz(size(uu,1)-1,size(uu,2)-1,size(level,1))=0;
N=size(level,1);
M=N*(N+1);
for level_z=2:size(level,1)
omg_xz(:,:,level_z)=omg(:,:,level_z)-level_z*(level_z+1)/M*(omg_xz(:,:,level_z)-omg_xz(:,:,size(level,1)));
end
3.2 计算的绘图程序和描述文件
%%%%%%%插值画图
X=60+2.5:160-2.5;
Y=10+2.5:60-2.5;
lonlon=double(lon(minlon+1:maxlon-1));
latlat=double(lat(minlat-1:-1:maxlat+1));
%修正前
for i=[200 500 700]
figure
level_z=find(level==i);
O=omg(:,:,level_z);
[X1 Y1 Z1]=griddata(lonlon,latlat,double(O(2:20,2:40)),X',Y,'cubic');%范围缩小一圈是起点是+
m_proj('miller','lon',[60 160],'lat',[10 60]);
m_grid('linestyle','none','box','fancy','tickdir','out');hold on
m_contourf(X1,Y1,Z1,30,'linestyle','none'); hold on
title([num2str(i),'hpa垂直速度分布图(修正前)'], 'Rotation', 0, 'FontSize', 16);%显示图标题
colormap(m_colmap('jet',10));%填充颜色搭配
colorbar('Location','eastoutside');%颜色栏
m_coast('color',[0 0 0]);%海岸线-黑色三原色
hold on;
saveas(gcf,[num2str(i),'hpa垂直速度分布图(修正前)','.jpg'])
close gcf %不加会重叠
end
%修正后
for i=[200 500 700]
figure
level_z=find(level==i);
O=omg_xz(:,:,level_z);
[X1 Y1 Z1]=griddata(lonlon,latlat,double(O(2:20,2:40)),X',Y,'cubic');%范围缩小一圈是起点是+
m_proj('miller','lon',[60 160],'lat',[10 60]);
m_grid('linestyle','none','box','fancy','tickdir','out');hold on
m_contourf(X1,Y1,Z1,30,'linestyle','none'); hold on
title([num2str(i),'hpa垂直速度分布图(修正后)'], 'Rotation', 0, 'FontSize', 16);%显示图标题
colormap(m_colmap('jet',10));%填充颜色搭配
colorbar('Location','eastoutside');%颜色栏
m_coast('color',[0 0 0]);%海岸线-黑色三原色
hold on;
saveas(gcf,[num2str(i),'hpa垂直速度分布图(修正后)','.jpg'])
close gcf %不加会重叠
end
3.3 绘制图形、统计表格和相关分析
3.4 存在的问题或遇到的问题或体会或小结