$\hspace{1cm}$本来是想着写一篇论文笔记，可是老板不在，实在是没有心情去写什么论文笔记了，突然想起来自己写过一篇饥荒游戏攻略，既然写了第一篇那就应该好好接下来写一个系列出来吧，由于最近也没有玩得太多，目前的存档在五十多天，对于许多刚开始玩饥荒的新手玩家来说第一个冬天是特别难熬的，主要因为冬天温度低，各种作物都不生长了，难以获得食物，其实冬天并不难过，将基地建在兔子洞附近，冬天就有源源不断的兔子肉可以吃，当然也可以养蜂蜜，吃蜂蜜也是极好的。
$\hspace{1cm}$闲话不多说，假设看这篇攻略的都是新手玩家，急需熬过第一个冬天，那么恭喜你，这个教程绝对是你过冬的神器。我们在开局选择秋天，白天可以选择长日，这样白天时间会略长一点，对于新手玩家来说要简单许多。好了我们角色还是选择威尔逊，刚开始的时候我们主要目标在于找到适合建基地的位置，当然针对不同的地图可能略有不同，首先收集基础材料，树枝、石头、干草、燧石等材料，就我个人的看法树枝、干草、木头各收集一组；石头和燧石多多益善，除了基础材料之外就是初级食材，初级食材是保命的东西必须要多多收集，主要是浆果和胡萝卜。一般在一边收集基础材料一边收集食物，然后要做的就是开图，开图的攻略有很多我这里也不会详细说明，总之就是到处跑，把整个地图分为几大块搞清楚，然后确定下面几个地方：1.矿区，矿区在开始的时候可以提供金子，这是建造基地必备的东西；2.沙漠地区，沙漠地区会有很多兔子洞，在这个地方建造基地是极好的，当然也可能有牛群，牛群也是一个好东西；3.找到红树林和里面的猪王，挖矿这个事情毕竟不能长久，所以找到猪王可以用肉跟猪王换金子，这个还是挺值得的；4.还有确定一些如森林，墓地，沼泽等地区，这些地方在刚开始的时候不建议进去深探，毕竟比较危险，整个地图会如图所示：

上面是我开对的地图，没有将所有位置都开发出来，不过找到了一些主要的位置，如草原，红树林，森林以及矿区等，这些位置找到了之后该有的资源如：黄金，石块等应该足够了，然后找到一个好地方建基地，还是这样一个原则，在兔子窝边上建造基地，这样能够方便过冬：

从上面两幅图我们可以看到我把基地建立在兔子窝边上，这样可以方便的捉兔子，然后在靠近河岸的地方，在这个地方猎狗不会从四面八方围过来，关于猎狗我在这里要提一句，刚开始猎狗就一两只比较少，带上猪人或者穿上木甲，用长矛就可以干死，到后来猎狗多了最好是种上食人花，这样猎狗就会被食人花打死完全不用担心。
$\hspace{1cm}$到这个时候应该在6-10天，事件还很充分，然后我们接下来的工作就是移植，我的意见是移植一组干草，一组树枝，一组浆果就足够了，这个工作大概要花费两天的时间，然后我们就要准备造一些必备的设置了，如干燥架，锅等，造这些东西估计也就两天时间，这样我们已经过了15天左右，然后我们去打死几只蜘蛛或者捉到兔子，带上猪人去杀蜘蛛，杀蜘蛛这样的事情是为了得到怪物肉，得到足够的怪物肉之后去猪村刷猪肉，杀死疯猪可以得到两个大肉和一个猪皮，所以杀疯猪比较划得来，一般来说喂猪吃四个怪物肉可以变疯猪，通过这种方法能够获取足够的猪肉和猪皮。得到猪肉和猪皮之后就可以回去建猪屋，这样以后利用猪人防御和杀猪都比较方便。
$\hspace{1cm}$我开的这个地图运气比较好，有干草彩蛋，但是附近有猪人守卫，不过这也难不倒我，我在边上种了一个蜘蛛卵，这样主任守卫和蜘蛛就会互殴，我只要默默的检肉就好了，然后在适当的时候把干草彩蛋收刮一下，爽到爆！！不过不是所有的地物都有彩蛋的，所以这个是可遇不可求的事情，至于怎么种蜘蛛卵，这个其实很简单，一个三级的蜘蛛巢穴，如果没有编程蜘蛛女王打爆他就有蜘蛛卵了。好了，到了这个时候应该过完一整个冬天是没有问题了。

MAP数学模型

令$z_{x,y}$和$g_{x,y}$分别为输入的辐亮度和传感器的响应，在位置$(x,y)$则传感器的成像模型可以表达为：
$$g_{x,y}=A_{x,y}z_{x,y}+B_{x,y}+n_{x,y}\hspace{1cm}(1)$$
则公式（1）的向量表示形式为：
$$\textbf{g}=\textbf{Az}+\textbf{B}+\textbf{n}\hspace{1cm}(2)$$
假设$\textbf{z}$为正确的影像，$\textbf{g}$为退化影像，则影像修复可以表示为：
$$\hat{x}=argmax_zp(\textbf{z}|\textbf{g})\hspace{1cm}(3)$$
则上式的贝叶斯公式表示为：
$$\hat{x}=argmax_z\frac{p(\textbf{g}|\textbf{z})}{p(\textbf{g})}\hspace{1cm}(4)$$
由于在$p(\textbf{z}|\textbf{g})$中$\textbf{g}$是独立的，因此公式（4）可以描述为：
$$\hat{x}=argmax_zp(\textbf{g}|\textbf{z})\hspace{1cm}(5)$$
公式（5）中等式右边第一项为：
$$p(\textbf{g}|\textbf{z})=\frac{1}{M_1}exp \{-\frac{1}{2}(\textbf{g-Az-B})^T\textbf{K}^{-1}\textbf{g-Az-B}\}\hspace{1cm}(6)$$
其中$\textbf{K}$为噪声的协方差矩阵，假设噪声在各个波段是独立的，则$\textbf{K}$为一个对角矩阵，则可以简化为：
$$p(\textbf{g}|\textbf{z})=\frac{1}{M_1}exp\{-\frac{1}{2}\lVert Q(\textbf{g}-\textbf{Az}-\textbf{B}) \rVert^2\}\hspace{1cm}(7)$$
其中$Q$为对角矩阵。
而公式(5)中等式右边第二项为：
$$p(\textbf{z})=\frac{1}{M_2}exp\{-\frac{1}{2\lambda}\Sigma_{x,y}\Sigma_{c\in C}\rho(d_c(z_{x,y})) \}\hspace{1cm}(8)$$
其中$M_2$是一个常数，$c$是一个影像类别数据集$C$中的一类，则数量$d_c(z_{x,y})$是像素$z_{x,y}$的空间测度，通常构成形式为一阶导或者是二阶导，$\rho(\cdot)$是Huber函数1：
$$\rho(i)=\left\{ \begin{array}{lr} i^2\quad |i|\le \mu \\ 2\mu|i|-\mu^2\quad |i|\gt\mu \end{array} \hspace{1cm}(9) \right.$$
对于$d_c(z_{x,y})$,在四个方向上计算二阶导，则有：
$$\begin{array}{lr} d^1_c(z_{x,y})=z_{x-1,y}-2z_{x,y}+z_{x+1,y}\quad(10)\\ d^2_c(z_{x,y})=z_{x,y-1}-2z_{x,y}+z_{x,y+1}\quad(11)\\ d^3_c(z_{x,y})=\frac{1}{\sqrt{2}}[z_{x-1,y-1}-2z_{x,y}+z_{x+1,y+1}]\quad(11)\\ d^4_c(z_{x,y})=\frac{1}{\sqrt{2}}[z_{x-1,y+1}-2z_{x,y}+z_{x+1,y-1}]\quad(11) \end{array}$$
将式（7）、（8）代入到公式（5）中有：
$$\hat{z}=argmin\lambda\lVert \textbf{Q}(\textbf{g}-\textbf{Az}-\textbf{B}) \rVert^2+\Sigma_{x,y}\Sigma_{c\in C}\rho(d_c(z_{x,y}))\hspace{1cm}(14)$$

假设和方法总纲

基于时序轮廓指数的方法通过匹配相似的时序光谱形状的像素来对位置像素值进行估计，此方法有三个主要的步骤：

1.时序光谱轮廓对准，对准的目的在于确保具有较高的轮廓相似度，且保证区别是由于土地覆盖响应差别引起的；
2.协同使用完整的和修复的轮廓，确保在进行脚哦正过程中有足够多的像素参与估计；
3.使用时序和空间信息解决突变点的模糊性的问题，采用KNN算法进行估计

（1）时序轮廓对准：
假设$x_i$和$y_i$为两个像素时序光谱，其中$i=1,...,t,$ $t$为时序光谱维数，假设两个光谱只在响应的强度上有差异，响应想去差异为$T$则在比较两个光谱相似性时需要减去强度的差异，比较$x_i$和$y_i-T$，且$T$的计算可以表示为：
$$T=\frac{1}{t}\Sigma^t_{t=1}(x_i-y_i)\,.\hspace{1cm}(1)$$
（2）使用完整的和未修复的像素时序轮廓：
其方法的有点可以描述为：
假设一个5时相的数据其中有一个时相数据发生了缺失表现为$11101$，其中$1$为可观察到数据的模式$0$为不可观察到数据，则对这样的缺失修复有（$11111\,,01111\,,11110$）等，我们将$11111$成为完整的轮廓，而在其他位置存在缺失的成为补充不完整轮廓。
（3）模糊性问题：
在对突变区域进行修复的过程中，如果只采用时间光谱轮廓信息则可能存在多种情况，而面对多种可能，求均值显然并不是最好的处理办法，因此考虑空间临近性，给空间临近性比较高的权值进行求解。

使用PBI方法进行缺失像素估计

PBI方法使用KNN的估计方法，KNN的估计方法优点在于操作简单，只使用最近邻的数据进行估计，且是非参数的方法，避免了复杂的参数求解的过程；而使用PBI方法进行缺失像素估计主要分为两个步骤：1.通过KNN方法选择与缺失像素光谱最相似的像素，对这些时序光谱进行分类；2.通过空间的临近性对不同类别的光谱给与不同的权重。

江湖风雨如晦，我自不动如山

侠客情似剑，文人笔如刀。
江湖六千里，有酒醉今朝。
九剑惊天下，武帝城墙高。
举世皆敌手，气运江湖老。
一剑两千甲，两袖青蛇摇。
风流曹官子，气散广陵道。
木剑出江湖，折剑有情长。
天下熙熙来，世人攘攘往。
立地入圣贤，情仇自此消。
江湖风雨近，庙堂心意高。
谁揽天下事，我自任逍遥。

金庸封闭古龙逝，江湖唯有英雄志

前一段时间看了孙晓的《英雄志》，不得不承认，这部小说当时后金庸时代的巅峰之作，文章构思精妙，结构紧凑，细节描写细致入微，同时又有着深刻的内涵。与燕垒生的《天行健》相比，在宏大场面上的描写要稍微逊色一些，然而在人物细节的描写和对武侠内涵的深刻阐述上《英雄志》堪称一流。
小说以明朝英宗为背景，描写了朝堂江湖的纷争，面对当世巨变（改朝换代），所有人都不能置身之外。不管是天下第一宁不凡还是剑神卓凌然；不管是权臣江充还是市井小民卢云。但是每个人都有着其追求和坚持的东西，这种追求和坚持才是侠之根本。
秦仲海那句此生不必跪人是多少人心中的理想！而卢云心中所求唯有正道而已，至于杨肃观，虽成修罗仍不改其道，那就是天道。弱肉强食这就是天道，从他对待吃羊的老虎的态度可以看出来。舍身求道，虽九死亦未悔。所以秦仲海最终战死沙场，杨肃观也难免一死，卢云舍生求道。但是心中所求就一定正确？卢云虽然心怀天下，然而却辜了一直深爱他的顾家小姐和银川公主。杨肃观为达目的可以利用所有人。
有时候不免会思考，到底什么是对什么是错，有的人也许生前风光无限死后却遗臭万年，有的人生前凄凄惨惨，最终落得青史留名，到底有什么值得？这不禁让我想到金庸的封笔大作《鹿鼎记》，韦小宝是一个不折不扣的小人物，没有郭靖黄蓉的侠气，没有令狐冲的豪放，也没有杨过的深情。这样一个人物如何能够称得上是侠，在江湖中，我们都是小人物，每个人都有自己的坚持也有自己的顾忌。有时候牺牲也许并不是伟大，偷生也不一定渺小。
英雄在世，难道就是宁叫我负天下人？难道牺牲身边爱我们的人去为了所谓的天下就是英雄?天行有常，连爱我们的人我们都可以舍弃，那拿什么去爱天下人，江湖太大，天下太广，君且待何如？

喜剧之王

本来想写一篇关于喜剧之王的观后感，可是又觉得怎么也写不出来，网上看了一篇觉得不错，跟大家分享一下，原文地址
寻找一份，能放飞梦想的坚持
文：米新磊
把所有的观察，怀疑，思索，挣扎，愤懑，不平都收束于幕后，沉淀一个老少咸宜的结论，那是喜剧。
一些重而痛的东西爬上脸庞，留下皱纹、沟壑、胡茬，而后一个转身，一副刀枪不入的表情抚平一切，那是喜剧演员。
尽可能放弃思考和敏锐。否则，从不老实归位的面部肌肉某个松弛的瞬间，透露出些微隔靴搔痒的真实，令人浮想联翩甚至悲从中来……
这是有一次在《看电影》扉页上看到的三段话。在这三段有点深刻的话语旁边，是周星驰和吉姆 凯瑞略显深沉的面容。这两位现当今东西方最出色的喜剧演员，演艺经历竟也出奇的相似。
都因为其独特的表演风格受到大众的喜爱。都因为其喜剧演员的身份在前期遭到主流奖项的排斥。而他们都因为自己默默的坚持最终拨开云天。
1998年，吉姆 凯瑞凭借《楚门的世界》一片展露了他隐藏在那张搞怪面孔下的精湛演技，夺得金球奖最佳男主角。
1999年，周星驰自编自导自演了《喜剧之王》。这部真实的喜剧展示了周星驰最让人难忘的真正演技。让无数影迷唏嘘不已感慨万千。
第一次看《喜剧之王》的时候还小，看完后第一感觉好像像鸡肋，食之无味，弃之可惜。现在想想原因，可能正是它不像当时眼中“正统”的周氏喜剧一样让我笑的癫狂吧。某些地方总是感觉笑的不痛不痒难以尽兴。
这便是我对它的第一印象。
第二次看的时候已经隐隐体会出其中的某些不同之处了，扫除了之前先入为主的不好的感觉，让我较为客观的从头到尾审视了一遍，然后得出了一个结论：其实，这是一部正剧。
第三次观影的时候我已经无可救药的陷入其中了。我也在笑，但觉得胸口很堵，有某些东西在心底涌动翻腾就是找不到出口。于是我发现我爱上了这部影片，以后每看一次都不会有丝毫的厌倦，都会有不同的感触，唯一不变的一点是，心中的暗涌依旧。
“努力∼奋斗∼” 尹天仇对着宽广的大海如是说。
“喂，阿姨，什么时候开工啊……”尹天仇对着被戛然挂断了无数次的电话如是说。
“其实呢，我是一个演员” 尹天仇一遍遍的对着每个人如是说。
我相信每个人都曾用一个梦想去鞭策过自己，刚开始仿佛什么都充满了希望，这是梦想的光辉点亮的希望。
可《肖申克的救赎》里的瑞德曾经严厉的警告安迪：“在这个地方，希望是个危险的东西。”
我明白他的意思，希冀绝对无法实现的希望，带来的痛苦会让把生命都吞嗤殆尽。
所以有了越来越多的人开始诅咒希望，他们曾经强烈希冀过的梦想但却变成了泡影。他们的诅咒给越来越多的摇摆不定着带来了不用辛苦坚持的借口，然后他们也加入了这个队伍……不能说他们是失败者，路是自己选的，走什么样的路是每个人的自由。可总会有这样的人，坚定的向梦想挪动，即便步伐小的让人绝望，却义无反顾的固执。比如，尹天仇。
有人说这是周星驰半自传体的作品，只有自身经历过才会有这样真切的感受。
他每天都去剧组跑龙套，即使阿姨对他下了驱逐令也仍然不放弃要找一个小小的角色。
他每天都翻开斯坦尼拉夫斯基的《演员的自我修养》，合上书本后的他显得格外平静。
他每次都会耐心的告诉骂他是“死跑龙套”的人说：“其实我是一个演员。”
他不会对谁发脾气不会对谁大吼大叫不会对谁有丝毫偏颇，好像周围的一切挫折都渺小的不足以关注。
可我可以清楚的看到他内心在挣扎，那种每天的失望带给他的不是麻木而是更大的失望。但灰色的事实却很快被一个坚定的信念还原成出了亮丽的色彩。“我是一个演员”，这个梦想是如此的绚丽如此强大能支持着他一直貌似平静的走下去。
但他很孤独，没有人去了解他丰富隽永的内心，虽然有梦想做伴，但却找不到倾诉的共鸣。一个人的旅途是残缺的，还好有飘飘的出现。
柳飘飘也是个曾经有梦的人，但现实刺痛的回忆让她把这些深深的埋到了最心底。她害怕回忆起来的那种焦灼，所以当第一次尹天仇问到她的初恋时她才会恼羞成怒转身离去。埋在心底的隐痛拿出来接受阳光的曝晒是需要勇气的，还好飘飘并没有像大多数人一样放弃梦想，只是把它藏了起来。当她遇到天仇后，她的梦想开始复苏了。他们不是因为痛苦而在一起寻找慰藉，他们是因为缘分而在一起放飞自己的梦想。
所以当尹天仇大声的对着飘飘喊：“我养你啊！”飘飘在车里无法抑止的泪流满面，这是梦想还原那一刻的悸动。
所以当尹天仇因为娟姐的器重得到男主角的角色时，他不可抑止的欣喜冲刷掉了以往所有的失落。这是梦想实现那一刻的激动。
如果梦想就这么实现，他们之前也就不会经历那么多的挫折了。我心里这样告诉自己，虽然我并不想这样的情景出现。
尹天仇失去了扮演男主角的机会，他死死的抓着剧本，想要紧紧的抓住一些东西。那是梦想的所有载体，咫尺之远，却眼看着要逝去。我终于还是看到了他失落的表情，那是一种根本无法再用平静来掩饰的失落，从来没想过放弃的梦想像外面的雨一样在凋零。
“你还是走吧！”剧组的人员最后还是忍不住说出了这句话。
我不知道周星驰在现实的奋斗过程中遇到过多少次这样的挫折，那种无可奈何的失落表情让我觉得心底都被掏空了一般。
同样的表情，同样的感觉也来自飘飘。看到和自己的偶像娟姐在一起的天仇，看到娟姐呵护关切的话语，她好像明白了些什么。勉强的微微一翘的嘴角，牵动着脸上因为倔强而受到的新伤。她想也许这就是现实吧，也许自己本来就不应该做梦吧！面对着驶远的车子，她自嘲一样的失落的说了一句：“就这样吧！”这句话说出来平静如斯，却承载了多少痛苦和压抑。
那一刻，我突然想，如果这不是一部电影，飘飘是否还会在几乎放弃的最后一刻，突然大喊着问一句：“你上次说养我，是真的吗”。如果不会的话，是不是真的“就这样了”呢？
于是，我也会想，如果这不是一部电影，尹天仇是否还会在那一次的失败之后最终得到了寻觅以久的成功呢？如果不会的话，他还能坚持的走多远呢？
毕竟，这不是现实，有人说现实不会出现这样的结果。
但是，这却无比真实，真实源自于最心底的向往和追寻。
人人心底都期待这样的结局，为什么要在转瞬间因为现实林林总总而败下阵来。
如果连去想的勇气都没有，梦想就真是一个悲哀的词了。
我决定心存感激，因为电影总是把最完整的人生和最美好的结局留给我们去回忆。我决定心存希望，因为尹天仇告诉我即使过程有多艰难也无所谓，携着希望的人总会找到幸福的出口。
在这个旅程中，总会有一个懂你的人，一直陪着你，直到放飞彼此的梦想为止。

TV模型数学理论

$\hspace{1cm}$TV模型与基于偏微分方程的主要用于修复较小的图像缺失区域和噪声，假设$u^0$是一个定义在二维平面$\Omega$上的光滑影像函数， $D$是一个待修复的区域，其区域直径为$d$，令$u^0|_d$为在区域$D$上$u^0$的限制，则影像的修复过程可以描述为获取$u_D$在$u^0|_d$上的拟合。
$\hspace{1cm}$如果区域$D$的直径$d$收缩到0，则可以采用线性插值的方法获取任意光滑的测试数据$u^0$：
$$\lVert u_D-d^0|_D\rVert_{\infty}=O(d^2)\hspace{1cm}(1)$$ 同样，影像修复过程可以使用k阶修复方式：
$$\lVert u_D-d^0|_D\rVert_{\infty}=O(d^{k+1})\hspace{1cm}(2)$$

二阶格林公式1影像修复

数学基础

假设$\Delta$表示拉普拉斯算子： $$\Delta u:=\frac{\delta^2u}{\delta x^2}+\frac{\delta^2u}{\delta y^2}\,.\hspace{1cm}(3)$$ 在区域$D$上的二阶格林公式为： $$\int_D(u\Delta v-v\Delta u)dxdy=\int_{\Gamma}(u\frac{\delta v}{\delta \textbf{n}}-v\frac{\delta u}{\delta\textbf{n}})ds\,,\hspace{1cm}(4)$$ 其中：$u$和$v$为在闭合环$D$上定义的二维函数，n为区域边界$\Gamma$上定义的朝外的法向量，$s$为边界的长度。
$\hspace{1cm}$假设$G(z_0,z)$是一个基于在$D$上基于珀松方程的格林函数，格林函数的意义为：$z_0$处所产生的脉冲对$z$处的场的作用，他是距离的函数，则对于任意一个$z_0=(x_0,y_0)\in D$，作为对应的点$z=(x,y)\in D$则$G(z_0,z)$求解为：
$$-\Delta G=\delta(z-z_0)\,,G|_{\Gamma}=0\,.\hspace{1cm}(5)$$ 则应用格林二阶函数（$u=u^0(z),v=-G(z_0,z)$）则讲$u,v$代入公式（4）中，可以得到：
$$u^0(z_0)=\int_\Gamma u^0(z(s))\frac{\delta(-G(z_0,z))}{\delta\textbf{n}}ds+\int_DG(z_0,z)(-\Delta u^0(z))dz\hspace{1cm}(6)$$ 其中$dz=dxdy\,.$公式(6)中的第一项表示为$u^h(z_0)$表示$f=u^0|_\Gamma$的光滑的扩展。

---

到公式（6）开始我就开始看看不懂了，简单的来说公式（4）就是格林公式，然后将$u,v$代入到公式（4）中得到公式（6），但是具体怎么得到的由于数学基础不够没法进行数学证明

---

$$dw_{z_0}=\frac{\delta(-G(z_0,z))}{\delta\textbf{n}}ds\hspace{1cm}(7)$$ 为边界$\Gamma$与源点$z_0$关联的和谐程度的度量，而非和谐程度可以表示为$u^a=u^0-u^h$且满足泊松方程2: $$\Delta u^a(z)=\Delta u^0(z)\,,z\in D\,,\qquad and\,\,u^a|_\Gamma=0\hspace{1cm}(8)$$
假设$d$为待修复影像区域$D$的直径$G(z_0,z)$是泊松分布的格林函数，则有以下条件成立：
$$\int_DG(z_0,z)dxdy\le\frac{d^2}{4}\,.\hspace{1cm}(9)$$
公式（9）的证明有兴趣的可以参看原文3 ，基于上述理论我们可以很简单的建立准确的基于二阶格林公式的影像修复方法：
（a）Linear inpainting via harmonic extension 假设我们只通过和谐扩展修复$u^0|_D$，则线性修复的表达为：
$$\lVert u^h-u^0|_D\rVert_\infty=O(d^2)\hspace{1cm}(10)$$
当$d\to0$时公式（10）成立，根据公式（6）均匀性修复误差正好就是非均匀性组份$u^a$。由于$u^0$为固定的光滑函数，因此存在常量$M$使得满足:
$$|\Delta u^0(z)|\le M\hspace{1cm}(11)$$
对于所有的$z\in D\,.$上式都成立，且对于$z_0\in D$根据公式（9）则有：
$$|u^a(z_0)|\le M\int_DG(z_0,z)dz\le \frac{Md^2}{4}\,.\hspace{1cm}(12)$$
（b）Cubic inpainting via Green’s formula 为了提高影像修复的精度，除了需要对整体结构进行修复外还需要对细节组份$u^a$进行修复。假设$u^\Delta_D$为线性修复的结果，然后我们修复$u^a|_D$通过$u^a_D$根据如下积分方程：
$$u^a_D(z_0)=\int_DG(z_0,z)(-u^\Delta_D(z))dz\,,\hspace{1cm}(13)$$
或者通过求解已知的泊松方程：
$$-\Delta u^a_D(z)=-u^\Delta_D(z)\,,\quad z\in D;\quad u^a_D|\Gamma=0\hspace{1cm}(14)$$
最后将这些新的细节信息加到均匀性修复过程中，可疑得到更加准确的修复结果$u_D$
$$u_D(z)=u^h(z)+u^a_D(z)\,.\hspace{1cm}(15)$$

局限性

实际上上述模型在大部分情况下是不适用的，因为：
1.影像通常情况下来说并不是平滑的函数，而是包含各种边界信息和非连续性的
2.影像在统计上也不能完全一致，由于存在影像噪声

全局变差模型的影像修复

模型数学基础

如上图所示，其中$D$为待修复区域，待修复区域具有分段光滑的边界$\Gamma$，$E$为待修复边界的外环。假设$u^0|_E$是被白噪声污染的区域，则基于变差修复模型为在扩展域找到变差函数$u$使得如下规则函数取得最小值：
$$R[u]:=\int_{E\bigcup D}r(|\Delta u|)dxdy\,,\hspace{1cm}(16)$$
在以下约束下：
$$\frac{1}{Area(E)\int_E|u-u^0|^2dxdy=\delta^2\hspace{1cm}(17)}$$
对式（16）（17）求解其欧拉拉格朗日方程为：
$$J_\lambda[u]=\int_{E\bigcup D}|\Delta u|dxdy+\frac{\lambda}{2}\int_E|u-u^0|^2dxdy\,,\hspace{1cm}(18)$$
其中$\lambda$为拉格朗日乘数子，则剩下的问题为求解这个欧拉拉格朗日方程。

$\hspace{1cm}$实在是不知道写什么题目了，随便打个题目凑凑数，这一段时间有很多事，又是女朋友离职去深证找工作，又是师兄和晓楠毕业，一连串的事情下来，自己都有些晕头转向了，一想到自己明年就要毕业了，到现在为止毕业论文什么的还没有着落，其实心里还是很着急。总是想着赶快做完了毕设好好找工作，可是到目前为止看了不少文献，可是总是感觉一团糟，不能将文章很好的串联起来，所以显得有些心烦意乱。 $\hspace{1cm}$也许每一个阶段都有着这一阶段的苦闷和烦恼吧，现在这一段时间让我想到了自己大三下那门通选课被挂了的时候，也跟现在差不多，总之就是很烦恼，总之就是很迷茫，感觉诸事不顺一样。其实相对来说自己还算是幸运的，毕竟在实验室一直学习也不用为生计发愁，目前文章也已经被录用，只等发表了。其实自己应该很感恩了，毕竟一切都是这么顺利。我还能奢望一些什么呢。老妈去韩国进修，总是时不时会晒一些照片，感觉过得挺舒服的，其实有信仰也不是一件坏事，毕竟人总是要有一些敬畏和感恩才好吧，如果把一切都归功于自己，从不感恩或者没有敬畏那才是可怕的吧。有时候觉得自己很奇怪，一方面觉得应该有敬畏和感恩，可是我却很怀疑上帝的存在，其实信仰这个东西是关乎心，在于人的吧。虽然我不相信上帝能主宰我的一切，可是我确实应该感谢上帝，应该学会感恩与敬畏。

$\hspace{1cm}$闲话说了不少，感觉还是槽点慢慢，就如同心中有一团火在烧，，让自己很难受，难道这也是一个周期性的症状么？在GitHub上一直在维护着的项目也停了下来，不知道应该怎样去维护了，主要是没有心思再去维护，不过还是应该再努努力吧。毕竟维护了这么久也不想就此荒废掉。刘艳去了深圳，或者有可能去广州，其实不管是广州还是深圳，总之在那边估计最好还是做开发的工作了，所以现在还是想学习一下怎么做web端的开发，现在也有一个项目在维护着，还是坚持做下去吧，毕竟多学一种技能总是不会有坏处的。
$\hspace{1cm}$我觉得自己会有一个毛病，那就是喜欢把所有的事情都放在心理，然后不停的去思考，去发酵，所以显得很沉重，压力很大，看来自己还是心不够宽，体不够胖呀，其实想想自己也不应该有这么大的压力，家里在市区有一栋房子，这就已经让我不应该为什么事情担忧了，不管在外面混得怎么样，总是有家可以回的不是么。还是应该放宽心，不要想的太多，只要能够好好的过自己的一生不就好了么！其实自己也有不少的愿望清单，相信也会一点点慢慢的实现吧。

引言

大部分纹理合成的方法可以通过使用何种纹理类型模型来进行分类，包括反应扩散模型，频率域方法，分形的理论与方法或者是统计和随机场方法。有一些使用混合模型的方法，包括确定性组份方法和随机组份方法。除了上述纹理之外还有扫描图像和手绘图像的纹理特征。本文着重于对随机纹理进行处理，目前对于纹理区分的理论主要基于以下假设：

当两种纹理在在频率和空间域的的频率分布是相同的，则这两种纹理难以区分出来，这个理论是建立在人类纹理感知的基础之上的。

影像金字塔

线性转换代表影像基函数的加权和，因此影像可以表示为： $$I(x,y)=\Sigma_iy_ig_i(x,y)\,,\hspace{1cm}(1)$$ 其中$I(x,y)$表示影像，$g_i(x,y)$为基函数，$y_i$为转换系数，转换系数是通过将信号投影到一系列投影函数中： $$y_i=\Sigma_{x,y}h_i(x,y)I(x,y)\,.\hspace{1cm}(2)$$ 比如傅立叶转换的基函数为空间频率变化的正弦函数，余弦函数。
影像金字塔也是一个典型的子带转换，影像金字塔的特征是基函数（投影函数）可以通过彼此的基函数的平移和膨胀得到（和膨胀的转换因子为$2^j$其中$j$为整数）。影像金字塔子带是通过卷积和降采样得到的，对于一系列的$j$降采样因子以$2^n$增长。

拉普拉斯金字塔

拉普拉斯金字塔的构造主要有两个步骤：1.缩减；2.扩张；缩减操作使用一个因子为2的降采样低通滤波器在各个唯独上进行处理，扩张操作使用一个因子为2的升采样因子，像素点之间的值用0填充。最常见的低通滤波核为：$\frac{1}{16}(1,4,6,4,1)\,.$一个完整的金字塔影像层包含两景影像：$l^0(低通滤波影像)$和$b^0$高通滤波影像，两景影像通过如下公式计算： $$l^0=Reduce(im)\\b^0=im-Expand(10)\,,$$ 其中$im$为原始的输入影像，原始影像可以通过$b^0$和$l^0$恢复：
$$recontructed-im=b^0-Expand(l^0)\,.\hspace{1cm}(3)$$ 下一层的金字塔也通过对$l^0$进行同样的操作获取，得到$l^{i+1}$和$b^{i+1}$

方向金字塔

有方向或者是长条型结构的纹理通过拉普拉斯金字塔不能够被探测到，因为拉普拉斯金字塔的基函数是径向对称的，为了合成各向异性的纹理，可以采用方向金字塔转换。方向金字塔和拉普拉斯金字塔一样，都是将影像分解为几个空间频带，另外对于每个空间频带，进一步划分为一系列的方向频带。

上图说明了使用方向金字塔进行纹理分析与合成的过程，左边的范式为分析部分，主要是构造方向金字塔的过程，右边是合成部分，通过纹理合成进行金字塔合成，在两个过程中的圆圈表示影像子带分解。方向金字塔首先通过一个径向对称频率响应的，高通滤波对应于四个方向。每一个连续层的金字塔都由上一层的低通滤波组成，$B_i$表示的是各个方向的带通滤波器。注意到金字塔影像像素远大于原始影像像素，但是只有低通滤波结果被降采样了，这样做能够在每个子带最小化混淆现象，避免混淆现象很重要，因为基于金字塔的纹理合成方法对每一个子带处理是独立的。

处理过程

直方图匹配

直方图匹配是对两景影像直方图进行匹配获取：（1）一幅影像的累积分布函数；（2）另一幅图像的逆累积分布函数。直方图匹配的处理过程为：

Match-histogram (im1,im2)
    im1-cdf = Make-cdf(im1)
    im2-cdf = Make-cdf(im2)
    inv-im2-cdf = Make-inverse-lookup-table(im2-cdf)
    Loop for each pixel do
        im1[pixel] =
            Lookup(inv-im2-cdf,
            Lookup(im1-cdf,im1[pixel]))

纹理匹配

纹理匹配的主要过程是对输入的噪声影像进行修改，使其与输入的纹理图像相似，1.匹配噪声图像的直方图和输入纹理的直方图；2.从修改后的噪声和纹理影像中构造金字塔；3.在两个金字塔数据循环，获取最佳的直方图匹配的金字塔子带；4.进行直方图匹配，生成初步的合成纹理。为了使得像素和金字塔直方图匹配都能够最好的匹配，进行迭代重新匹配影像直方图和子带金字塔直方图。

Match-texture(noise,texture)
    Match-Histogram (noise,texture)
    analysis-pyr = Make-Pyramid (texture)
    Loop for several iterations do
        synthesis-pyr = Make-Pyramid (noise)
    Loop for a-band in subbands of analysis-pyr
         for s-band in subbands of synthesis-pyr
         do
         Match-Histogram (s-band,a-band)
    noise = Collapse-Pyramid (synthesis-pyr)
    Match-Histogram (noise,texture)

目前文章1并没有对其收敛性进行证明，但是文章给出经验在处理过程中总是能够收敛，由于滤波器并不是完美的，迭代次数过多会引入重构误差。
除了以上处理以外还包括辩解处理，颜色处理，纹理混合等操作。

结论

局限性

文章2提出的算法还是有着其局限性：1.输入的纹理应该是一个均匀的纹理，如果纹理不均匀，则合成的结果会出现污斑现象；2.算法对于准周期性纹理和随机镶嵌的纹理并不适用；3.对于局部存在方向性，但是不同影像部分方向不同的纹理并不适用；

结论

文章提出了一种纹理合成的方法，其主要的有点在于本文提出的方法十分简单，在处理过程中并不需要理解复杂的纹理结构。

引言

$\hspace{1cm}$原则上，使用纹理合成的方法能够对数字影像进行修复，特别是在一些纹理结构比较强的区域，通过纹理合成的方法进行修复可疑获取比较好的影像修复效果，然而如果待修复的区域包含额外的颜色或者是强度梯度信息，则采用纹理合成的方法就不能达到理想的效果了。影像修复的方法可以作为纹理合成方法的补充，影像修复方法一般是通过待修复区域周围像素信息传递对待修复区域进行修复的，与基于纹理合成的方法不同，影像修复的方法能够正确处理影像强度信息，但是对于纹理细节并不能够进行很好的恢复。

修复方法

整个修复方法的处理流程为：

Fig2
$\hspace{1cm}$从Fig2的处理流程我们可以看出整个处理过程为：首先对影像进行分解得到高频和低频的分量，对高频和低频的分量分别进行处理。低频分量主要包含影像的强度信息，包含较少的细节信息，因此对低频的修复着重修复影像强度，而基于PDE1的影像修复方法对于纹理细节比较少的区域能够较好的保留其强度信息，而对于高频分量，将其通过高斯金字塔分解为n+1层，假设每一层为$H_i\,(i=0,...,n)\,.$从第$H_n$层开始采用多分辨率纹理合成方法2对每一层的待修复区域进行修复。将修复后的低频成分和高频成分结合得到修正后的影像。

频率分解

$\hspace{1cm}$处理过程中第一步就是要获取影像$I$分解后得到高频和低频分量。首先通过离散DCT3变换将原始影像变换为一系列的子带，然后我们选取前$k$个子带（前$k$个子带包含了主要信息），然后通过前$k$个子带进行DCT逆变换，通过逆变换后就可以得到低频影像$L_k$。而相应的高频影像则通过原始影像减去低频影像得到：$H_k:=I-L_k\,.$在这个过程中低频组份选取的$k$值确定了低频组份的模糊程度，$k$值越大则低频组份中的细节信息越丰富，$k$值越小则说明低频组份越模糊，实际上$k$值的选取对影像修复效果有比较大的影响。在校正过程中我们希望能够尽量大的去除低频成分中的细节信息只保留强度信息。而我们面临的问题在于如何确定$k$值的问题，对于这个问题文章给出了一个解决方案4。首先获取高频成分后计算高频成分$H_k$的自相关矩阵5$A_k$：$A_k:=DCT^{-1}(DCT(H_k)\cdot DCT(H_k))\,.$如果输入影像不是一个矩形，则非矩形部分通过0值补全，求得$A_k$之后剩下的工作为求取$A_k$的标准差，取$k$值使得$A_k$的标准差小于0.001，则最小的$k$值为选取的值，文章认为此时能够取得比较好的校正效果。

高斯影像金字塔

$\hspace{1cm}$通过上述过程分别得到高频和低频分量之后剩下的步骤为对高频分量进行高斯金字塔分解，实际上高斯金字塔分解是一个很简单的过程：
$$\begin{eqnarray*} && \omega(u,v)=\widehat{\omega}(u)\,\widehat{\omega}(v)\\ && \widehat{\omega}(0)=0.4,\widehat{\omega}(\pm1)=0.25,\widehat{\omega}(\pm1)=0.05 \end{eqnarray*}$$
高斯金字塔6必须要求影像大小满足$(p2^N+1)\times(q2^N+1)$像素个数且（$p\,,q\,,N\in\mathbb N$）在这样大小的影像下可以保证影像高斯金字塔具有$N+1$层。在此情况下我们将影像分解为$n$层且$n\ll N+1\,.$整个过程可以杯描述为：在纹理合成的第$i$层像素信息可以从更高层$(i+1,...,n)$的核中获取。在实际处理过程中进行三层金字塔分解就能获取比较好的效果。

纹理合成

纹理合成的难点有三点7:1）如何定义纹理单元和纹理内容；2）如何重构纹理概率分布；3）如何在二维合成过程中进行线性化。单像素合成：假设$I$是一个通过纹理样本影像合成的影像$I_{smp}\subset I_{real}$,其中$I_{real}$是真实的无限的纹理。假设$p\in I$是一个像素，$w(p)\subset I$是一个以$p$为中心的一个模版，定义$d(w_1,w_2)$为两个块之间的距离测度，假设除$p$以外$w(p)$中其他像素值都是已知的，在基于MRF的纹理合成模型为：
$$\Omega(p)= { w^`\in I_real:d(w^`,w(p))=0}$$
通过模板$w(p)$在整个影像上做匹配，找到最佳匹配点，采用距离测度则为找到一个匹配区域使得距离最小，而距离测度描述为：$d=d_{SSD}*G$，其中$G$为二维高斯核。找到最佳匹配点之后可以通过纹理合成算法对点$p$进行纹理修复。

处理结果

$\hspace{1cm}$从上图结果可以看出文章算法具有比较好的校正效果，整个校正过程有两个参数起重要作用：1.DCT变换后子带个数的选取，子带个数的选取决定了低频组份$L_k$；2.对高频成分进行高斯金字塔分解。实际上通过标准差小于0.001的范围选取的$k$值能够获取比较好的处理结果。这样对于子带数目$k$的选取就可以自动化进行。对高频组份的金字塔分解选取合适的分解层数是一个比较困难的问题，通常在处理中使用两层或三层金字塔。

$\hspace{1cm}$Bandelet处理的过程：总的来说Bandelet处理过程是一个自底向上的处理过程；
其主过程为：

1.首先对影像进行小波变换，得到小波变换后的影像，小波变换可以通过小波行变换进行一维处理；
2.完成小波变换之后进行四叉树分解，在分解的过程中通过四叉树分成最小的正方形格网，对于每个格网计算其最佳几何流方向；
3.进行自底向上的过程的最佳方向的处理，这一个处理主要是对步骤2中根据四叉树分割结果进行合并，合并的原则为，分别计算合并之后影像块的拉格朗日系数和四个小块拉格朗日系数的和，比较合并后的拉格朗日系数和四个小块拉格朗日系数之和的大小，如果合并之后的拉格朗日系数小于四个小块的拉格朗日系数之和，则说明可以进行合并，否则不进行合并，如果进行合并，则层级系数+1，重复此过程直到无法进行合并。

主过程中步骤2的处理步骤为：

1.对每一小块影像进行角度采样，设每一景影像的大小为：$L*L$，则将圆周角$[0,\pi)$等角度离散化为$L^2-1$个，则角度的可能取值为： $$\theta=\frac{k\pi}{L^2-1}\,,k=0,1,...,L^2-2$$
2.构造一个$L*L$大小的网格，即构造和$L*L$大小的数组，对于每一个角度$\theta$计算在网格点的正交投影误差 $$t=-sin(\theta)*i+cos(\theta)*j$$ 然后对$t$按照从小到大进行排序，然后对排序后的数据进行一维小波变换；
3.对于小波变换后的数据进行量化，量化的主要方法为： $$Q(x)=\left\{\begin{array}{lr}0,\,|x|\le T\\sgn(x)(q+1/2)T,qT\leq|x|<(q+1)T\end{array}\right.$$ 量化后得到量化后的小波变换系数；
4.首先获取几何流比特数$R_g=\lceil log_2(length(\theta)) \rceil$，大于阈值的系数编码比特数$R_b$，$R_b$等于大于阈值的小波系数个数乘以7，最后计算拉格朗日系数如下： $$L(f_{\theta},R)=||f_{\theta}-\tilde{f_{\theta}})||+\lambda T^2(R_g+R_b)$$ 其中$\lambda$为拉格朗日乘子，取$\lambda=3/28$，$\tilde{f_{\theta}}$为量化之后的小波信号，$f_{\theta}$为量化之前的小波信号，这样得到拉格朗日系数；
5.获取拉格朗日系数之后找到最小的拉格朗日系数，然后将最小的拉格朗日系数对应的角度作为子块对应的角度，然后沿着此角度方向，通过步骤2的步骤进行排序并计算小波变换系数。

总的来说Bandelet变换主要是通过定义一个最佳几何流方向，然后沿着最佳几何流方向进行小波变换，在边缘或几何流变化剧烈的影像块能够有比较细的分割，而对于几何流变化不那么剧烈的影像块进行合并。这样保证了最大限度的保留影像的纹理特征和影像压缩的目的。

DCT——离散余弦变换

DCT变换数学基础

1D序列的离散余弦变换为：
$$C(u)=\alpha(u)\Sigma_{x=0}^{N-1}f(x)cos[\frac{\pi(2x+1)u}{2N}]\hspace{1cm} (1)$$
其中$u=0,1,2,...,N-1\,.$ $f(x)$为原始数据值，$C(u)$为转换后频率值，则对于(1)式存在逆变换
$$f(x)=\Sigma_{u=0}^{N-1}\alpha(u)C(u)cos[\frac{\pi(2x+1)u}{2N}]\hspace{1cm} (2)$$
其中$x=0,1,2,...,N-1\,.$在式(1)和式(2)中$\alpha(u)$定义为：
$$\begin{equation} \alpha(u)= \left\{ \begin{array}{lr} &\sqrt{\frac{1}{N}} \\ &\sqrt{\frac{2}{N}} \end{array} \right. \hspace{1cm}(3) \end{equation}$$
2D离散余弦变换公式为：
$$C(u,v)=\alpha(u)\alpha(v)\Sigma_{x=0}^{N-1}\Sigma_{y=0}^{N-1}f(x,y)cos[\frac{\pi(2x+1)u}{2N}]cos[\frac{\pi(2y+1)v}{2N}]\hspace{1cm}(4)$$
对于公式(4)有$u,v=0,1,3,...,N-1$并且$\alpha(u)$和$\alpha(v)$都在公式(3)中定义了。则2D离散的DCT变换的逆变换为：
$$f(x,y)=\Sigma_{u=0}^{N-1}\Sigma_{v=0}^{N-1}\alpha(u)\alpha(v)C(u,v)cos[\frac{\pi(2x+1)u}{2N}]cos[\frac{\pi(2y+1)v}{2N}]\hspace{1cm}(5)$$

3D余弦变换：
3D余弦变换的主要是对序列影像进行的变换，3D
余弦变换的公式同2D，1D余弦变换相同，首先对各个影像单独进行余弦变换，然后对于每一个影像像元的序列组成的一维数据进行余弦变换，从而得到3D余弦变换的结果。

DCT变换的属性

1.去相关性：主要优势为移除邻域像素之间的冗余。
2.数据能量压缩：
3.可分离性：
公式(4)可以被描述为：
$$C(u,v)=\alpha(u)\alpha(v)\Sigma_{x=0}^{N-1}cos[\frac{\pi(2x+1)u}{2N}]\Sigma_{y=0}^{N-1}f(x,y)cos[\frac{\pi(2x+1)v}{2N}]\hspace{1cm}(6)$$
可分离性的主要优势在于能够通过1D的处理分别进行行列处理。
4.对称性：对称性表现在进行行转换和进行列转换都是独立的，一个分离的对称变换可以被描述为：
$$T=AfA\hspace{1cm}(7)$$
其中$A$为对称转换矩阵，且其中的元素$a(i,j)$的值为：
$$a(i,j)=\alpha(j)\Sigma_{j=0}^{N-1}cos[\frac{\pi(2j+1)i}{2N}]$$
这样，在计算过程中可以事先计算出矩阵$A$然后通过变换矩阵进行处理
5.正交性：正交性指的是DCT变换的基矩阵是正交的，整个逆转换的过程通过公式(6)可得：
$$f=A^-1TA^-1$$
由于$A$是正交矩阵，因此其转置和逆是相同的，在逆变换的过程中可以直接通过求转置快速处理。
在了解余弦变换的过程中顺便完成了余弦变换的代码

FMM影像修复方法

数学描述

$\hspace{1cm}$以上为FMM方法的示意图，图中$p$为待修复的点，$\delta\Omega$为待修复区域的边界$\Omega$为待修复的区域，对于一个待修复点$p$周围的的邻域小区域$B_\epsilon(p)$，区域大小为$\epsilon$，如Figure2图（a）所示，则未知点$p$的像素值应该由已知的区域$B_\epsilon(p)$得到，首先考虑灰度影像，彩色影像是灰度影像的扩展。如果区域$\epsilon$足够小的话，我们认为使用一个一阶导的近似可以进行修复，进行修复的公式为：
$$I_{q}(p)=I(q)+\Delta I(q)(p-q)\,.(1)$$
$\hspace{1cm}$上述公式表明在位置$p$的校正结果与$q$点位置和$q$点位置和$p$点位置距离差值的一阶近似，也就是$q$点的1阶近似（在区域比较小的时候比较适用）。接下来我们修复$p$点实用区域$B_\epsilon(p)$内所有像素进行处理，对每一个像素都有一个权值，则$p$的修复后的像素值为：
$$I(p)=\frac{\Sigma_{q\in{B_\epsilon(p)}}\omega(p,q)[I(q)+\Delta I()q(p-q)]}{\Sigma_{q\in{B_\epsilon(p)}}\omega(p,q)}\,.(2)$$
其中$\omega(p,q)$为权函数，权函数的设计既要考虑到$p$点信息的传递也要考虑到在$B_\epsilon(p)$区域影像的细节。

FMM影像修复

$\hspace{1cm}$上一章节描述了影像修复的基本原理和数学基础，则对整个区域的$\Omega$的修复方法为，对处于$\delta\Omega$边界上的点根据公式（2）进行修复，然后进行迭代，直到整个区域都修复完成，在这里有一个修复顺序的问题，在处理过程中使用待修复边界点与原始边界点的距离进行判断，距离原始边界越越近的点首先进行修复，这是因为距离边界越近的点修复的可靠性越大。对于上述操作可以采用FMM方法进行处理，简单的说FMM方法就是求解 $Eikonal$方程1

$$|\Delta T|=1\hspace{1cm}on \hspace{3pt}\Omega,\hspace{1cm}with\,T=0\,on\,\delta\Omega\,.(3)$$
其中方程（3）中对$T$的求解就是$\Omega$中像素距离初始边界$\delta\Omega$的距离，整个求解过程的伪代码为：

$\delta\Omega_i\,$=boundary of region to inpaint
$\delta\Omega\,$= $\,\delta\Omega_i$
while($\delta\Omega\,$not empty)
{
$p\,$=$\,$pixel of $\delta\Omega\,$closest to $\delta\Omega_i\,$
inpaint $\,p$ using Eqn (2)
advance $\delta\Omega\,$into $\Omega$
}

FMM方法的主要优势在于此方法明确平衡了区分已知和未知区域窄带，确定了下一个需要进行修复的像素，这里的像素窄带表示的就是我们的边界$\delta\Omega$，则对于每一个像素除了储存像素的灰度值之外还储存一个标志变量，这个标志变量包含三个值

• BAND：像素值属于窄带，其值经过更新。
• KNOWN：在边界$\delta\Omega$之外的点，也就是在已知区域的点，其T和灰度值I是已知的。
• INSIDE：在边界的像素值点，在边界区域中，T和I的值目前是未知的。

则FMM算法的初始值和传播相位如下：首先设置T在边界$\delta\Omega$上或边界外值为0，在待修复边界内设置极大值，对所有的影像都按上述条件设置初始值

while (NarrowBand not empty)
{
extract P(i,j) = head(NarrowBand); /* STEP 1 */
f(i,j) = KNOWN;
for (k,l) in (i1,j),(i,j1),(i+1,j),(i,j+1)
if (f(k,l)!=KNOWN)
{
if (f(k,l)==INSIDE)
{
f(k,l)=BAND; /* STEP 2 */
inpaint(k,l); /* STEP 3 */
}
T (k,l) = min(solve(k1,l,k,l1), /* STEP 4 */
solve(k+1,l,k,l1),
solve(k1,l,k,l+1),
solve(k+1,l,k,l+1));
insert(k,l) in NarrowBand; /* STEP 5 */
}
float solve(int i1,int j1,int i2,int j2)
{
float sol = 1.0e6;
if (f(i1,j1)==KNOWN)
if (f(i2,j2)==KNOWN)
{
float r = sqrt(2(T(i1,j1)T(i2,j2))*(T(i1,j1)T(i2,j2)));
float s = (T(i1,j1)+T(i2,j2)r)/2;
if (s>=T(i1,j1) && s>=T(i2,j2)) sol = s;
else
{ s += r; if (s>=T(i1,j1) && s>=T(i2,j2)) sol = s; }
}
else sol = 1+T(i1,j1));
else if (f(i2,j2)==KNOWN) sol = 1+T(i1,j2));
return sol;
}

对于单个像素的修复伪代码如下：

void inpaint(int i,int j)
{
for (all (k,l) in Bε(i,j) such that f(k,l)!=OUTSIDE)
{
r = vector from (i,j) to (k,l);
dir = r * gradT(i,j)/length(r);
dst = 1/(length(r)*length(r));
lev = 1/(1+fabs(T(k,l)T(i,j)));
w = dir*dst*lev;
if (f(k+1,l)!=OUTSIDE && f(k1,l)!=OUTSIDE &&
f(k,l+1)!=OUTSIDE && f(k,l1)!=OUTSIDE)
gradI = (I(k+1,l)I(k1,l),I(k,l+1)I(k,l1));
Ia += w * (I(k,l) + gradI * r);
s += w;
}
I(i,j) = Ia/s;
}

FMM算法细节

首先使用FMM方法计算初始待修复边界外的区域，获取待修复边界外部距离。由于我们只需要距离边界点小于$\epsilon$的点，因此只在边界外小于$\epsilon$的区域内计算T，这样我们就提高了计算速度，另外我们对边界内的点运行FMM进行处理，因此，整个影像上的T可以表现为：
$$\begin{equation} T(p)= \left\{ \begin{array}{lr} &T_{in}(p)\, & if \,p\in\Omega \\ -&T_{out(p)}\,&if\,p\not\in\Omega\, \end{array} \right. \,.(5) \end{equation}$$
然后通过$3\times 3$滤波算子进行滤波，然后计算$\Delta T$通过中心差分。

以上为文修复的效果，从上图可以看出，与BSCB方法相比，文章提出的方法在细节上表现并没有比BSCB方法好，但是文章的方法处理效率要远高于BSCB方法2。

局限性

对于以上采用微分校正的方法只能校正比较小的，或者线性的缺失区域，缺失范围较大的区域，校正方法会造成过度平滑，并不能很好进行处理，因此对于比较大的缺失需要进一步处理，但是采用微分校正的方法可疑获取校正趋势3。

结合纹理合成的多分辨率影像修复方法

引言

$\hspace{1cm}$原则上，使用纹理合成的方法能够对数字影像进行修复，特别是在一些纹理结构比较强的区域，通过纹理合成的方法进行修复可疑获取比较好的影像修复效果，然而如果待修复的区域包含额外的颜色或者是强度梯度信息，则采用纹理合成的方法就不能达到理想的效果了。影像修复的方法可以作为纹理合成方法的补充，影像修复方法一般是通过待修复区域周围像素信息传递对待修复区域进行修复的，与基于纹理合成的方法不同，影像修复的方法能够正确处理影像强度信息，但是对于纹理细节并不能够进行很好的恢复。

修复方法

整个修复方法的处理流程为：

fig 1

$\hspace{1cm}$从Fig2的处理流程我们可以看出整个处理过程为：首先对影像进行分解得到高频和低频的分量，对高频和低频的分量分别进行处理。低频分量主要包含影像的强度信息，包含较少的细节信息，因此对低频的修复着重修复影像强度，而基于PDE1的影像修复方法对于纹理细节比较少的区域能够较好的保留其强度信息，而对于高频分量，将其通过高斯金字塔分解为n+1层，假设每一层为$H_i\,(i=0,...,n)\,.$从第$H_n$层开始采用多分辨率纹理合成方法2对每一层的待修复区域进行修复。将修复后的低频成分和高频成分结合得到修正后的影像。

频率分解

$\hspace{1cm}$处理过程中第一步就是要获取影像$I$分解后得到高频和低频分量。首先通过离散DCT3变换将原始影像变换为一系列的子带，然后我们选取前$k$个子带（前$k$个子带包含了主要信息），然后通过前$k$个子带进行DCT逆变换，通过逆变换后就可以得到低频影像$L_k$。而相应的高频影像则通过原始影像减去低频影像得到：$H_k:=I-L_k\,.$在这个过程中低频组份选取的$k$值确定了低频组份的模糊程度，$k$值越大则低频组份中的细节信息越丰富，$k$值越小则说明低频组份越模糊，实际上$k$值的选取对影像修复效果有比较大的影响。在校正过程中我们希望能够尽量大的去除低频成分中的细节信息只保留强度信息。而我们面临的问题在于如何确定$k$值的问题，对于这个问题文章给出了一个解决方案4。首先获取高频成分后计算高频成分$H_k$的自相关矩阵5$A_k$：$A_k:=DCT^{-1}(DCT(H_k)\cdot DCT(H_k))\,.$如果输入影像不是一个矩形，则非矩形部分通过0值补全，求得$A_k$之后剩下的工作为求取$A_k$的标准差，取$k$值使得$A_k$的标准差小于0.001，则最小的$k$值为选取的值，文章认为此时能够取得比较好的校正效果。

高斯影像金字塔

$\hspace{1cm}$通过上述过程分别得到高频和低频分量之后剩下的步骤为对高频分量进行高斯金字塔分解，实际上高斯金字塔分解是一个很简单的过程：
$$\begin{eqnarray*} && \omega(u,v)=\widehat{\omega}(u)\,\widehat{\omega}(v)\\ && \widehat{\omega}(0)=0.4,\widehat{\omega}(\pm1)=0.25,\widehat{\omega}(\pm1)=0.05 \end{eqnarray*}$$
高斯金字塔6必须要求影像大小满足$(p2^N+1)\times(q2^N+1)$像素个数且（$p\,,q\,,N\in\mathbb N$）在这样大小的影像下可以保证影像高斯金字塔具有$N+1$层。在此情况下我们将影像分解为$n$层且$n\ll N+1\,.$整个过程可以杯描述为：在纹理合成的第$i$层像素信息可以从更高层$(i+1,...,n)$的核中获取。在实际处理过程中进行三层金字塔分解就能获取比较好的效果。

纹理合成

纹理合成的难点有三点7：1）如何定义纹理单元和纹理内容；2）如何重构纹理概率分布；3）如何在二维合成过程中进行线性化。单像素合成：假设$I$是一个通过纹理样本影像合成的影像$I_{smp}\subset I_{real}$,其中$I_{real}$是真实的无限的纹理。假设$p\in I$是一个像素，$w(p)\subset I$是一个以$p$为中心的一个模版，定义$d(w_1,w_2)$为两个块之间的距离测度，假设除$p$以外$w(p)$中其他像素值都是已知的，在基于MRF的纹理合成模型为：
$$\Omega(p)= { w^`\in I_real : d(w^`,w(p))=0}$$
通过模板$w(p)$在整个影像上做匹配，找到最佳匹配点，采用距离测度则为找到一个匹配区域使得距离最小，而距离测度描述为：$d=d_{SSD}*G$，其中$G$为二维高斯核。找到最佳匹配点之后可以通过纹理合成算法对点$p$进行纹理修复。

处理结果

$\hspace{1cm}$从上图结果可以看出文章算法具有比较好的校正效果，整个校正过程有两个参数起重要作用：1.DCT变换后子带个数的选取，子带个数的选取决定了低频组份$L_k$；2.对高频成分进行高斯金字塔分解。实际上通过标准差小于0.001的范围选取的$k$值能够获取比较好的处理结果。这样对于子带数目$k$的选取就可以自动化进行。对高频组份的金字塔分解选取合适的分解层数是一个比较困难的问题，通常在处理中使用两层或三层金字塔。

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

$\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一直重复以上步骤直到求解结果稳定