基于 Google Earth Engine 构建人工林语义分割数据集:从图源下载到 500×500 样本裁剪
一、写在前面
最近在做人工林遥感语义分割任务时,我构建了一个面向高分辨率人工林分割的数据集。论文中的数据集覆盖亚洲、欧洲、北美和非洲四个区域,其中北美子集使用 NAIP 航空影像,并包含原始影像、NDVI、EVI 和像素级标注;其他区域主要来自高分辨率卫星影像和航空底图场景。大场景影像统一裁剪为 500×500 像素的小样本,最终形成四个区域共 4728 组 image-mask 样本。
数据集地址:
https://huggingface.co/datasets/jkhsdfja/ArtificialForest-HRSegDataset
这篇文章主要记录一个可复现的数据构建流程:如何使用 Google Earth Engine(GEE)寻找遥感图源、下载影像、计算植被指数,并进一步裁剪成语义分割数据集。
二、整体流程
整个流程可以分成 6 步:
- 在 GEE 中确定研究区域 ROI
- 选择合适的遥感影像源
- 进行云掩膜、影像合成和裁剪
- 计算 NDVI、EVI 等植被指数
- 导出 GeoTIFF 影像
- 本地裁剪为 500×500 图像块,并制作 mask 标签
需要注意的是,GEE 可以直接导出其数据目录中的公开影像,例如 Sentinel-2、Landsat、NAIP 等。
如果使用 Google Earth 网页端底图作为人工目视参考,一般不能直接从 GEE 导出 Google Earth 底图本身,更规范的做法是使用 GEE 中公开可下载的数据源,并在 QGIS、ArcGIS 或 LabelMe 中完成标注。
三、在 GEE 中获取图源
1. 北美区域:使用 NAIP 高分辨率影像
NAIP 是美国农业部提供的航空影像,分辨率较高,通常包含 RGB 和近红外波段,适合做植被分割、林地提取和人工林边界分析。
在 GEE 中可以这样加载 NAIP:
// ===============================
// 1. 定义研究区域
// ===============================
// 可以在 GEE 左侧 Geometry Imports 中手动绘制 ROI
var roi = geometry;
// ===============================
// 2. 加载 NAIP 影像
// ===============================
var naip = ee.ImageCollection("USDA/NAIP/DOQQ")
.filterBounds(roi)
.filterDate("2020-01-01", "2023-12-31");
// 多景影像镶嵌
var image = naip.mosaic().clip(roi);
// NAIP 常见波段:R, G, B, N
var rgb = image.select(["R", "G", "B"]);
var nir = image.select("N");
// ===============================
// 3. 计算 NDVI
// NDVI = (NIR - Red) / (NIR + Red)
// ===============================
var ndvi = image.normalizedDifference(["N", "R"]).rename("NDVI");
// ===============================
// 4. 计算 EVI
// EVI = 2.5 * (NIR - Red) / (NIR + 6*Red - 7.5*Blue + 1)
// ===============================
var evi = image.expression(
"2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))",
{
"NIR": image.select("N"),
"RED": image.select("R"),
"BLUE": image.select("B")
}
).rename("EVI");
// ===============================
// 5. 合并 RGB + NDVI + EVI
// ===============================
var outImage = rgb.addBands(ndvi).addBands(evi);
// 可视化
Map.centerObject(roi, 14);
Map.addLayer(rgb, {min: 0, max: 255}, "NAIP RGB");
Map.addLayer(ndvi, {min: -1, max: 1, palette: ["brown", "white", "green"]}, "NDVI");
Map.addLayer(evi, {min: -1, max: 1, palette: ["brown", "white", "blue"]}, "EVI");
// ===============================
// 6. 导出到 Google Drive
// ===============================
Export.image.toDrive({
image: outImage,
description: "NAIP_RGB_NDVI_EVI",
folder: "GEE_Plantation_Dataset",
fileNamePrefix: "naip_rgb_ndvi_evi",
region: roi,
scale: 1,
maxPixels: 1e13,
fileFormat: "GeoTIFF"
});
如果研究区在美国,优先推荐 NAIP,因为它包含近红外波段,可以直接计算 NDVI 和 EVI。
四、非北美区域:使用 Sentinel-2 获取可复现图源
如果研究区不在美国,NAIP 无法覆盖,可以使用 Sentinel-2。
Sentinel-2 分辨率虽然没有 NAIP 那么高,但它公开、稳定、可复现,并且具有红光、绿光、蓝光、近红外等波段,可以计算 NDVI 和 EVI。
// ===============================
// 1. 定义研究区域
// ===============================
var roi = geometry;
// ===============================
// 2. Sentinel-2 云掩膜函数
// ===============================
function maskS2Clouds(image) {
var qa = image.select("QA60");
var cloudBitMask = 1 << 10;
var cirrusBitMask = 1 << 11;
var mask = qa.bitwiseAnd(cloudBitMask).eq(0)
.and(qa.bitwiseAnd(cirrusBitMask).eq(0));
return image.updateMask(mask).divide(10000);
}
// ===============================
// 3. 加载 Sentinel-2 SR 数据
// ===============================
var s2 = ee.ImageCollection("COPERNICUS/S2_SR_HARMONIZED")
.filterBounds(roi)
.filterDate("2023-01-01", "2023-12-31")
.filter(ee.Filter.lt("CLOUDY_PIXEL_PERCENTAGE", 20))
.map(maskS2Clouds);
// 中值合成,降低云和噪声影响
var image = s2.median().clip(roi);
// Sentinel-2 波段
// B2: Blue, B3: Green, B4: Red, B8: NIR
var rgb = image.select(["B4", "B3", "B2"]).rename(["R", "G", "B"]);
// ===============================
// 4. 计算 NDVI
// ===============================
var ndvi = image.normalizedDifference(["B8", "B4"]).rename("NDVI");
// ===============================
// 5. 计算 EVI
// ===============================
var evi = image.expression(
"2.5 * ((NIR - RED) / (NIR + 6 * RED - 7.5 * BLUE + 1))",
{
"NIR": image.select("B8"),
"RED": image.select("B4"),
"BLUE": image.select("B2")
}
).rename("EVI");
// ===============================
// 6. 合并导出
// ===============================
var outImage = rgb.addBands(ndvi).addBands(evi);
Map.centerObject(roi, 12);
Map.addLayer(rgb, {min: 0, max: 0.3}, "Sentinel-2 RGB");
Map.addLayer(ndvi, {min: -1, max: 1, palette: ["brown", "white", "green"]}, "NDVI");
Map.addLayer(evi, {min: -1, max: 1, palette: ["brown", "white", "blue"]}, "EVI");
Export.image.toDrive({
image: outImage,
description: "S2_RGB_NDVI_EVI",
folder: "GEE_Plantation_Dataset",
fileNamePrefix: "s2_rgb_ndvi_evi",
region: roi,
scale: 10,
maxPixels: 1e13,
fileFormat: "GeoTIFF"
});
这里导出的 GeoTIFF 包含 5 个通道:
R
G
B
NDVI
EVI
后续可以根据任务需求只使用 RGB,也可以将 NDVI、EVI 作为额外模态输入模型。
五、为什么要计算 NDVI 和 EVI?
在人工林分割任务中,仅使用 RGB 影像有时会遇到几个问题:
1. 人工林和天然林纹理相似
2. 草地、灌木和人工林光谱相近
3. 阴影、裸土和建筑物容易造成误分割
4. 不同地区、不同卫星影像之间存在亮度差异
NDVI 和 EVI 可以提供更直接的植被信息:
NDVI:反映植被绿度和长势
EVI:对高覆盖植被、土壤背景和大气影响更稳定
在我的论文中,模型采用了光学图像与 NDVI、EVI 对齐融合的方式,用于增强人工林区域的表征能力。
六、导出大图后,本地裁剪为 500×500 样本
论文中的大场景影像约为 5000×5000 像素,并统一裁剪成 500×500 像素的小图块,用于语义分割网络训练。
下面给出一个本地裁剪 GeoTIFF 的 Python 脚本。1. 安装依赖
pip install rasterio numpy tqdm
2. 裁剪影像为 500×500
import os
import rasterio
from rasterio.windows import Window
from tqdm import tqdm
import numpy as np
def crop_geotiff_to_tiles(
image_path,
out_dir,
tile_size=500,
stride=500,
min_valid_ratio=0.6
):
os.makedirs(out_dir, exist_ok=True)
with rasterio.open(image_path) as src:
width = src.width
height = src.height
profile = src.profile.copy()
tile_id = 0
for y in tqdm(range(0, height - tile_size + 1, stride)):
for x in range(0, width - tile_size + 1, stride):
window = Window(x, y, tile_size, tile_size)
tile = src.read(window=window)
# 判断有效像素比例,避免导出大量空白图块
valid_mask = np.any(tile != 0, axis=0)
valid_ratio = valid_mask.sum() / (tile_size * tile_size)
if valid_ratio < min_valid_ratio:
continue
tile_profile = profile.copy()
tile_profile.update({
"height": tile_size,
"width": tile_size,
"transform": src.window_transform(window)
})
out_path = os.path.join(out_dir, f"image_{tile_id:05d}.tif")
with rasterio.open(out_path, "w", **tile_profile) as dst:
dst.write(tile)
tile_id += 1
print(f"裁剪完成,共生成 {tile_id} 个图块。")
if __name__ == "__main__":
image_path = "naip_rgb_ndvi_evi.tif"
out_dir = "tiles/images"
crop_geotiff_to_tiles(
image_path=image_path,
out_dir=out_dir,
tile_size=500,
stride=500
)
如果希望增加样本数量,可以把 stride 改小,例如:
stride=250
这样会产生重叠裁剪样本。
七、如何制作语义分割 mask?
GEE 下载的是影像,不会自动生成精确的人工林 mask。
语义分割数据集中的 mask 一般需要人工标注或者半自动标注。
常见方法有三种:
1. 在 QGIS / ArcGIS 中绘制人工林矢量边界
2. 使用 LabelMe 标注 polygon
3. 使用已有林地矢量数据作为初始标签,再人工修正
如果你已经在 QGIS 中绘制了人工林 polygon,可以导出为 shapefile 或 geojson,然后使用 Python 栅格化成 mask。1. 安装依赖
pip install geopandas rasterio shapely numpy
2. 矢量标注转 mask
import rasterio
import geopandas as gpd
from rasterio.features import rasterize
import numpy as np
def vector_to_mask(image_path, vector_path, out_mask_path):
"""
image_path: 参考影像路径
vector_path: 人工林 polygon 标注文件,例如 shp 或 geojson
out_mask_path: 输出 mask tif
"""
with rasterio.open(image_path) as src:
transform = src.transform
height = src.height
width = src.width
crs = src.crs
profile = src.profile.copy()
gdf = gpd.read_file(vector_path)
# 保证矢量和影像坐标系一致
if gdf.crs != crs:
gdf = gdf.to_crs(crs)
shapes = [(geom, 1) for geom in gdf.geometry if geom is not None]
mask = rasterize(
shapes=shapes,
out_shape=(height, width),
transform=transform,
fill=0,
dtype=np.uint8
)
profile.update({
"count": 1,
"dtype": "uint8"
})
with rasterio.open(out_mask_path, "w", **profile) as dst:
dst.write(mask, 1)
print("mask 已保存:", out_mask_path)
if __name__ == "__main__":
vector_to_mask(
image_path="naip_rgb_ndvi_evi.tif",
vector_path="plantation_label.geojson",
out_mask_path="plantation_mask.tif"
)
生成的 mask 中:
0:背景
1:人工林
八、同步裁剪 image 和 mask
语义分割任务要求 image 和 mask 严格一一对应。
因此,裁剪时必须使用相同的窗口。
import os
import rasterio
from rasterio.windows import Window
import numpy as np
from tqdm import tqdm
def crop_image_and_mask(
image_path,
mask_path,
image_out_dir,
mask_out_dir,
tile_size=500,
stride=500,
min_valid_ratio=0.6
):
os.makedirs(image_out_dir, exist_ok=True)
os.makedirs(mask_out_dir, exist_ok=True)
with rasterio.open(image_path) as img_src, rasterio.open(mask_path) as mask_src:
width = img_src.width
height = img_src.height
img_profile = img_src.profile.copy()
mask_profile = mask_src.profile.copy()
tile_id = 0
for y in tqdm(range(0, height - tile_size + 1, stride)):
for x in range(0, width - tile_size + 1, stride):
window = Window(x, y, tile_size, tile_size)
image_tile = img_src.read(window=window)
mask_tile = mask_src.read(1, window=window)
valid_mask = np.any(image_tile != 0, axis=0)
valid_ratio = valid_mask.sum() / (tile_size * tile_size)
if valid_ratio < min_valid_ratio:
continue
img_tile_profile = img_profile.copy()
img_tile_profile.update({
"height": tile_size,
"width": tile_size,
"transform": img_src.window_transform(window)
})
mask_tile_profile = mask_profile.copy()
mask_tile_profile.update({
"height": tile_size,
"width": tile_size,
"transform": mask_src.window_transform(window),
"count": 1,
"dtype": "uint8"
})
image_out_path = os.path.join(image_out_dir, f"image_{tile_id:05d}.tif")
mask_out_path = os.path.join(mask_out_dir, f"mask_{tile_id:05d}.tif")
with rasterio.open(image_out_path, "w", **img_tile_profile) as dst:
dst.write(image_tile)
with rasterio.open(mask_out_path, "w", **mask_tile_profile) as dst:
dst.write(mask_tile.astype(np.uint8), 1)
tile_id += 1
print(f"同步裁剪完成,共生成 {tile_id} 对 image-mask 样本。")
if __name__ == "__main__":
crop_image_and_mask(
image_path="naip_rgb_ndvi_evi.tif",
mask_path="plantation_mask.tif",
image_out_dir="dataset/images",
mask_out_dir="dataset/masks",
tile_size=500,
stride=500
)
裁剪完成后,可以得到类似结构:
dataset/
├── images/
│ ├── image_00000.tif
│ ├── image_00001.tif
│ └── ...
└── masks/
├── mask_00000.tif
├── mask_00001.tif
└── ...
九、将 GeoTIFF 转成 PNG/JPG 训练格式
有些语义分割代码习惯读取 PNG/JPG。
如果只使用 RGB 三通道,可以把 GeoTIFF 转为 PNG。
import os
import rasterio
import numpy as np
from PIL import Image
from tqdm import tqdm
def tif_to_png(tif_dir, png_dir):
os.makedirs(png_dir, exist_ok=True)
for name in tqdm(os.listdir(tif_dir)):
if not name.endswith(".tif"):
continue
tif_path = os.path.join(tif_dir, name)
with rasterio.open(tif_path) as src:
img = src.read([1, 2, 3]) # 读取 RGB 三通道
img = np.transpose(img, (1, 2, 0))
# 如果是 0-1 浮点数据,转成 0-255
if img.max() <= 1.0:
img = img * 255
img = np.clip(img, 0, 255).astype(np.uint8)
out_name = name.replace(".tif", ".png")
out_path = os.path.join(png_dir, out_name)
Image.fromarray(img).save(out_path)
print("转换完成。")
if __name__ == "__main__":
tif_to_png(
tif_dir="dataset/images",
png_dir="dataset/images_png"
)
mask 也可以转为单通道 PNG:
import os
import rasterio
import numpy as np
from PIL import Image
from tqdm import tqdm
def mask_tif_to_png(mask_tif_dir, mask_png_dir):
os.makedirs(mask_png_dir, exist_ok=True)
for name in tqdm(os.listdir(mask_tif_dir)):
if not name.endswith(".tif"):
continue
tif_path = os.path.join(mask_tif_dir, name)
with rasterio.open(tif_path) as src:
mask = src.read(1)
mask = mask.astype(np.uint8)
out_name = name.replace(".tif", ".png")
out_path = os.path.join(mask_png_dir, out_name)
Image.fromarray(mask).save(out_path)
print("mask 转换完成。")
if __name__ == "__main__":
mask_tif_to_png(
mask_tif_dir="dataset/masks",
mask_png_dir="dataset/masks_png"
)
十、划分训练集、验证集和测试集
import os
import random
import shutil
def split_dataset(
image_dir,
mask_dir,
out_dir,
train_ratio=0.7,
val_ratio=0.2,
seed=42
):
random.seed(seed)
names = sorted(os.listdir(image_dir))
names = [n for n in names if n.endswith((".png", ".jpg", ".tif"))]
random.shuffle(names)
total = len(names)
train_num = int(total * train_ratio)
val_num = int(total * val_ratio)
train_names = names[:train_num]
val_names = names[train_num:train_num + val_num]
test_names = names[train_num + val_num:]
splits = {
"train": train_names,
"val": val_names,
"test": test_names
}
for split, split_names in splits.items():
os.makedirs(os.path.join(out_dir, split, "images"), exist_ok=True)
os.makedirs(os.path.join(out_dir, split, "masks"), exist_ok=True)
for name in split_names:
image_path = os.path.join(image_dir, name)
mask_name = name.replace("image", "mask")
mask_path = os.path.join(mask_dir, mask_name)
shutil.copy(image_path, os.path.join(out_dir, split, "images", name))
shutil.copy(mask_path, os.path.join(out_dir, split, "masks", mask_name))
print("数据集划分完成:")
print("train:", len(train_names))
print("val:", len(val_names))
print("test:", len(test_names))
if __name__ == "__main__":
split_dataset(
image_dir="dataset/images_png",
mask_dir="dataset/masks_png",
out_dir="ArtificialForest-HRSegDataset"
)
十一、下载我的公开数据集
如果只是想直接使用整理好的人工林语义分割数据集,可以从 Hugging Face 下载:
https://huggingface.co/datasets/jkhsdfja/ArtificialForest-HRSegDataset
方式一:使用 git clone
git lfs install
git clone https://huggingface.co/datasets/jkhsdfja/ArtificialForest-HRSegDataset
方式二:使用 huggingface_hub
pip install huggingface_hub
from huggingface_hub import snapshot_download
snapshot_download(
repo_id="jkhsdfja/ArtificialForest-HRSegDataset",
repo_type="dataset",
local_dir="./ArtificialForest-HRSegDataset"
)
十二、数据集构建时需要注意的问题
1. 图像和 mask 必须严格对齐
遥感语义分割中,影像和标签只要错位几个像素,模型训练就会明显变差。
因此,裁剪 image 和 mask 时必须使用同一个窗口。
2. mask 不能使用双线性插值
如果对 mask resize,必须使用最近邻插值:
Image.NEAREST
或者:
cv2.INTER_NEAREST
mask 是类别索引图,不是普通灰度图,不能对类别编号做平均。
3. 尽量保留复杂场景
如果只保留边界清晰、背景简单的样本,模型在真实场景中容易失效。
论文数据集中保留了云影、山地、林地与非林地过渡区、建筑物和裸土干扰等复杂情况,用于提高模型泛化能力。
4. 植被指数要和原图保持同尺度
如果 RGB 图像是 500×500,NDVI 和 EVI 也要裁剪成同样大小。
否则在多模态融合时会出现空间错位。
十三、总结
本文记录了一个从 GEE 获取遥感图源并构建人工林语义分割数据集的完整流程。核心步骤包括:
1. 在 GEE 中选择研究区域
2. 下载 NAIP 或 Sentinel-2 影像
3. 计算 NDVI、EVI 等植被指数
4. 导出 GeoTIFF
5. 本地裁剪为 500×500 图像块
6. 制作人工林 mask
7. 划分训练集、验证集和测试集
对于人工林语义分割任务,RGB 图像可以提供纹理和空间结构信息,NDVI 和 EVI 可以提供植被生理信息,两者结合更有利于区分人工林、天然林、草地、裸土和建筑物等复杂背景。
公开数据集地址:
https://huggingface.co/datasets/jkhsdfja/ArtificialForest-HRSegDataset
后续我也会继续整理人工林语义分割模型训练、mask 可视化、mIoU 计算和多模态输入构建等相关内容。
AtomGit 是由开放原子开源基金会联合 CSDN 等生态伙伴共同推出的新一代开源与人工智能协作平台。平台坚持“开放、中立、公益”的理念,把代码托管、模型共享、数据集托管、智能体开发体验和算力服务整合在一起,为开发者提供从开发、训练到部署的一站式体验。
更多推荐



所有评论(0)