本发明属于低信噪比下的硅单晶生长图像检测技术领域，涉及一种基于PSO算法的自适应随机共振的硅单晶生长图像的去噪方法。
背景技术：
硅单晶被广泛应用于制造光伏电池及集成电路。硅单晶的制备方法通常是将多晶硅在高温下融化，再用直拉法或区熔法从熔体中生长出棒状的硅单晶。
在直拉法晶体生长过程中，晶体直径是需要检测并控制的重要宏观参数。在晶体生长过程中，引晶-放肩阶段中需要晶体直径由籽晶向目标直径过渡；在等径阶段需要测量晶体直径，从而控制晶体稳定生长；在收尾阶段需要逐步减小晶体直径，直到晶体直径缩小为一点并与硅液面脱离。因此，精确测量晶体直径在整个晶体生长过程中非常重要。晶体直径测量通常是利用图像处理的方法通过检测晶体与硅液的固液交界面—弯月面位置或大小实现的。然而，在晶体成像过程中，会不可避免的受到环境中的噪声的干扰，导致晶体图像被污染，从而影响晶体图像弯月面的检测，也就影响了晶体直径的测量。当晶体图像中含有的噪声量较多时，会导致无法检测弯月面，即无法检测晶体的直径。因此，有必要对低信噪比下的晶体图像做有效的去噪处理，提高图像的质量，使得弯月面得以检测，进而准确检测晶体直径。
基于随机共振理论的图像去噪是指在由非线性系统、图像信号、噪声三者组成的系统中，非线性系统、图像信号、噪声可以发生协同作用，从而将部分噪声能量转化为有用的图像信号能量，增强图像信号，从而提高图像质量。传统的图像处理去噪算法都是将噪声视为有害的，通过分析噪声的特征，以滤除噪声为目标，通过去除噪声来提高图像的质量。当图像噪声强度较大时，传统方法在滤除噪声的同时将不可避免地滤除图像中的有用信息，而随机共振的图像去噪方法在去除噪声的同时，依然可以保留图像中的特征，提高图像的质量。
技术实现要素：
本发明的目的是提供一种基于PSO自适应随机共振的硅单晶生长图像的去噪方法，提高低信噪比下的硅单晶生长图像质量，进而提高后续硅单晶弯月面检测的精确度。
为实现上述目的本发明采用以下技术方案：
低信噪比下的硅单晶生长图像的自适应随机共振去噪方法，包括如下步骤：
步骤1：将硅单晶生长图像做灰度映射与降维扫描得到一维信号，其中，灰度映射是将图像的灰度值范围线性映射至[0,1]之间，降维扫描是将图像按照逐行扫描，得到一维信号；
步骤2：建立双稳态系统模型，并将步骤1中得到的一维信号作为双稳态系统输入信号，经过数值求解后得到一维输出信号；
步骤3：对步骤2中得到的一维输出信号按行逆扫描成图像，并进行灰度映射，使其对比度线性拉伸，灰度值均匀分布在[0,1]之间，提高图像的对比度；
步骤4：通过粒子群优化算法(PSO)调整双稳态系统的参数a和b，选用无参考图像的评价指标—Donoho噪声标准差作为PSO算法的适应度函数，当PSO没有达到最大迭代次数，则返回执行步骤2；计算不同参数a和b下对应的双稳态系统输出图像的Donoho噪声标准差，当PSO达到最大迭代次数时，则执行步骤5；
步骤5：寻找步骤4中得到的所有Donoho噪声标准差中的最小值，该最小值对应的参数a和b的值即为最优参数解，该参数下的随机共振图像，即为最终的输出图像。
作为本发明进一步的方案，步骤2中建立的双稳态随机共振模型如下：
其中，u(t)为一维输入信号，u(t)＝s(t)+η(t)，其中s(t)为输入信号中不包含噪声的有用信号，η(t)为具有高斯分布的白噪声，统计平均值＜η(t)＞＝0，自相关函数为＜η(t)η(t')＞＝2Dδ(t-t')，t'为t的时间延时，D为噪声强度；
U(x)是非线性双稳态势函数，其表达式为：
其中a，b为非线性双稳系统的结构参数，且有a＞0、b＞0，因此，双稳态随机共振模型表示如下：
当噪声驱动下的双稳态势阱间跃迁率与平均等待时间匹配时系统发生随机共振，双稳态势阱间跃迁率称为Kramers跃迁率rK，它取决于噪声分布和强度，定义为：
跃迁速率与参数a，b有关，当a，b取得适当值时，噪声和有用信号在两个势阱中做同频率运动，此时即发生随机共振，噪声能量转移到有用信号上，从而增强有用信号，达到去噪的目的；
对图像做双稳态随机共振处理时，图像信号f(x,y)经过步骤1中的降维扫描和灰度映射转换成了一维信号u(n)，该信号中包含有用的图像信号和噪声，当双稳态系统参数a，b取得适当值时，该信号经过双稳态系统会发生随机共振，噪声的能量转移到有用的图像信号上，从而实现图像增强检测。
针对双稳态随机共振系统方程，采用四阶龙格-库塔算法(Runge-Kutta)对其求解：
其中x(n)为系统输出第n次采样值；u(n)＝s(n)+η(n)为输入信号与噪声之和的第n次采样值，在前面得到的图像的一维信号u中包含有图像信息和噪声信息；h为采样步长，在实际中为采样间隔。
作为本案发明进一步的方案，其中步骤4中选用粒子群的参数优化算法，对双稳态系统参数a和b进行优化，具体步骤如下：
首先初始化双稳态系统参数a和b，即粒子的最大位置Xmax，再初始化粒子群的各项参数，具体包括：寻优维数dimension，粒子个数popsize_num，学习因子c1、c2，最大飞行速度Vmax，算法的最大迭代次数Imax；
计算PSO中第i代的第k个粒子的适应度函数值fitnessi,k时，选用面向噪声失真的无参考图像指标—改进的Donoho噪声标准差作为PSO的适应度函数，Donoho噪声标准差按如下方法计算：
对观测图像进行正交小波变换，记Gj(u,v)为第j尺度上的小波系数，由小波变换的线性特征，得到：
Gj(u,v)＝Fj(u,v)+Nj(u,v) (6)；
其中，Fj(u,v)和Nj(u,v)分别为原始图像和噪声在第j尺度上的小波系数，Nj仍然是零均值的高斯噪声，而Fj通常遵从广义高斯分布并且均值为零，记Fj的方差为σ2F,j,Nj的方差为σ2N,j,可由Gj对原始图像小波系数Fj进行估计，即：
Fj(u,v)＝λGj(u,v) (7)；
式中λ为取定的系数，按照最小均方误差准则，可得λ的最佳取值为：
根据小波系数与滤波器之间的关系，得到：
σN,j2＝(||Lj-1||σn)2 (9)；
式中，σn为空域中的噪声标准差；||Lj-1||为滤波器的范数，即：
对于确定的子带，||Lj-1||已知，σn可由Donoho公式粗略估计得到:
s＝median(D)/0.6745 (11)；
式中，D为含噪图像在高频对角子带HH1中的系数的绝对值；median表示取中值,这样σ2N,j可求，由此σ2F,j可估计为：
σ2F,j＝σ2G,j-σ2N,j (12)；
其中，σ2G,j的表达式为：
M和N分别为相应子带的行数和列数，这样，取j＝1时，通过式(7)～(13)就可以得到图像中噪声的方差估计值，该方差估计值即为粒子的适应度函数值，即：
fitnessi,k＝σ2i,k (14)；
在每一个粒子的fitnessi,k中寻找自身的极小值Pi(t)＝(Pi1,Pi2,…,Pid)，然后在整个种群中的时刻t搜索得到全局最优解Pg(t)＝(Pg1,Pg2,…,Pgd)，当粒子找到了过程中发现的两个最优解，则同时按式(15)和式(16)更新种群中所有粒子的飞行速度Vi(t)和位置Xi(t)：
Vi(t+1)＝w×Vi(t)+c1×r1×(Pi(t)-Xi(t))+c2×r2×(Pg(t)-Xi(t) (15)；
Xi(t+1)＝Xi(t)+Vi(t+1) (16)；
式中w为惯性权重，w较大则可较好的进行全局搜索，w较小则可较好的进行局部搜索，采用了线性递减w方法来选取合适的w；c1＝c2＝2为学习因子；r1,r2为两个随机数，均匀分布在区间[0,1]；
w的计算公式如下：
当PSO算法达到最大迭代次数Imax时，算法终止。
与现有技术相比，本发明具有以下优点：本发明利用随机共振具有将噪声能量转移到信号上从而增强有用信号的特性，将含噪声的硅单晶生长图像经过随机共振处理，可有效去除噪声，提高图像的质量；利用粒子群算法(PSO)优良的寻优特性，对双稳态系统的参数进行优化，实现参数的自适应选择，快速达到随机共振的状态，大大提高处理的效率，克服了系统参数选择难的问题；选用面向噪声失真的无参考量化指标—改进的Donoho噪声标准差，作为PSO算法的适应度函数，也作为图像去噪的评价指标，符合实际应用情况中无法获得原始图像的特点，提高了实用价值。
附图说明
图1为本发明方法中算法的总体流程图；
图2为本发明方法中双稳态势阱图；
图3为本发明中双稳态系统发生随机共振时，粒子在两个势阱中的跃迁示意图；
图4为基于随机共振原理的图像处理示意图；
图5为本发明中改进的Donoho噪声标准差随高斯白噪声强度变化图；
图6为放肩阶段和等径阶段的硅单晶图像分别经过自适应随机共振处理后的结果图；
图7为本发明实施例中对等径阶段含噪硅单晶图像进行自适应随机共振去噪和其他三种图像处理算法去噪的结果图；
图8为本发明实施例中对放肩阶段含噪硅单晶图像直接检测弯月面，和经过自适应随机共振去噪后检测弯月面的结果图；
图9为本发明实施例中对等径阶段含噪硅单晶图像直接检测弯月面，和经过自适应随机共振去噪后检测弯月面的结果图。
具体实施方式
下面结合附图和具体实施例对本发明作进一步详细阐述。
本发明针对8英寸硅单晶生产炉中的晶体生长图像进行去噪处理，其算法的总体示意图如图1所示，含噪的硅单晶图像经过灰度映射和降维扫描后得到相应的一维信号；将此一维信号作为双稳态系统的输入信号；采用粒子群寻优算法(PSO)，以改进的Donoho噪声标准差作为PSO的适应度函数(也是图像去噪效果的评价指标)，寻找双稳态系统的最佳参数，该参数下得到的随机共振图像作为最终的输出图像。
低信噪比下的硅单晶生长图像的自适应随机共振去噪方法，具体按照以下步骤实施：
步骤1：将硅单晶生长图像做灰度映射与降维扫描得到一维信号。其中，灰度映射是将图像的灰度值范围由[0,255]线性映射至[0,1]之间，降维扫描是将图像按照逐行扫描，得到一维信号；
步骤2：建立双稳态系统模型，并将步骤1中得到的一维信号作为双稳态系统输入信号，经过数值求解后得到一维输出信号；
步骤2具体实施过程：
建立双稳态随机共振模型如下：
其中，u(t)为降维后的一维图像输入信号，u(t)＝s(t)+η(t)，其中s(t)为输入信号中不包含噪声的有用信号，η(t)为具有高斯分布的白噪声，统计平均值＜η(t)＞＝0，自相关函数为＜η(t)η(t')＞＝2Dδ(t-t')，t'为t的时间延时，D为噪声强度。
U(x)是非线性双稳态势函数，其表达式为：
其中a，b为非线性双稳系统的结构参数，且有a＞0、b＞0。因此，双稳态随机共振模型表示如下：
双稳态势函数的曲线如图2所示，粒子在两个势阱内运动的示意图，如图3所示。
当噪声驱动下的双稳态势阱间跃迁率与平均等待时间匹配时系统发生随机共振。双稳态势阱间跃迁率称为Kramers跃迁率rK，它取决于噪声分布和强度，定义为：
可以看到，跃迁速率与参数a，b有关，当a，b取得适当值时，噪声和有用信号在两个势阱中做同频率运动，此时即发生随机共振，噪声能量转移到有用信号上，从而增强有用信号，达到去噪的目的。
对硅单晶生长图像做双稳态随机共振处理时，图像信号f(x,y)经过步骤1中的降维扫描和灰度映射转换成了一维信号u(n)，该信号中含有图像中的有用信号和噪声，当双稳态系统参数a，b取得适当值时，该信号经过双稳态系统会发生随机共振，噪声的能量转移到有用的图像信号上，从而增强有用信号。基于随机共振原理的图像处理示意图如图4所示。
针对双稳态随机共振系统方程，采用四阶龙格-库塔算法(Runge-Kutta)对其求解：
其中x(n)为系统输出第n次采样值；u(n)＝s(n)+η(n)为输入信号与噪声之和的第n次采样值，在前面得到的图像的一维信号u中包含有图像信息和噪声信息；h为采样步长，在实际中为采样间隔。
步骤3：对步骤2中得到的一维输出信号按行逆扫描成图像，并进行灰度映射，使其对比度线性拉伸，灰度值均匀分布在[0,1]之间，提高图像的对比度；
步骤4：通过粒子群寻优算法(PSO)调整双稳态系统的参数a和b，选用无参考图像的评价指标—Donoho噪声标准差作为PSO算法中的适应度函数，当PSO没有达到最大迭代次数，则返回执行步骤2，计算不同参数a和b下对应的双稳态系统输出图像的Donoho噪声标准差，当PSO达到最大迭代次数时，则执行步骤5；
步骤4的具体实施过程：
选用粒子群的参数优化算法(PSO)，对双稳态系统参数a和b进行优化，具体步骤如下：
首先初始化双稳态系统参数a和b，即粒子的最大位置Xmax，再初始化粒子群的各项参数，具体包括：寻优维数dimension，粒子个数popsize_num，学习因子c1、c2，最大飞行速度Vmax，算法的最大迭代次数Imax。
计算PSO中第i代的第k个粒子的适应度函数值fitnessi,k时，选用面向噪声失真的无参考图像指标：改进的Donoho噪声标准差作为PSO的适应度函数，Donoho噪声标准差按如下方法计算：
对观测图像进行正交小波变换，记Gj(u,v)为第j尺度上的小波系数。由小波变换的线性特征，可得：
Gj(u,v)＝Fj(u,v)+Nj(u,v) (6)；
其中，Fj(u,v)和Nj(u,v)分别为原始图像和噪声在第j尺度上的小波系数，Nj仍然是零均值的高斯噪声，而Fj通常遵从广义高斯分布并且均值为零。记Fj的方差为σ2F,j,Nj的方差为σ2N,j,可由Gj对原始图像小波系数Fj进行估计。即：
Fj(u,v)＝λGj(u,v) (7)；
式中λ为取定的系数，按照最小均方误差准则，可得λ的最佳取值为：
根据小波系数与滤波器之间的关系，可得：
σN,j2＝(||Lj-1||σn)2 (9)；
式中，σn为空域中的噪声标准差；||Lj-1||为滤波器的范数，即：
对于确定的子带，||Lj-1||已知，σn可由Donoho公式粗略估计得到:
s＝median(D)/0.6745 (11)；
式中，D为含噪图像在高频对角子带HH1中的系数的绝对值。median表示取中值,这样σ2N,j可求，由此σ2F,j可估计为：
σ2F,j＝σ2G,j-σ2N,j (12)；
其中，σ2G,j的表达式为：
M和N分别为相应子带的行数和列数。这样，取j＝1时，通过上式就可以得到图像中噪声的方差估计了。该方差估计即为粒子的适应度函数值，即：
fitnessi,k＝σ2i,k (14)；
在每一个粒子的fitnessi,k中寻找自身的极小值Pi(t)＝(Pi1,Pi2,…,Pid)，再在整个种群中在时刻t搜索得到的全局最优解Pg(t)＝(Pg1,Pg2,…,Pgd)。当粒子找到了过程中发现的两个最优解，则同时按式(15)和式(16)更新种群中所有粒子的飞行速度Vi(t)和位置Xi(t)：
Vi(t+1)＝w×Vi(t)+c1×r1×(Pi(t)-Xi(t))+c2×r2×(Pg(t)-Xi(t) (15)；
Xi(t+1)＝Xi(t)+Vi(t+1) (16)；
式中w为惯性权重，w较大则可较好的进行全局搜索，w较小则可较好的进行局部搜索，采用了线性递减w方法来选取合适的w；c1＝c2＝2为学习因子；r1,r2为两个随机数，均匀分布在区间[0,1]。
w的计算公式如下：
当PSO算法达到最大迭代次数Imax时，算法终止。
步骤5：找到步骤4中得到的所有Donoho噪声标准差中的最小值，此时对应的参数a和b的值即为最优参数解，输出该参数下的随机共振图像，即为最终的输出图像。
为了验证上述基于随机共振的硅单晶生长图像的去噪效果，针对8英寸硅单晶生产炉中拍摄的硅单晶图像，选取等径阶段和放肩阶段的硅单晶图像作为对象进行处理。首先把图像转化成对应的灰度图像，然后作为算法的输入图像进行灰度映射和降维扫描处理；建立双稳态随机共振模型，设置双稳态系统的初始参数a＝1，b＝1；设置PSO寻优算法中的粒子数量popsize_num＝100个，寻优维数dimension＝2，学习因子c1＝c2＝2，粒子的最大位置Xmax＝20，最大飞行速度Vmax＝20，算法的最大迭代次数Imax＝100；另外，将本发明方法与其他几种传统的去噪算法进行了对比。
图5为改进后的Donoho噪声标准差随高斯白噪声强度变化图，可以看到，改进的Donoho噪声标准差与噪声强度D近似成正比例关系，因此用它来评价图像的噪声量是完全可行的。
图6为放肩阶段和等径阶段的晶体图像分别经过双稳态随机共振处理后的结果图，放肩阶段的硅单晶含噪图经过PSO自适应寻优参数a＝1.2354，b＝3.5720,等径阶段的硅单晶含噪图经过PSO自适应寻优参数：a＝2.3301，b＝7.6276。
图7为含噪的硅单晶生长图像经过自适应随机共振去噪处理的结果，并与其它几种图像去噪算法：低通滤波、维纳滤波、中值滤波进行了对比。从图7中可以看出，自适应随机共振去噪后的晶体图的亮度明显增加。
表1为本发明实施例中对等径阶段含噪硅单晶图像进行自适应随机共振和其他三种图像处理算法去噪后图像的Donoho噪声标准差；
表1
如表1所示含噪的硅单晶生长图像经过自适应随机共振去噪处理的结果，并与其它几种图像去噪算法处理后的图像的Donoho噪声标准差的对比。Donoho噪声标准差表示的是图像中的含噪量，通常该指标的数值越小代表图像中噪声含量越小。由表1可以看出，原含噪图像的Donoho噪声标准差为70.1304，经过表中的几种去噪算法处理后，图像的Donoho噪声标准差降低到10左右，但自适应随机共振方法的标准差值最低，具有较好的去噪优势。
图8为放肩阶段硅单晶的弯月面和区域边缘的检测图，图a为原始放肩阶段含噪的硅单晶图像；图b、图c为原始含噪硅单晶图像直接检测弯月面和区域边缘的结果图；图d为经过自适应随机共振去噪处理后的图像；图e、图f为图d的弯月面和区域边缘检测图。
图9为等径阶段硅单晶的弯月面和区域边缘的检测图，图a为原始等径阶段含噪的硅单晶图像；图b、图c为原始含噪硅单晶图像直接检测弯月面和区域边缘的结果图；图d为经过自适应随机共振去噪处理后的图像；图e、图f为图d的弯月面和区域边缘检测图。从图中可以看到，经过自适应随机共振去噪处理后的图像，检测出来的硅单晶弯月面清晰度有明显的提高，有利于通过后续方法精确检测硅单晶的直径，为保证等径阶段硅单晶稳定生长奠定基础。
以上所述为本发明较佳实施例，对于本领域的普通技术人员而言，根据本发明的教导，在不脱离本发明的原理与精神的情况下，对实施方式所进行的改变、修改、替换和变型仍落入本发明的保护范围之内。