中国科学院大气物理研究所研究员冯锦明课题组Matlab-MetToolbox工具箱简介

dapmood · 发表于 2021-2-24 17:31

EDA365欢迎您登录！

您需要登录才可以下载或查看，没有帐号？注册

x

MetToolbox工具箱简介：

1 风矢量图：

clc;clear;close all
[x,y] = meshgrid(70:10:140,15:10:55);
u = cos(x).*y*10;
v = sin(x).*y*10;
figure(1);
% % m_proj('Equidistant Cylindrical','long',[70 140],'lat',[15 55]);
% m_proj('lambert','lon',[70,140],'lat',[15,55]);
% ma=shaperead('F:/RMeteoInfo/data/map/bou2_4l.shp');
% % m_line( [ma(：）.X], [ma(：）.Y],'color',[0.5,0.5,0.5]);%绘制范围内的地图
% % m_line([ma(：）.X],[ma(：）.Y],'color','k');%绘制范围内的地图
% m_plot([ma(：）.X],[ma(：）.Y],'color','k');%绘制范围内的地图
% hold on
h=quiver(x,y,u,v,0.5,'r-');
hl=quiverlegend(h,10);
%设置x轴范围和刻度
set(gca,'XLim',[70 140]);%X轴的数据显示范围
set(gca,'XTick',[70:10:140]);%设置要显示坐标刻度
set(gca,'XTickLabel',[70:10:140]);%给坐标加标签
%设置y轴范围和刻度
set(gca,'YLim',[15 55]);%X轴的数据显示范围
set(gca,'YTick',[15:10:55]);%设置要显示坐标刻度
set(gca,'YTickLabel',[15:10:55]);%给坐标加标签

2 高版本m_map工具箱指南针&比例尺：

clc;clear;close all
m_proj('miller','lat',[-77 77]);
m_coast('patch',[.7 1 .7],'edgecolor','none');
m_grid('box','fancy','linestyle','-','gridcolor','w','backcolor',[.2 .65 1]);
m_ruler([.4 .9],.3,'tickdir','out','ticklen',[.007 .007]);
% 需要更新m_map
% m_northarrow(76,48,0.8,'type',4,'aspect',1.5);
m_northarrow(-125.4,49.7,.2,'type',2);
% northarrow(76,48,0.8,'type',4,'aspect',1.5);
title('Great Circle Routes','fontsize',14,'fontweight','bold');
set(gcf,'color','w'); % Need to do this otherwise 'print' turns the lakes black

3 利用Map Toolbox不同投影全球地形的绘制

clc;clear;close all
load etopo
load topo
figure(1);
ctopo=topo(90:160,50:160);
clegend=[1 70 50];
worldmap(topo,topolegend);
meshm(topo,topolegend);
dEMCmap(topo)
axis off
setm(gca,'Origin',[0 105 0])
ptopo = topo; ptopo(topo<0) = 0;
ptopo = 100*ptopo/(6371*1000);
figure(2);
axesm('MapProjection','ortho','Origin',[0 0 15])
meshm(topo,topolegend,size(topo),topo/10000);
shading interp
demcmap(topo)
lightm(23.5,80,'Color',[1 1 0.2])
material([.6 0.8 0]);
lighting phong
axis off
daspectm('km',20000)
figure(3);
axesm('MapProjection','globe','Origin',[0 105 0])
meshm(topo,topolegend,size(topo),topo/100000);
shading interp
% [lat,lon]=meshgrat(topo,topolegend);
% suRFm(lat,lon,topo,ptopo)
demcmap(topo)
lightm(23.5,20,'Color','y')
material([.6 1 0]);
lighting phong
axis off
view([200 0])
z=topo;
figure(4);
surf(z)
shading interp
demcmap(z)
view([0 90])
axis off
light('position',[1 180 1000],'color','w')
material([0.6 1 0]);
lighting phong

4 利用uniteprv合并东北三省与内蒙四个区域

clc;clear;close all
load prvpatch;
figure
hold on
wnprv={'内蒙','黑龙江','吉林','辽宁'};
nw=1;
for i=1:38
if any(cellfun(@(x) ~isempty(strmatch(x,s(i).tag)),wnprv))
      disp([nw i])
      if nw==1
         long=s(i).long;
         lat=s(i).lat;
      else
         [long lat]=uniteprv(long,lat,s(i).long,s(i).lat);
      end
      nw=nw+1;
end
plot(s(i).long,s(i).lat)
end
% [long lat]=pactpt(long,lat,0.1);
patch(long,lat,'m');
mapnorm;

5 填充带条纹的等值线

clc;clear;close all
z=peaks;
x=1:size(z,2);
y=1:size(z,1);
cv=-8:2:8;
cmap=cusjet(10);
ifirst=0;
ndelta=1;
iswbg=true;
bardir='horiz';
% csty={'r.-',linewidth,markersize,linecolor,markeredgecolor,markerfacecolor};
% tex={ctex,nnl,ldir,csty}
csty={'r.-'};
tex={-4,[100 80],1,{'r:'}; 2,[100 80],3,{'b-'}};
figure
[C,h,hb]=freecontourf2(x,y,z,cv,cmap,ifirst,ndelta,iswbg,bardir,tex);

6 空间插值&白化

clc;clear;close all
load('mdata/chnmon756.mat');
idate=find(stat.date>=19881200&stat.date<=19981100);
x=stat.long;
y=stat.lat;
z=stat.prec(idate,：);
z=mean(z)';%mm/day
lo=[72:1:137];
la=[15:1:55];
[clo,cla]=chngon('g');
[lo2 la2]=meshgrid(lo,la);
in=inpolygon(lo2,la2,clo,cla);
cv=[0.1 1:5];
% cv=[0:5];
% cmap=nclcmap(142,numel(cv)+1);
load('colorbar-mat/rainbow.mat');
cmap=rainbow;
zi=cus_griddata(x,y,z,lo,la,'v4');
figure
hold on
% [C,h,hb]=freecontourf(lo,la,zi,cv,cmap,true);
[C,h,hb]=freecontourf(lo,la,zi,cv);
set(h,'edgecolor','none')
set(h,'linew',2)
chnmask2(lo,la,'g')
hc=mapnorm(lo([1 end]),la([1 end]),false,false,true,true,true,0,[10 5],1,'g',2);
grid off
set(gcf,'pos',[89 200 520 445])
set(gca,'pos',[0.0913 0.1504 0.8747 0.7730]);
ssaxes(4,0.14,[3 4 1.6 1.1 1.35 0.98],'nanhai2',2);
set(gca,'linew',2);
set(gca,'fontsize',14)

7 WRF模拟嵌套区域绘制

clc;clear;close all
si=[1, 107, 94, 47];
sj=[1, 47, 37, 122];
nx=[179 169 151 145]; %nested grid number should be 3n+1
ny=[161 199 157 151];
dxy=[30000 10000 3333.33 3333.33];
grat=[1 1/3 1/3 1/3];
pid=[1 1 2 2];
% param=[nx(1) ny(1) 100 35 30 60 36000];
% param=[nx(1) ny(1) 105 35 (nx(1)+1)/2 (ny(1)+1)/2 30 60 115 dxy([1 1])];
param=[nx(1) ny(1) 102 36.5 (nx(1)+1)/2 (ny(1)+1)/2 30 60 102 dxy([1 1])];
[x y]=wrf_nestxy(si,sj,nx,ny,grat,pid);
nb=10;
nbb=nb*[1 1/3];
figure,hold on
maplam(param,false,false,true,true,true,0,[10 5],3,'h',2,':');
for i=2:length(pid)
[xx yy]=meshgrid(x{i},y{i});
xx2=[xx(1,：) xx(:,end)' xx(end,end:-1:1) xx(end:-1:1,1)'];
yy2=[yy(1,：) yy(:,end)' yy(end,end:-1:1) yy(end:-1:1,1)'];
plot(xx2,yy2,'m','linew',2)
end
[xx yy]=meshgrid(x{1},y{1});
xx2=[xx(1,1)+nbb(1) xx(1,end)'-nbb(1) xx(1,end)-nbb(1) xx(1,1)+nbb(1) xx(1,1)+nbb(1)];
yy2=[yy(1,1)+nbb(1) yy(1,1)'+nbb(1) yy(end,1)-nbb(1) yy(end,1)-nbb(1) yy(1,1)+nbb(1)];
plot(xx2,yy2,'b','linew',2)
load prvpatch
[clong clat]=polyjoin({s(1:34).long},{s(1:34).lat});
[cx cy]=deg2grd(clong,clat,param);
plot(cx,cy,'k')

8 Analysis of power spectrum

clc;clear;close all
% Analysis of power spectrum
% load enso_sst.mat
zsim=xlsread('data2.xlsx');
zi=zsim(:,1);
zi=zi-mean(zi); % must remove mean.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%fft
Fs = 1;             % Sampling frequency
T = 1/Fs;                   % Sample time
L = length(zi);             % Length of signal
t = (0：L-1)*T;             % Time vector
y=zi;
figure,hold on
% NFFT = 2^nextpow2(L); % Next power of 2 from length of y
NFFT = L; % Next power of 2 from length of y
Y = fft(y,NFFT);
P=Y.*conj(Y)/NFFT; %power spectrum density
P=smooth(P,7);
f=Fs/2*(0:NFFT/2-1)/(NFFT/2);
P=P(1:NFFT/2);
df=1/NFFT;
svar=(P(1)+P(end))*df/2+sum(P(2:end-1))*df;
svar0=var(y);
P=P/svar*svar0;
% Plot single-sided amplitude spectrum.
plot(f,P,'r')
title('Single-Sided Amplitude Spectrum of y(t)')
xlabel('Frequency (cycles/year)')
ylabel('Variance')

9 Taylor图

clc;clear;close all
zsim=xlsread('data2.xlsx');
% Example for ploting taylor diagram
lss={ 'o','v','s','^','p','o','v','s','^','p','o','v','s','^','p','o','v','s','^','p','o'};
lcc=[1.0000       0       0
0 0.8000       0
1.0000       0       0
0 0.9000 0.9000
1.0000       0 1.0000
1.0000 0.8000       0
1.0000       0       0
0 0.8000       0
0       0 1.0000
0 0.9000 0.9000
1.0000       0 1.0000
1.0000 0.8000       0
0       0 1.0000
0 0.8000       0
0       0 1.0000
0 0.9000 0.9000
1.0000       0 1.0000
1.0000 0.8000       0
1.0000       0       0
0 0.8000       0
0       0 1.0000];
zobs=zsim(:,1);
zsim=zsim(:,1:end);
stag=cellstr(num2str((1:size(zsim,2))'));
stdd=std(zsim)/std(zobs);
ccorr=corr(zobs,zsim);
z1=zobs-mean(zobs);
z2=zsim-mean(zsim);
rmse=rms(z2-z1)/std(zobs);
rmse2=sqrt(stdd.*stdd-2*stdd.*ccorr+1);
[hp,hp_std,hp_rmse,ht_corr,ht_rmse,hl]=TaylorDiag(zobs,zsim,lcc,lss,stag);
set(hp,'markersize',8);
set(ht_rmse,'color','k');

SsaaM7 · 发表于 2021-2-24 18:37

不愧是中国科学院大气物理研究所，Matlab-MetToolbox工具箱简介介绍的真的太棒

帐号		自动登录	找回密码
密码			注册

中国科学院大气物理研究所研究员冯锦明课题组Matlab-MetToolbox工具箱简介

EDA365欢迎您登录！

推荐内容 /1