CN114494306B - 边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法 - Google Patents

Abstract

本发明公开一种边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法，首先，将待处理的甲骨拓片图像进行对比度增强、形态学操作和连通域几何属性分析，提取出甲骨卜辞文字；其次，利用边缘梯度的自相关矩阵对文字边缘的局部曲率进行分析，把其划分为平滑区域和笔划拐点区域；然后，对于平滑区域的文字边缘轮廓，采用线性插值方法对其进行平滑重构，而对于拐点区域的文字边缘轮廓，则利用局部梯度协方差和最小均方误差法计算最优的线性插值权重，实现边缘定向插值重构。

Description

边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法

技术领域

本发明涉及数字图像处理和古文字信息处理的交叉领域，尤其是一种能保持甲骨卜辞字形和笔划的整体风格，修复质量高、处理速度快、鲁棒性强、字形轮廓平滑度和边缘特征保真度好，有效克服字迹断裂、突起状毛刺、边缘锯齿效应和边缘特征退化的边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法。

背景技术

甲骨卜辞文字是镌刻于兽骨或龟甲上的一种古文字，其内容主要以记载商代王族的占卜记录为主，反映了距今3600多年前的殷商时期历史。作为可考证的最早文字***，甲骨卜辞文字对中国乃至世界文明溯源均有重要的研究意义。然而，一方面，甲骨文字数量众多，字形繁杂，直接书写于龟甲和兽骨，其书写材料不够光滑，本身又存在齿缝、兆纹、盾纹、刻痕等纹理，加之深埋地下数千年受到严重腐蚀，字形往往模糊不清，并且常受到点状噪声、片状斑纹和固有纹理的强烈干扰；另一方面，甲骨文在挖掘出土时受到损坏，字形边缘轮廓往往出现断裂和残缺等情况。这就导致人工释读和计算机文字识别工作困难重重。

甲骨文字修复的根本目的是借助计算机图像处理技术，在甲骨拓片图像上自动定位甲骨卜辞文字的边缘轮廓位置，并对断裂、毛刺、污迹、粗细不均、锯齿效应明显的边缘轮廓进行复原，进而为刻辞文本的自动识别等环节提供较高质量的字形图像或字形图像集合。它对于保证后续的文字分割和字形复原质量、提高字形特征提取和自动识别精度具有基础且重要的意义，亦是计算机辅助甲骨文字形研究、碎片缀合、文意释读等研究工作的基础任务之一。

目前，典型的图像修复算法大体分为三类，即基于偏微分方程的图像修复方法、基于纹理合成的图像修复方法、基于稀疏表示的图像修复方法。首先，基于偏微分方程的图像修复算法的基本原理是利用信息扩散原理，将待修复区域的周边信息沿着等照度线方向将迭代到破损区域，从外向内扩散式地得到修复信息，使破损区域越来越小直至修复完成。例如，Shen等人提出了一种整体变分算法；Chan等人提出了一种曲率驱动扩散算法等。该类方法适用于处理划痕等小尺度破损图像，但是不适合修复较大受损区域图像，会造成其边缘模糊。其次，基于纹理合成的图像修复方法从视觉认知规律的相似性出发，在未受损区域中选取与待修复块相似的图像块或纹理结构，并利用纹理合成技术对待修复区域进行复原。例如，Hu等人提出一种用于唐卡图像的修复与预测方法；张东等人提出一种小波变换与纹理合成相结合的图像修复方法。该方法有利于保证图像纹理结构的相似性和连续性，不仅可用来修复受损区域较小的图像，也被广泛应用于破损区域较大的图像修复中，其不足之处在于要求待修复图像的纹理特征具备一定规律性。最后，基于稀疏表示的图像修复方法将信号投影到一组由原子组成的完备/过完备字典上，再利用不同原子的线性组合来表示原始图像信号的稀疏特征，从而选择少量原子就能获得原始信号的非线性近似，实现待修复图像的重构。例如，Elad等人提出了一种稀疏形态成分分解模型；Sahoo等人提出了一种基于局部图像块的稀疏近似模型；李清泉等人结合人机交互和稀疏表示提出了一种壁画修复方法。该类方法比较适合提取和修复结构繁杂的图像特征信息，如不规则的边缘、纹理与角度等，对平滑区域和结构简单区域的修复效果却不够理想，其应用范围有限。总体来看，基于偏微分方程的图像修复方法无法处理甲骨拓片图像中的齿缝、盾纹和片状斑纹；基于纹理合成的图像修复方法无法抑制甲骨拓片图像中的点状噪声，且由于甲骨文字分布的规律性不强，其修复效果不理想；基于稀疏表示的图像修复方法尚不能抵抗甲骨固有纹理的干扰，会在修复过程中引起更多伪边缘。同时，由于甲骨文字在凿刻过程中会受到坚硬的龟甲和兽骨纹理的影响，导致文字边缘大多有突起状毛刺和明显的锯齿效应，这是上述方法均不能修复的。因此，上述三类图像修复方法尚无法对甲骨卜辞文字进行有效处理，而且它们均需要多次迭代，其计算效率较低，时间复杂度高。

针对古文字字形的修复，张显俊等人提出了一种基于傅里叶描述子的书法碑刻文字字形修复方法，采用8-邻域搜索法对填充后的字形图像进行轮廓提取，得到连通的字形外轮廓，再利用基于傅立叶描述子的特征提取和边缘平滑技术对字形轮廓进行重构，获得了平滑后的碑刻文字轮廓。马小虎等人对甲骨文的字形轮廓进行特征点提取，并使用贝塞尔曲线拟合甲骨文字的轮廓，能够把图像形式的甲骨文转化为曲线轮廓描述的甲骨文，但是无法实现残损字形轮廓的修复。考虑到甲骨文字边缘具有明显的锯齿效应，其轮廓震荡剧烈，顾绍通提出了一种基于自适应阈值和分形几何的甲骨拓片字形复原方法，通过分析甲骨拓片的噪声特点以及甲骨字形边缘的分形特征，计算甲骨拓片图像边缘的分形维数并对甲骨文字的锯齿状边缘进行压缩变换，实现了文字边缘的自适应平滑。杨亦鸣等人引进泊松分布来刻画甲骨拓片图像的连通区域面积，提出了一种基于泊松分布和分形几何的甲骨拓片字形复原方法，通过比较每个连通区域与拓片上所有连通区域的数学期望的差值来识别拓片上的噪声区域，再根据甲骨拓片字形边缘的分形维数和字形轮廓线上当前点与左右相邻点形成的向量夹角，对甲骨拓片的文字区域进行判别，从而达到平滑甲骨字形边缘的目的。为了更好地保持被平滑了的文字边缘的分形特征，顾绍通进一步构造了一种分形插值函数，对提取到的甲骨文字特征点进行分形插值与逼近，进而在保留甲骨字形整体风格的同时，对文字边缘进行平滑。不过，由于甲骨拓片中固有纹理和片状斑纹的分形维数与真正的甲骨文字往往存在一定相似性，上述两种方法仅对甲骨拓片图像上的点状噪声具有很好的抵抗能力，却对其固有纹理和片状斑纹的干扰不够鲁棒。而且，在不断逼近原始字形的过程中，分形插值操作对于迭代次数非常敏感，若迭代参数选取不当，反而会加剧字形退化甚至改变文字形状。于是，顾绍通进一步提出了一种基于数学形态学的甲骨拓片字形复原方法，分别沿着水平方向、竖直方向以及斜线方向，对字形边缘进行基于形态学膨胀和腐蚀操作，从而平滑甲骨文字边缘的凹凸形态。虽然该方法对突起状毛刺的复原表现出一定效果，但是形态学操作在平滑边缘的同时却忽略了原有文字轮廓和曲率的先验约束，不可避免地导致边缘特征出现退化，容易消除曲率较大的文字拐点。

综上所述，目前尚没有一种可同时抵抗龟甲和兽骨上存在的点状噪声、片状斑纹和固有纹理干扰的，修复质量高、处理速度快、鲁棒性强、字形轮廓平滑度和边缘特征保真度好，尤其是能有效克服字迹断裂、突起状毛刺、边缘锯齿效应、边缘特征退化的具备自适应能力的甲骨卜辞文字轮廓自动修复方法。

发明内容

本发明是为了解决现有技术所存在的上述技术问题，提供一种能保持甲骨卜辞字形和笔划的整体风格，修复质量高、处理速度快、鲁棒性强、字形轮廓平滑度和边缘特征保真度好，尤其是能够有效克服字迹断裂、突起状毛刺、边缘锯齿效应和边缘特征退化的边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法。

本发明的技术解决方案是：一种边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法，其特征在于按照如下步骤进行：

步骤1.输入待处理的甲骨卜辞拓片图像I_ori；

步骤2.利用形态学开运算估计I_ori的背景亮度，得到背景图像I_background；

步骤3.从I_ori中减去I_background，得到前景图像I_foreground；

步骤4.采用stretchlim方法为图像I_foreground计算出自适应灰度变换的最佳下界low和最佳上界high，进而将I_foreground进行对比度拉伸运算，使其亮度值的范围处于[low,high]区间，得到对比度增强后的甲骨拓片图像I_enhanced；

步骤5.利用最大类间方差法计算全局最佳分割阈值T_seg，并利用阈值T_seg对I_enhanced进行粗分割，进而将I_enhanced转换为二值图像I_bw，其中，前景文字区域的灰度为白色，背景区域的灰度为黑色；

步骤6.对二值图像I_bw中的每个8-连通区域进行标记，令8-连通区域的总数为N_con；

步骤7.测量二值图像I_bw中每个8-连通区域的属性；

步骤7.1计算I_bw中第i个8-连通区域的面积A_i、边缘坐标集合C_i、像素坐标集合P_i，所述i为正整数且1≤i≤N_con；

步骤7.2利用第i个8-连通区域的边缘坐标集合C_i和像素坐标集合P_i计算外接矩形框，得到其左上角坐标和右下角坐标/>并显示该连通区域的外接矩形框；

步骤8.计算所有8-连通区域的最大面积A_max、最小面积A_min、平均面积A_avg和面积分布的标准差A_σ，并绘制所有8-连通区域的面积分布直方图；

步骤9.判断二值图像I_bw的噪声区域并进行填充，得到去噪后的图像I′_bw，令i←1，并令N_text←0，所述N_text表示I′_bw中包含的文字区域数量；

步骤9.1对于图像I_bw的第i个8-连通区域，若其面积A_i小于A_avg-A_σ，则将当前的8-连通区域判定为点状噪声区域，并将其填充为背景；

步骤9.2对于图像I_bw的第i个8-连通区域，若其面积A_i大于A_avg+A_σ，则将当前的8-连通区域判定为片状斑纹噪声区域，并将其填充为背景；

步骤9.3对于图像I_bw的第i个8-连通区域，若其外接矩形框的高宽比或/>则将当前的8-连通区域判定为齿缝、兆纹、盾纹和刻痕等固有的非文字区域并将其填充为背景，否则将当前的8-连通区域判定为文字区域，并令N_text←N_text+1，所述ρ_max和ρ_min表示预设的常量；

步骤9.4令i←i+1，若i≤N_con，则返回步骤9.1，否则令j←1并转入步骤10；

步骤10.从图像I′_bw中取出第j个文字区域，对应第j个甲骨卜辞单字；

步骤11.从图像I′_bw中取出第j个文字区域的外接矩形框所包含的像素集合，进而建立一幅子图像I_j，令其高度为h_j像素，宽度为w_j像素；

步骤12.利用Canny算子对I_j进行边缘检测，得到边缘图像I_edge中；

步骤13.利用沿着水平方向的差分算子(-2,-1,0,1,2)计算I_edge的水平梯度I_x，并利用沿着竖直方向的差分算子(-2,-1,0,1,2)^T计算I_edge的竖直梯度I_y，所述上标T表示转置运算；

步骤14.根据公式(1)，计算水平梯度I_x的Hadamard乘积，得到

其中，“*”表示Hadamard乘积的运算符；

步骤15.根据公式(2)，计算竖直梯度I_y的Hadamard乘积，得到

步骤16.根据公式(3)，计算水平梯度I_x与竖直梯度I_y的Hadamard乘积，得到

步骤17.建立高斯平滑算子G，其尺寸为h_G×h_G，并利用该算子对 和/>进行卷积滤波，所述h_G表示预设的常量；

步骤18.计算子图像I_j中每个像素的特征值；

步骤18.1令c←1，l←1；

步骤18.2根据公式(4)，建立一个大小为2×2阶的自相关矩阵M，用于表示子图像I_j中位于第c行、第l列的像素的梯度自相关；

所述表示/>中位于第c行、第l列的元素值，/>表示/>中位于第c行、第l列的元素值，/>表示/>中位于第c行、第l列的元素值；

步骤18.3根据公式(5)，计算自相关矩阵M的特征值R(c,l)；

R(c,l)＝det(M)-κ·trace²(M) (5)

所述det(M)表示矩阵M所对应的行列式的值，trace(M)表示矩阵M的迹，κ表示预设的加权常数；

步骤18.4令l←l+1，若l≤w_j，则返回步骤18.2，否则，转入步骤18.5；

步骤18.5令c←c+1，l←1，若c≤h_j，则返回步骤18.2，否则，转入步骤18.6；

步骤18.6根据公式(6)，计算全部特征值的最大值R_max(c,l)，并令β←0.1×R_max(c,l)；

所述max{·}表示取最大值函数，β表示用于判断特征点的预设阈值；

步骤19.在子图像I_j中寻找特征角点；

步骤19.1令c←1，l←1，并初始化一幅边缘特征点图像I_corner，其高度为h_j像素，宽度为w_j像素；

步骤19.2若R(c,l)≥β，则将子图像I_j中位于第c行、第l列的像素判定为特征角点并令I_corner(c,l)←1，否则令I_corner(c,l)←0，所述I_corner(c,l)表示边缘特征点图像I_corner中位于第c行、第l列的像素值；

步骤19.3令l←l+1，若l≤w_j，则返回步骤19.2，否则，转入步骤19.4；

步骤19.4令c←c+1，l←1，若c≤h_j，则返回步骤19.2，否则，转入步骤20；

步骤20.初始化一幅图像其高度为2×h_j像素，宽度为2×w_j像素，用于存储修复后的第j个甲骨卜辞单字图像；

步骤21.利用线性插值方法对图像的全部像素值进行平滑重构，从而保持原始甲骨卜辞文字的主要字形及其轮廓信息；

步骤21.1根据公式(7)，对图像中偶数行、偶数列的全部像素点进行平滑插值重构；

所述表示图像/>中位于第2m行、第2n列的像素值，I_j(m,n)表示图像I_j中位于第m行、第n列的像素值，且m∈{1,2,3,…,h_j}，n∈{1,2,3,…,w_j}；

步骤21.2根据公式(8)，对图像中奇数行、偶数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m-1行、第2n列的像素值，表示图像/>中位于第2m-2行、第2n列的像素值，γ表示预设的线性加权系数；

步骤21.3根据公式(9)，对图像中偶数行、奇数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m行、第2n-1列的像素值，表示图像/>中位于第2m行、第2n-2列的像素值；

步骤21.4根据公式(10)，对图像中奇数行、奇数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m-1行、第2n-1列的像素值，表示图像/>中位于第2m-1行、第2n-2列的像素值；

步骤22.将边缘特征点图像I_corner视为拐点区域，进而根据I_j与之间的几何对偶一致性，利用I_j的梯度协方差系数计算出最优的线性插值权重，对/>进行梯度协方差引导的边缘定向插值；

步骤22.1令c←1，l←1；

步骤22.2若I_corner(c,l)＝1，则转入步骤22.3，否则转入步骤22.9；

步骤22.3在I_j中以第c行、第l列的像素为中心建立一个大小为W×W的窗口，并将该窗口中的全体像素值按照光栅扫描顺序保存到一个列向量y中，得到所述W为一个预设常数，/>表示I_j中位于第/>行、第/>列的像素值，I_j(c,l)表示I_j中位于第c行、第l列的像素值，/>表示I_j中位于第行、第/>列的像素值；

步骤22.4建立一个阶数为4×W²的矩阵C，对于列向量y的第K个元素I_j(c′,l′)，将I_j(c′-1,l′-1)、I_j(c′+1,l′-1)、I_j(c′-1,l′+1)、I_j(c′+1,l′+1)保存到矩阵C的第K列，所述1≤K≤W²，I_j(c′,l′)表示I_j中位于第c′行、第l′列的像素值，I_j(c′-1,l′-1)表示I_j中位于第c′-1行、第l′-1列的像素值，I_j(c′+1,l′-1)表示I_j中位于第c′+1行、第l′-1列的像素值，I_j(c′-1,l′+1)表示I_j中位于第c′-1行、第l′+1列的像素值，I_j(c′+1,l′+1)表示I_j中位于第c′+1行、第l′+1列的像素值；

步骤22.5根据公式(11)，计算基于最小均方误差的最优线性插值权重向量α；

α＝(CC^T)^-1×(Cy) (11)

步骤22.6根据公式(12)，对进行梯度协方差引导的边缘定向插值；

所述表示图像/>中位于第2c+1行、第2l+1列的像素值，表示图像/>中位于第2c行、第2l列的像素值，/>表示图像/>中位于第2c+2行、第2l列的像素值，/>表示图像/>中位于第2c+2行、第2l+2列的像素值，/>表示图像/>中位于第2c行、第2l+2列的像素值，α₀表示向量α的第1个分量，α₁表示向量α的第2个分量，α₂表示向量α的第3个分量，α₃表示向量α的第4个分量；

步骤22.7根据公式(13)，对进行梯度协方差引导的边缘定向插值；

所述表示图像/>中位于第2c+1行、第2l列的像素值，表示图像/>中位于第2c+1行、第2l-1列的像素值；

步骤22.8根据公式(14)，对进行梯度协方差引导的边缘定向插值；

所述表示图像/>中位于第2c行、第2l+1列的像素值，表示图像/>中位于第2c-1行、第2l+1列的像素值；

步骤22.9令l←l+1，若l≤w_j，则返回步骤22.2，否则，转入步骤22.10；

步骤22.10令c←c+1，l←1，若c≤h_j，则返回步骤22.2，否则，得到第j个甲骨卜辞单字图像的轮廓修复结果并转入步骤23；

步骤23.令j←j+1，若j≤N_text，则返回步骤10，否则，输出全部甲骨卜辞单字图像的轮廓修复结果。

本发明与现有技术相比，具有如下优点：第一，对甲骨拓片图像进行对比度增强、形态学操作和连通域几何属性分析等预处理，使甲骨卜辞文字与点状噪声、片状斑纹以及齿缝、兆纹、盾纹、刻痕等固有纹理相互分离，有利于后续步骤更准确地提取出卜辞文字的字迹特征；第二，利用边缘梯度的自相关矩阵及其特征值对甲骨文字边缘点的局部曲率进行分析，从而将文字边缘划分为平滑区域和笔划拐点区域，有利于检测出甲骨文字的突起状毛刺、边缘锯齿；第三，对于平滑区域的边缘轮廓，采用线性插值方法对其进行平滑重构，从而保持原始甲骨卜辞文字的主要字形及其轮廓信息。对于拐点区域的边缘轮廓，则根据低分辨率图像与高分辨率图像之间的几何对偶一致性，利用局部梯度协方差和最小均方误差法计算出最优的线性插值权重，实现了梯度协方差引导的边缘定向插值重构，进而在保持笔划拐点和边缘特征不变的情况下，以主要边缘的方向来抑制和修正突起状毛刺、边缘锯齿，并修复小尺度断裂的字迹；第四，无需手工处理甲骨拓片图像，可以批量自动处理甲骨文字图像，在时间效率上明显优于人工修复方法。因此，本发明公开的边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法为甲骨文研究者，特别是为基于计算机视觉的甲骨文字识别研究者提供了一种提高文字质量和识别率的有效途径，具有保持甲骨卜辞字形和笔划的整体风格，修复质量高、处理速度快、鲁棒性强、字形轮廓平滑度和边缘特征保真度好的特点，能够有效克服字迹断裂、突起状毛刺、边缘锯齿效应和边缘特征退化。

附图说明

图1是本发明对《甲骨文合集》00032号甲骨拓片图像的文字修复结果图。

图2是本发明对《甲骨文合集》00102号甲骨拓片图像的文字修复结果图。

图3是本发明对《甲骨文合集》00137号甲骨拓片图像的文字修复结果图。

具体实施方式

本发明提供一种边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法，按照如下步骤进行：

步骤1.输入待处理的甲骨卜辞拓片图像I_ori；

步骤2.利用形态学开运算估计I_ori的背景亮度，得到背景图像I_background；

步骤3.从I_ori中减去I_background，得到前景图像I_foreground；

步骤4.采用stretchlim方法为图像I_foreground计算出自适应灰度变换的最佳下界low和最佳上界high，进而将I_foreground进行对比度拉伸运算，使其亮度值的范围处于[low,high]区间，得到对比度增强后的甲骨拓片图像I_enhanced；

步骤5.利用最大类间方差法计算全局最佳分割阈值T_seg，并利用阈值T_seg对I_enhanced进行粗分割，进而将I_enhanced转换为二值图像I_bw，其中，前景文字区域的灰度为白色，背景区域的灰度为黑色；

步骤6.对二值图像I_bw中的每个8-连通区域进行标记，令8-连通区域的总数为N_con；

步骤7.测量二值图像I_bw中每个8-连通区域的属性；

步骤7.1计算I_bw中第i个8-连通区域的面积A_i、边缘坐标集合C_i、像素坐标集合P_i，所述i为正整数且1≤i≤N_con；

步骤7.2利用第i个8-连通区域的边缘坐标集合C_i和像素坐标集合P_i计算外接矩形框，得到其左上角坐标和右下角坐标/>并显示该连通区域的外接矩形框；

步骤8.计算所有8-连通区域的最大面积A_max、最小面积A_min、平均面积A_avg和面积分布的标准差A_σ，并绘制所有8-连通区域的面积分布直方图；

步骤9.判断二值图像I_bw的噪声区域并进行填充，得到去噪后的图像I′_bw，令i←1，并令N_text←0，所述N_text表示I′_bw中包含的文字区域数量；

步骤9.1对于图像I_bw的第i个8-连通区域，若其面积A_i小于A_avg-A_σ，则将当前的8-连通区域判定为点状噪声区域，并将其填充为背景；

步骤9.2对于图像I_bw的第i个8-连通区域，若其面积A_i大于A_avg+A_σ，则将当前的8-连通区域判定为片状斑纹噪声区域，并将其填充为背景；

步骤9.3对于图像I_bw的第i个8-连通区域，若其外接矩形框的高宽比或/>则将当前的8-连通区域判定为齿缝、兆纹、盾纹和刻痕等固有的非文字区域并将其填充为背景，否则将当前的8-连通区域判定为文字区域，并令N_text←N_text+1，所述ρ_max和ρ_min表示预设的常量，本实施例中，令ρ_max＝0.65，ρ_min＝0.25；

步骤9.4令i←i+1，若i≤N_con，则返回步骤9.1，否则令j←1并转入步骤10；

步骤10.从图像I′_bw中取出第j个文字区域，对应第j个甲骨卜辞单字；

步骤11.从图像I′_bw中取出第j个文字区域的外接矩形框所包含的像素集合，进而建立一幅子图像I_j，令其高度为h_j像素，宽度为w_j像素；

步骤12.利用Canny算子对I_j进行边缘检测，得到边缘图像I_edge中；

步骤13.利用沿着水平方向的差分算子(-2,-1,0,1,2)计算I_edge的水平梯度I_x，并利用沿着竖直方向的差分算子(-2,-1,0,1,2)^T计算I_edge的竖直梯度I_y，所述上标T表示转置运算；

步骤14.根据公式(1)，计算水平梯度I_x的Hadamard乘积，得到

其中，“*”表示Hadamard乘积的运算符；

步骤15.根据公式(2)，计算竖直梯度I_y的Hadamard乘积，得到

步骤16.根据公式(3)，计算水平梯度I_x与竖直梯度I_y的Hadamard乘积，得到

步骤17.建立高斯平滑算子G，其尺寸为h_G×h_G，并利用该算子对 和/>进行卷积滤波，所述h_G表示预设的常量，本实施例中，令h_G＝7；

步骤18.计算子图像I_j中每个像素的特征值；

步骤18.1令c←1，l←1；

步骤18.2根据公式(4)，建立一个大小为2×2阶的自相关矩阵M，用于表示子图像I_j中位于第c行、第l列的像素的梯度自相关；

所述表示/>中位于第c行、第l列的元素值，/>表示/>中位于第c行、第l列的元素值，/>表示/>中位于第c行、第l列的元素值；

步骤18.3根据公式(5)，计算自相关矩阵M的特征值R(c,l)；

R(c,l)＝det(M)-κ·trace²(M) (5)

所述det(M)表示矩阵M所对应的行列式的值，trace(M)表示矩阵M的迹，κ表示预设的加权常数，本实施例中，令κ＝0.05；

步骤18.4令l←l+1，若l≤w_j，则返回步骤18.2，否则，转入步骤18.5；

步骤18.5令c←c+1，l←1，若c≤h_j，则返回步骤18.2，否则，转入步骤18.6；

步骤18.6根据公式(6)，计算全部特征值的最大值R_max(c,l)，并令β←0.1×R_max(c,l)；

所述max{·}表示取最大值函数，β表示用于判断特征点的阈值；

步骤19.在子图像I_j中寻找特征角点；

步骤19.1令c←1，l←1，并初始化一幅边缘特征点图像I_corner，其高度为h_j像素，宽度为w_j像素；

步骤19.2若R(c,l)≥β，则将子图像I_j中位于第c行、第l列的像素判定为特征角点并令I_corner(c,l)←1，否则令I_corner(c,l)←0，所述I_corner(c,l)表示边缘特征点图像I_corner中位于第c行、第l列的像素值；

步骤19.3令l←l+1，若l≤w_j，则返回步骤19.2，否则，转入步骤19.4；

步骤19.4令c←c+1，l←1，若c≤h_j，则返回步骤19.2，否则，转入步骤20；

步骤20.初始化一幅图像其高度为2×h_j像素，宽度为2×w_j像素，用于存储修复后的第j个甲骨卜辞单字图像；

步骤21.利用线性插值方法对图像的全部像素值进行平滑重构，从而保持原始甲骨卜辞文字的主要字形及其轮廓信息；

步骤21.1根据公式(7)，对图像中偶数行、偶数列的全部像素点进行平滑插值重构；

所述表示图像/>中位于第2m行、第2n列的像素值，I_j(m,n)表示图像I_j中位于第m行、第n列的像素值，且m∈{1,2,3,…,h_j}，n∈{1,2,3,…,w_j}；

步骤21.2根据公式(8)，对图像中奇数行、偶数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m-1行、第2n列的像素值，表示图像/>中位于第2m-2行、第2n列的像素值，γ表示预设的线性加权系数，本实施例中，令γ＝0.5；

步骤21.3根据公式(9)，对图像中偶数行、奇数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m行、第2n-1列的像素值，表示图像/>中位于第2m行、第2n-2列的像素值；

步骤21.4根据公式(10)，对图像中奇数行、奇数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m-1行、第2n-1列的像素值，表示图像/>中位于第2m-1行、第2n-2列的像素值；

步骤22.将边缘特征点图像I_corner视为拐点区域，进而根据I_j与之间的几何对偶一致性，利用I_j的梯度协方差系数计算出最优的线性插值权重，对/>进行梯度协方差引导的边缘定向插值；

步骤22.1令c←1，l←1；

步骤22.2若I_corner(c,l)＝1，则转入步骤22.3，否则转入步骤22.9；

步骤22.3在I_j中以第c行、第l列的像素为中心建立一个大小为W×W的窗口，并将该窗口中的全体像素值按照光栅扫描顺序保存到一个列向量y中，得到所述W为一个预设常数，/>表示I_j中位于第/>行、第/>列的像素值，I_j(c,l)表示I_j中位于第c行、第l列的像素值，/>表示I_j中位于第行、第/>列的像素值，本实施例中，令W＝5；

步骤22.4建立一个阶数为4×W²的矩阵C，对于列向量y的第K个元素I_j(c′,l′)，将I_j(c′-1,l′-1)、I_j(c′+1,l′-1)、I_j(c′-1,l′+1)、I_j(c′+1,l′+1)保存到矩阵C的第K列，所述1≤K≤W²，I_j(c′,l′)表示I_j中位于第c′行、第l′列的像素值，I_j(c′-1,l′-1)表示I_j中位于第c′-1行、第l′-1列的像素值，I_j(c′+1,l′-1)表示I_j中位于第c′+1行、第l′-1列的像素值，I_j(c′-1,l′+1)表示I_j中位于第c′-1行、第l′+1列的像素值，I_j(c′+1,l′+1)表示I_j中位于第c′+1行、第l′+1列的像素值；

步骤22.5根据公式(11)，计算基于最小均方误差的最优线性插值权重向量α；

α＝(CC^T)^-1×(Cy) (11)

步骤22.6根据公式(12)，对进行梯度协方差引导的边缘定向插值；

所述表示图像/>中位于第2c+1行、第2l+1列的像素值，表示图像/>中位于第2c行、第2l列的像素值，/>表示图像/>中位于第2c+2行、第2l列的像素值，/>表示图像/>中位于第2c+2行、第2l+2列的像素值，/>表示图像/>中位于第2c行、第2l+2列的像素值，α₀表示向量α的第1个分量，α₁表示向量α的第2个分量，α₂表示向量α的第3个分量，α₃表示向量α的第4个分量；

步骤22.7根据公式(13)，对进行梯度协方差引导的边缘定向插值；

/>

所述表示图像/>中位于第2c+1行、第2l列的像素值，表示图像/>中位于第2c+1行、第2l-1列的像素值；

步骤22.8根据公式(14)，对进行梯度协方差引导的边缘定向插值；

所述表示图像/>中位于第2c行、第2l+1列的像素值，表示图像/>中位于第2c-1行、第2l+1列的像素值；

步骤22.9令l←l+1，若l≤w_j，则返回步骤22.2，否则，转入步骤22.10；

步骤22.10令c←c+1，l←1，若c≤h_j，则返回步骤22.2，否则，得到第j个甲骨卜辞单字图像的轮廓修复结果并转入步骤23；

步骤23.令j←j+1，若j≤N_text，则返回步骤10，否则，输出全部甲骨卜辞单字图像的轮廓修复结果。

从《甲骨文合集》中选取第00032号、00102号、00137号甲骨拓片图像进行实验，采用本发明进行甲骨卜辞文字修复的结果如图1-图3所示，其中，图1(a)、图2(a)和图3(a)为原始甲骨拓片图像，图1(b)、图2(b)和图3(b)为采用本发明将背景和噪声区域填充后的图像，图1(c)、图2(c)、图2(e)和图3(c)分别为图1(b)、图2(b)、图3(b)中用白色矩形框标出的文字的原始字迹，图1(d)、图2(d)、图2(f)和图3(d)分别为图1(b)、图2(b)、图3(b)中用白色矩形框标出的卜辞文字采用本发明修复后的边缘轮廓。

从图1(b)、图2(b)和图3(b)可见，本发明能够在保留卜辞文字区域的基础上，有效地去除龟甲和兽骨上的点状噪声、片状斑纹和固有纹理。而且，从图1(d)可见，本发明修复了甲骨卜辞文字左上侧发生断裂的相连笔划和右下侧孔洞的多余豁口；从图2(d)和图2(f)可见，本发明有效缓解了甲骨卜辞文字的边缘锯齿，并突出了图2(c)上方的圆形孔洞，避免其出现边缘特征退化而影响文字识别精度；从图3(d)可见，本发明消除了文字顶端的突起状毛刺以及倾斜笔划的边缘锯齿。因此，本发明完整地保持了甲骨卜辞文字的字形特点，抑制了字迹断裂、突起状毛刺、边缘锯齿效应和边缘特征退化等现象，实现了文字边缘轮廓的平滑处理，进而更加有利于提高后续的甲骨文字自动识别的精度。

Claims (1)

1.一种边缘梯度协方差引导的甲骨卜辞文字轮廓修复方法，其特征在于按照如下步骤进行：

步骤1.输入待处理的甲骨卜辞拓片图像I_ori；

步骤2.利用形态学开运算估计I_ori的背景亮度，得到背景图像I_background；

步骤3.从I_ori中减去I_background，得到前景图像I_foreground；

步骤4.采用stretchlim方法为图像I_foreground计算出自适应灰度变换的最佳下界low和最佳上界high，进而将I_foreground进行对比度拉伸运算，使其亮度值的范围处于[low,high]区间，得到对比度增强后的甲骨拓片图像I_enhanced；

步骤5.利用最大类间方差法计算全局最佳分割阈值T_seg，并利用阈值T_seg对I_enhanced进行粗分割，进而将I_enhanced转换为二值图像I_bw，其中，前景文字区域的灰度为白色，背景区域的灰度为黑色；

步骤6.对二值图像I_bw中的每个8-连通区域进行标记，令8-连通区域的总数为N_con；

步骤7.测量二值图像I_bw中每个8-连通区域的属性；

步骤7.1计算I_bw中第i个8-连通区域的面积A_i、边缘坐标集合C_i、像素坐标集合P_i，所述i为正整数且1≤i≤N_con；

步骤7.2利用第i个8-连通区域的边缘坐标集合C_i和像素坐标集合P_i计算外接矩形框，得到其左上角坐标和右下角坐标/>并显示该连通区域的外接矩形框；

步骤8.计算所有8-连通区域的最大面积A_max、最小面积A_min、平均面积A_avg和面积分布的标准差A_σ，并绘制所有8-连通区域的面积分布直方图；

步骤9.判断二值图像I_bw的噪声区域并进行填充，得到去噪后的图像I′_bw，令i←1，并令N_text←0，所述N_text表示I′_bw中包含的文字区域数量；

步骤9.1对于图像I_bw的第i个8-连通区域，若其面积A_i小于A_avg-A_σ，则将当前的8-连通区域判定为点状噪声区域，并将其填充为背景；

步骤9.2对于图像I_bw的第i个8-连通区域，若其面积A_i大于A_avg+A_σ，则将当前的8-连通区域判定为片状斑纹噪声区域，并将其填充为背景；

步骤9.3对于图像I_bw的第i个8-连通区域，若其外接矩形框的高宽比或/>则将当前的8-连通区域判定为齿缝、兆纹、盾纹和刻痕等固有的非文字区域并将其填充为背景，否则将当前的8-连通区域判定为文字区域，并令N_text←N_text+1，所述ρ_max和ρ_min表示预设的常量；

步骤9.4令i←i+1，若i≤N_con，则返回步骤9.1，否则令j←1并转入步骤10；

步骤10.从图像I′_bw中取出第j个文字区域，对应第j个甲骨卜辞单字；

步骤11.从图像I′_bw中取出第j个文字区域的外接矩形框所包含的像素集合，进而建立一幅子图像I_j，令其高度为h_j像素，宽度为w_j像素；

步骤12.利用Canny算子对I_j进行边缘检测，得到边缘图像I_edge中；

步骤13.利用沿着水平方向的差分算子(-2,-1,0,1,2)计算I_edge的水平梯度I_x，并利用沿着竖直方向的差分算子(-2,-1,0,1,2)^T计算I_edge的竖直梯度I_y，所述上标T表示转置运算；

步骤14.根据公式(1)，计算水平梯度I_x的Hadamard乘积，得到

其中，“*”表示Hadamard乘积的运算符；

步骤15.根据公式(2)，计算竖直梯度I_y的Hadamard乘积，得到

步骤16.根据公式(3)，计算水平梯度I_x与竖直梯度I_y的Hadamard乘积，得到

步骤17.建立高斯平滑算子G，其尺寸为h_G×h_G，并利用该算子对和/>进行卷积滤波，所述h_G表示预设的常量；

步骤18.计算子图像I_j中每个像素的特征值；

步骤18.1令c←1，l←1；

步骤18.2根据公式(4)，建立一个大小为2×2阶的自相关矩阵M，用于表示子图像I_j中位于第c行、第l列的像素的梯度自相关；

所述表示/>中位于第c行、第l列的元素值，/>表示/>中位于第c行、第l列的元素值，/>表示/>中位于第c行、第l列的元素值；

步骤18.3根据公式(5)，计算自相关矩阵M的特征值R(c,l)；

R(c,l)＝det(M)-κ·trace²(M) (5)

所述det(M)表示矩阵M所对应的行列式的值，trace(M)表示矩阵M的迹，κ表示预设的加权常数；

步骤18.4令l←l+1，若l≤w_j，则返回步骤18.2，否则，转入步骤18.5；

步骤18.5令c←c+1，l←1，若c≤h_j，则返回步骤18.2，否则，转入步骤18.6；

步骤18.6根据公式(6)，计算全部特征值的最大值R_max(c,l)，并令β←0.1×R_max(c,l)；

所述max{·}表示取最大值函数，β表示用于判断特征点的预设阈值；

步骤19.在子图像I_j中寻找特征角点；

步骤19.1令c←1，l←1，并初始化一幅边缘特征点图像I_corner，其高度为h_j像素，宽度为w_j像素；

步骤19.2若R(c,l)≥β，则将子图像I_j中位于第c行、第l列的像素判定为特征角点并令I_corner(c,l)←1，否则令I_corner(c,l)←0，所述I_corner(c,l)表示边缘特征点图像I_corner中位于第c行、第l列的像素值；

步骤19.3令l←l+1，若l≤w_j，则返回步骤19.2，否则，转入步骤19.4；

步骤19.4令c←c+1，l←1，若c≤h_j，则返回步骤19.2，否则，转入步骤20；

步骤20.初始化一幅图像其高度为2×h_j像素，宽度为2×w_j像素，用于存储修复后的第j个甲骨卜辞单字图像；

步骤21.利用线性插值方法对图像的全部像素值进行平滑重构，从而保持原始甲骨卜辞文字的主要字形及其轮廓信息；

步骤21.1根据公式(7)，对图像中偶数行、偶数列的全部像素点进行平滑插值重构；

所述表示图像/>中位于第2m行、第2n列的像素值，I_j(m,n)表示图像I_j中位于第m行、第n列的像素值，且m∈{1,2,3,…,h_j}，n∈{1,2,3,…,w_j}；

步骤21.2根据公式(8)，对图像中奇数行、偶数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m-1行、第2n列的像素值，/>表示图像/>中位于第2m-2行、第2n列的像素值，γ表示预设的线性加权系数；

步骤21.3根据公式(9)，对图像中偶数行、奇数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m行、第2n-1列的像素值，/>表示图像/>中位于第2m行、第2n-2列的像素值；

步骤21.4根据公式(10)，对图像中奇数行、奇数列的全部像素点进行平滑插值重构，边界像素采用奇对称延拓原则进行重构；

所述表示图像/>中位于第2m-1行、第2n-1列的像素值，表示图像/>中位于第2m-1行、第2n-2列的像素值；

步骤22.将边缘特征点图像I_corner视为拐点区域，进而根据I_j与之间的几何对偶一致性，利用I_j的梯度协方差系数计算出最优的线性插值权重，对/>进行梯度协方差引导的边缘定向插值；

步骤22.1令c←1，l←1；

步骤22.2若I_corner(c,l)＝1，则转入步骤22.3，否则转入步骤22.9；

步骤22.3在I_j中以第c行、第l列的像素为中心建立一个大小为W×W的窗口，并将该窗口中的全体像素值按照光栅扫描顺序保存到一个列向量y中，得到所述W为一个预设常数，/>表示I_j中位于第/>行、第/>列的像素值，I_j(c,l)表示I_j中位于第c行、第l列的像素值，/>表示I_j中位于第行、第/>列的像素值；

步骤22.4建立一个阶数为4×W²的矩阵C，对于列向量y的第K个元素I_j(c′,l′)，将I_j(c′-1,l′-1)、I_j(c′+1,l′-1)、I_j(c′-1,l′+1)、I_j(c′+1,l′+1)保存到矩阵C的第K列，所述1≤K≤W²，I_j(c′,l′)表示I_j中位于第c′行、第l′列的像素值，I_j(c′-1,l′-1)表示I_j中位于第c′-1行、第l′-1列的像素值，I_j(c′+1,l′-1)表示I_j中位于第c′+1行、第l′-1列的像素值，I_j(c′-1,l′+1)表示I_j中位于第c′-1行、第l′+1列的像素值，I_j(c′+1,l′+1)表示I_j中位于第c′+1行、第l′+1列的像素值；

步骤22.5根据公式(11)，计算基于最小均方误差的最优线性插值权重向量α；

α＝(CC^T)^-1×(Cy) (11)

步骤22.6根据公式(12)，对进行梯度协方差引导的边缘定向插值；

所述表示图像/>中位于第2c+1行、第2l+1列的像素值，/>表示图像/>中位于第2c行、第2l列的像素值，/>表示图像/>中位于第2c+2行、第2l列的像素值，/>表示图像/>中位于第2c+2行、第2l+2列的像素值，表示图像/>中位于第2c行、第2l+2列的像素值，α₀表示向量α的第1个分量，α₁表示向量α的第2个分量，α₂表示向量α的第3个分量，α₃表示向量α的第4个分量；

步骤22.7根据公式(13)，对进行梯度协方差引导的边缘定向插值；

所述表示图像/>中位于第2c+1行、第2l列的像素值，/>表示图像/>中位于第2c+1行、第2l-1列的像素值；

步骤22.8根据公式(14)，对进行梯度协方差引导的边缘定向插值；

所述表示图像/>中位于第2c行、第2l+1列的像素值，/>表示图像/>中位于第2c-1行、第2l+1列的像素值；

步骤22.9令l←l+1，若l≤w_j，则返回步骤22.2，否则，转入步骤22.10；

步骤22.10令c←c+1，l←1，若c≤h_j，则返回步骤22.2，否则，得到第j个甲骨卜辞单字图像的轮廓修复结果并转入步骤23；

步骤23.令j←j+1，若j≤N_text，则返回步骤10，否则，输出全部甲骨卜辞单字图像的轮廓修复结果。