Matlab处理气象数据（十七）季节变化

最新推荐文章于 2022-11-03 21:09:53 发布
最新推荐文章于 2022-11-03 21:09:53 发布
阅读量5.6k 收藏 60
点赞数 10
分类专栏： Matlab 气象 文章标签： matlab
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/ymyz1229/article/details/106062525
版权
气象 同时被 2 个专栏收录
23 篇文章 63 订阅
订阅专栏
Matlab
18 篇文章 70 订阅
订阅专栏

17.1 各省冬季和夏季的DTR趋势

以江西为例，其他省份只需将下面代码中的“江西”替换为其他省份的名称即可。

%%
clc,clear
load('E:\数据\201507\各省MASK\江西');%MASK是各省掩膜
load('E:\数据\201501\Tmax1');%导入NCEP数据最高温度 72*128*396
load('E:\数据\201501\Tmax2');%导入观测数据最高温度 72*128*396
load('E:\数据\201501\Tmin1');%导入NCEP数据最低温度 72*128*396
load('E:\数据\201501\Tmin2');%导入观测数据最低温度 72*128*396
DTR1=Tmax1-Tmin1;
DTR2=Tmax2-Tmin2;

%冬季DTR温度 72*128*32
for i=1:32
	wDTR1(:,:,i)=squeeze(nanmean(DTR1(:,:,(12+12*(i-1)):(14+12*(i-1))),3));
	wDTR2(:,:,i)=squeeze(nanmean(DTR2(:,:,(12+12*(i-1)):(14+12*(i-1))),3));
end

%全国夏季DTR温度 72*128*33
for i=1:33
	smDTR1(:,:,i)=squeeze(nanmean(DTR1(:,:,(6+12*(i-1)):(7+12*(i-1))),3));
	smDTR2(:,:,i)=squeeze(nanmean(DTR2(:,:,(6+12*(i-1)):(7+12*(i-1))),3));
end

clear Tmax1 Tmax2 Tmin1 Tmin2 DTR1 DTR2

%%

%冬季

%对NCEP数据
T_mean=squeeze(nanmean(wDTR1,3));
for i=1:32
	wDTR1(:,:,i)=squeeze(wDTR1(:,:,i))-T_mean;
end
T12=wDTR1;%删除边界区域外的值
for i=1:72
    for j=1:128
		if isnan(MASK(i,j))==1
			T12(i,j,:)=nan;
		end
    end
end

T12=reshape(T12,[72*128 32]);
T12=nanmean(T12,1);
NCEP=T12;
p1=polyfit(1:32,T12,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target1=polyval(squeeze(p1(:)),1:32);
clear wDTR1 T12 T_mean  

%对观测数据
T_mean=squeeze(nanmean(wDTR2,3));
for i=1:32
	wDTR2(:,:,i)=squeeze(wDTR2(:,:,i))-T_mean;
end

clear i j T_mean  

T22=wDTR2;%删除边界区域外的值
for i=1:72
    for j=1:128
		if isnan(MASK(i,j))==1
			T22(i,j,:)=nan;
		end
    end
end

T22=reshape(T22,[72*128 32]);
T22=nanmean(T22,1);
Obs=T22;
p2=polyfit(1:32,T22,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target2=polyval(squeeze(p2(:)),1:32);
clear j i wDTR2 T22 

%求P值
t=1:32;
t=t';
t=[ones(32,1),t];
NCEP = NCEP';
[b1,bint1,r1,rint1,stats1]= regress(NCEP,t,0.05);
Obs = Obs';
[b2,bint2,r2,rint2,stats2]= regress(Obs,t,0.05);

%作图
clf
subplot(1,2,1)
plot(NCEP,'r','linewidth',2)
hold on
plot(Obs,'b','linewidth',2)
plot(target1,'r','linewidth',2)%趋势线
plot(target2,'b','linewidth',2)
legend({'NCEP','Observed'},'Location','Northwest')
set(gca,'xtick',[2 7 12 17 22 27 32]); 
set(gca,'xticklabel',{'1980','1985','1990','1995','2000','2005','2010'}); 
axis([ -inf inf -2 2]);
h=xlabel('Year');
set(h, 'FontSize',8,'FontWeight','Bold')
h=ylabel('Temp. Anom. (\circC)');
set(h, 'FontSize',8,'FontWeight','Bold')
set(gca, 'FontSize',8,'FontWeight','Bold','tickdir','out')
xlim([1 33])%x轴范围锁定
box off %去掉外框

%添加斜率方程、R2、P值
h=text(9,-1.4,['y(NCEP)=' num2str(roundn(b1(1,:),-3)) '+' num2str(roundn(b1(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(9,-1.6,['R^2(NCEP)='  num2str(roundn(stats1(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(9,-1.8,['P(NCEP)='  num2str(roundn(stats1(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(22,-1.4,['y(Obs)=' num2str(roundn(b2(1,:),-3)) '+' num2str(roundn(b2(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(22,-1.6,['R^2(Obs)='  num2str(roundn(stats2(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(22,-1.8,['P(Obs)='  num2str(roundn(stats2(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(20,1.8,['Jiangxi(winter)']) %添加省份名称
set(h, 'FontSize',15)

%%

%夏季

%对NCEP数据
T_mean=squeeze(nanmean(smDTR1,3));
for i=1:33
smDTR1(:,:,i)=squeeze(smDTR1(:,:,i))-T_mean;
end
T12=smDTR1;%删除边界区域外的值
for i=1:72
    for j=1:128
		if isnan(MASK(i,j))==1
			T12(i,j,:)=nan;
		end
    end
end

T12=reshape(T12,[72*128 33]);
T12=nanmean(T12,1);
NCEP=T12;
p1=polyfit(1:33,T12,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target1=polyval(squeeze(p1(:)),1:33);

clear smDTR1 T12 T_mean 

%对观测数据
T_mean=squeeze(nanmean(smDTR2,3));
for i=1:33
	smDTR2(:,:,i)=squeeze(smDTR2(:,:,i))-T_mean;
end
clear i j T_mean 
T22=smDTR2;%删除边界区域外的值
for i=1:72
    for j=1:128
		if isnan(MASK(i,j))==1
			T22(i,j,:)=nan;
		end
    end
end

T22=reshape(T22,[72*128 33]);
T22=nanmean(T22,1);
Obs=T22;
p2=polyfit(1:33,T22,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target2=polyval(squeeze(p2(:)),1:33);
clear j i smDTR2 T22 

%求P值
t=1:33;
t=t';
t=[ones(33,1),t];
NCEP = NCEP';
[b1,bint1,r1,rint1,stats1]= regress(NCEP,t,0.05);
Obs = Obs';
[b2,bint2,r2,rint2,stats2]= regress(Obs,t,0.05);

%作图
subplot(1,2,2)
plot(NCEP,'r','linewidth',2)
hold on
plot(Obs,'b','linewidth',2)
plot(target1,'r','linewidth',2)%趋势线
plot(target2,'b','linewidth',2)
legend({'NCEP','Observed'},'Location','Northwest')
set(gca,'xtick',[2 7 12 17 22 27 32]); 
set(gca,'xticklabel',{'1980','1985','1990','1995','2000','2005','2010'}); 
axis([ -inf inf -2 2]);
h=xlabel('Year');
set(h, 'FontSize',8,'FontWeight','Bold')
h=ylabel('Temp. Anom. (\circC)');
set(h, 'FontSize',8,'FontWeight','Bold')
set(gca, 'FontSize',8,'FontWeight','Bold','tickdir','out')
xlim([1 33])%x轴范围锁定
box off %去掉外框

%添加斜率方程、R2、P值
h=text(9,-1.4,['y(NCEP)=' num2str(roundn(b1(1,:),-3)) '+' num2str(roundn(b1(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(9,-1.6,['R^2(NCEP)='  num2str(roundn(stats1(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(9,-1.8,['P(NCEP)='  num2str(roundn(stats1(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(22,-1.4,['y(Obs)=' num2str(roundn(b2(1,:),-3)) '+' num2str(roundn(b2(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(22,-1.6,['R^2(Obs)='  num2str(roundn(stats2(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(22,-1.8,['P(Obs)='  num2str(roundn(stats2(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(18,1.8,['Jiangxi(summer)']) %添加省份名称
set(h, 'FontSize',15)

得到结果如下图所示：
在这里插入图片描述

江西省冬季和夏季的DTR趋势

17.2 全国平均温度的四季变化趋势

clc;clear
load('E:\数据\201507\MASK');%MASK是城市区域掩膜
load('E:\数据\201501\T1');%导入NCEP数据平均温度 72*128*420
load('E:\数据\201501\T2');%导入观测数据平均温度 72*128*420
 
%全国春季年平均温度 72*128*35
for i=1:35
    sT1(:,:,i)=squeeze(nanmean(T1(:,:,(3+12*(i-1)):(5+12*(i-1))),3));
    sT2(:,:,i)=squeeze(nanmean(T2(:,:,(3+12*(i-1)):(5+12*(i-1))),3));
end
 
%%
%对NCEP数据
T_mean=squeeze(nanmean(sT1,3));
for i=1:35
	sT1(:,:,i)=squeeze(sT1(:,:,i))-T_mean;
end
clear i j T_mean 
 
T12=sT1;%删除边界区域外的值
for i=1:72
    for j=1:128
		if isnan(MASK(i,j))==1
			T12(i,j,:)=nan;
		end
    end
end
 
T12=reshape(T12,[72*128 35]);
T12=nanmean(T12,1);
NCEP=T12;
p1=polyfit(1:35,T12,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target1=polyval(squeeze(p1(:)),1:35);
clear j i T12 
 
 
%对观测数据
T_mean=squeeze(nanmean(sT2,3));
for i=1:35
	sT2(:,:,i)=squeeze(sT2(:,:,i))-T_mean;
end
clear i j T_mean
 
T22=sT2;%删除边界区域外的值
for i=1:72
    for j=1:128
		if isnan(MASK(i,j))==1
			T22(i,j,:)=nan;
		end
    end
end
 
T22=reshape(T22,[72*128 35]);
T22=nanmean(T22,1);
Obs=T22;
p2=polyfit(1:35,T22,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target2=polyval(squeeze(p2(:)),1:35);
clear j i T22 
 
 
%求P值
t=1:35;
t=t';
t=[ones(35,1),t];
NCEP = NCEP';
[b1,bint1,r1,rint1,stats1]= regress(NCEP,t,0.05);
Obs = Obs';
[b2,bint2,r2,rint2,stats2]= regress(Obs,t,0.05);
 
 
%作图
clf
subplot(1,2,1)
plot(NCEP,'r','linewidth',2)
hold on
plot(Obs,'b','linewidth',2)
plot(target1,'r','linewidth',2)%趋势线
plot(target2,'b','linewidth',2)
 
legend({'NCEP','Observed'},'Location','Northwest')
set(gca,'xtick',[2 7 12 17 22 27 32]); 
set(gca,'xticklabel',{'1980','1985','1990','1995','2000','2005','2010'}); 
axis([ -inf inf -2 2]);
h=xlabel('Year');
set(h, 'FontSize',8,'FontWeight','Bold')
h=ylabel('Temp. Anom. (\circC)');
set(h, 'FontSize',8,'FontWeight','Bold')
set(gca, 'FontSize',8,'FontWeight','Bold','tickdir','out')
xlim([1 35])%x轴范围锁定为1~35
box off %去掉外框
 
%添加斜率方程、R2、P值
h=text(10,-1.5,['y(NCEP)=' num2str(roundn(b1(1,:),-3)) '+' num2str(roundn(b1(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(10,-1.7,['R^2(NCEP)='  num2str(roundn(stats1(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(10,-1.9,['P(NCEP)='  num2str(roundn(stats1(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(23,-1.5,['y(Obs)=' num2str(roundn(b2(1,:),-3)) '+' num2str(roundn(b2(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(23,-1.7,['R^2(Obs)='  num2str(roundn(stats2(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(23,-1.9,['P(Obs)='  num2str(roundn(stats2(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(26,1.8,['spring(a)']) %添加名称
set(h, 'FontSize',15)
%可以用 delete(h) 来撤销句柄
 
%%
%对NCEP数据
T_mean=squeeze(nanmean(sT1,3));
for i=1:35
	sT1(:,:,i)=squeeze(sT1(:,:,i))-T_mean;
end
clear i j T_mean 
 
T12=sT1;%删除边界区域外的值
for i=1:72
    for j=1:128
		if isnan(MASK(i,j))==1
			T12(i,j,:)=nan;
		end
    end
end
 
%抠非城市的
for i=1:72
    for j=1:128
		if isnan(sT1(i,j,:))==1
			sT1(i,j,:)=0;
		end
	end
end
for i=1:72
    for j=1:128
		if isnan(T12(i,j,:))==1
			T12(i,j,:)=0;
		end
	end
end
T13=sT1-T12;
for i=1:72
    for j=1:128
		if T13(i,j,:)==0
			T13(i,j,:)=nan;
		end
	end
end
 
T13=reshape(T13,[72*128 35]);
T13=nanmean(T13,1);
NCEP=T13;
p1=polyfit(1:35,T13,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target1=polyval(squeeze(p1(:)),1:35);
clear j i sT1 T13 
 
 
%对观测数据
T_mean=squeeze(nanmean(sT2,3));
for i=1:35
	sT2(:,:,i)=squeeze(sT2(:,:,i))-T_mean;
end
clear i j T_mean
 
T22=sT2;%删除边界区域外的值
for i=1:72
    for j=1:128
		if isnan(MASK(i,j))==1
			T22(i,j,:)=nan;
		end
    end
end
 
%抠非城市的
for i=1:72
    for j=1:128
		if isnan(sT2(i,j,:))==1
			sT2(i,j,:)=0;
		end
	end
end
for i=1:72
    for j=1:128
		if isnan(T22(i,j,:))==1
			T22(i,j,:)=0;
		end
	end
end
T23=sT2-T22;
for i=1:72
    for j=1:128
		if T23(i,j,:)==0
			T23(i,j,:)=nan;
		end
	end
end
 
T23=reshape(T23,[72*128 35]);
T23=nanmean(T23,1);
Obs=T23;
p2=polyfit(1:35,T23,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target2=polyval(squeeze(p2(:)),1:35);
clear j i sT2 T23 
 
 
%求P值
t=1:35;
t=t';
t=[ones(35,1),t];
NCEP = NCEP';
[b1,bint1,r1,rint1,stats1]= regress(NCEP,t,0.05);
Obs = Obs';
[b2,bint2,r2,rint2,stats2]= regress(Obs,t,0.05);
 
 
%作图
subplot(1,2,2)
 
plot(NCEP,'r','linewidth',2)
hold on
plot(Obs,'b','linewidth',2)
plot(target1,'r','linewidth',2)%趋势线
plot(target2,'b','linewidth',2)
 
legend({'NCEP','Observed'},'Location','Northwest')
set(gca,'xtick',[2 7 12 17 22 27 32]); 
set(gca,'xticklabel',{'1980','1985','1990','1995','2000','2005','2010'}); 
axis([ -inf inf -2 2]);
h=xlabel('Year');
set(h, 'FontSize',8,'FontWeight','Bold')
h=ylabel('Temp. Anom. (\circC)');
set(h, 'FontSize',8,'FontWeight','Bold')
set(gca, 'FontSize',8,'FontWeight','Bold','tickdir','out')
xlim([1 35])%x轴范围锁定为1~35
box off %去掉外框
 
%添加斜率方程、R2、P值
h=text(10,-1.5,['y(NCEP)=' num2str(roundn(b1(1,:),-3)) '+' num2str(roundn(b1(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(10,-1.7,['R^2(NCEP)='  num2str(roundn(stats1(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(10,-1.9,['P(NCEP)='  num2str(roundn(stats1(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(23,-1.5,['y(Obs)=' num2str(roundn(b2(1,:),-3)) '+' num2str(roundn(b2(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(23,-1.7,['R^2(Obs)='  num2str(roundn(stats2(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(23,-1.9,['P(Obs)='  num2str(roundn(stats2(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(26,1.8,['spring(b)']) %添加名称
set(h, 'FontSize',15)
%可以用 delete(h) 来撤销句柄

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

全国城市区域(a)和非城市区域(b)平均温度的四季变化趋势

类似地,可以得到全国城市区域和非城市区域最高温度和最低温度的四季变化趋势，此处不再赘述。

17.3 全国DTR的四季变化趋势

clc;clear
load('E:\数据\201507\MASK');%MASK是城市区域掩膜
load('E:\数据\201501\Tmax1');%导入NCEP数据最高温度 72*128*396
load('E:\数据\201501\Tmax2');%导入观测数据最高温度 72*128*396
load('E:\数据\201501\Tmin1');%导入NCEP数据最低温度 72*128*396
load('E:\数据\201501\Tmin2');%导入观测数据最低温度 72*128*396
 
DTR1=Tmax1-Tmin1;
DTR2=Tmax2-Tmin2;
 
%全国春季DTR温度 72*128*33
for i=1:33
	sDTR1(:,:,i)=squeeze(nanmean(DTR1(:,:,(3+12*(i-1)):(5+12*(i-1))),3));
	sDTR2(:,:,i)=squeeze(nanmean(DTR2(:,:,(3+12*(i-1)):(5+12*(i-1))),3));
end
 
 
%%
%对NCEP数据
T_mean=squeeze(nanmean(sDTR1,3));
for i=1:33
	sDTR1(:,:,i)=squeeze(sDTR1(:,:,i))-T_mean;
end
clear i j T_mean
 
T12=sDTR1;%删除边界区域外的值
for i=1:72
	for j=1:128
		if isnan(MASK(i,j))==1
			T12(i,j,:)=nan;
		end
	end
end
 
T12=reshape(T12,[72*128 33]);
T12=nanmean(T12,1);
NCEP=T12;
p1=polyfit(1:33,T12,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target1=polyval(squeeze(p1(:)),1:33);
clear j i T12
 
 
%对观测数据
T_mean=squeeze(nanmean(sDTR2,3));
for i=1:33
	sDTR2(:,:,i)=squeeze(sDTR2(:,:,i))-T_mean;
end
clear i j T_mean
 
T22=sDTR2;%删除边界区域外的值
for i=1:72
	for j=1:128
		if isnan(MASK(i,j))==1
			T22(i,j,:)=nan;
		end
	end
end
 
T22=reshape(T22,[72*128 33]);
T22=nanmean(T22,1);
Obs=T22;
p2=polyfit(1:33,T22,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target2=polyval(squeeze(p2(:)),1:33);
clear j i T22
 
%求P值
t=1:33;
t=t';
t=[ones(33,1),t];
NCEP = NCEP';
[b1,bint1,r1,rint1,stats1]= regress(NCEP,t,0.05);
Obs = Obs';
[b2,bint2,r2,rint2,stats2]= regress(Obs,t,0.05);
 
%作图
clf
subplot(1,2,1)
plot(NCEP,'r','linewidth',2)
hold on
plot(Obs,'b','linewidth',2)
plot(target1,'r','linewidth',2)%趋势线
plot(target2,'b','linewidth',2)
legend({'NCEP','Observed'},'Location','Northwest')
set(gca,'xtick',[2 7 12 17 22 27 32]);
set(gca,'xticklabel',{'1980','1985','1990','1995','2000','2005','2010'});
axis([ -inf inf -2 2]);
h=xlabel('Year');
set(h, 'FontSize',8,'FontWeight','Bold')
h=ylabel('Temp. Anom. (\circC)');
set(h, 'FontSize',8,'FontWeight','Bold')
set(gca, 'FontSize',8,'FontWeight','Bold','tickdir','out')
xlim([1 33])%x轴范围锁定为1~33
box off %去掉外框

%添加斜率方程、R2、P值
h=text(8,-1.5,['y(NCEP)=' num2str(roundn(b1(1,:),-3)) '+' num2str(roundn(b1(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(8,-1.7,['R^2(NCEP)='  num2str(roundn(stats1(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(8,-1.9,['P(NCEP)='  num2str(roundn(stats1(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(21,-1.5,['y(Obs)=' num2str(roundn(b2(1,:),-3)) '+' num2str(roundn(b2(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(21,-1.7,['R^2(Obs)='  num2str(roundn(stats2(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(21,-1.9,['P(Obs)='  num2str(roundn(stats2(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(26,1.8,['spring(a)']) %添加名称
set(h, 'FontSize',15)
%可以用 delete(h) 来撤销句柄
 
 
%%
%对NCEP数据
T_mean=squeeze(nanmean(sDTR1,3));
for i=1:33
	sDTR1(:,:,i)=squeeze(sDTR1(:,:,i))-T_mean;
end
clear i j T_mean
 
T12=sDTR1;%删除边界区域外的值
for i=1:72
	for j=1:128
		if isnan(MASK(i,j))==1
			T12(i,j,:)=nan;
		end
	end
end
 
%抠非城市的
for i=1:72
	for j=1:128
		if isnan(sDTR1(i,j,:))==1
			sDTR1(i,j,:)=0;
		end
	end
end
for i=1:72
	for j=1:128
		if isnan(T12(i,j,:))==1
			T12(i,j,:)=0;
		end
	end
end
T13=sDTR1-T12;
for i=1:72
	for j=1:128
		if T13(i,j,:)==0
			T13(i,j,:)=nan;
		end
	end
end
 
T13=reshape(T13,[72*128 33]);
T13=nanmean(T13,1);
NCEP=T13;
p1=polyfit(1:33,T13,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target1=polyval(squeeze(p1(:)),1:33);
clear j i sDTR1 T13
 
 
%对观测数据
T_mean=squeeze(nanmean(sDTR2,3));
for i=1:33
	sDTR2(:,:,i)=squeeze(sDTR2(:,:,i))-T_mean;
end
clear i j T_mean
 
T22=sDTR2;%删除边界区域外的值
for i=1:72
	for j=1:128
		if isnan(MASK(i,j))==1
			T22(i,j,:)=nan;
		end
	end
end
 
%抠非城市的
for i=1:72
	for j=1:128
		if isnan(sDTR2(i,j,:))==1
			sDTR2(i,j,:)=0;
		end
	end
end
for i=1:72
	for j=1:128
		if isnan(T22(i,j,:))==1
			T22(i,j,:)=0;
		end
	end
end
T23=sDTR2-T22;
for i=1:72
	for j=1:128
		if T23(i,j,:)==0
			T23(i,j,:)=nan;
		end
	end
end
 
T23=reshape(T23,[72*128 33]);
T23=nanmean(T23,1);
Obs=T23;
p2=polyfit(1:33,T23,1);  % p输出两个值第一个是a  第二个是b   y=ax+b
target2=polyval(squeeze(p2(:)),1:33);
clear j i sDTR2 T23
 
%求P值
t=1:33;
t=t';
t=[ones(33,1),t];
NCEP = NCEP';
[b1,bint1,r1,rint1,stats1]= regress(NCEP,t,0.05);
Obs = Obs';
[b2,bint2,r2,rint2,stats2]= regress(Obs,t,0.05);
 
%作图
subplot(1,2,2)
plot(NCEP,'r','linewidth',2)
hold on
plot(Obs,'b','linewidth',2)
plot(target1,'r','linewidth',2)%趋势线
plot(target2,'b','linewidth',2)
legend({'NCEP','Observed'},'Location','Northwest')
set(gca,'xtick',[2 7 12 17 22 27 32]);
set(gca,'xticklabel',{'1980','1985','1990','1995','2000','2005','2010'});
axis([ -inf inf -2 2]);
h=xlabel('Year');
set(h, 'FontSize',8,'FontWeight','Bold')
h=ylabel('Temp. Anom. (\circC)');
set(h, 'FontSize',8,'FontWeight','Bold')
set(gca, 'FontSize',8,'FontWeight','Bold','tickdir','out')
xlim([1 33])%x轴范围锁定为1~33
box off %去掉外框
 
%添加斜率方程、R2、P值
h=text(8,-1.5,['y(NCEP)=' num2str(roundn(b1(1,:),-3)) '+' num2str(roundn(b1(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(8,-1.7,['R^2(NCEP)='  num2str(roundn(stats1(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(8,-1.9,['P(NCEP)='  num2str(roundn(stats1(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(21,-1.5,['y(Obs)=' num2str(roundn(b2(1,:),-3)) '+' num2str(roundn(b2(2,:),-3)) '*x' '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(21,-1.7,['R^2(Obs)='  num2str(roundn(stats2(:,1),-2)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(21,-1.9,['P(Obs)='  num2str(roundn(stats2(:,3),-4)) '']);
set(h, 'FontSize',10,'FontWeight','Bold')
h=text(26,1.8,['spring(b)']) %添加名称
set(h, 'FontSize',15)
%可以用 delete(h) 来撤销句柄

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

全国城市区域(a)和非城市区域(b)DTR的四季变化趋势

相关链接：
Matlab处理气象数据——目录

岩鸣杨子
关注
  • 10
    点赞
  • 踩
  • 60
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
使用matlab处理地震数据
10-15
里面包含了matlab代码以及测试的地震数据。以及Wail A.Mousa]Processing of seismic reflection data Using matlab书籍PDF版，包含了偏移反褶积动校正一系列流程。
Matlab数据可视化分析气象笔记.md
最新发布
09-04
本文首先以气象数据分析为例,介绍了数据加载、统计计算和绘图可视化的过程,给出了绘制温湿度变化折线图的代码。然后以温度时间序列分析为例,讲解了绘制时间序列图、进行季节性分解、以及使用时间序列工具进行预测的...
MODIS16天数据月合成/季合成
ronvicki的博客
11-11 9057
1.MRT重投影、转换格式 重投影后命名为MOD13Q1.A2019129.250m_16_days_NDVI.tif 2.重命名 5月单独一期，即序号为129的一期。 首先去掉前面的MOD13Q1. 再将001.250替换为0101 将129.250替换为0501 得到A20200501_16_days_NDVI.tif 3.按月合成 ...
matlab_一种贝叶斯算法，用于检测变化点并将时间序列分解为趋势、季节性和突变
07-06
一种快速、通用的贝叶斯模型平均算法，用于将时间序列或1D序列数据分解为单个组件，例如突变、趋势和周期/季节变化，如赵等人（2019）所述。BEAST适用于变化点检测（即断点或结构中断）、非线性趋势分析、时间序列分解和时间序列分割。查看使用BEAST的选定研究列表。 BEAST是用C/C++实现的，但可以从R和Matlab中访问 在matlab run eval(webread('http://b.link/beast',weboptions('cert',''))) In R, run install.packages("Rbeast")
季节优化算法（Seasons optimization algorithm，SOA）附matlab代码
qq_59747472的博客
11-03 380
This paper introduces a new stochastic bio-inspired optimization algorithm, denoted as seasons optimization (SO) algorithm.This algorithm is inspired by the growth cycle of trees in diferent seasons of a year. It is an iterative and population-based algori
让用户输入一个月份，判断这个月是哪个季节？
weixin_46680533的博客
01-14 3199
让用户输入一个月份，判断这个月是哪个季节？假定3到 4月是春季，5到8月是夏季，9到10是秋季，11、12、1、2 月是冬季 a = int (input(“请输入一个月份：”)) if a == 3 or a == 4: print("%s月是春季"%a) elif 5 <= a <=8: print("%s月是夏季"%a) elif 9 <= a <=10: print("%s月是秋季"%a) else: print("%s月是冬季"%a) ...
气象基础知识matlab,气象类专业Matlab课程教学探索与思考
weixin_39843738的博客
03-17 469
气象类专业Matlab课程教学探索与思考2019-05-30 15:00:33次浏览寇蕾蕾【摘 要】气象类专业学生需要处理来自地面站、探空站、气象雷达、气象卫星等各种观测手段的不同类型和格式的资料，Matlab具有强大的数值计算和可视化显示功能，使得处理和分析各种气象资料变得方便和快捷。本文介绍了气象类专业Matlab课程教学工作的相关问题，从教学目标的制定、教学内容的设计、教学方法的选择对气象类...
Matlab处理气象数据（一）读取NCEP数据
设计の诗
05-10 1万+
数据获取 NCEP数据air.2m.mon.mean.nc，下载地址 读取 使用matlab命令读取NCEP数据。 常用的读取NC头文件的命令是 nc_dump 和 nc_dump 。 执行nc_dump命令： nc_dump('E:\数据\空间数据\NECP\air.2m.mon.mean.nc ') 得到如下结果： NetCDF-3 Classic E:\数据\空间数据\NECP\air.2m.mon.mean.nc { dimensions: lon = 192 ; lat = 94 ; le
wind_data_assessment.zip_matlab气象数据
07-14
该压缩包"wind_data_assessment.zip_matlab气象数据"是一个包含MATLAB代码和数据库的集合，主要用于处理和分析气象塔的风力数据。这个程序旨在帮助用户评估风能潜力，通过对历史风速和风向数据的处理，进行数据质量...
人工智能-项目实践-回归分析-MATLAB气象数据回归预测模型
05-21
在本项目实践中，我们将深入探讨如何使用人工智能技术，特别是回归分析方法，来建立一个MATLAB气象数据回归预测模型。回归分析是一种统计学工具，用于研究变量间的关系，并预测未来趋势。在这个项目中，我们将利用...
混合层厚度计算
04-14
计算全球大洋的内区和西边界混合层厚度及流函数
ra_mld_海洋_混合层深度_matlab图像处理_MLD是什么计量_
10-01
计算海洋混合层和障碍层深度并通过M_Map工具包画出全球海洋MLD和BL的分布图
基于MATLAB实现3种气象数据的读取和绘图 .pdf
05-07
利用MATLAB在Window环境下实现MICAPS、GRIB和NetCDF三种气象数 据的读取，以及物理量图绘制的方法，可为MATLAB在气象领域的应用提供一些参考
MATLAB进行季节调整步骤
04-18
MATLAB进行季节调整步骤
matlab气象学：气象因子的小波功率谱与小波分析图象
04-24
小波分析是一种强大的工具，它能够揭示数据中的局部特征和时间变化，这在处理气象数据时特别有用。MATLAB作为一种广泛使用的数学计算和编程环境，是进行这种复杂分析的理想选择。本项目将详细介绍如何利用MATLAB进行...
【气象水文】Python和Matlab处理VIIRS卫星数据（代码+数据）.zip
05-20
《Python与Matlab在气象水文领域处理VIIRS卫星数据详解》 在现代气象水文研究中，遥感数据的应用日益广泛，其中VIIRS（Visible Infrared Imaging Radiometer Suite）卫星数据因其高分辨率和全球覆盖特性，成为重要...
Matlab处理气象数据（四）观测数据的预处理和计算
设计の诗
05-10 1万+
观测数据为中国气象数据网上获得的中国地面温度月值0.5°×0.5°格点数据集，时间范围是1961年1月至2013年12月。这套数据为txt格式，包含头文件。头文件信息为： NCOLS 128 NROWS 72 XLLCORNER 72.000000000000 YLLCORNER 18.000000000000 CELLSIZE 0.50000000000000 NODATA_VALUE -99.00000 从中可以得知，有128
毕业设计 - 天气数据分析
HUXINY的博客
10-08 1万+
文章目录1 前言2 项目简介3 开始分析3.1 海洋对当地气候的影响3.2 导入数据集3.3 温度数据分析3.4 湿度数据分析3.5 风向频率玫瑰图3.6 计算风速均值的分布情况4 最后-毕设帮助 1 前言 Hi，大家好，这里是丹成学长，今天向大家介绍 一个项目 基于GRU的 电影评论情感分析 大家可用于 毕业设计 毕设帮助，开题指导，技术解答 ????746876041 2 项目简介 本实例将对意大利北部沿海地区的气象数据进行分析与可视化。首先会运用 Python 中 matplotlib 库对数据进行
Matlab处理气象数据——目录
热门推荐
设计の诗
05-10 1万+
整理了一下研究生期间做的一些气象数据处理方面的资料，讲研究过程和matlab代码一并贴出来，以便存档和查阅。 论文原文见知网 目录 1.读取NCEP数据 2.调用NCEP数据，进行预处理 3.NCEP数据的简单计算 4.观测数据的预处理和计算 5.NCEP数据的剪裁、插值和两套数据的对比 6.数据的面积加权 7.分段趋势线的做法 8.数据的逐点变化趋势 8.1 分别对NCEP和观测数据进行逐点趋势分析 8.2 平均、最高、最低温度差值的趋势 8.3 退耕还林前后的温度差值对比 9
matlab处理气象数据txt
07-30
在使用MATLAB处理气象数据的过程中，可以使用以下步骤进行操作。 首先，需要导入气象数据txt文件。可以使用MATLAB中的文本导入工具，如readtable函数，逐行读取txt文件中的数据，并存储为一个数据表格。在导入数据时，可以设置适当的参数，例如指定分隔符、数据类型和列名等。 其次，可以对数据进行预处理。根据具体需求，可以进行数据筛选、缺失值处理、异常值检测和处理，以及对数据进行平滑、插值或滤波等操作。预处理数据可以提高后续分析的准确性和可靠性。 然后，可以进行数据分析和统计。根据具体问题，可以采用MATLAB中的各种统计工具和函数，如均值、方差、相关性分析、频率分布分析等。还可以进行时序分析，例如时间序列分析、周期性分析和趋势分析等，以了解气象数据的特点和规律。 最后，可以通过数据可视化的方式呈现分析结果。MATLAB提供了丰富的绘图和可视化工具，如绘制折线图、散点图、柱状图、等高线图、热力图等。可以根据数据的特点和需求，选择适当的图形进行展示和解释。 总之，MATLAB是一个强大的工具，可以方便地处理和分析气象数据txt文件。通过合理的数据处理和分析，可以从中获取有用的信息，为气象研究和应用提供支持。

“相关推荐”对你有帮助么？

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值