基于飞桨实现高光谱影像和全色影像融合

项目背景

高光谱影像因其高光谱分辨率（具有几百个光谱波段）所提供的丰富光谱信息在许多实际应用中都大放光彩，如图像分类、异常检测、变化检测和定量农业等领域。然而由于信噪比和卫星重访周期等原因，高光谱影像的空间分辨率很低，如Modis卫星影像的空间分辨率为500m，使得影像的空间细节严重丢失，极大地限制了其在应用中的范围和精度。

高光谱全色融合是指融合具有高光谱、低空间分辨率的高光谱影像、以及高空间分辨率的单波段全色影像，来得到具有高光谱、高空间分辨率的影像，这是提升高光谱空间分辨率的一种有效的方式。

然而现有的融合方法主要关注于小比例融合任务（即高光谱和全色影像的空间分辨率比值），缺乏实际应用性。如资源一号02D卫星获取的高光谱影像具有166个波段，空间分辨率为30m，而全色影像空间分辨率为2.5m，二者分辨率之比达到12。针对这类大比例分辨率，目前还没有很好的研究成果。

因此本文基于飞桨框架首次聚焦于大比例融合任务（比例为16），并针对融合问题的病态性（即从单波段全色影像预测多波段高光谱影像的反射率），本文提出了一种基于高光谱投影丰度空间的融合网络。此外，根据全色强度和丰度特征之间的线性和非线性的关系，还设计了细节注入的子网络来将全色细节注入到低维的丰度空间而非原始高维空间，有效地缓解了该融合问题的病态性，提升了融合影像的质量。

技术方案

基于全色特征和丰度特征之间的线性和非线性关系，本文构建了PDIN子网络来注入全色细节。

1.线性关系

全色和丰度特征之间的线性关系的推导如上所示，由第四个公式可知，全色强度为丰度的线性组合，而全色强度本身为高光谱的线性组合，因此将全色特征注入到丰度空间和注入到高光谱空间具有等价性。

2.非线性关系

图1：丰度特征和全色强度之间的模拟散点分布图

该非线性关系通过模拟实验得到。首先从USGS光谱库中选择光谱曲线来模拟纯净端元，然后随机生成丰度像元数据，该丰度计算得到的标准差即为图1横轴。随机得到的丰度乘以端元即为构建的高光谱数据，然后乘以worldview2的光谱响应函数生成全色数据。图1即计算得到的丰度标准差和全色强度的散点分布图。由于该分布形状和鱼很相似，故本文称该分布为‘鱼分布’。

运行下方代码即可生成图1分布图。

1import numpy as np2import matplotlib.pyplot as plt3srf_pan = np.load('/home/aistudio/data/data124050/srf_pan.npy')    # 光谱响应函数45filename = '/home/aistudio/work/file/envi_plot3.txt'        # 选择的各种地物纯净光谱6curve_data =np.loadtxt(filename, skiprows=13)[:, 1:]78pixel_num = 50009abun = np.random.rand(pixel_num, curve_data.shape[1])
10abun_sum = np.sum(abun, axis=1)
11abun = abun / np.expand_dims(abun_sum, axis=1)
12
13intensity = np.dot(srf_pan, np.dot(curve_data[:, :], np.transpose(abun))[:srf_pan.shape[0], :]) / np.sum(srf_pan)
14x_var = np.std(abun, axis=1)    # 每个像元端元丰度的标准差
15
16plt.scatter(x_var, intensity, color='b', s=10)
17plt.xlabel('Abundance STD')
18plt.ylabel('PAN Intensity')
19plt.show()

3.细节注入子网络（PDIN）

图1所示的分布为正确的分布，但实际在网络前向传播过程中，一些模块如上采样模块生成的数据分布并不正确，我们的想法是通过调整数据的分布来变相的注入全色细节。该调整过程设想如下图：

图2：假设的数据分布变换的过程

如图2所示，我们假设绿点是不正确的点（其像元丰度存在误差），而蓝点为符合‘鱼分布’的正确的点。我们希望将绿点移动到蓝点区域以符合‘鱼分布’。由于纵轴是对应的全色强度，是已有的数据，故其能起到监督或约束的作用。即在绿点移动过程中，对应丰度像元所解码得到的‘全色强度’应不偏离真实全色强度。

绿点的移动通过改变丰度像元标准差可以完成，将像元乘以一个权重即可。权重的意思是从图2可以看到，不同的绿点和红线（绿蓝边界）的距离不同，即需要移动的距离也不同，因此丰度像元乘以的权重应自适应每个绿点。由于篇幅有限，具体的推导过程请查看源论文。

图3：细节注入子网络（PDIN）

下面是基于飞桨构建的PDIN类函数：

1# 利用线性关系和非线性关系构建的融合网络2class pan_abunstd(nn.Layer):3    def __init__(self, endmember_num=30):4        super(pan_abunstd, self).__init__()5        self.multiply_weight = simple_net_res(1, endmember_num, kernelsize=3) # 乘性权重，改变方差6        self.plus_weight = simple_net_res(1, endmember_num, kernelsize=3)  # 加性权重，改变均值78    def forward(self, w1, w2, b1, b2, pan, abun, dim_band=1):9        update_weight = F.sigmoid(F.relu(pan - w1 * paddle.std(abun, axis=1, keepdim=True) - b1))     # 应该更新的位置权重
10        update_weight2 = F.sigmoid(F.relu(w2 * paddle.std(abun, axis=1, keepdim=True) + b2 - pan))
11        update_weight = update_weight + update_weight2
12        result0 = self.multiply_weight(pan) * update_weight * abun + self.plus_weight(pan) * update_weight
13        return result0 + abun

其中forward函数里的w1, w2, b1, b2在主类函数里通过paddle.create_parameter()函数定义。

4.总体网络框架

图4：融合网络框架Pgnet

本文提出的融合网络分为四部分，如图4所示。

Part 1为编码部分，即将输入高光谱影像投影到低维丰度空间中，采用了Conv+BN+LReLU的网络结构进行投影。

Part 2为上采样部分，由于本文关注于大比例融合任务，故采用了两步上采样来逐渐提升影像分辨率，绿色方块为上采样block，包含bicubic内插和单卷积操作。紫色方块Double-net为双卷积操作，用于匹配丰度的分辨率和提取特征。

图5：SAB block

Part 3为注意力增强部分，采用了联级的像素级注意力机制来进一步增强特征，如图5所示。并且SAB block中加入了PDIN来进一步注入全色空间细节。

Part 4为解码部分，即将丰度特征还原为高光谱影像。值得注意的是，这部分的卷积参数可以看成线性混合模型中的端元。

实验结果

本文在JiaXing，Chikusei，XiongAn 和Real dataset四个数据集上验证了所提出的Pgnet的融合性能。此外，为了验证网络在小比例融合任务上的泛化性能，我们进行了ratio 4 和 8的实验。消融实验也被设计来验证各个网络模块的必要性和分布变化规律的正确性，即在PDIN中，通过改变数据分布来注入空间细节并不只是设想，实验结果证实了模拟实验所提出的假设。

数据生成代码，第一步：生成低分辨率高光谱和全色数据；第二步：切割成训练patch：

1def generate_data(path, path_srf):    # 第一步2    ratio_hs = 163    # spectral response function of worldview24    srf_pan = np.expand_dims(np.expand_dims(np.load(path_srf), axis=-1), axis=-1)56    noise_mean = 0.07    noise_var = 0.00018    dowmsample16 = Downsampler(ratio_hs)    # 高斯核降采样9    original_msi = np.float32(np.load(path))
10    band, row, col = original_msi.shape
11original_msi = original_msi[:, :(row - row%(ratio_hs*4)), :(col - col%(ratio_hs*4))]
12
13    # ratio 16
14    temp_blur = dowmsample16(paddle.to_tensor(np.expand_dims(original_msi, axis=0), dtype='float32'))
15    print(temp_blur.shape)
16    temp_blur = np.squeeze(temp_blur.numpy())
17    _, rows, cols = temp_blur.shape
18    blur_data = []
19    for i3 in range(temp_blur.shape[0]):    # 加入噪声
20        blur_data.append(np.expand_dims(temp_blur[i3, :, :] +
21                                        np.random.normal(noise_mean, noise_var ** 0.5, [rows, cols]), axis=0))
22    blur_data = np.concatenate(blur_data, axis=0)
23    print('blur_data.shape:' + str(blur_data.shape))
24
25    # simulated pan image
26    temp_pan = np.expand_dims(np.sum(original_msi * srf_pan, axis=0) / np.sum(srf_pan), axis=0)
27    print('temp_pan.shape:' + str(temp_pan.shape))
28
29    return original_msi, temp_blur, temp_pan
30
31def crop_data(hrhs, lrhs, pan, ratio=16, train_ratio=0.8):        # 第二步
32    training_size = 64  # training patch size
33    testing_size = 256  # testing patch size
34    idx = int(lrhs.shape[2] * train_ratio)    # 切割的坐标索引值
35
36    '''产生训练和测试数据'''
37    train_hs_image, train_hrpan_image, train_label = \
38        Crop_traindata(lrhs[:, :, :idx],
39                        pan[:, :, :idx*ratio],
40                        hrhs[:, :, :idx*ratio],
41                        ratio=ratio, size=training_size,
42                        test=False)
43    test_hs_image, test_hrpan_image, test_label = \
44        Crop_traindata(lrhs[:, :, idx:],
45                        pan[:, :, idx*ratio:],
46                        hrhs[:, :, idx*ratio:],
47                        ratio=ratio, size=testing_size,
48                        test=True)
49    return train_hs_image, train_hrpan_image, train_label, test_hs_image, test_hrpan_image, test_label

训练代码：

1#  定义数据和模型2dataset0 = Mydata(train_hs_image, train_hrpan_image, train_label)        #  训练数据3train_loader = paddle.io.DataLoader(dataset0, num_workers=0, batch_size=opt.batch_size,4                                 shuffle=True, drop_last=True)56dataset1 = Mydata(test_hs_image, test_hrpan_image, test_label)        #  测试数据7test_loader = paddle.io.DataLoader(dataset1, num_workers=0, batch_size=opt.test_batch_size,8                       shuffle=False, drop_last=False)9
10model = Pg_net(band=opt.in_nc, endmember_num=opt.endmember, ratio=ratio)
11
12L2_loss = nn.loss.MSELoss()        #  loss 函数
13samloss = SAMLoss()
14
15scheduler = optim.lr.StepDecay(opt.lr, opt.step, gamma=opt.momentum, verbose=False)    #  学习率调整
16optimizer = optim.Adam(learning_rate=scheduler, parameters=model.parameters())
17for epoch in range(opt.num_epochs):
18    time0 = time.time()
19    loss_total = 0.0
20
21    scheduler.step(epoch)
22    model.train()
23    for i, (images_hs, images_pan, labels) in enumerate(train_loader()):
24        result = model(images_hs, images_pan)
25
26    loss_l2 = L2_loss(result, labels)
27        loss_sam = samloss(result, labels)
28        loss = loss_l2 + 0.01*loss_sam
29
30        loss.backward()
31        optimizer.step()
32        optimizer.clear_gradients()
33        loss_total += loss.numpy()
34
35    if ((epoch+1) % 10) == 0:
36        print('epoch %d of %d, using time: %.2f , loss of train: %.4f' %
37            (epoch + 1, opt.num_epochs, time.time() - time0, loss_total))
38
39# 测试结果输出：
40model.eval()
41image_all = []
42with paddle.no_grad():
43    for (images_hs, images_pan, _) in test_loader:
44        outputs_temp = model(images_hs, images_pan)
45        image_all.append(outputs_temp)
46    a = paddle.concat(image_all, axis=0)
47return a

JiaXing 实验结果：

图6：各个方法的融合结果在JiaXing数据集上。其中第一列为真实影像，后续每列分别为GSA, SFIM, Wavelet, MTF_GLP_HPM, CNMF, Hyperpnn2, HSpeSet2, HyperKite方法的结果，本文方法在最后一列 (详细方法请看论文)。其中奇数行为融合高光谱影像，偶数行是和真值的差值。

表1：定量实验结果在JiaXing数据集上。

从上述实验结果可以看出，相比于其他方法，本文方法展现出了更好的光谱和空间保真度。包括多个指标的定量结果也表明我们的方法取得了最好的性能。由于空间关系，其他三个数据集的结果并未展现，但我们的方法均取得了最好的结果。

数据分布变换假设的验证

为了验证提出的PDIN对于数据分布变换的效果，我们对测试影像的中间结果进行了分析，结果如下：

图 7：中间结果变化过程。其中每行代表不同的数据集，以第一行JiaXing数据为例，前三个中间结果散点图（横轴为丰度标准差，纵轴为全色强度）分别代表第二个上采样block后、第二个PDIN后和Part 3的末端。第四到六幅影像为对应解码得到的高光谱影像，最后一个为真实影像。

如图7所示，从第一列到第二列，即经过了一次PDIN的细节注入后，数据分布变得些许像模拟数据的分布了（鱼分布），对应的高光谱影像质量也有所提升（第五列）。在经过了带有PDIN的注意力机制后，即第三列，数据分布变成了和模拟数据几乎相等的分布，对应的高光谱质量也和真实影像非常接近（第六列）。该图所展示的中间结果的变换验证了提出的PDIN能够通过纠正数据分布来注入细节和提升影像质量，并且也验证了‘鱼分布’的正确性。

结论

本文关注于大比例的高光谱和全色融合任务，使用飞桨框架构建了基于投影丰度子空间的逐级上采样融合网络。利用丰度和全色特征之间的线性和非线性关系构建了全色细节注入子网络PDIN，有效地将全色细节注入到了丰度特征中，缓解了融合问题的病态性。定性和定量的实验结果验证了所提出的融合网络的有效性和高效性。对于中间结果的可视化也验证了PDIN能够通过纠正数据分布来注入全色细节和提升影像质量。该轻量型的网络设计（只有49100个参数）和所取得的实验结果使得高光谱和全色融合的实际应用迈出了一步。

欢迎fork项目，一起交流学习：

基于飞桨实现高光谱影像和全色影像融合相关推荐

高光谱数据集_基于飞桨实现高光谱反演：通过遥感数据获取土壤某物质含量
[飞桨开发者说]马云飞,枣庄学院地理信息科学专业在读,计算机视觉技术爱好者,研究方向为神经网络在遥感领域的应用等. 高光谱反演是什么? 高光谱反演是使用遥感卫星拍摄的高光谱数据以及实地采样化验的某物质 ...
基于PaddleSeg的眼底血管分割——使用飞桨助力医学影像分析
基于PaddleSeg的眼底血管分割--使用飞桨助力医学影像分析一.项目背景研究表明,各类眼科疾病以及心脑血管疾病会对视网膜血管造成形变.出血等不同程度的影响.随着生活水平的提高,这类疾病的发病率 ...
基于飞桨复现图像分类模型TNT，实现肺炎CT分类
本项目介绍了TNT图像分类模型,讲述了如何使用飞桨一步步构建TNT模型网络结构,并尝试在新冠肺炎CT数据集上进行分类.由于作者水平有限,若有不当之处欢迎批评指正. TNT模型介绍 TNT模型全称是Tr ...
基于“飞桨”的深度学习智能车
本文作者:吴东昱,北京钢铁侠科技深度学习算法工程师前言我在观察历届智能车竞赛以及教学实验中发现,采用传统视觉算法的视觉智能车只能在特定赛道中行驶,一旦赛道环境改变,必须修改大量的代码才能运行. ...
基于飞桨 DeepLabV3+实现人体肾组织图像中肾小球识别
基于飞桨 DeepLabV3+实现人体肾组织图像中肾小球识别一.项目背景估计显示,地球上有超过70亿人,银河系中有3000亿颗恒星.相比之下,成年人体内含有37万亿个细胞.确定这些细胞之间的功能和 ...
基于飞桨图像分类套件PaddleClas的柠檬分类竞赛实战
前情提要通过之前教程中的学习,相信大家对于如何搭建一个分类网络已经清晰了.那么我们不禁会想,有没有更快速的尝试模型及技巧的方法呢?因为我们在上一次课程中使用的代码都需要自己进行开发,自己写需要很 ...
强烈推荐 | 基于飞桨的五大目标检测模型实战详解
机器视觉领域的核心问题之一就是目标检测(object detection),它的任务是找出图像当中所有感兴趣的目标(物体),确定其位置和大小.对于人类来说,目标检测是一个非常简单的任务.然而,计算机能 ...
基于飞桨的稻米加工品质监测系统
基于飞桨的稻米加工品质监测系统本项目拟解决目前农业中普遍存在的稻米过度加工问题,通过开发基于PaddlePaddle的密集粘连小目标的图像分割与目标检测算法,量化各环节稻米的加工指标,进而实现产线工 ...
基于飞桨复现语义分割网络HRNet，实现瓷砖缺陷检测
点击左上方蓝字关注我们 [飞桨开发者说]路星奎,沈阳化工大学信息工程学院研究生在读,PPDE飞桨开发者技术专家,研究方向为图像分类.目标检测.图像分割等内容简介本项目讲述了HRNet网络结构,并尝 ...
搭建基于飞桨的OCR工具库，总模型仅8.6M的超轻量级中文OCR，单模型支持中英文数字组合识别、竖排文本识别、长文本识别的PaddleOCR
介绍基于飞桨的OCR工具库,包含总模型仅8.6M的超轻量级中文OCR,单模型支持中英文数字组合识别.竖排文本识别.长文本识别.同时支持多种文本检测.文本识别的训练算法. 相关链接 PaddleOCR ...

基于飞桨实现高光谱影像和全色影像融合

基于飞桨实现高光谱影像和全色影像融合相关推荐

最新文章

热门文章