手上刚好整理过一批 NSCAT 3 级风场数据,翻了翻当时的处理脚本和踩坑记录,干脆把这些东西梳理成一篇完整的说明。无论你是刚接触散射计数据、准备做海面风场分析,还是想在论文里引用这套历史数据,这篇内容应该能帮你少走不少弯路。
1. 这套数据是什么来头:NSCAT 任务背景与数据定位
NSCAT(NASA Scatterometer,美国宇航局散射计)搭载在1996年8月发射的日本先进地球观测卫星(ADEOS-I)上,是当时全球第一颗业务化运行的Ku波段散射计。它通过测量海面粗糙度引起的雷达后向散射系数变化,反演海面10米高度处的风矢量,覆盖全球无冰海洋面积约78%。1997年6月,ADEOS-I卫星因太阳能帆板故障导致整星失效,NSCAT实际有效观测时间只有大约9个多月。虽然寿命短,但它为后续QuikSCAT、ASCAT等散射计提供了关键的技术验证,也是1990年代中期唯一能够提供高分辨率全球海洋风场的星载传感器。
JPL(喷气推进实验室)基于NSCAT L2B轨道风矢量产品,制作了这套3级每日网格化浏览图像。它本质上是一种快速可视化产品,把不规则分布的轨道风场重采样到等经纬度网格,生成每日全球风场分布图,输出为图像和NetCDF数据两种形式。因为NSCAT的观测窗口比较特殊,而且L2B数据使用门槛较高,很多利用NSCAT研究热带气旋、海洋环流、海气相互作用的论文都会直接采用这套L3产品。
你需要明确一个关键点:NSCAT L3浏览图像的定位是"快速查看"和"统计分析基础"。它的空间分辨率是1°×1°,时间分辨率是每日一双(升轨和降轨分别合成),适合分析天气尺度以上的风场特征,不适合做中尺度锋面或台风内部精细结构的个例研究。后面我会详细解释为什么网格化之后分辨率会变粗,以及哪些应用场景下你需要回到L2B。
数据覆盖时间从1996年9月15日到1997年6月29日,每日一个文件,全球范围覆盖。下载途径一般是NASA物理海洋学分布式数据档案中心(PO.DAAC),现在通过Earthdata Search界面检索"NS CAT Level 3"即可找到。老用户在1998-2005年间常用的FTP路径已经发生了变化,HTTPS和OPeNDAP接口现在是主要的数据访问方式。
2. 文件里到底有什么:NetCDF结构、变量含义与元数据解读
这套3级产品最常用的格式是NetCDF,文件命名大致遵循s-nscat-l3-daily-YYYYMMDD-...nc的形式。打开文件后,你首先会看到维度定义和几个核心变量。对于不熟悉NetCDF的读者,可以把它理解成一种自描述的数组容器,变量名、单位、缺失值都写在文件内部,比纯二进制和ASCII要友好得多。
核心变量大致如下:
| 变量名 | 维度 | 单位 | 含义 |
|---|---|---|---|
wind_speed | (lat, lon) | m/s | 海面10米高度风速 |
wind_direction | (lat, lon) | degree | 风的来向,气象学角度 |
u_wind/v_wind | (lat, lon) | m/s | 风矢量的纬向和经向分量 |
n_measurements | (lat, lon) | count | 参与平均的观测次数 |
time | 1 | days since 1996-01-01 | 参考时间 |
latitude/longitude | 1D | degree | 网格中心坐标 |
quality_flag | (lat, lon) | - | 质量标记 |
rain_flag | (lat, lon) | - | 降雨影响标记 |
其中u_wind和v_wind是最便于直接使用的变量,因为风场分析经常涉及散度、涡度、通量计算,直接用分量比风速风向转换方便得多。wind_direction是气象学中的来向定义,0度表示北风(风从北往南吹),90度表示东风。做风矢量图的时候要注意,在数学坐标系中绘制时,需要做u=-wind_speed*sin(dir*pi/180)或者直接用已有的u_wind、v_wind,不要自己拿速度和方向去转,否则方向和速度对应不上很容易出错。
n_measurements这个变量很多人会忽略,但它非常有用。它表示某网格内有多少个有效观测值参与了平均,数值越高,该格点的风场可信度越高。如果某个网格的n_measurements等于0,说明没有升轨或降轨数据覆盖。需要注意的是,NSCAT的一个轨道覆盖带宽度约600公里,相邻轨道之间在中低纬度存在缝隙,所以不是每天每个格点都有数据。在做逐日图时,这些n_measurements=0的区域要设置成NaN,不能当作缺测风场处理,更不能把0当成风速值来画图。
元数据部分记录了数据生成的版本号、输入L2B产品版本、网格化算法简介、以及JPL数据处理团队的联系方式。早期版本(如V2、V3)与后续重处理版本(如V4)之间,风速和风向存在细微的系统性差异。建议你在论文方法部分明确写出你使用的是哪个版本,因为不同版本的风速均值差异最多可达0.3-0.5 m/s,虽然看起来小,但在气候统计分析里足以改变趋势和显著性结论。
文件内部还包含了一些辅助质量信息,包括每个网格使用的观测时间窗、平均后向散射系数等,但这些字段在不同发布版本中差异较大,需要查阅当时的数据说明文档(README)来确定具体含义。早期JPL发布的数据说明文档往往以PDF形式附带在数据目录中,比NetCDF全局属性更详细。
3. 读取与预处理:从NetCDF到可分析数据的关键步骤
拿到NetCDF文件后,最直接的处理方式是使用Python的xarray库。不仅因为xarray对NetCDF的维度、变量、坐标处理非常优雅,还因为它内置了resample和groupby方法,便于做多日合成和气候态计算。下面是一个典型的读取和预处理流程,也是我实际处理这套数据时的标准操作。
3.1 读取并生成掩码
import xarray as xr import numpy as np ds = xr.open_dataset('s-nscat-l3-daily-19961001-...nc') lat = ds['latitude'].values lon = ds['longitude'].values # 读取U/V分量 u = ds['u_wind'].squeeze() v = ds['v_wind'].squeeze() speed = np.sqrt(u**2 + v**2) # 质量与降雨标记 qual = ds['quality_flag'].squeeze() rain = ds['rain_flag'].squeeze() nmeas = ds['n_measurements'].squeeze() # 有效数据掩码:质量标记合格、无降雨、观测数大于0 mask = (qual == 0) & (rain == 0) & (nmeas > 0) u = u.where(mask) v = v.where(mask) speed = speed.where(mask)这里有几个重要的细节需要说明。第一,quality_flag==0通常表示质量合格,非0值表示数据受各种因素污染,包括陆地污染、海冰误判、极端风速不可信等。如果你忽略标记直接使用所有数据,沿海岸线和高纬度地区会出现很多高风速的奇异点。第二,rain_flag对于Ku波段散射计非常重要。Ku波段(14 GHz)信号对雨滴敏感,降雨会显著改变海面粗糙度,导致反演风速严重偏高(尤其是大雨条件下,风速可能被高估5-10 m/s)。rain_flag标记了受降雨影响的观测,在气候态统计分析中应该剔除。第三,nmeas>0这个条件看上去简单,却是最容易翻车的地方。有些处理工具在网格化填充时会把空值填为0,如果你不设这个条件,直接把所有风速为0的区域当成"静风",在全球图上就会看到轨道缝隙处出现大片0风速带,这完全不符合实际海况。
3.2 建立陆地/海冰掩码
虽然NSCAT L3产品已经标注了大部分陆地像元,但近岸网格受陆地回波污染的情况仍会出现。一种常见做法是使用NSCAT L3自带的land_mask或sea_ice_flag变量(如果有),另一种是使用外部海岸线数据生成掩码。
# 如果文件内没有掩码变量,可以用外部海岸线生成 # 这里使用Natural Earth的1:110m海岸线示意 import geopandas as gpd from shapely.geometry import Point world = gpd.read_file('ne_110m_land.shp') lon_2d, lat_2d = np.meshgrid(lon, lat) mask_land = np.zeros(lon_2d.shape, dtype=bool) for i in range(lat_2d.shape[0]): for j in range(lat_2d.shape[1]): pt = Point(lon_2d[i,j], lat_2d[i,j]) mask_land[i,j] = world.contains(pt)不过这种双重循环在1°网格上跑全球(180×360)倒还好,但如果以后处理0.25°数据,就会慢得让人崩溃。更好的做法是用regionmask库,或者用cartopy的feature接口先生成岸线栅格。例如:
import cartopy.crs as ccrs import cartopy.feature as cfeature from cartopy.io import shapereader # 用cartopy在1°网格上生成陆地掩码 import xarray as xr land_110 = cfeature.NaturalEarthFeature( 'physical', 'land', '110m', edgecolor='face', facecolor='none') lon_b, lat_b = np.meshgrid(lon, lat) land_mask = np.zeros(lon_b.shape, dtype=bool) # 这里示意用经纬度范围初筛,实际可用shapely的contains # regionmask库一行即可完成,推荐 import regionmask mask_region = regionmask.defined_regions.natural_earth_v5_0_0.land_110() land_mask = mask_region.mask(lon_2d, lat_2d).values land_mask = ~np.isnan(land_mask) # True表示陆地我的经验是:对于NSCAT这种1°浏览产品,用文件自带的land_flag(如果有)加上一个简单的近岸缓冲(比如距离海岸线0.5°以内的格点全部剔除)就够用了,不需要精确到像元的陆地掩码。原因是L3本身的网格化分辨率有限,海岸附近的风场反演误差远大于分辨率损失,你在海岸线附近看到的"强风"很可能是陆地回波泄漏导致的假信号。
3.3 升轨与降轨数据的区分使用
NSCAT的每日L3文件通常包含升轨(ascending)和降轨(descending)的分别合成变量,或者通过orbit_flag标记不同轨道方向。为什么需要区分?因为升轨和降轨的观测地方时不同。对NSCAT来说,ADEOS卫星的轨道设计使升轨大约在地方时10:30左右经过赤道,降轨大约在22:30左右。海面风场具有明显的日变化(尤其是海陆风影响区域和热带对流活动区),如果把两个时刻的数据混在一起做日均,会引入至少1 m/s量级的内噪声。
在1996-1997年的NSCAT L3版本中,某些时段文件里升轨和降轨变量是分开命名的,比如u_wind_asc、u_wind_desc。在做日平均风场时,我通常先分别检查两个方向风速的差异,如果差异超过2 m/s的区域比例大于10%,说明当天天气系统活跃,日均值需要谨慎解读。实际做气候态时,建议升轨和降轨分别计算平均值,再对两个平均场做算术平均,这样等效于对日变化做了最简单的订正。
4. 画风场图的实际操作:一张可用的风矢量图是怎样生成的
对于海洋风场来说,一张好的浏览图像应当同时表达风速大小和方向。JPL原始浏览图像一般以颜色填充表示风速,叠加箭头表示风矢量方向。下面以Python的matplotlib和cartopy为例,给出一个可以复制使用的绘图流程,并说明几个容易忽略的细节。
4.1 风场填色图的绘制逻辑
import matplotlib.pyplot as plt import cartopy.crs as ccrs import cartopy.feature as cfeature # 以某一天数据为例(示意,假设已经从xarray中取出u、v、speed) fig = plt.figure(figsize=(14, 6)) ax = plt.axes(projection=ccrs.PlateCarree()) ax.set_global() # 绘制风速填色 cf = ax.pcolormesh(lon, lat, speed, cmap='Spectral_r', vmin=0, vmax=20, shading='auto') # 叠加海岸线 ax.add_feature(cfeature.LAND, facecolor='lightgray') ax.add_feature(cfeature.COASTLINE, linewidth=0.5) # 叠加风矢量箭头(直接绘制U/V分量) # 为了清晰,每3个格点画一个箭头 q = ax.quiver(lon[::3], lat[::3], u[::3, ::3], v[::3, ::3], scale=400, width=0.002, alpha=0.7, transform=ccrs.PlateCarree()) cbar = fig.colorbar(cf, ax=ax, orientation='horizontal', pad=0.05, shrink=0.8) cbar.set_label('Wind speed (m/s)') ax.set_title('NSCAT L3 Daily Wind Vector 1996-10-01') plt.show()这段代码基本复现了JPL浏览图像的核心内容。使用Spectral_r的原因是这个colormap在0-20 m/s范围内的视觉层次分明,低风速偏蓝、高风速偏红,符合海洋气象学社区的习惯。如果你需要更接近JPL原始图像的配色,可以用ncl_default或viridis以外的一些气象专用colormap,比如cmcolor或cmocean的balance。
4.2 出图时容易被忽略的坐标陷阱
在绘制这套L3数据时,有一个非常常见的坑:文件里的longitude坐标可能是0到360,也可能是-180到180,不同版本不统一。如果你不做诊断直接画,从0到360的经度范围在PlateCarree投影上会导致地图中心位置偏移,图形看起来像被"拉断"了。建议读过文件后先检查lon.min()和lon.max(),如果大于180,利用xarray的assign_coords方法转换为-180到180,或者直接对lon进行((lon + 180) % 360) - 180变换。
另一个坑是关于quiver的投影参数。在cartopy中,当你用ax.quiver绘制U/V分量时,变换需要设置transform=ccrs.PlateCarree(),否则风矢量在非等经纬度投影下会指向错误的方位。NSCAT本身是全球覆盖的,我最常用的投影是PlateCarree(适合快速查看全球图)和Orthographic(适合极地视角)。如果做区域图,比如热带太平洋或北大西洋,建议用Mercator或LambertConformal,这时风矢量方向会随纬度变化,必须依赖cartopy正确处理矢量旋转,不要让箭头"飞"出图框。
4.3 多久采样一个箭头合适?
很多初学者一上来就把每个格点的风矢量都画上去,结果整个图形变成一片密密麻麻的黑色箭头,什么都看不清。1°网格的全球图有180×360个格点,每个格点上都画一个箭头显然不现实。我的经验是:全球图每4-5个格点取一个箭头,区域图每2个格点取一个箭头。取样的方式不要用np.arange简单隔点选取,而要确认取样后仍然覆盖了关键系统(比如气旋中心、急流带),必要时可以先用scipy.ndimage.uniform_filter对u/v做轻微平滑,再降采样。这样能保证箭头方向不是噪声。
对NSCAT L3来说,原始数据本身已经是日均产品,网格内的风场已经被平滑过,再做平滑对趋势分析影响不大,但要注意:不要对wind_speed填色图做过多的平滑,因为浏览图像的一个核心功能是呈现空间变率,如果你为了美观把风速场平滑到看不出锋面结构,那就失去意义了。
5. 版本坑与质量标记:为什么同一区域数值和别人对不上
使用NSCAT L3数据,最让人头疼的问题是你得到的数值可能和别人论文里的数值存在系统偏差。这不一定是操作错误,更可能是使用了不同的数据版本和不同的质量筛选策略。整理一下常见的原因,方便你排查。
5.1 网格化算法差异
NSCAT L3历史上经历过多次重处理。最早的V1版本把轨道数据简单平均到网格,V2引入了距离加权插值,V3加入了海冰和降雨标记,V4修正了风和海况状态的边界处理。不同版本的同一日平均风速,在赤道辐合带和高纬度西风带可以相差0.5-1 m/s。如果论文里没有写明确版本号,很难直接对比。
JPL在1998年发布的第一版官方文档中把L3产品描述为"为用户提供快速浏览能力",这让很多人误以为它是研究级产品。但从后续评估来看,L3在热带地区的风速偏差和均方根误差确实比同时段的L2B轨道产品要大。原因是网格化过程中的平滑效应会损失部分风速极值和梯度信息,而风场反演里最需要关心的恰恰是热带气旋、冷空气爆发这类极端事件。对于气候平均意义上的风速研究,L3完全够用;但对于台风个例和强对流系统的风场结构分析,还是应该回到L2B。
5.2 降雨标记是"推荐使用"还是"必须使用"?
我在前文已经提到rain_flag的重要性,这里再展开说说。Ku波段散射计发展早期,雨的影响曾经是一个令人头疼的问题。雨滴对雷达波的散射和衰减同时起作用,导致σ0显著升高。NSCAT的降雨标记是根据辐射计和散射计联合估计得到的,它不仅标记了大雨,也标记了中等降雨影响。数据文档中建议所有应用都使用降雨标记来排除数据,但实际操作中,很多人为了保留更多样本,会直接忽略这个标记。
我的建议是:做气候态统计时,必须剔除rain_flag=True的像元;在个例分析中,如果要保留,至少要在图上标出降雨影响区域,并合理评估这部分数据的可靠性。快速估算显示,热带地区的降雨影响像元比例可达10-20%,如果你不剔除,全球平均风速会被抬高约0.3-0.4 m/s。这个偏量看似小,但在研究年际变化和趋势时足以导致结论反转。
5.3 质量标记的不同位域含义
NSCAT L3的quality_flag是一个位掩码字段,不同数值组合对应不同的质量问题。有些版本用0表示好、1表示坏,有些版本则定义了一系列位。如果你只判断数值是否等于0,有可能会把"可用但质量较差"的像元全部剔除,也可能保留了"看起来很好但实际是陆地污染"的像元。最稳妥的方法是在数据集中找到质量标记的定义说明(通常在全局属性flag_meanings中),用位运算检查每一个位:
# 假设某版本flag第0位表示陆地污染,第1位表示海冰,第2位表示降雨 bad_land = (qual & 1) != 0 bad_ice = (qual & 2) != 0 bad_rain = (qual & 4) != 0对于NSCAT L3的多数版本,直接判断qual == 0是可行的简化操作,但在正式分析中,建议解析每一位的质量信息,至少把陆地污染和海冰误判单独提取出来,看看它们在空间上是否集中在某些区域。这样可以避免把系统误差当成信号。
6. 从浏览图像到实际研究:NSCAT L3的典型应用与统计预处理
不管你是做海洋环流、海气通量、还是极地科学,拿到NSCAT L3以后最常做的操作首先是验证数据合理性,然后按月或按季合成气候态,最后再做具体分析。这里分享两个层面的经验。
6.1 数据合理性验证
在信任任何数据集之前,先做基本的合理性检查。对L3风场,我会依次检查这几项:
- 全球风速平均值应大致在6-8 m/s之间(取决于季节和是否剔除降雨区域),如果全球日均风速低于5 m/s或者高于10 m/s,说明可能没有正确过滤陆地/冰面像元,或者大雨污染极其严重。
- 画一张任意一天的全球风速图,确认低纬度存在明显的信风带结构,中纬度西风带风速较强,极地东风带弱且破碎。
- 检查北极和南极海冰边缘区的风速,那里经常出现虚假的高风速(海冰表面回波和海面不同,散射计反演算法在海冰区域失效)。用MODIS或NSIDC海冰密集度数据做掩膜,能有效剥离这部分虚假信号。
- 和浮标数据做点对点比较时,要注意时空匹配窗口。NSCAT的观测瞬时值和浮标小时平均在风速上存在约1 m/s的差异,风向差异约20度(尤其是在低风速条件下)。散射计反演在4-6 m/s以下的低风速段风向误差很大,因为此时海面粗糙度主要由风浪而不是涌浪决定,方向信号弱。
6.2 月平均合成与缺测补插
NSCAT数据只有9个多月,跨了1996年9月到1997年6月,四季还算齐全,但没有覆盖完整年度周期。要做准气候态,可以把同一月份的若干天数据做平均,比如1996年10月1-31日和1996年10月的日平均再平均。但因为轨道空隙的存在,某些格点一个月内可能只有一半的天数有有效观测,直接平均会造成偏差。建议在做逐月平均前,先统计每个格点的有效天数,有效天数少于10天的格点直接标记为缺测。
缺测区域的补插是另一个大坑。全球L3网格上的空值带(轨道缝隙)是有规律分布的,中低纬度地区每天都存在。如果你对全球风场做EOF分析或计算散度,这些空值会导致结果出现条带状虚假信号。常见做法是使用最优插值(OI)或克里金插值,但我不建议在逐日尺度上做复杂的空间补插。原因很简单,NSCAT数据的轨道缝隙内可能确实存在中小尺度天气系统,而你用周边插值补出来的风场,本质上是一种平滑猜测,在后续物理量计算中可能放大噪声。
如果你的目标是一个月的平均风场,更好的做法是先把每天的U/V网格数据放到一个三维数组中,对每个格点的时间序列求平均,同时保留有效天数。不同区域的有效天数差异很大(极地轨道重叠多,有效天数多;赤道附近轨道相邻间距大,有效天数少),直接用平均场绘图的可靠性在不同纬度带差别很大。出图时把有效天数也画出来,能让读者对可信度有直观感受。
6.3 NSCAT与后续散射计数据的时间衔接
NSCAT失效后,直到1999年QuikSCAT发射,中间约两年的全球散射计风场存在空缺。如果你想做更长时序的风场变化分析,必然要把NSCAT和QuikSCAT衔接起来。这里有个明显的系统偏差问题:NSCAT和QuikSCAT虽然都是Ku波段散射计,但天线配置、入射角范围、反演算法版本不同,两者在相同海域的风速存在约0.3-0.5 m/s的偏差,风向偏差小于5度。直接拼接使用会导致1996-1997年和1999-2009年的风速序列出现跳变。
解决方法是使用重叠期数据做统计校准。NSCAT和QuikSCAT实际没有重叠观测时段,因此无法直接做同平台交叉定标。常用的校准桥梁是NDBC浮标和TAO/TRITON浮标阵列。你可以在两个时间段内,分别把卫星数据和浮标数据做回归,然后通过浮标这个"公共参考"把NSCAT调整到QuikSCAT的基准。这个两步校正的过程会累积一定误差,但比直接拼接要可信得多。
此外,NSCAT的一个独特价值在于它能同时提供VV和HH两种极化的背向散射系数。这在后续的QuikSCAT(只有VV极化)和ASCAT(VV极化为主)上是没有的。不过L3浏览产品并不直接包含极化信息,如果你需要研究极化比模型,应该直接使用NSCAT L1.7或L2A级数据。在L3尺度上,NSCAT更多被用来提供海面风场的大尺度分布特征。
7. 从产品文档到实测效果的几个结论
说回到JPL原始浏览图像本身。JPL在发布L3浏览图像时,以彩色图形式直观展示每日全球风场。这种图像适合新闻稿、科普和快速检查,但正式研究不能把JPG图片作为数据源。原因是JPEG压缩会损失色标精度,每一个色标对应的风速区间在0.5-1 m/s量级,同时图像本身不携带精确的空间坐标和元数据。正确做法是直接下载NetCDF数据,自己在本地绘制。
另一个容易被忽视的问题是:PODAAC网站上的NSCAT L3产品在经过系统迁移后,部分旧版文件的命名和变量结构发生了细微变化。建议下载文件后立刻用ncdump -h查看变量列表,确认你要用的变量确实存在。例如,早期版本的wind_speed在文件里可能叫ws,后来统一改成wind_speed。如果直接套用旧脚本,很容易报KeyError。
最后想提醒一点:NSCAT这套浏览图像虽然"只"有9个月数据,但它捕捉到了1997-1998年厄尔尼诺事件发展初期热带太平洋风场异常的演变过程。结合1997年后期的浮标和再分析资料,NSCAT风场数据仍然是研究这次强厄尔尼诺事件中大气海洋响应不可多得的观测资料。如果你正在研究那个时段的海气过程,花点时间把这套L3数据处理干净,往往能得到一些再分析资料里被平滑掉的真实信号。