CNN深度学习Landsat遥感影像地物分类Python源码包实战

发布时间:2026/10/11 9:50:07
CNN深度学习Landsat遥感影像地物分类Python源码包实战 简介这是一份面向遥感与深度学习入门者的CNN地物分类完整项目基于PyTorch与Python实现Landsat影像处理全流程适合高校学生、科研人员及从事智慧城市、农业监测等领域的开发者参考。资源包含数据预处理脚本、模型训练与预测代码、样本切片生成工具并附带TIF影像样例、模型权重文件及XML/TFW辅助参数可支撑毕业设计、课程作业或项目初期演示。压缩包共10个文件涵盖3个Python脚本、2个TIFF图像、2个XML元数据、1个H5权重模型、1个TFW坐标参数及1份说明文档整体大小约14.88MB结构紧凑、便于直接上手。目前已有88人学习使用。通过该资源读者可快速理解深度卷积网络在遥感影像地物分类中的完整工作流掌握从Landsat数据切片、训练分类模型到预测新影像的关键步骤适合希望在真实地理数据上开展实验并快速复现结果的学习者。1. CNN深度学习遥感影像地物分类这个Landsat Python源码包到底能干什么直接把一段Landsat原始影像丢给CAD软件去解译能把人盯到眼瞎。这个源码包做的事情很简单把遥感影像地物分类这件事拆成一条标准的深度学习流水线从Landsat多波段数据切块、训练一个七分类CNN模型到对整幅影像做推理预测并输出带地理坐标的分类结果图。我用它跑通了一次完整的原始tif进、分类tif出流程sequence可复现不是那种缺胳膊少腿的占位工程。这个包适合两类人。一类是做遥感或地信相关课程设计、毕业设计的学生需要一个能跑通、能讲清楚原理的完整项目另一类是刚接触遥感影像深度学习分类的从业者想知道Landsat数据预处理和推理落地的坑到底在哪。注意一个关键点这个包是Keras/TensorFlow体系的.h5权重文件不是PyTorch的.pth下载前先搞清楚你的环境能不能跑TF后面会详细说。2. 把Landsat大影像喂给CNN之前createImageChips.py切片原理与实操2.1 为什么要切片而不是把整幅tif直接塞进CNN遥感影像和日常图片分类最不一样的地方在于尺寸。一张Landsat 8 OLI影像覆盖范围动辄上万×上万像元全色波段分辨率15米多光谱30米一张Scene数据量从几百MB到几个GB不等。CNN的输入尺寸是固定的比如这个项目里的3×3卷积堆叠结构输入往往是64×64或128×128大小的图像块你想把整幅影像直接输入网络显存和计算量都会直接爆炸完全没有可操作性。所以第一步永远是切块遥感深度学习领域叫Image Chips或者Tile。1_createImageChips.py做的事情就是把大tif按固定尺寸切成小方块并为每一块生成对应的标签。这里要理解一个隐含逻辑监督分类需要带标签的训练数据Landsat原始影像本身没有标签你需要一份已经标注好的分类结果图通常是人工目视解译或者基于已有土地覆盖数据生成的label.tif脚本通过读取影像标签的对应关系来自动生成训练样本标签图里每个像元的整数值就代表了地物类别比如0代表水体、1代表植被、2代表建设用地、3代表农田、4代表裸地等具体看项目的类别定义文件。# 1_createImageChips.py 核心逻辑关键代码节选 import rasterio import numpy as np import os from tqdm import tqdm def create_chips(image_path, label_path, save_dir, chip_size64, step32): with rasterio.open(image_path) as src_img, rasterio.open(label_path) as src_lab: img_data src_img.read().transpose(1, 2, 0) # 转为 HWC 格式 lab_data src_lab.read(1) # 标签通常是单波段 height, width img_data.shape[:2] # 按窗口滑动切割step chip_size 时为重叠切块可以增广样本 for row in range(0, height - chip_size 1, step): for col in range(0, width - chip_size 1, step): img_chip img_data[row:rowchip_size, col:colchip_size, :] lab_chip lab_data[row:rowchip_size, col:colchip_size] # 去掉全为nodata的无效块避免把边界值学进去 if np.all(lab_chip -1): continue np.save(os.path.join(save_dir, fchip_{row}_{col}_img.npy), img_chip) np.save(os.path.join(save_dir, fchip_{row}_{col}_lab.npy), lab_chip)这个脚本的运行逻辑并不复杂用rasterio分别读入原始影像和标签图然后按滑窗方式切块切片尺寸和步长是两个最关键的参数。chip_size决定CNN输入的固定尺寸要根据你的网络结构来定如果网络里是全连接层输入尺寸必须固定一致这个项目的3×3卷积CNN输入尺寸在模型里已经写死了改尺寸就得改网络结构。step参数控制切块之间的重叠程度stepchip_size时是紧邻切割没有重叠stepchip_size/2时生成约两倍的数据量相当于免费的数据增广。重叠切的代价是训练样本之间有相关性可能轻微过拟合但样本量不足时重叠是有效的补救手段。我一般会先检查一下输出的chip数量和类别分布。常见做法是统计每个类别像元占比如果某几类样本很少训练出来的模型大概率会把稀有类别直接忽略掉。备选的切块策略包括类别均衡采样单独为稀有类别的区域加密切块后面避坑部分细说。2.2 Landsat七个波段是什么以及波段怎么选这个包的训练权重叫CNN_7class_3by3.h5这里的7class指最终要分七类地物有的项目里7class也可能指7个波段输入具体看模型第一层的input_shape。Landsat数据通常使用的波段组合里Landsat 8的OLI传感器有9个波段常用的有蓝B2、绿B3、红B4、近红外B5、短波红外SWIR1/B6和SWIR2/B7。如果创建训练数据时选了7个波段大概率是B2到B7加一个热红外或者全色波段你需要打开example.tif和模型脚本确认一下输入通道数到底是多少。波段选择直接影响地物分类精度不同波段对地物的响应差别很大。水体在近红外波段吸收很强反射率很低植被在近红外和高反射区形成红边特征建设用地在短波红外和可见光上反射都比较高。所以你的CNN输入波段组合里有没有近红外对水体、植被分类精度影响非常大。如果脚本里写的是3波段RGB输入分类精度会肉眼可见地下降。我建议你先用rasterio直接读一下tif的元数据和波段数量确认输入通道数匹配再看下一步。# 查看tif波段数和元数据确认通道匹配 python -c import rasterio with rasterio.open(example.tif) as src: print(src.count, src.height, src.width, src.transform) print(src.crs) 这段bash命令用来快速检查影像通道数和空间参考如果你下载的Landsat数据波段顺序和原项目不一致就得在切块之前先做波段重排。常见做法是用GDAL的gdal_translate或者rasterio的read(index)按索引重排波段。务必在切块前完成波段顺序的校正不然训练出来的模型和预测时用的波段分布不一致分类结果会一塌糊涂。波段选择还影响数据的存储量。7波段float32的tif和3波段uint8的tif数据量相差接近十倍切块后的npy文件占用磁盘空间更大。我一般会把原始影像反射率转成uint8再存npy以节省内存损失一点精度但对CNN训练影响不大。3. 训练一个七分类CNN2_trainModel.py的网络结构与参数解析3.1 小卷积核堆叠为什么适合遥感影像地物分类这个项目的权重名带3by3字样说明网络用的卷积核主要是3×3尺寸。3×3卷积是现在几乎所有CNN的基础组件两个3×3卷积堆叠起来的感受野等于一个5×5卷积但参数量更少非线性更强。遥感地物分类和ImageNet识别不一样的地方在于地物边界是连续且渐变的一个像元可能是几种地物的混合3×3这种小卷积核能在保证局部细节的前提下逐层扩大感受野非常适合像素级别的分类任务。七分类任务相对简单不需要太过深厚的网络。一个经典的轻量CNN结构就够了三层3×3卷积加池化展平后接全连接层做分类。Landsat数据本身是大范围成像但每块chip是局部场景网络能提取的特征包括光谱特征各波段的灰度组合关系和纹理特征局部空间结构这些特征通过3×3卷积逐层抽象最后在全连接层组合成分类决策。这里不需要预训练权重遥感影像和自然图像的特征分布差异大直接用随机初始化训练反而更省事。# 2_trainModel.py 核心模型定义简化节选对应CNN_7class_3by3.h5 import tensorflow as tf from tensorflow.keras import layers, models def build_cnn(input_shape(64, 64, 7), num_classes7): inputs tf.keras.Input(shapeinput_shape) # 第一个卷积块3x3卷积通道升到32下采样一次 x layers.Conv2D(32, (3, 3), activationrelu, paddingsame)(inputs) x layers.Conv2D(32, (3, 3), activationrelu, paddingsame)(x) x layers.MaxPooling2D((2, 2))(x) # 第二个卷积块通道升到64再下采样 x layers.Conv2D(64, (3, 3), activationrelu, paddingsame)(x) x layers.Conv2D(64, (3, 3), activationrelu, paddingsame)(x) x layers.MaxPooling2D((2, 2))(x) # 第三个卷积块通道升到128 x layers.Conv2D(128, (3, 3), activationrelu, paddingsame)(x) x layers.Conv2D(128, (3, 3), activationrelu, paddingsame)(x) x layers.GlobalAveragePooling2D()(x) # 用全局平均池化替代Flatten减少参数量 # 分类头七类地物对应7个输出 outputs layers.Dense(num_classes, activationsoftmax)(x) model models.Model(inputs, outputs) return model model build_cnn() model.compile( optimizertf.keras.optimizers.Adam(learning_rate1e-4), losssparse_categorical_crossentropy, metrics[accuracy] )这个网络结构有三个值得注意的地方。第一用GlobalAveragePooling2D而不是Flatten加全连接层为了压缩参数数量防止在样本量不大的情况下过拟合同时让网络对输入尺寸有一定鲁棒性当然预测前还是要reshape到固定尺寸。第二padding全都用same这样卷积不会缩小特征图尺寸边界信息保留得更充分。第三loss用的是sparse_categorical_crossentropy而不是categorical_crossentropy两者的区别在于标签是一维整数索引还是one-hot编码如果你的训练标签是npy里直接存的整数类别号就用sparse版本少一步预处理。学习率的设置上1e-4对于这种小数据集是安全值数据集很小的时候学习率太高会在损失曲面震荡。如果要调参我一般先跑5个epoch看loss下降趋势如果loss纹丝不动就把学习率提到3e-4或1e-3如果loss先降后弹上去那就是过拟合了调低学习率的同时加一点Dropout。3.2 数据集划分和训练流程验证集intended用途切块脚本生成的是海量的npy文件训练的时候不能一次性全load进内存需要写一个数据生成器DataGenerator边读边训练。常见做法是继承tensorflow.keras.utils.Sequence在每个epoch打乱样本顺序按照batch_size读取chunk样本。验证集从所有样本里按比例随机抽或者按区块抽——按区块抽更符合遥感实际情况因为空间相邻的样本高度相似随机抽取验证集会高估模型的真实泛化能力。# DataGenerator 简洁实现节选 import numpy as np from tensorflow.keras.utils import Sequence class ChipGenerator(Sequence): def __init__(self, img_paths, lab_paths, batch_size32): self.img_paths img_paths self.lab_paths lab_paths self.batch_size batch_size self.indexes np.arange(len(img_paths)) def __len__(self): return int(np.ceil(len(self.img_paths) / self.batch_size)) def __getitem__(self, idx): idxs self.indexes[idx * self.batch_size : (idx 1) * self.batch_size] imgs np.stack([np.load(self.img_paths[i]) for i in idxs]) labs np.stack([np.load(self.lab_paths[i]).ravel() for i in idxs]) # 每个chip取众数作为该chip的标签避免混合类别干扰 from scipy import stats label np.array([stats.mode(l).mode for l in labs]) return imgs, label这里有个细节值得注意每个64×64的chip内可能有多种地物怎么定义这个chip的标签三种常见做法取像元众数、取中心像元类别、或者做逐像元分类语义分割。这个项目的h5权重最后预测的是整幅图的逐像元分类结果所以训练时更合理的做法其实是按像元做分割训练但很多人图省事用众数做标签这会导致边界区域预测不准。如果你发现预测结果里地物边界锯齿感极强且有小块错误斑块八成是标签生成方式太粗暴。我一般会在训练脚本里加一个选项按中心像元标签训练边界精度会有可见提升。训练时长和显存占用也要有预期。CPU上训练64×64×7的输入、batch_size32、跑50个epoch几万样本可能需要十几个小时GPU上可以缩短到一小时内。训练结束后模型会保存为CNN_7class_3by3.h5这个h5包含了完整网络结构和权重预测的时候直接load_model就能用不需要重新定义网络。4. 从训练好的模型到整幅分类图3_predictNewData.py推理全流程4.1 滑窗推理与拼接为什么你不能直接resize整幅影像去预测拿到训练好的h5权重接下来要对一幅全新的Landsat影像做分类预测输出结果就是new_class.tif。脚本的预测流程和切块是镜像关系用同样的窗口尺寸从左到右从上到下扫描整幅影像每个窗口经过CNN得到一个类别预测再按窗口位置把预测结果拼回大图。这里面最关键的参数是重叠率——推理时的滑窗步长不一定要等于切块的步长。# 3_predictNewData.py 滑窗预测核心逻辑节选 import numpy as np import rasterio from tensorflow.keras.models import load_model model load_model(CNN_7class_3by3.h5) with rasterio.open(example.tif) as src: img src.read().transpose(1, 2, 0).astype(np.float32) profile src.profile.copy() height, width img.shape[:2] chip_size 64 step 32 # 重叠推理步长小于chip_size减少边界伪影 pred_map np.zeros((height, width), dtypenp.uint8) count_map np.zeros((height, width), dtypenp.float32) for row in range(0, height - chip_size 1, step): for col in range(0, width - chip_size 1, step): chip img[row:rowchip_size, col:colchip_size, :] # 归一化方式必须与训练时严格一致这是最常见的翻车点 chip (chip - chip_mean) / chip_std chip_input np.expand_dims(chip, axis0) pred model.predict(chip_input, verbose0)[0] # 得到7个类别的概率 pred_cls np.argmax(pred, axis-1).astype(np.uint8) # 投票合成重叠区域累加最后取多数类别 pred_map[row:rowchip_size, col:colchip_size] pred_cls count_map[row:rowchip_size, col:colchip_size] 1.0 # 重叠区域取平均多数投票避免边缘拼接痕迹 pred_map (pred_map / np.maximum(count_map, 1.0)).astype(np.uint8) # 以原始影像的坐标信息写tif确保new_class.tif能和原图叠加 with rasterio.open(new_class.tif, w, **profile) as dst: dst.write(pred_map, 1)预测脚本最容易出错的地方集中在两个环节。一是归一化参数必须和训练时完全一致训练时用的mean和std如果在训练脚本里写过预测时就老老实实用同一个值不能偷懒用当前影像自己算均值否则光谱分布整体偏移分类精度断崖式下跌。二是为什么推理用重叠窗口边缘的预测通常比中心区域噪声更大重叠推理让每个像元被多个窗口预测到取多数投票得到平滑结果。这个技巧能有效消除分类图上的棋盘格效应。输出的new_class.tif旁边还带着new_class.tfw、new_class.tif.aux.xml、new_class.tif.xml几个文件。.tfw是World File文件保存了仿射变换六参数左上角坐标、像元分辨率、旋转系数让分类图能和原始影像在GIS软件里精确叠加.aux.xml是GDAL自动生成的辅助文件记录统计信息和色彩映射。这些配套文件说明脚本在写tif时正确地复制了原始影像的地理参考信息如果你拿到的分类结果能直接在ArcGIS/QGIS里和Landsat原图叠在一起说明坐标链路是通的。4.2 理解输出结果七类地物的分类图怎么读new_class.tif里每个像元的值是0到6的整数对应七类地物。渲染的时候最好给它配一套颜色映射水体用蓝色0、植被用绿色1、建设用地用红色2或灰色、农田用黄色、裸地用棕色这样一眼能看出分类的分布趋势。如果直接在ArcGIS里打开new_class.tif看到的是灰度图或彩色混杂的伪彩色图可以右键属性→符号系统→唯一值指定类别字段和色带。有的项目还会额外生成一个类别名称与颜色的CSV或QGIS样式文件README里应该写明每个整数代表什么地物类别——如果没写那就去看训练时用的label图对应的类别表。我处理这个包时先看的是数据说明文档里有没有类别映射表这个东西比模型权重还重要因为很多分类结果要写进报告里。一个需要留意的行为是CNN逐像元分类的结果会呈现椒盐噪声就是孤立的小块类别区域在连续地物中出现常见于建设用地和裸地这种光谱相近的类别。这不是bug是逐像元分类的「黑匣子」特性。后续处理一般会用众数滤波或条件随机场CRF做后处理平滑消除孤立的分类噪声。5. 避坑与常见问题排查训练、预测、文件兼容性三个环节的6条血泪经验5.1 坐标信息丢失分类图和原图叠加不上现象生成的new_class.tif单独打开正常但拖进QGIS和原图叠加时跑到地球另一边或者完全不显示。原因预测脚本在写入结果tif时只复制了profile的一部分transform仿射变换没带出来或者原始影像本身就带了地理坐标但中间有人用numpy存npy再读回来时把坐标信息丢在第一步切块之前。解决用rasterio.open读取原始影像后保留src.profile写入时直接**profile传参这就够了。我遇到过更隐蔽的情况profile里有nodata参数没设置导致背景区域被当成0类参与后续统计。所以写完tif后务必用rasterio重新读一次打印transform和crs确认落在正确的坐标范围。5.2 预测结果全是某一类比如整幅图都变成植被现象模型评估时精度还不错但预测新影像时所有像元都输出同一个类别。原因训练样本和预测影像的光谱范围不一致。典型场景是训练数据来自某个季节的Landsat比如夏季植被茂盛预测用的是另一个季节的影像比如冬季落叶光谱分布差异超过模型泛化能力。另一个可能是归一化参数写死了某个地区的均值方差换一景影像后分布偏移。解决做预测之前先统计待预测影像各波段的均值和标准差和训练数据的归一化参数对比差异超过30%就要考虑做直方图匹配或者把训练数据重新采样一批不同时相的影像。做工程项目时训练集覆盖不同时相/不同区域是铁律。5.3 Keras版本不一致导致h5模型加载失败现象load_model(CNN_7class_3by3.h5)报错提示Unknown layer或者Unable to load尤其是老代码用了tensorflow.keras里的特殊层。原因h5文件里保存的网络结构包含层配置信息新版Keras移除了某些旧层名称。常见于用keras和tensorflow.keras混用的情况或者TF 1.x的模型在TF 2.x下加载。解决优先用custom_objects从keras或tensorflow.keras导入对应层来加载from tensorflow.keras.models import load_model import tensorflow.keras.layers as layers # 如果模型用了自定义层或旧层名按需传入custom_objects model load_model(CNN_7class_3by3.h5, custom_objects{GlobalAveragePooling2D: layers.GlobalAveragePooling2D})如果这一步还是报错一个可行的替代方案是忽略模型架构手动在Python里重新定义相同的网络结构然后只加载权重model.load_weights(CNN_7class_3by3.h5)注意这时候h5里得有model这个key否则就要打开h5文件查看结构。h5文件本质是HDF5格式可以查看内部结构没耐心的可以直接跑一遍模型定义脚本看能否成功加载。这类问题在老的公开项目里出现频率很高遇到不要慌先确认TF版本和Keras版本在纯环境里装指定的老版本TF反而最快。5.4 训练时显存溢出OOM但数据集本身不现象切块阶段产出几千个npy文件后开始训练加载几十个样本就把显存撑爆程序直接被杀。原因Sequence的__getitem__一次性读了整个batch的所有npy进内存如果chip文件很大或者batch_size设得过高显存崩溃。更隐性的原因是从磁盘读npy再转成numpy数组的过程本身也很耗时卡在IO上。解决先把batch_size降到16或8试跑一个epoch确认能过再往上加。另外npy文件如果存的是float64需要转成float32能省一半内存。如果还想更省可以在__getitem__里只从磁盘读取当前batch的索引对应的文件路径不要用全量加载。用tf.data.Dataset.from_generator配合num_parallel_calls也能提升IO效率。5.5 标签类别和训练时的类别含义对不上现象预测结果里2类是水体但训练标签里2类是建设用地整个类别是错位的。原因切块和训练的数据集来自不同的标签文件两份标签图的类别编码不一致。有的是0-based有的是1-based有的把255当背景。解决切块完成后、训练之前先打印label的独立值分布np.unique(lab_data)确认整数类别范围和你预期一致。然后跑一次预测用example.tif里已知的地物类型比如明显的水体区域做目视验证再看预测结果是不是正确。这是整个流程里最便宜的一次检查能省一整天调参时间。5.6 类别不平衡稀有类别精度极低现象训练出来的模型总体精度80%但看一眼混淆矩阵农田类的recall只有20%大部分农田被分成了裸地或草地。原因Landsat影像里各类地物面积天然不平衡如果农田占比不到5%切块得到的样本也极少模型倾向把所有样本判成多数类。类别不平衡在遥感分类里是常态不是个例。解决最简单的方式是在DataGenerator里加类别权重model.fit的class_weight参数可以直接传权重字典。另一个更激进的做法是切块时对稀有类别的区域做重叠加密采样比如步长减半再切一次。数据增强上对稀有类别做随机旋转90度、翻转可以低成本增加样本多样性。如果这些做完还是不行就要考虑用Focal Loss替换交叉熵损失专门针对难分类样本。我自己的经验是先试class_weight见效最快且不动网络结构再试Focal Loss能进一步压坪稀有类别的错误率最后才考虑采样策略调整因为那要重新生成npy文件时间成本高。6. 进阶玩法用example.tif快速验证整套流程并做分类后处理与精度检查拿到这个源码包我建议你按从后往前的顺序验证先跑3_predictNewData.py输入是example.tif输出new_class.tif对照README或已有的new_class.tif如果包里带了预测结果看看两个分类图是否一致。如果差异很大就要检查训练出来的模型权重是否和示例结果来自同一套训练数据。example.tif的存在非常关键相当于给了你一条基准答案验证完流程正确之后再替换成你自己的Landsat影像。在验证完流程通顺后进阶的重点放在分类后处理上。CNN直接输出的new_class.tif会带有椒盐噪声和细碎图斑做两步后处理就能显著提升成图质量。第一步用众数滤波scipy.ndimage.median_filter配合size3或5去掉孤立点第二步做多数滤波或主/次要分析scipy.ndimage的binary_opening或者用GDAL的gdal_sieve去掉小图斑。# 后处理示例先用中值滤波去噪声再剔除小图斑 from scipy.ndimage import median_filter import numpy as np # 对分类结果逐类别处理避免类别边界被滤波模糊 filtered median_filter(pred_map, size5, modenearest) # 剔除面积小于阈值比如50个像元的孤立图斑 from skimage import morphology cleaned morphology.remove_small_objects( filtered.astype(bool), min_size50 ).astype(np.uint8) * filtered最终检查精度的常规做法是拿分类结果和已有的高分辨率影像或Google Earth影像做目视对比选几个典型区域水体边界、建设用地与裸地交接处放大看判断分类是否正确。定量评估的话需要一份测试区域的真实标签Ground Truth用sklearn.metrics里的confusion_matrix和classification_report算总体精度、Kappa系数、每类精确率和召回率。你下载这个源码包后至少应该做一次完整的预测→后处理→精度报告流程才算真正把这个项目用起来。我自己的习惯是每次拿到这种遥感分类项目第一件事不是看模型结构而是先跑预测脚本把example.tif跑一遍确认输出分类图分布合理然后手动设一个小的验证区域用影像自带的波段数据通过NDVI等指数估算水体、植被、裸地的分布再对比分类结果图的对应位置。从那以后我每次做类似的Landsat分类任务都强制走一遍先预测后验证再训练的顺序能少掉很多后期的返工。这个包的完整代码链路做得比较规整跑通之后你替换成自己的数据会省很多事希望帮到你。本文还有配套的精品资源点击获取