遥感图像PCA融合的并行算法研究与实现

行宽， 也即每一维波段数据的长度。具体可分为计
算每个波段的平均值和计算每两个波段间的协方
１５４
微电子学与计算机
２００６ 年第 ２３ 卷第 １０ 期
差两步。 （ ２） 计算协方差矩阵的特征向量 ｖ 及特征值 ｄ （ ３） 将特征向量矩阵 ｖ 乘以 ＭＳ 数据矩阵得到
主成分矩阵 Ｐ， 记 Ｐ 中对应于最大特征值的一行 Ｐ１ 为第一主成分。
算， 得到最终的协方差矩阵。这一步骤涉及到 ｒｏｏｔ
进程和其他进程间的通信有 ３ 次。
（ ２） 计算协方差矩阵的特征向量 ｖ 及 特 征 值
ｄ。因为这一步骤的计算量较小， 所以由主进程完
成， 并将计算结果广播到其他非主进程结点。
（ ３） 在 每 个 结 点 上 进 行 Ｒ、Ｇ、Ｂ 向 Ｐ１、Ｐ２、Ｐ３ 的 转换， 记 Ｐ 中对应于最大特征值的一行 Ｐ１ 为第一主 成分。
在 遥 感 图 像 融 合 技 术 中 ， 主 成 分 分 析 ＰＣＡ （Ｐｒｉｎｃｉｐａｌ－ Ｃｏｍｐｏｎｅｎｔ－ Ａｎａｌｙｓｉｓ） 是一种经典的融合 方法， 然而由于实现运算量大， 计算复杂度高， 随着 所获取的遥感图像数据量的不断增大， 该方法无法 满足一些时效性要求， 从而影响了 ＰＣＡ 融合方法的 广泛应用。因此， 随着并行化技术的发展， 研究高效 实用的 ＰＣＡ 融合并行算法具有较高的理论与实用 价值。本文针对已有的 ＰＣＡ 融合串行算法， 提出了 一种新的 ＰＣＡ 融合并行算法 Ｐ－ ＰＣＡ（Ｐａｒａｌｌｅｌ Ｐｒｉｎ－
（ ４） 对每个结点， 将 Ｐ 矩阵的第一主成分用经
过拉伸后的全色图像 ｐａｎ 的数据替换。由于 ｐａｎ 图
像 拉 伸 需 要 用 到 ｐａｎｍａｘ， ｐａｎｍｉｎ， 以 及 Ｐ１ｍａｘ， Ｐ１ｍｉｎ， 如 何 计算出全局的最大最小值有多种方法。我们采用的
方法是： 各个结点先计算本结点上的最大最小值，
ＰＣＡ， 并对该并行算法进行了通信优化， 在机群系统上进行了实现。针对 ＩＫＯＮＯＳ 图像进行实验的结果表明该算法
可获得良好的并行加速比， 并行效率较高， 为遥感图像的实际应用提供了有益的指导和借鉴。
关键词： 遥感图像融合， ＰＣＡ， 并行算法
中图分类号： ＴＰ３１
文献标识码： Ａ
文章编号： １０００－ ７１８０（ ２００６） １０－ ０１５３－ ０３
ｅｎｄｎｕｍ 为该结点所读取图像的结束象素的标号。然
பைடு நூலகம்
后将局部值传给 ｒｏｏｔ 进程， 再由 ｒｏｏｔ 进程算出整幅
图像的 Ｒ、Ｇ、Ｂ 后广播给其它结点。其它结点按照
ｅｎｄｎｕｍ
! 算式（ ２） 计算， 得到一个结点局部的
（Ｘｉ－ Ｘ ）（Ｙｉ－
ｉ＝ｂｅｇｉｎｎｕｍ
Ｙ ）值， 再次传入主结点， 由主结点进行加法及剩余计
２００６ 年第 ２３ 卷第 １０ 期
微电子学与计算机
１５３
遥感图像 ＰＣＡ 融合的并行算法研究与实现
胡 冰 周海芳 王攀峰 刘衡竹
（ 国防科技大学 计算机学院， 湖南 长沙 ４１００７３）
摘 要： 文章针对已有的遥感图像 ＰＣＡ 融合串行算法， 提出了一种基于数据并行的新的 ＰＣＡ 融合并行算法 Ｐ－
（Ｘｉ－ Ｘ ）（Ｙｉ－ Ｙ ）
ｉ＝１
ｉ＝（Ｎ－ １）Ｓ／Ｎ＋１
ｉ＝（Ｎ－ １）Ｓ／Ｎ＋１
（ｓ－ １）
（ ２）
每个结点先分别计算该结点上图像的局部
ｅｎｄｎｕｍ
! ｘｉ， 其中 ｘ 表示 Ｒ、Ｇ、Ｂ 三 个 波 段 中 的 一 个 ，
ｉ＝ｂｅｇｉｎｎｕｍ
ｂｅｇｉｎｎｕｍ 为该结点所读取图像的起始象素的标号，
Ａｂｓｔｒ ａｃｔ： Ｔｈｉｓ ｐａｐｅｒ ｄｉｓｃｕｓｓｅｓ ａ ｎｅｗ ｐａｒａｌｌｅｌ ａｌｇｏｒｉｔｈｍ ｏｆ ＰＣＡ ｉｍａｇｅ ｆｕｓｉｏｎ Ｐ－ ＰＣＡ， ｂａｓｅｄ ｏｎ ｔｈｅ ｔｒａｄｉｔｉｏｎａｌ ＰＣＡ ｉｍ－ ａｇｅ ｆｕｓｉｏｎ ａｌｇｏｒｉｔｈｍ， ａｎｄ ｍｏｄｉｆｉｅｓ ｓｏｍｅ ｄｅｔａｉｌｓ ｔｏ ｉｍｐｒｏｖｅ ｔｈｅ ｅｆｆｉｃｉｅｎｃｙ． Ｗｅ ａｌｓｏ ｔｅｓｔ ｔｈｅ ｐａｒａｌｌｅｌ ａｌｇｏｒｉｔｈｍ ｗｉｔｈ ｔｈｅ ＭＳ ｉｍａｇｅ ａｎｄ Ｐａｎ ｉｍａｇｅ ｏｆ ｔｈｅ ＩＫＯＮＯＳ ｏｎ ｔｈｅ ｐｃ ｃｌｕｓｔｅｒｓ ｓｙｓｔｅｍ ｗｉｔｈ ｇｏｏｄ ｒｅｓｕｌｔｓ， ｗｈｉｃｈ ｃａｎ ｂｅ ａ ｖａｌｕａｂｌｅ ｒｅｆｅｒｅｎｃｅ ｆｏｒ ｐｒａｃｔｉｃａｌ ａｐｐｌｉｃａｔｉｏｎｓ ｉｎ ｒｅｍｏｔｅ ｓｅｎｓｉｎｇ． Ｋｅｙ ｗｏｒ ｄｓ： Ｒｅｍｏｔｅ ｓｅｎｓｉｎｇ ｉｍａｇｅ ｆｕｓｉｏｎ， ＰＣＡ， Ｐａｒａｌｌｅｌ ａｌｇｏｒｉｔｈｍ
２ ＰＣＡ 融合的串行算法
主成分分 析 在 数 学 上 称 为 Ｋ－ Ｌ 变 换［２］， 是 在 统
计 特 征 基 础 上 进 行 的 一 种 多 维 （多 波 段 ）正 交 线 性 变
换。实际操作是将原来的各个变量（这些变量中部分
有相关关系）重新组合， 组合后的新变量是互不相关
的。此方法用于对全色图像和多光谱图像数据进行
其中， !表示!， !表示 ! ， !表示 ! 。
１
ｉ＝１
２
ｉ＝１＋ｓ／Ｎ
ｎ
ｉ＝（Ｎ－ １）Ｓ／Ｎ＋１
由上式可推得， 对于每个结点来说， 只需要计算∑
ＸＹ、∑Ｘ、∑Ｙ 三项的值， 再传递给主进程， 主进程便
能够由此计算出协方差矩阵。即主进程与其他进程
间的通信由原有的 ３ 次减少到一次即可， 比优化前
融合时， 具有显著的优势［３］。为了更好的理解本文提
出的并行算法， 首先简要介绍 ＰＣＡ 串行算法的主要
步骤。
（ １） 将 ＭＳ 图像作为矩阵输入， 计算输入矩阵
的协方差矩阵， 计算公式为：
ｓ
!（Ｘｉ－ Ｘ ）（Ｙｉ－ Ｙ ）
σｘｙ＝ＣＯＶ（Ｘ，Ｙ）＝ ｉ ＝ １ Ｓ－ １
（ １）
Ｘ、Ｙ 分别表示不同的波段， ｓ 表示 输 入 矩 阵 的
（ ４） 将 Ｐ１ 用 经 过 拉 伸 后 的 全 色 图 像 ｐａｎ 的 数 据进行替换，得到新的 Ｐ 矩阵。
（ ５） ｖ 矩阵的转置 ｖ’乘以 Ｐ 矩阵， 反变换回 ＲＧＢ 坐标系统。
分析上述串行算法， 可以看出在第 ３ 步和第 ５ 步中， 存在着大量的矩阵乘法， 是串行算法的计算 瓶颈。而在第 １ 步中， 由于需要多次扫描整个全局 的数据进行计算， 计算量更是相当大。
３ ＰＣＡ 融合并行算法 Ｐ－ ＰＣＡ 及其优化 为了提高算法效率， 满足实时性要求， 我们设
计了可扩展性能好、执行效率高的 ＰＣＡ 融合并行算 法 Ｐ－ ＰＣＡ（ Ｐａｒａｌｌｅｌ Ｐｒｉｎｃｉｐａｌ－ Ｃｏｍｐｏｎｅｎｔ－ Ａｎａｌｙｓｉｓ） 。 ３．１ 数据划分
常用的并行编程模型主要有： 任务并行， 流水 线并行和数据并行三类。由于 ＰＣＡ 融合方法对各个 象素单元的操作基本相同， 因此我们选择数据并行 模型。对数据集的划分采用的是垂直分割方法， 将 图像以水平线分割到每个结点， 保证每个进程间的 交互尽量少。设结点数为 Ｎ， 图像的总行数值为 ｍ， 若 ｍ 能整除 Ｎ， 则每个结点分配到的图像行数为 ｍ／ Ｎ； 若不能整除， 则非 ｒｏｏｔ 进程的结点分配到的图像 行数为 ｍ／Ｎ， 剩余的行数分配给 ｒｏｏｔ 结点进程， 为 ｍ－ ｍ／Ｎ×（ Ｎ－ １） 。这种数据划分简单易行， 且满足负 载均衡的要求。 ３．２ 并行算法流程
减少了进程间的通信量， 减少了无谓的通信等待时
表 １ 各种图像的执行时间（ 单位 ｓ）
进程数／图 １０２４×１０２４ ２０４８×２０４８ ４０９６×４０９６ ８１９２×８１９２
像大小
１
０．４２５８ １．６７２３ ６．７３６０ ２６．１６３１
２
０．２２２６ ０．８６５２ ３．４７３９ １３．４５９３
４
方差串行计算式 （ １） 所得到的最直接最简易的变
形， 在式（ ２） 中计算协方差矩阵时需要很多全局数
据， 因而要求各进程间进行多次数据通信来交换数
据， 但是频繁的数据通信会严重的影响并行效率。
为使并行程序的执行效率更高， 我们针对第 １ 步进
行通信优化处理， 经过数学变形后得到算式（ ４） ， 它
０．１１５５ ０．４２６５ １．７０８９ ６．６１３５
８
０．０５９４ ０．２１５１ ０．８５３８ ３．２８１０
１６
０．０３５２ ０．１２５８ ０．４９０６ １．８２６６
间， 加快了执行速度。 对于整个过程而言， 进程间的通信次数由原有
１ 引言 全 色 （ Ｐａｎ） 图 像 和 多 光 谱 （ ＭＳ） 图 像 是 最 常 见
的 两 种 遥 感 图 像 。 ＭＳ 图 像 包 含 了 地 面 目 标 丰 富 的光谱信息， 但它的空间分辨率不高； 另一方面， Ｐａｎ 图 像 尽 管 具 有 非 常 丰 富 的 空 间 细 节 信 息 ， 却 不包含光谱信息。实际应用中往往需要能同时具 备较高空间分辨率和光谱分辨率的遥感图像。多 源遥感图像的融合技术是解决这一问题的有效 途径， 能产生具有高空间分辨率的多光谱图像， 对 于 从 图 像 中 提 取 更 丰 富 、更 有 用 信 息 ， 它 是 一 种 既 有 效 又 经 济 的 方 法［１］。
ｉ＝１ ｉ＝１
＝
１ ｓ－ １
&
)
’
*
! ! ! （’
’ ’
Ｘｉ Ｙｉ ＋
Ｘｉ Ｙｉ ＋…＋
Ｘｉ Ｙｉ ）－
* * *
’ ’
１
２
ｎ
* *
’
*
’
*
! ! ! ! ! ! ’
１’
（ ’
ｓ’’
(
１
Ｘｉ ＋
２
Ｘｉ ＋…＋