企业官网建设流程全解析

简介：《MODIS数据综合处理软件 V1.0》是一套面向遥感、GIS及生态环境研究人员的桌面工具，集成MODIS陆面产品批量处理、多核并行计算与可视化分析功能，支持NDVI、EVI、ET/PET、LST、LAI、GPP/NPP等常用产品，适用于科研数据预处理、资源环境监测与教学演示等场景。软件附详细使用手册，安装包与手册打包为一个zip压缩包，共2000个文件；文件类型以JavaScript、JSON、HTML等前端程序为主，兼有Python脚本、PDF文档、CSS样式、XML配置及Markdown笔记，能覆盖软件运行、用户查阅与二次开发等多种用途，压缩包整体176.28MB，目录结构清晰，解压即可部署。已有140人学习下载。借助该软件，用户无需编写复杂代码即可完成多文件批量导入、一键处理与结果出图；多核并行机制可显著缩短大规模数据计算时间，可视化工具则便于快速把握地表参数的时空变化趋势，显著提升MODIS数据处理与制图效率，是相关领域数据分析流程中一款实用型辅助软件。

1. 用 MODIS 数据综合处理软件之前，先回答“能不能直接用”这个问题

把“MODIS数据综合处理软件 V1.0”这串字拆开看，它不像某个商业软件那么神秘——它想解决的，是每个用过 MODIS 数据的人都撞过的墙：从 NASA 下载下来的 HDF 文件，用 GIS 打开一片黑，或者看到一个波段根本不知道是什么投影。我最早拿到 MOD11A1 地表温度产品时，直接把它当普通栅格丢进 ArcGIS，出来的图数值上万分，我一度以为卫星坏了，后来才意识到是没做综合处理。所谓“综合处理软件 V1.0”，在我的理解里就是一套把下载、拼接、重投影、裁剪、质量控制和批量导出的流程固化成脚本和参数模板的工具箱。它面向的人很具体：做地表温度、植被指数、水体或土地利用变化研究，又不想每次都在 GUI 里点半小时鼠标的从业者和研究生。这篇文章就按这个思路，从数据格式认知讲到落地代码，再讲清楚热词里那个灵魂拷问——MODIS 下载的地表温度数据到底能不能直接用。答案是：能，但绝对不能直接画图。

2. 先认清 MODIS 的“脾气”：HDF 格式、产品家族与正弦投影

2.1 HDF 文件结构和子数据集：为什么打开一个文件有十几个图层

MODIS 的 L2/L3 标准产品是 HDF4 格式，后缀 .hdf。HDF 本身是一个容器，里面可以装多个科学数据集（SDS），专业说法叫子数据集。你用普通软件打开一个 MOD11A1.hdf，看到的是“栅格”和“属性表”，但其实里面装了两类东西：一类是实际的地理数据，例如白天地表温度（LST_Day_1km）、夜间地表温度、质量控制波段（QC_Day）、观测时间、视角等；另一类是元数据，包括投影信息、尺度因子、有效值范围、填充值。这两个东西混在一个文件里，是新手最容易困惑的地方。

用 GDAL 的命令行工具看一眼就知道里面有什么：

gdalinfo MOD11A1.A2020001.h23v04.006.2020001030000.hdf

输出会列出类似这样的子数据集：

Subdatasets: SUBDATASET_1_NAME=HDF4_EOS:EOS_GRID:MOD11A1.A2020001.h23v04.006.2020001030000.hdf:MOD_Grid_Daily_1km_LST:LST_Day_1km SUBDATASET_2_NAME=HDF4_EOS:EOS_GRID:MOD11A1.A2020001.h23v04.006.2020001030000.hdf:MOD_Grid_Daily_1km_LST:QC_Day

这里每一行是一个独立子数据集。LST_Day_1km 是数值型地表温度灰度值，QC_Day 是 8 位无符号整数质量控制位域。我的习惯是先把这两层拆出来做预处理，而不是一上来就放到拼接工具里——因为 QC 在你做重采样和裁剪时会被插值污染，一旦污染，后续质量筛选全是错的。拆层的代码后面会说，这里先建立概念：所谓“处理一个 MODIS 文件”，第一步永远是“从 HDF 里把这个子数据集正确地抠出来”。

2.2 MODIS 产品家族：你下载的到底是什么东西

MODIS 数据有几十种产品，命名规则是 MOD（Terra）或 MYD（Aqua）开头，后接三位产品号。常见的一些做地表覆盖和气候研究会用到的产品，我整理如下：

区分这些产品很重要，因为“综合处理”的流程会根据产品不同而改变。比如 MOD13A1 本身是 16 天合成，处理时不需要再做时间平均；而 MOD11A1 是瞬时值，一天里 Terra 和 Aqua 各过境两次，处理时如果要把白天和晚上分开，就得同时读文件名里的时间字段和产品内部的时间 SDS。我在处理地表温度时，一般会把日产品和 8 天合成产品分开两条流水线，因为前者适合做极端温度事件分析，后者适合做气候态平均。

除了产品号，文件名里还有一个关键信息：分片编号。MODIS 标准产品采用正弦投影的分片网格（tile grid），整个地球被横竖切成若干块，文件名里 h23v04 这样的编号就代表横向第 23 片、纵向第 4 片。中国的陆地基本覆盖在 h23v04、h24v04、h25v04、h26v04、h27v04 等几片上，跨片研究区必须把所有覆盖到的片子都下载下来再做拼接，不能只下中间一片。

2.3 正弦投影是怎么回事：为什么 WGS84 坐标的人容易栽在这里

MODIS 标准产品的默认投影是正弦投影（Sinusoidal），GDAL 能识别，但 ArcGIS 的老版本和一些在线平台不一定认，打开后经常出现整个图层定位到海里或者影像拉伸变形的现象。正弦投影本质上是一种等积伪圆柱投影，赤道附近变形小，高纬度变形大。它在全球尺度做统计时面积保持得好，但到了区域尺度，想要和矢量边界、气象站点或者其他卫星资料叠加，就必须重投影到你自己的目标坐标系。

常见做法有两种：一是转成 WGS84 地理坐标系（EPSG:4326），适合做全球或大区域制图；二是转成 UTM（根据研究区经度选带号），适合做局部定量分析。我的经验是：做地表温度时间序列，用 WGS84 就够了，因为你不涉及面积算量；做森林覆盖或蒸散发涉及面积统计的，必须转 UTM 或者转 Albers 等积投影。这一步看着简单，但如果在拼接前不做，后面裁剪出来的数据就会出现边界错位。

所以，综合处理软件 V1.0 的第一版核心应该包含四件事：HDF 子数据集提取、拼接、重投影、裁剪。这四个步骤的顺序不能乱。先拼接再重投影，比先重投影再拼接稳定得多，因为原始 tile 之间是严格无缝的，一旦各自重投影，两个相邻 tile 的边界会因为重采样产生细缝和重叠，后面拼起来要处理羽化和接边线，属于给自己找罪受。

3. 用 GDAL 搭一条预处理流水线：从 HDF 抠数据到拼接与裁剪

3.1 下载与目录组织：给 V1.0 定一个不会被自己搞乱的目录结构

我一般在 NASA Earthdata Search 上下载数据，选定时间范围和 tile 号后，一次把覆盖研究区的所有 tile 都勾选下载。下载后先做文件归档，否则处理到一半你会发现某个日期的文件缺失，蹲在电脑前抓瞎。推荐的目录结构是：

MODIS/ raw/ # 原始 HDF 文件，按 tile 和时间子目录 lst_day/ # 拆分出的 LST 子数据集 lst_night/ qc/ mosaic/ # 拼接后的 GeoTIFF clipped/ # 按研究区裁剪后的结果 log/ # 处理日志

这个结构看着简单，但能避免一个很常见的翻车：不同时间的文件混在一起，跑了批量脚本后发现 2020 年第 200 天和 2021 年第 200 天拼在了一起。我的做法是把 HDF 文件名里的日期提取出来，写进输出文件名，格式统一为LST_YYYYDOY.tif，比如LST_2020200.tif，这样后面做时间序列排序时，按文件名排序就是按时间排序。

3.2 抠子数据集：gdal_translate 提取单波段

这一步是流水线的第一环。用 gdal_translate 把 HDF 里的特定子数据集拉出来，转成 GeoTIFF。这里要特别叮嘱：子数据集全名很长，每次手敲都会疯掉，我都是先跑一遍 gdalinfo 把名字复制出来，或者用脚本动态获取。

gdal_translate -of GTiff \ "HDF4_EOS:EOS_GRID:MOD11A1.A2020001.h23v04.006.2020001030000.hdf:MOD_Grid_Daily_1km_LST:LST_Day_1km" \ ./lst_day/LST_day_2020001.tif

这条命令做的事是：打开 HDF 文件定位到 LST_Day_1km 子数据集，写成单个 GeoTIFF。输出文件保留了原始投影信息和地理变换参数。参数说明：-of GTiff指定输出格式；坐标和投影信息默认继承源文件，所以此时如果直接打开，看到的还是正弦投影。我一般在这一步不做任何重采样，只为把数据从 HDF 容器里救出来，避免每次都要处理 HDF 嵌套结构。

3.3 按日期拼接：gdalwarp 还是 gdal_merge.py

拼接有很多工具，但处理 1km 分辨率的地表温度产品，我推荐用 gdalwarp，因为它能同时完成镶嵌和重投影，一步到位。gdal_merge.py 也可以拼，但它只做简单拼接，不处理投影不一致的问题，两片不同日期或不同轨道的影像只要投影参数有一丁点差异，拼完就会出现明显的接缝。

gdalwarp -t_srs EPSG:4326 -r bilinear -overwrite \ ./lst_day/LST_day_2020001_tile1.tif \ ./lst_day/LST_day_2020001_tile2.tif \ ./mosaic/LST_day_2020200_wgs84.tif

这里把覆盖研究区的两个 tile 拼接并重投影到 WGS84。-t_srs EPSG:4326是目标坐标系，-r bilinear是重采样方法。地表温度这种连续变量，用双线性插值比最近邻更平滑；但如果你在处理分类产品或者需要严格保留原始像元值的场景，比如火点检测，要改用-r near。重采样方法选错是新手常踩的坑，后面避坑章节还会展开。

再看一遍这条命令的逻辑：输入是两个 tile 的白天 LST GeoTIFF，输出是拼接好的 WGS84 GeoTIFF。gdalwarp 会自动读取每个文件的覆盖范围，找到重叠区，按默认策略做叠加。对于 MOD11A1 这种每天多轨数据，有时同一天相邻轨道的重叠区数值不一致，我的做法是拼完后用 QC 做掩膜后处理，而不是在拼接时加羽化参数。羽化会让真正的高温异常被抹平，得不偿失。

3.4 按矢量边界裁剪：把数据切到研究区

拼接完成了，如果研究区只是一个城市或者一个流域，下一步是裁剪。最常见的方式是用一个 shapefile 边界做掩膜裁剪：

gdalwarp -cutline ./shp/study_area.shp -crop_to_cutline \ -t_srs EPSG:4326 -dstnodata -9999 \ ./mosaic/LST_day_2020200_wgs84.tif \ ./clipped/LST_day_2020200_study.tif

-cutline指定边界矢量，-crop_to_cutline让输出范围严格跟随边界范围。这里有一个细节：边界矢量的坐标系最好和目标栅格一致，如果不一致，gdalwarp 会自动做坐标变换，但你要确保 shapefile 里定义了正确的投影信息，否则裁剪结果可能在研究区外或者完全空白。

到这一步，一个 MODIS 文件已经从 HDF 变成了一个带正确投影、按边界裁剪好的 GeoTIFF。但这才完成“综合处理”的一半，因为还没有做质量控制。而缺少质量控制的地表温度数据，就是热词里问到的“能不能直接用”的根源。

4. 地表温度数据到底能不能直接用：质量控制、尺度因子与无效值

4.1 灰度值 ≠ 真实温度：0.02 这个尺度因子是怎么用的

先直接回答热词那个问题：modis 下载地表温度数据可以直接用吗？答案是“不可以直接画图”，原因有两个。第一，LST 产品里存储的是灰度值（DN 值），不是物理温度。MOD11A1 的白天地表温度波段，灰度值乘以 0.02 才得到开尔文温度，再减去 273.15 才等于摄氏度。第二，影像里不全是有效观测，有云遮挡、有填充值、有质量不佳的像元。如果你不处理这两件事，做出来的地表温度图会同时包含大量异常低值（云顶温度被当成地表温度）和诡异的高值。

读取 LST 灰度值的脚本我一般这样写：

from osgeo import gdal import numpy as np # 打开 HDF 中的 LST 子数据集 lst_ds = gdal.Open( "HDF4_EOS:EOS_GRID:MOD11A1.A2020001.h23v04.006.2020001030000.hdf:" "MOD_Grid_Daily_1km_LST:LST_Day_1km" ) lst_dn = lst_ds.ReadAsArray().astype(np.float32) # 尺度因子 0.02，单位是开尔文 lst_k = lst_dn * 0.02 lst_c = lst_k - 273.15

这里lst_dn是原始灰度值，lst_c是摄氏温度。有一个易错点：填充值 0 乘以 0.02 等于 0，0 再减 273.15 成了 -273.15，这个值会严重拉低后续统计的最低值。所以必须先剔除填充值再算物理量。正确的顺序是：先找填充值和无效值，再乘尺度因子，最后才做后续运算。

4.2 QC 波段：质量不是用眼睛看的，是拿 bit 算出来的

MOD11A1 质量控制波段（QC_Day）是 8 位无符号整数，它不是一个“0 到 255 的等级分数”，而是每两位 bit 记录一种质量信息的位域。最常见的做法是提取 bit0 和 bit1 判断综合质量等级：0 表示质量好，1 表示质量一般，2 表示云遮挡，3 表示云阴影。很多人拿到 QC 后直接做阈值筛选，比如保留 QC < 64，这是完全错误的，因为 QC=64 的二进制是 01000000，bit0-1 是 0（质量好），却被误杀了。

正确的 QC 筛选方式是用位运算：

from osgeo import gdal import numpy as np qc_ds = gdal.Open( "HDF4_EOS:EOS_GRID:MOD11A1.A2020001.h23v04.006.2020001030000.hdf:" "MOD_Grid_Daily_1km_LST:QC_Day" ) qc = qc_ds.ReadAsArray() # bit0-1: 0=良好, 1=一般, 2=云, 3=云阴影 quality = qc & 0b11 valid_mask = (quality == 0) | (quality == 1)

取qc & 0b11位与运算，得到低两位的十进制值。这样好质量的像元，即使它的高比特位有其他信息（比如云检测标记、日/夜标记），也不会被误过滤。这一点是地表温度处理里的分水岭，理解它之后，你再去读其他 MODIS 产品的 QA 文档，会发现套路都一样。

4.3 拼起来的完整处理：温度、质量、掩膜一次搞定

把上面几段组合成一个完整的最小流程：

from osgeo import gdal import numpy as np product_path = ( "HDF4_EOS:EOS_GRID:MOD11A1.A2020001.h23v04.006.2020001030000.hdf:" "MOD_Grid_Daily_1km_LST:" ) lst_dn = gdal.Open(product_path + "LST_Day_1km").ReadAsArray().astype(np.float32) qc = gdal.Open(product_path + "QC_Day").ReadAsArray() # 第一步：根据 QC 生成有效像元掩膜 good_pixels = ((qc & 0b11) == 0) | ((qc & 0b11) == 1) # 第二步：剔除填充值（DN 值为 0 是无效观测） good_pixels &= (lst_dn != 0) # 第三步：尺度因子换算 + 无效值抑制 lst_k = np.where(good_pixels, lst_dn * 0.02, np.nan) lst_c = lst_k - 273.15 # 统计去看一眼，确认数值范围合理 print("温度范围（摄氏度）：", np.nanmin(lst_c), "~", np.nanmax(lst_c))

过程的逻辑是：先算质量控制得到哪些像元可信，再剔掉填充值，最后才做单位换算。顺序不能错，如果先换算再筛选，填充值 0 会污染温度统计。np.where在这里把无效像元直接置为np.nan，后续做任何统计运算都不会再被它干扰。

到这一步，地表温度数据才算“能用了”。数据科学家经常说 Garbage in garbage out，在 MODIS 这里，不过滤 QC 就是 Garbage in garbage out 的说明书级案例。用直方图检查一下处理后的lst_c：正常的 LST 白昼应该在 -20 到 60 摄氏度之间，如果出现大量负值，往往是云遮挡没滤干净；如果出现大量超过 60 摄氏度的像元，先检查是不是裸土或者火点，再检查尺度因子是不是漏乘了。

5. 综合处理中的 4 个常见坑：现象、原因和排查顺序

5.1 拼接后的影像出现明显的十字接缝或暗带

现象：用 gdalwarp 拼完两个 tile 后，影像中间出现一条清晰的分界线，一边亮一边暗，或者在重叠区出现几何错位。

原因：这种情况绝大多数不是算法问题，而是输入的两片影像日期不一致。MODIS 逐日产品每天由多次过境拼接而成，云覆盖不同，地表温度也会随时段变化。你下数据时选了同一天，但下载到的是不同轨道来源的 tile，它们虽然日期相同，但过境时刻差了一个多小时，正午地表温度变化很快，拼在一起自然有色差。

解决：先查两个 tile 的元数据，确认过境时间（文件名里有精确到秒的时间字段）。处理地表温度时，同一天的 Terra 白天数据如果过境时间相差超过 30 分钟，我一般直接放弃拼接，改用 MOD11A2 的 8 天合成产品，或者把研究区限定在单一过境覆盖范围内。拼接前额外检查一遍文件时间，比拼完再调试快得多。

5.2 QC 波段用阈值筛选后，好的像元全被删了

现象：用qc < 64之类的条件筛完，影像变成大片空洞，连夏天晴空万里的时候都是碎的。

原因：QC 是按位域存储的，不是数值越高质量越差。MOD11A1 的 QC_Day 为 8 位，高 bit 位存放云检测标志、邻近云等信息，数值大小和质量好坏没有单调关系。我见过有人用qc == 0筛数据，结果影像只剩零星十几个像元，因为大部分像元的 bit2-3 或 bit4-5 并不为 0，但它们本身不代表质量问题。

解决：永远用位与运算。qc & 0b11取低两位，再判断是否为 0 或 1。如果要做更严格的质量控制，再去看产品文档里 bit2-3 和 bit4-5 的定义，不要用十进制阈值。

5.3 温度图像数值在几万甚至十几万，明显不是温度

现象：做完处理，输出影像的像元值在 10000 到 30000 之间，画出来像是外星地貌。

原因：忘了乘尺度因子，直接把 DN 值当物理量使用了。MOD11A1 的 LST 灰度值通常在 7500 到 14000 之间，乘 0.02 后才是 150 到 280 开尔文。这是热词“能不能直接用”里最典型的错误。

解决：乘尺度因子是硬性步骤，写进任何脚本的前三行。我的习惯是把scale = 0.02初始化在脚本顶部，并加注释说明它是 MOD11A1 的固定属性。处理完输出后，立刻用np.nanmin和np.nanmax打印数值范围，范围不对就停下来查，不要等到画图才暴露。

5.4 按边界裁剪后影像全部空白或者偏移了几百公里

现象：用-cutline裁剪完，输出的 GeoTIFF 打开是一片黑，或者数据在研究区边界外，像是一块被平移过的影像。

原因：两个坐标参考不一致。最常见的是栅格还是正弦投影，而 shapefile 是 WGS84 经纬度，虽然 gdalwarp 会自动做投影变换，但如果 shapefile 缺少.prj文件，GDAL 不知道它的坐标系，就会按经纬度数值直接当坐标用，结果裁剪范围跑到非洲或者大海。

解决：先跑gdalinfo shapefile.shp确认矢量有投影信息，若没有.prj，在 GIS 软件里给 shapefile 手动指定 WGS84。裁剪前，先gdalinfo看一眼栅格范围和矢量范围是否大致吻合，在同一个经纬度体系下，两者的 extent 应该存在重叠。这一步只要养成习惯，几乎可以杜绝裁剪空白的坑。

5.5 时间序列里某一天的温度断崖式偏低，其他天都正常

现象：做 2020 年逐日 LST 曲线，大部分日期在 20 到 35 度之间，某一天突然掉到 -5 度，画出来的曲线像心电图。

原因：那天的研究区恰好被云覆盖，而且 QC 筛选条件太宽（保留了 quality==1 或 2 的数据）。MOD11A1 的 QC 位域里，云掩膜信息在高 bit 位，仅用qc & 0b11会保留“一般质量”的像元，而这些像元里很多是云边缘。

解决：做时间序列时，质量控制收紧到quality == 0再加一道云掩膜检查。具体做法是读 MOD11A1 附属的云掩膜波段（如果产品提供），或者直接用 MOD35 云产品做二次过滤。我处理长时间序列时还有个习惯：对每个日期计算云覆盖比例，云覆盖超过 40% 的日期直接标记为缺测，不参与平均和统计。宁可少一天数据，也不要一个假值混进序列里。

6. 把 V1.0 变成自己的批处理骨架：脚本模板与验证方法

到这一步，你已经有了完整的单日处理链路。把它固化成批处理骨架，才算得上“软件”。我通常会把处理流程抽象成一个 Python 脚本，输入是 HDF 文件列表，输出是裁剪后的温标 GeoTIFF 加一张 QC 概览图。骨架结构如下：

from osgeo import gdal import numpy as np import glob def process_one_day(hdf_path, output_dir, cutline_shp): # 1. 列出 HDF 内所有子数据集 ds = gdal.Open(hdf_path) subdatasets = ds.GetSubDatasets() # 2. 按名称定位 LST 和 QC 层 def find_sds(keyword): for name, desc in subdatasets: if keyword in name: return name return None lst_name = find_sds("LST_Day_1km") qc_name = find_sds("QC_Day") # 3. 读取数据并做 QC 筛选与尺度换算 lst_dn = gdal.Open(lst_name).ReadAsArray().astype(np.float32) qc = gdal.Open(qc_name).ReadAsArray() good = ((qc & 0b11) == 0) & (lst_dn != 0) lst_c = np.where(good, lst_dn * 0.02 - 273.15, np.nan) # 4. 输出 GeoTIFF（原始投影，后面统一交给 gdalwarp 处理） driver = gdal.GetDriverByName("GTiff") out_path = f"{output_dir}/LST_C_{lst_name.split(':')[-1]}.tif" out_ds = driver.Create(out_path, lst_c.shape[1], lst_c.shape[0], 1, gdal.GDT_Float32) out_ds.GetRasterBand(1).WriteArray(lst_c) out_ds.SetGeoTransform(ds.GetGeoTransform()) out_ds.SetProjection(ds.GetProjection()) out_ds.FlushCache() return out_path # 批量入口：遍历某个日期的所有 tile for hdf in sorted(glob.glob("./raw/MOD11A1.A2020200.*.hdf")): process_one_day(hdf, "./lst_day", "./shp/study_area.shp")

这段脚本把每天多个 tile 从 HDF 中拆出并换算成摄氏温度，输出文件仍保留正弦投影。后续再用 gdalwarp 统一做投影转换和批量拼接。脚本里的find_sds函数解决了一个很实际的问题：HDF 子数据集全名在不同版本的 MODIS 产品里可能有细微差异，按关键字定位比硬编码全名更稳。

骨架搭好后，验证环节不能省。我的验证习惯是用气象站点的实测气温数据做参考，但注意：LST 是地表温度，不是 1.5 米气温，两者差值白天可能达到 5 到 15 摄氏度，夜间差值小很多。比较时主要看两点：一是 LST 和站点气温的相关系数是否显著；二是高温期、低温期出现的日期是否对齐。如果你处理的是裸土区域，LST 白天波动幅度会比气温大得多，不要因为这个就觉得数据处理错了。

还有一个我花过不少时间才养成的习惯：每次批量处理完，导出一张随意日期的 LST 分布图和对应 QC 掩膜图，用眼睛扫一遍。MODIS 处理里的很多问题——投影错、尺度因子漏乘、云没滤净——在这种抽查下会在十分钟内暴露，比事后发现数据不可用再重跑强得多。V1.0 的代码不追求一次写对，追求的是出问题时能快速定位到是哪个环节。毕竟这行当的常态就是：前面省下的检查时间，后面都会以通宵重算的方式还回去。希望这篇笔记能帮你少踩几个我已经踩过的坑。

本文还有配套的精品资源，点击获取