本质上是数据插值,把粗分辨率数据插值成细网格
记录一个例子:
B3.m程序:
clear
clc
%% 读取数据
Lon_Data=ncread('2020-Jan-(1-31)-(0-24)-(150-450)-RH.nc','longitude');
Lat_Data=ncread('2020-Jan-(1-31)-(0-24)-(150-450)-RH.nc','latitude');
RH_Data=ncread('2020-Jan-(1-31)-(0-24)-(150-450)-RH.nc','r');
level=ncread('2020-Jan-(1-31)-(0-24)-(150-450)-RH.nc','level');
%% 计算250,300hPa平均湿度,6点和12点平均
r250=squeeze(RH_Data(:,:,find(level==250),:));
r300=squeeze(RH_Data(:,:,find(level==300),:));
r250_mean=nanmean(r250,3); %平均函数
r300_mean=nanmean(r300,3);
%% 读.shp 地区边界
map=shaperead('广西.shp');
bou_x=[map(:).X];%提取经度
bou_y=[map(:).Y];%提取纬度
%% 掩摸,抠数据
latlim_GX = [21 26.5];
lonlim_GX = [104 112.5];
lat_GX=21:0.02:26.5;
lon_GX=104:0.02:112.5;
M=length(lat_GX);
N=length(lon_GX);
R = georasterref('RasterSize',[M N],'ColumnsStartFrom','south','LatitudeLimits',latlim_GX,'LongitudeLimits',lonlim_GX);
Z = zeros(R.RasterSize);
Z = vec2mtx(bou_y,bou_x,Z,R,'filled');
Z(Z>1)=nan;
Z(Z==1)=0;
%% 湿度插值
[X,Y]=meshgrid(lon_GX,lat_GX);
[plon,plat]=meshgrid(Lon_Data,Lat_Data);
r_250=griddata(double(plon),double(plat),r250_mean',X,Y);
r_250(isnan(Z))=nan;
r_300=griddata(double(plon),double(plat),r300_mean',X,Y);
r_300(isnan(Z))=nan;
%% 272.6hpa高度插值
l=[250 300];
r_272=nan(276,426);
for m=1:276
for n=1:426
y=[r_250(m,n),r_300(m,n)];
r_272(m,n)= interp1(l,y,272.6); %线性插值
clear y
end
end
%% 绘图
subplot(2,2,1)
contourf(X,Y,r_272,100,'LineStyle','none');
hold on;
plot(bou_x,bou_y,'r','LineWidth',1.5);
hold off;
box on
hc = colorbar;
set(gca,'Clim',[16 26]);
set(hc,'linewi',1)
set(hc,'YTick',[16:1:26])
set(get(hc,'Title'),'string','(%)','fontsize',16,'fontweight','Bold');
set(gcf,'color','white')
set(gca,'Xlim',[104 113]);
set(gca,'Ylim',[21 27]);
set(gca,'XMinorTick','on','TickDir','out');
set(gca,'fontweight','Bold','fontsize',16,'linewi',1)
xlabel('经度','FontSize',16);
ylabel('纬度','Fontsize',16);
set(gca,'fontweight','Bold','fontsize',15,'linewidth',2)
title('272.6 hPa','fontweight','Bold','fontsize',18)
效果图: