南信大气象统计实习5，合成分析

耿瑞满

%有问题需要加QQ：3588430252
%m_map包以及othercolor需要提前下载
%elnino.m判断是否为厄尔尼诺年（粗略估计）
%合成分析
%对厄尔尼诺年判断
%12为X方向，18为Y方向，516为时间T维数
clc;clear;
fid = fopen('sstpx.grd');
[sst_jp,count]=fread(fid,inf,'float32');
sst_jp=reshape(sst_jp,12,18,516);


%% 画图,对-999处转换成NaN,并进行griddata插值
gcf=figure(1);
set(figure(1),'color','white')
addpath 'D:\Codes_Projects\MATLAB_codes\othercolor'
addpath 'C:\Program Files\MATLAB\R2023a\toolbox\m_map'
addpath 'D:\Codes_Projects\MATLAB_codes\Internship'
addpath 'D:\Codes_Projects\MATLAB_codes\practice_and_trial\fill\area_fill'

%转化-999成Nan，因为缺测值较多，插值过后仍然有缺测，将NaN改为0
lon=120:10:290;
lat=-27.5:5:27.5;
[lon0_grid,lat0_grid]=meshgrid(lon,lat);
sst_jp(sst_jp==-999)=0;

%griddata插值
[elnino_year,sig]=elnino(sst_jp);
%exalple_year=elnino_year(1);

exalple_year=2;
start_end=12*(exalple_year-1)+1:12*(exalple_year-1)+12;
lon_grid = linspace(min(lon), max(lon), 171);  
lat_grid = linspace(min(lat), max(lat), 56);  
[lon_grid, lat_grid] = meshgrid(lon_grid, lat_grid );
sst_grid= griddata(lon0_grid,lat0_grid,mean(sst_jp(:,:,start_end),3),lon_grid, lat_grid,'natural');


%色标显示，不画等值线
m_proj('Equidistant Cylindrical','long',[min(lon) max(lon)],'lat',[min(lat) max(lat)]);%选中范围
[c,h]=m_contourf(lon_grid,lat_grid,sst_grid,'ShowText','on');
level=[-3.5,-3,-2.5,-2,-1.5,-1,-0.5,0,0.5,1,1.5,2];
clim([-3.5 2]);
cbar=colorbar('southoutside','ticklength',0);
set(cbar,'Ticks',-3.5:0.5:2,'TickLabels',level);
set(h,'LineColor','none');


%加上色标单位 t=strcat(t,'mm/year');
t=get(cbar,'YTickLabel');
set(cbar,'Ticks',-3.5:0.5:2,'TickLabels',level,'YTickLabel',t);
colormap(othercolor('BuDOr_12'));

%绘制范围内的地图
ma=shaperead('D:\Maps\world\World_countries.shp');
m_line([ma(:).X],[ma(:).Y],'color','k');
m_grid('ytick',-25:5:25,'xtick',120:10:290,'tickdir','out','linest','none','fontname','Times','fontsize',7,'linewidth',1.5);
set(get(cbar,'ylabel'),'string','Sea Level anomality/(mm/year)','fontsize',12,'fontweight','bold');

%确定厄尔尼诺年份
cer_date=datetime(1948,1,1) + calyears(exalple_year);
title(sprintf('Year of %d SST',year(cer_date)),'FontSize',12);
print(gcf,sprintf('year=%dSST.jpg',exalple_year),'-r400','-djpeg')