Matlab海洋绘图：海面热通量分布——以南海为例 - 海洋地质热流测量

图片赏析

1 H3 ^/ b s3 r. b0 t% U* B0 B8 ~

前言

- {" L7 l/ d$ S, l1 g& l. p

+ l9 h5 w1 W/ y0 h. H" s

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

数据来源：

% G* L( T; f8 p2 n

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

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

% a7 Z+ j( F! n5 j; H) l: i

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

; E. J- h+ y3 L( B" } a

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

7 i1 @ R, I2 o0 S0 M. f

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

空间精度：0.25°

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

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

" m5 D6 |7 B' k& w4 Y4 D

源代码

0 R. [2 D( @9 L7 z* b! S

%南海1998年12月31日某时的海面热通量数据及分布图 %% ( h# R! k% J' m9 } j6 A! ^6 \ %第一步：数据的展示 tic;Path="E:\OneDrive - Inc\Share_science\NC文件\";%文件夹地址 6 Z6 i" E. o, f0 Q" S3 Z9 @ @ filename="SEAFLUX-OSB-CDR_V02R00_FLUX_D19981231_C20160824.nc";%文件名 2 K% C6 K3 l) r1 f* E; M8 { abs_path=fullfile(Path,filename);%合并，绝对路径文件 clearvars -except abs_path;%清除除了abs_path以外所有变量 & Q1 w, x( w W ncdisp(abs_path);%读取所选文件 , u6 p" P$ U# f8 M8 C3 a- p whos;toc; & N8 J; D. t6 y# g %% 3 D* U$ q7 d2 n. h. [9 ?7 ? %第二步：数据读取 tic;Lat_ini=ncread(abs_path,lat); ( y3 V9 L" N/ w* J8 K" L! ~ Lon_ini=ncread(abs_path,lon); Time_ini=ncread(abs_path,time); 0 ^6 T3 f( \) l; s, b; n Lflx_ini=ncread(abs_path,surface_upward_latent_heat_flux); Sflx_ini=ncread(abs_path,surface_upward_sensible_heat_flux); 4 E' _4 d8 T7 h4 a9 ]8 T# b7 A toc;clear abs_path; %% %第三步：确定边界/截取数据 0 \4 x* y% G+ k9 W- ]6 ] %空间边界：南海范围北纬2°至北纬25° 东经104°至东经124° tic;Boundary=[104,124;2,25;]%修改范围[经度;纬度],注意顺序 " k T' c# e; N- N# P Spatial_Resolution=abs(Lat_ini(1)-Lat_ini(2));Sp_Re=Spatial_Resolution;%自动确定精度 * l z. C. l) V9 J& R: F- J Lon_start_No=find(Lon_ini-(Boundary(1))<=Sp_Re&Lon_ini-(Boundary(1))>=0);%找到离104°最近格点的编号（行号） ( e; q/ k, W# c7 b/ G- l' \7 G Lon_end_No=find(Lon_ini-(Boundary(3))<=Sp_Re&Lon_ini-(Boundary(3))>=0);%找到离124°最近格点的编号（行号） Lat_start_No=find(Lat_ini-(Boundary(2))<=Sp_Re&Lat_ini-(Boundary(2))>=0);%找到离2°最近格点的编号（行号） & G0 h& W: W3 k- a- ^( o t. j Lat_end_No=find(Lat_ini-(Boundary(4))<=Sp_Re&Lat_ini-(Boundary(4))>=0);%找到离25°最近格点的编号（行号） Lon=Lon_ini(Lon_start_No:Lon_end_No); Lat=Lat_ini(Lat_start_No:Lat_end_No);%北半球可能要颠倒一下编号，因为是从北极往赤道计数的 clear Lat_ini;clear Lon_ini;%回收内存 Lat=double(Lat);Lon=double(Lon);%双精度化 % @0 j! y5 ] U latlim=[min(Lat),max(Lat)];lonlim=[min(Lon),max(Lon)]; Time=9%输入时间，注意三小时的倍数 Time_No=Time/3+1; %截取数据并降维度 2 G b1 T5 K8 Q, D) l3 `; m7 L' k Lflx=squeeze(Lflx_ini(Lon_start_No:Lon_end_No,Lat_start_No:Lat_end_No,Time_No));;%北半球有时可能要颠倒一下编号，因为是从北极往赤道计数的 + t' ^) l M$ k+ y/ g8 E Sflx=squeeze(Sflx_ini(Lon_start_No:Lon_end_No,Lat_start_No:Lat_end_No,Time_No));;%北半球有时可能要颠倒一下编号，因为是从北极往赤道计数的 clear Sp_Re;toc;%清除不再使用的变量 %% %第四步：绘图 4 z3 s2 t2 i1 c" S, u9 s. t tic;figure; m_proj(mercator,lat,latlim,lon,lonlim); - j; o9 p$ z6 N( p& v [longrid,latgrid]=meshgrid(Lon,Lat); & l0 `9 e* J, _1 B' C m_pcolor(longrid,latgrid,Sflx); %pcolor或contourf tic;m_gshhs_f(patch,k);toc;%添加海岸线：c<l<h<f % Y" n* h, \# z% K% e/ g1 z+ q m_grid(box,fancy,tickdir,in);%网格化 %brighten(.5);%亮化 1 Z/ E/ S: S% w, B3 \ colormap(hsv); h = colorbar(h);%色标 h.Label.String = W·m^{-2} h.Location = eastoutside;%色标位置 title([1998年12月31日,num2str(Time),:00,南海海面感热通量],... FontName,黑体,fontsize,12,FontWeight,"bold") %标题 & i& s) m0 [; c; }; _, v3 T9 ?. i xlabel(Longitude);ylabel(Latitude) ;set(gca,fontsize,9); toc; + M3 y5 N; T B

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

+ V0 `, X' m* d