基于条带波修复的影像云阴影去除

$\hspace{1cm}$对于遥感影像来说，去除云阴影并对对缺失数据进行填充是一个比较困难的问题，本文1提出了一个从遥感影像云或云阴影中修复数据的方法。本文提出的方法基于带波滤波转换和多分辨率几何分组方法，该方法主要有两个步骤：1.影像不同区域的几何流弯曲是通过带通转换和多尺度分组方式获取的。第一步完成后能够获取比较好的结构信息，在云遮挡区域，然后进行第二步纹理填充，在云遮挡区域通过区域中几何弯曲流信息扩散进行合成。

带波变换和多分辨率分组

$\hspace{1cm}$在影像边界处的影像规则是各向异性的规则，因此尽管影像在轮廓处是不连续的，影像在平行与边界弯曲的方向的切线方向是可区分的。通过带波变换2采用各向异性的规则构建正交向量，正交向量在函数具有最大平整度的方向延长。

带波变换

$\hspace{1cm}$Bandelet变换是一种基于边缘的图像表示方法，能自适应地跟踪图像的几何正则方向。构造Bandelet变换的中心思想是定义图像中的几何特征为矢量场，矢量长表示了图像空间结构的灰度值变化的局部正则方向，带波基并不是预先确定的，而是通过优化最终的应用结果来自适应的选择具体的基。对于几何正则图像，在局部范围内几何流是平行的。

Bandelet变换算法流程：
（a）将哈尔小波系数矩阵分解成边长为8的小方块,将小方块内的系数加载到图形处理器的共享存储器中，对于每个小方块，在（0，π)区间内，等间隔地抽取32个角度值作为几何流的方向，用符号INF表示无几何流方向；
（b)采用投影值公式，计算每个小方块内哈尔小波系数的坐标在抽取的32个几何流方向上的正交投影值；
（c）将正交投影值按从小到大排序，得到一个长为64的排序索引,对于无几何流方向，排序索引为数字1～64；
（d）将小方块内的哈尔小波系数，按照排序索引进行重新排序，得到一个一维信号，对一维信号进行一维哈尔小波变换，得到小波变换后的一维信号，对小波变换后的一维信号进行量化处理，得到量化后的一维信号，对量化后的一维信号进行一维哈尔小波逆变换得到重构信号；
（e）计算估计误差，选择使估计误差最小的方向作为小方块的最佳几何流方向；
（f）采用投影值公式，计算每个小方块内哈尔小波系数的坐标在最佳几何流方向上的正交投影值，将正交投影值按从小到大排序，得到一个长为64的排序索引，若最佳几何流方向是无几何流方向，排序索引为数字1～64,，使用排序索引重新排序小方块内的哈尔小波系数，得到一维信号，对一维信号，进行一维哈尔小波变换，得到哈尔小波变换后的一维信号，使用排序索引，对哈尔小波变换后的一维信号进行重新排序，得到Bandelet变换系数矩阵；
（g）将Bandelet变换系数矩阵转换为整型值；
（h）将Bandelet变换系数矩阵从显存传回内存；

$\hspace{1cm}$Bandelet先定义了一种能表征图像局部正则方向的几何矢量线；再对图像的支撑区间$S$进行二进剖分$S=\cup_i\Omega_i$，当剖分足够细时，每一个剖分区间$\Omega_i$中最多只包含图像的一条轮廓线(边缘)。在所有不包含轮廓线的局部区域$\Omega_i$，图像灰度值的变化是一致正则的，因此，在这些区域内不定义几何矢量线的方向。而对于包含轮廓线的局部区域，几何正则的方向就是轮廓的切线方向。根据局部几何正则方向，在全局最优的约束下，计算区域$\Omega_i$上矢量场$\tau(x_1,x_2)$的矢量线，再沿矢量线将定义在$\Omega_i$的区间小波进行Bandelet化(bandeletization)以生成Bandelet基，以能够充分利用图像本身的局部几何正则性。Bandelet化过程实际上是沿着矢量线进行小波变换的过程，就是所谓的弯曲小波变换，于是所有的剖分区域$\Omega_i$上的Bandelet的几何构成了一组$L^2(S)$上的标准正交基3。
$\hspace{1cm}$在几何方向上，假设图像方程是规则的，则对于这种几何正则图像，图像沿平行于边缘线的方向的变化是正则的，而垂直于变换线的方向则变化剧烈。图像中的的几何流一般指定义在图像的支撑区域上的向量场$\tau(x_1,x_2)$,用于描述每一点正则变化的方向。Bandelet基是Bandelet变换的核心，为了构建bandelet基，应该首先建立四叉树结构，以便于寻找图像的几何方向。一般先对图像进行二维双正交离散小波变换，对变换后的多尺度分解系数进行处理假设原始图像为$f(x,y)\,,$则$f(x,y)$的二维小波变换为：

$$\left\{ \begin{array}{lr} &<f,\psi_{j,m}^`>=f*\psi_{j,m}^`(m_12^j,m_22^j) \hspace{1cm} (j\in z,m1 m2\in z^2,s\in{H,V,D})\\ &\varphi_{j,m}(x)=\varphi_{j,m1}(x_1)\varphi_{j,m2}(x_2) \\ &\psi_{j,m}^V(x)=\psi_{j,m1}(x_1)\varphi_{j,m2}(x_2)\\ &\psi_{j,m}^D(x)=\psi_{j,m1}(x_1)\psi_{j,m2}(x_2)\\ &\psi_{j,m}^D(x)=\varphi_{j,m1}(x_1)\psi_{j,m2}(x_2)\, \end{array} \hspace{1cm}(1) \right.$$

其中$\varphi_{j,m}$表示粗略逼近$\psi_{j,m}^H(x)\,,\psi_{j,m}^V(x),\psi_{j,m}^D(x)\,,$分别是水平、垂直和对角方向子带的高频信息，而尺度函数$\psi(t)$和小波函数$\varphi(t)$的伸缩和平移量分别为：
$$\begin{array}{lr} \psi_{j,m}(t)=2^{-j/2}\psi(2^{-i}t-k)\\ \varphi(t)=2^{-j/2}\varphi(2^{-j}t-k) \end{array} \hspace{1cm}(2)$$
对于变换所得的多尺度图像$f_j$,除了低频部分外，每个尺度下有三个方向的变换，进一步用二进四叉树剖分的方法对每个方向的变换系数进行处理。二进分割时先将子带分成四个小子带，然后对每个小子带在下一层分割中继续分割为四个更小的子带，直到达到预设的最小尺度$J_{min}$。获取分割子带后需要计算每个子带的最佳几何流，计算最佳几何流的流程为：

(1)对角度进行采样：设要处理的子带的大小$L*L$，则将圆周角$[0,\pi)$等间隔离散为$L^2-1$个，则可能的取值为： $$\theta=\frac{k\pi}{L^2-1}\,,k=0,0,...,L^2-2\hspace{1cm}(3)$$ 对于无几何流的情形，则标记为$Inf$，则一共可以获取$L^2$个采样角度：$(\theta,Inf)$。
(2)曲波变换：设子带大小$L*L$，则可以构造一个$L*L$大小的网格，利用公式(4)计算在采样角度上逐个网格点的正交投影误差： $$t=-sin(theta)*x(i)+cos(theta)*y(j)\hspace{1cm}(4)$$ 式中$x(i)\,,y(j)$分别是所构造的网格点的坐标，$theta$是(3)中所得的采样角度。网格点根据误差值大小，对网格点进行从小到大的排列，根据网格点的位置，找出对应位置的子带系数，对子带系数再进行一次一维小波变换，以实现曲波变换。
(3)几何流计算：采用Lagrangian系数计算几何流，设用$R_b$表示进行量化后，Bandelet系数编码所需的比特数，量化公式$Q(x)$为：
$$Q(x)= \left\{ \begin{array}{lr} 0,\,|x|\le T\\ sgn(x)(q+1/2)T,qT\leq|x|<(q+1)T \end{array} \hspace{1cm}(5) \right.$$
其中，$T$是设定的量化阈值，一般是根据经验事先设定，$R_g$代表几何流编码的比特数（也就是方向个数），$q=\lfloor(|x|/T)\rfloor$。拉格朗日系数为： $$L(f_{\theta},R)=|||f_{\theta}-\tilde{f_\theta}||+\lambda T^2R_g+R_b\hspace{1cm}(6)$$
上式中$\lambda$为拉格朗日橙子，$\tilde{f_\theta}$表示重构出的一维信号，是步骤(2)中所得的一维小波系数采用量化公式(5)量化后的结果。得到最佳几何流方向后剩下的步骤为对Bandelet块内存储的小波

带波算子正交基

$\hspace{1cm}$带通近似4可以通过计算正交阿尔伯特正交基的多项式阈值来获取。而阿尔伯特转换可以认为是一个适应于不规则采样的一个弯曲小波变换。
$$b_{j,l,n}^k(x)=\Sigma_pa_{l,n}[p]\psi_{j,p}^k(x)\,.\hspace{1cm}(7)$$
$\hspace{1cm}a_{l,n}[p]$是阿尔伯特转换的相关系数，$a_{l,n}[p]$取决于局部几几何流，因为局部几何流定义了弯曲采样的位置

在位置$n$和分辨率为$2^j$的情况下，带波函数被定义了。指数$k$表示小波方向。阿尔伯特转换引入了一个新的分辨率因子$2^l>2^j$，这个分辨率因子定义了带波函数的延展率，带波$b_{j,l,n}^k(x)$继承了小波$\psi_j^k(x)\,.$的规则。

分割几何流

正交带波族取决于对每一个分辨率$2^j$和方向$k$定义的局部适应流。此种平行流通过对曲线的积分进行特征化。为了通过多项式流你和几何模型，我们需要将在S中的一系列小波系数进行分割。

影像修复

设置$x$是在关联几何场中任意对于$x\in\mathbb{R^n}$，我们构建了n维向量空间，叫做正切空间，表示为$T_x(\mathbb{R^n})$，关联几何场是一个一一映射，对于每一个点$x\in\mathbb{R^n}$都有一个向量$v\in T_x(\mathbb{R^n})$。假设一个点开始于$x_0\in\mathbb{R^n}$在时间$t=0$处，根据关联场中描述对的速度方向移动，这样其轨迹可以被描述为：$C:[0,\varpropto)\to\mathbb{R^n}$，而轨迹表示微分方程的曲线积分，其初始条件为$C(0)=x_0$，方程形式为： $$\frac{dC}{dt}(t)=f(C(t))\hspace{1cm}(8)$$

对于方程(8)在任意时刻，其积分形式可以描述为： $$C(t)=x_0+\int_0^tf(C(s))ds\hspace{1cm}(9)$$
求解趋势线问题转换为求解微分方程的问题，积分曲线通过固定的速传播到云遮挡区域。为了保证积分曲线的连续性，我们求解以下方程的最小值: $$l=\int|C'(t)|dt\hspace{1cm}(10)$$ 方程（10）中$C^`(t)$是$C$的方向导数，对每一个$\theta$的方向通过带波转换获取：
$$C'(t)=(x'(t)cos(\theta),y'(t)sin(\theta))\hspace{1cm}(11)$$
整个实验过程：
1.定义修复区域和修复区域边界
2.在边界处修复填补影像信息
2.1在分辨率$2^j$的影像上计算二维正交小波变换
2.2计算分割后的小波系数
$$m_{even}=argmin_{m\in\Omega_{even}}\Sigma_{|n|<P}|a[m-n]-a[m_{odd}-n]|^2\hspace{1cm}(12)\,.$$
2.3根据相关场定义相邻像素对的权重均值和差异
2.4计算边界点的积分流曲线
2.5修复方程为：
$$C(t_n+1)=C(t_n)+\Delta tC^`(t_n)\hspace{1cm}(13)$$
2.6一直重复以上步骤直到求解结果稳定