Matlab海洋绘图：海面热通量分布——以南海为例

) O# V9 Q, D2 r; I7 y7 E2 K4 P3 }

图片赏析

# U! e. k. r% X/ M0 G: F5 A; T

. K$ Z9 x6 f; }0 @; ~; y

5 P: l& B. D" L& [; S7 s

前言

' o& D0 V- r! p1 J/ q% A

南海的经纬度：南海范围北纬2°至北纬25° 东经104°至东经124°

, C! T4 a% j. V8 y" n4 X" S+ y

数据来源：

National Climatic Data Center​www.ncdc.noaa.gov/cdr/atmospheric/ocean-heat-fluxes

1 C9 g8 r$ X! \% T7 D

数据链接提供者：广东海洋大学海洋与气象学院海洋科学系纪同学

数据来源机构：NOAA美国国家海洋与大气管理局-NCDC环境信息国家中心

产品名称：Ocean Heat Fluxes海洋热通量

时间精度：每日3-hourly记录

空间精度：0.25°

注：本程序需调用UBC英属哥伦比亚大学Rich Pawlowicz 开发的程序包m_map；

亦需调用gshhs高精度海岸线数据。

源代码

%南海1998年12月31日某时的海面热通量数据及分布图 %% 2 H. p' e4 v9 F5 x- r4 x %第一步：数据的展示 0 M1 }! a4 N5 f" o% h6 F8 l) V) S tic;Path="E:\OneDrive - Inc\Share_science\NC文件\";%文件夹地址 E$ w C3 g7 [4 n+ f6 e) [' @ filename="SEAFLUX-OSB-CDR_V02R00_FLUX_D19981231_C20160824.nc";%文件名 abs_path=fullfile(Path,filename);%合并，绝对路径文件 ( Z7 Q# Y+ c- `$ {; d# ]6 ] clearvars -except abs_path;%清除除了abs_path以外所有变量 ; f$ r) F8 n$ \, Z ncdisp(abs_path);%读取所选文件 1 K# D; o$ B; M: C; e6 x whos;toc; & I9 `5 e1 L% t7 h' v [+ J %% %第二步：数据读取 . U! i: u+ w% d tic;Lat_ini=ncread(abs_path,lat); Lon_ini=ncread(abs_path,lon); Time_ini=ncread(abs_path,time); Lflx_ini=ncread(abs_path,surface_upward_latent_heat_flux); Sflx_ini=ncread(abs_path,surface_upward_sensible_heat_flux); toc;clear abs_path; %% T$ X6 B2 e+ T, L/ F- _. \ %第三步：确定边界/截取数据 %空间边界：南海范围北纬2°至北纬25° 东经104°至东经124° # B! c4 @9 J! }, F0 k& } tic;Boundary=[104,124;2,25;]%修改范围[经度;纬度],注意顺序 , s8 A( f4 Q( V) e, [# c/ \ Spatial_Resolution=abs(Lat_ini(1)-Lat_ini(2));Sp_Re=Spatial_Resolution;%自动确定精度 Lon_start_No=find(Lon_ini-(Boundary(1))<=Sp_Re&Lon_ini-(Boundary(1))>=0);%找到离104°最近格点的编号（行号） " @. _& W$ R% U- i9 h& `5 @ Lon_end_No=find(Lon_ini-(Boundary(3))<=Sp_Re&Lon_ini-(Boundary(3))>=0);%找到离124°最近格点的编号（行号） 3 b0 T e$ Y0 q Lat_start_No=find(Lat_ini-(Boundary(2))<=Sp_Re&Lat_ini-(Boundary(2))>=0);%找到离2°最近格点的编号（行号） Lat_end_No=find(Lat_ini-(Boundary(4))<=Sp_Re&Lat_ini-(Boundary(4))>=0);%找到离25°最近格点的编号（行号） / i$ z' T/ m" i. U# k Lon=Lon_ini(Lon_start_No:Lon_end_No); Lat=Lat_ini(Lat_start_No:Lat_end_No);%北半球可能要颠倒一下编号，因为是从北极往赤道计数的 clear Lat_ini;clear Lon_ini;%回收内存 6 D9 u% a4 e( l- Y$ H; E. z Lat=double(Lat);Lon=double(Lon);%双精度化 . _6 t+ B9 ?5 _. `: e, _, n latlim=[min(Lat),max(Lat)];lonlim=[min(Lon),max(Lon)]; " N9 |9 L# y% @ Q% | Time=9%输入时间，注意三小时的倍数 - ^( {( V# ]# K; B1 [. K n. Z Time_No=Time/3+1; ! K2 a4 ?" s$ z) _% Y; q5 t6 f" S %截取数据并降维度 Lflx=squeeze(Lflx_ini(Lon_start_No:Lon_end_No,Lat_start_No:Lat_end_No,Time_No));;%北半球有时可能要颠倒一下编号，因为是从北极往赤道计数的 ; e8 _2 i& Y! S Sflx=squeeze(Sflx_ini(Lon_start_No:Lon_end_No,Lat_start_No:Lat_end_No,Time_No));;%北半球有时可能要颠倒一下编号，因为是从北极往赤道计数的 clear Sp_Re;toc;%清除不再使用的变量 %% %第四步：绘图 tic;figure; ! q' ~5 ?" f1 c m_proj(mercator,lat,latlim,lon,lonlim); [longrid,latgrid]=meshgrid(Lon,Lat); 9 B: g \+ C m" w9 v m_pcolor(longrid,latgrid,Sflx); %pcolor或contourf ; n' H! c7 n8 e$ Q1 P. y$ M tic;m_gshhs_f(patch,k);toc;%添加海岸线：c<l<h<f m_grid(box,fancy,tickdir,in);%网格化 %brighten(.5);%亮化 colormap(hsv); 2 q4 u$ G9 O% \' E3 v h = colorbar(h);%色标 h.Label.String = W·m^{-2} h.Location = eastoutside;%色标位置 title([1998年12月31日,num2str(Time),:00,南海海面感热通量],... FontName,黑体,fontsize,12,FontWeight,"bold") %标题 + n$ Q5 p6 S7 y F! V xlabel(Longitude);ylabel(Latitude) ;set(gca,fontsize,9); . X9 }# W+ b% A: p9 c- Q; _ toc; 9 Q1 P2 i# f/ Q$ P6 N6 K; a

: m! E& c4 g8 j: N: U

附上另外俩残次品绘图结果，数据来源：IOCADS

+ g1 Q4 [8 A2 M

& _. Q" r" `! [! u* Q$ _