CN104504710A - X射线光栅相位衬度成像中摩尔条纹的识别方法及装置 - Google Patents

Abstract

本发明公开了一种X射线光栅相位衬度成像中的摩尔条纹的识别方法及装置，所述方法包括：步骤1：对待识别的摩尔条纹图像进行光照不均匀修正；步骤2：对经过光照不均匀修正后的摩尔条纹图像进行滤波；步骤3：对滤波后的摩尔条纹图像进行二值化，获得二值化图像；步骤4：对所述二值化图像进行细化，以提取二值化图像中摩尔条纹的中心线；步骤5：根据所提取的摩尔条纹的中心线识别摩尔条纹的精确位置。本发明提供的方法实现了摩尔条纹角度和方向计算的自动化，避免的仪器使用者的主观因素导致的误差，并且仪器调节速度大大提高，精确度也得到了很好的保证。

Description

X射线光栅相位衬度成像中摩尔条纹的识别方法及装置

技术领域

本发明涉及一种条纹图像的角度与周期的自动检测方法，特别的，涉及X射线光栅相位衬度成像装置对准过程中的莫尔条纹的角度和周期的自动识别方法。

背景技术

自1895年德国科学家伦琴首次发现X射线至今的100多年里，X射线因其极强的穿透本领，被广泛用于物体的成像领域。传统的X射线成像方法主要是基于物体对X射线的吸收，最终获得的成像衬度的好坏很大程度上取决于物体内部各部分对X射线吸收特性差异的大小。在医学领域，由于人体软组织部分对X射线的吸收很少，这就意味着X射线吸收成像方法运用于人体的软组织病变诊断有很大的局限性。

自上世纪九十年代开始，随着第三代同步辐射装置的发展，硬X射线相位衬度成像技术应运而生。目前已经有多种X射线相位衬度成像技术得以发展。其机理简单地说，就是利用了X射线穿透物体后其相位发生的移动进行成像。相比较吸收衬度成像，相位衬度成像的优势在于，同样剂量的X射线穿透软组织，相位移动产生的变化比射线强度吸收产生的变化大得多，因此所得到的X射线图像衬度将会得到很大的提高，参见参考文献[1]。

X射线光栅步进相衬成像方法是目前发展较为成熟的一种X射线相衬成像方法，由于其可利用通用X射线机产生的多色、非相干性光进行成像，目前被广泛采用，参见参考文献[2]。现在普遍采用的X射线光栅相位衬度成像方法，是Pfeiffer F等人于2006年首次提出的，该方法采用三块不同功能的光栅，实现了在通用X射线机上完成相位衬度成像。实验中，光源光栅主要作用是将普通X射线光源分割成一系列互不相干的线光源。物体样本放置于相位光栅之前，单个X射线线光源是部分相干的，可以与相位光栅产生泰伯效应，最后通过放置于探测器之前的分析光栅，获取相位变化信息，参见参考文献[3]。

采用该方法获取相位衬度图像的一个关键环节是相位光栅和分析光栅的对准，对准的精度对获取的图像质量有明显的影响。对准方法为：先调整分析光栅，使其光栅线水平，之后通过调整相位光栅，使得相位光栅的泰伯自成像与分析光栅恰好完全重合。判定对准是否完成是依据相位光栅的泰伯自成像与分析光栅形成的莫尔条纹角度和周期判定的。依据摩尔条纹的相关知识可知，当摩尔条纹垂直并且周期相等并趋于无限大时就可以近似认为对准完成。而目前X射线光栅相衬成像实验中对准时摩尔条纹的角度与周期皆依据肉眼判定，其弊端为带有个人主观意识，精度无法保证并且速度较慢，参见参考文献[4]。

参考文献：

[1]Chapman L D，Tomlinson W C，Johnston R E，Washburn D，Pisano E，Gmur N，Zhong Z，Menk R，Arfelli F，Sayers D 1997phys.med.biol.42 2015

[2]Atsushi MOMOSE，Recent Advances in X-ray Phase Imaging，Japanese Journal of Applied Physics，Vol.44，No.9A，2005，pp.6355-6367

[3]Franz Pfeiffer，TimmWeitkamp，Oliver Bunk，Christian David，Phaseretrieval and differentialphase-contrast imaging with low-brillianceX-raysources，nature physics VOL 2 APRIL 2006

[4]PavloBaturin，Mark Shafer，Optimization of grating-basedphase-contrast imaging setup，Medical Imaging 2014：Physics of MedicalImaging，Vol.9033，90334

发明内容

为了实现X射线光栅相衬成像对准过程中摩尔条纹角度和周期的准确快速测算，从而提高仪器对准精度并获得高质量的相衬图像。

本发明提出了一种X射线光栅相位衬度成像中的摩尔条纹的识别方法，其包括：

步骤1：对待识别的摩尔条纹图像进行光照不均匀修正；

步骤2：对经过光照不均匀修正后的摩尔条纹图像进行滤波；

步骤3：对滤波后的摩尔条纹图像进行二值化，获得二值化图像；

步骤4：对所述二值化图像进行细化，以提取二值化图像中各摩尔条纹的初始位置信息；

步骤5：根据所提取的各摩尔条纹的初始位置信息识别各摩尔条纹的精确位置。

本发明还提出了一种X射线光栅相位衬度成像中的摩尔条纹的识别装置，其包括：

修正模块：对待识别的摩尔条纹图像进行光照不均匀修正；

滤波模块：对经过光照不均匀修正后的摩尔条纹图像进行滤波；

二值化模块：对滤波后的摩尔条纹图像进行二值化，获得二值化图像；

细化模块：对所述二值化图像进行细化，以提取二值化图像中各摩尔条纹的初始位置信息；

识别模块：根据所提取的各摩尔条纹的初始位置信息识别各摩尔条纹的精确位置。

与现有技术相比，本发明提供的方案实现了摩尔条纹角度和方向计算的自动化，避免的仪器使用者的主观因素导致的误差，并且仪器调节速度大大提高，精确度也得到了很好的保证。

附图说明

图1为X射线光栅步进相位衬度成像***的构成示意图；

图2为两块光栅的几种典型相对位置和相应的莫尔条纹示意图；

图3为本发明中X射线光栅相位衬度成像中的摩尔条纹图像的识别方法的流程图；

图4(a)-(e)为实现本发明X射线光栅相位衬度成像中的摩尔条纹图像的识别方法的软件界面和处理步骤示意图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明作进一步的详细说明。

摩尔条纹产生原理如下：由泰伯效应可知相位光栅在泰伯距离产生光栅自成像，我们将分析光栅放置于特定的泰伯距离，二者的条纹相交形成摩尔条纹，并由放置于分析光栅后面的CCD探测器接收获得原始图像。

图1示出了X射线光栅步进相位衬度成像***的构成。如图1所示，***从右至左依次为X射线光源，光源光栅，试验样品，相位光栅，分析光栅和CCD探测器。

图2示出了分析光栅和相位光栅的几种典型相对位置和相应的莫尔条纹。如图2所示，第一行各图左侧为分析光栅示意图，右侧为相位光栅的泰伯自成像，第二行为他们叠加后形成的莫尔条纹。

如图3所示，本发明公开了一种X射线光栅相位衬度成像中的摩尔条纹图像的识别方法，其包括：

步骤1：获取条纹图像，并对条纹图像进行光照不均匀修正；其中，所述条纹图像可以为摩尔条纹，也可以为一般的条纹图像；

实际实验中的X射线光源为点光源，因此获取的条纹图像必然存在光照不均匀的状况。不均匀的光照会使图像质量下降，显示效果变差更重要的是会影响后续摩尔条纹计算的精度。因此对获取的原始图像，本发明最先要进行光照不均匀修正。

针对点光源的特征可知，图像中心亮度最大，各个像素点亮度大小为该像素点到光照中心像素点距离的函数，往四周亮度逐渐衰减，符合点光源照射度的距离平方反比定律与点光源照射度的距离平方反比余弦定律。针对该特征，本发明提出了一种简单可行的修正方法。

点光源照射度的距离平方反比定律，即假设点光源的光强为I_θ，点光源到探测器照射中心的距离为I，则照射中心光照度为而点光源照射度的距离平方反比余弦定律，即假设点光源照射位置与照射面不垂直，仍假设点光源的光强为I_θ，点光源到探测器平面的距离为I，点光源与照射点法线方向夹角为θ，则照射位置光照度为据此假定获取的图像灰度值矩阵为A＝(a_ij)_m×n，光照中心位置像素点为a_pq，有这里的k为光照度到图像灰度值的转换系数，则对任意一位置像素点灰度值为

$a_{\mathrm{ij}}=k\frac{I_{\mathrm{\θ}}}{I^{\′2}}\cos \mathrm{\θ}=k\frac{I_{\mathrm{\θ}}}{{(I/\cos \mathrm{\θ})}^2}\mathrm{COS\θ}=k\frac{I_{\mathrm{\θ}}}{I^2}\·{\cos }^3\mathrm{\θ}$

其中I′为光源到待修正的当前像素点(i，j)的距离，因此我们可以得到当前像素点(i，j)的灰度值a_ij＝a_pq·cos³θ。据此，对于获取的图像上像素点原始灰度值a_ij，我们修正其灰度值为a′_ij＝a_ij/cos³θ，其中而I可以计算测得，设待修正的当前像素点(i，j)到光照中心的距离为i，j为原始图像上当前像素点(i，j)的像素坐标，p，q为光照中心位置像素点(p，q)的像素坐标，即 $a_{\mathrm{ij}}^{\′}=a_{\mathrm{ij}}\·{(I^2/(I^2+{(i-p)}^2+{(j-q)}^2))}^{-\frac{3}{2}}.$

故计算中指定光照中心坐标(p，q)，测量光源到传感器距离I，对图像上的像素点(i，j)便可完成光照不均匀修正。

步骤2：对经过光照不均匀修正的摩尔条纹图像进行滤波。

由于实际环境条件的限制，CCD本身特性等因素影响，获得的原始图像不免存在噪声干扰。CCD图像常见的噪声包括椒盐噪音、脉冲噪音、高斯噪音等，另外对本发明而言，光栅线的影像也是噪声信息。噪声信息会对后续图像处理结果有明显影响，因此需要对图像进行滤波处理。

常用的图像滤波方法包括空间域滤波和频域滤波。空间域滤波是指直接对图像的像素灰度进行处理，依据每个像素点的灰度值特征进行滤波。常见的空间域滤波方法包括直方图均衡化，中值滤波法，均值滤波法等。空间域滤波方法简单直观，但是滤波效果往往不够理想，因此本发明中优先采用频域滤波法。

频域滤波方法通常做法为：先将图像进行快速傅里叶变换，选取频域滤波函数进行滤波，将滤波后的图像进行傅里叶反变换，得到滤波后的图像，从而滤除噪声，提高图像质量。

当相位光栅和分析光栅发生干涉时会形成莫尔条纹，假设在X，Y平面内，坐标点表示为(x，y)，两块光栅的周期分别为d₁和d₂，令第一块光栅即相位光栅或分析光栅的栅线与Y轴平行，第二块光栅即分析光栅或相位光栅的栅线与Y轴顺时针成θ角，假定两光栅间隙为零，第一块光栅的透过函数的傅里叶变换为f₁(T)，第二块光栅的透过函数的傅里叶变换为f₂(T)，它们在空间形成的干涉场光强分布为：

$\begin{array}{c}f\left(T\right)=f_1\left(T\right)\·f_2\left(T\right)=\frac{a_{01}{a}_{02}}{4}+\frac{a_{02}}{2}\×\underset{n=1}{\overset{\∞}{\mathrm{\Σ}}}{a}_n\cos \left(n\frac{2\mathrm{\π}}{d_1}x\right)+\frac{a_{01}}{2}\×\underset{m=1}{\overset{\∞}{\mathrm{\Σ}}}{a}_m\cos \left(m\frac{2\mathrm{\π}}{d_1}\right(x\cos \mathrm{\θ}\\ -y\sin \mathrm{\θ}\left)\right)+\left\{\underset{n=1}{\overset{\∞}{\mathrm{\Σ}}}{a}_n\cos (n-\frac{2\mathrm{\π}}{d_1}x)\right\}\×\left\{\underset{m=1}{\overset{\∞}{\mathrm{\Σ}}}{a}_m\cos \left(m\frac{2\mathrm{\π}}{d_2}(x\cos \mathrm{\θ}-y\sin \mathrm{\θ})\right)\right\}\end{array}$

其中，x，y为摩尔条纹上任意点的坐标，T为傅里叶变换周期；a₀₁，a_n为第一块光栅透过函数的傅里叶变换系数，a₀₂，a_m为第二块光栅透过函数傅里叶变换系数，分析可知，上式等式右边第一项不含相位因子，其代表背景光，第二项含有频率成分其包含了第一块光栅即相位光栅的结构信息，第三项含有频率成分其包含了第二块光栅即分析光栅的结构信息，而第四项含有两光栅的和频和差频成分，属于摩尔条纹信息。因此第一项和第四项为我们需要保留的有用信息，第二项和第三项需要滤除。当然，实际得到的图像包含着更复杂的噪声频率成分，在滤波过程中都应该加以滤除。

根据以上分析的图像频率信息特征，本发明采用了一种低通与带通结合的组合滤波器，低通用于滤出背景成分，带通用以滤出莫尔条纹信息。滤波器种类，本发明中可以选择理想滤波器，指数滤波器，巴特沃夫滤波器等，其中理想滤波器形式较为简单并且效果达到要求，故为本发明优先采用的方法。即滤波函数为H(u，v)，有

$H(u,v)=\left\{\begin{array}{c}1,u^2+v^2\≤r_1\\ 1,{(u-u_0)}^2+{(v-v_0)}^2\≤r_2\\ 0,\mathrm{others}\end{array}\right.$

其中r₁为低通滤波器滤波半径，r₂为带通滤波器滤波半径，u，v为滤波平面上点的坐标，u₀，v₀为带通滤波中心坐标。

故假设滤波前图像信息傅里叶变换为F(u，v)，滤波后图像傅里叶变换应该为G(u，v)＝F(u，v)H(u，v)，对其进行傅里叶反变换，得到滤波后的图像。

步骤3：对进行了滤波后的原始图像进行二值化。

二值化顾名思义就是将图像分成感兴趣区域和不感兴趣区域两部分。本发明中，图像亮纹为感兴趣区域，二值化之后将所述感兴趣区域用1表示，而其他区域用0来表示。基于图像分割理论的图像二值化方法已经很成熟，现在通常采用的方法为阈值化分割，阈值化分割算法主要有两个步骤：确定需要的分割阈值和将分割阈值与像素值比较以划分像素。阈值化分割方法包括局部阈值法，全局阈值法，手动阈值法等，其中确定合适的阈值是图像分割的关键，而阈值的提取方法也是多种多样。不同的阈值分割方法适合于不同特征的图像，本发明需要处理的图像信息相对简单，并且经过预处理，图像质量也比较好，因此大部分阈值分割方法都适用。

步骤4：对进行了二值化后图像的感兴趣区域进行细化，以提取图像摩尔条纹的初始位置信息。

经过上述处理，图像信息已经大大简化，但是本发明的目标是获得摩尔条纹的精确位置，因此需要进一步提取莫尔条纹的位置信息，即图像细化。图像细化就是提取图像的主干信息，以莫尔条纹提取为例，就是获取摩尔条纹的中心点坐标和角度。

步骤5：根据所提取的摩尔条纹的初始位置信息识别摩尔条纹的精确位置。

经过细化处理提取的条纹初始位置信息已经大致满足精度需要，但是为了进一步提高精度，本发明提出了一种基于灰度方差权重的迭代计算方法，以识别摩尔条纹的精确位置。该方法具体为：

步骤51：获得当前待识别摩尔条纹的初始位置信息，设其中心点为(x₀，y₀)，长度为l，角度为θ₀；

步骤52：设定迭代计算次数，计算偏置角度α和划线数量n，以(x₀，y₀)为中心，l为长度，向左右侧各倾斜α划n条直线；

步骤53：设各直线上像素点灰度值的方差依次为σ₁，σ₂，σ₃，……，σ_2n。同时由之前描述可知，各直线对应的直线偏角β_i也不难计算，如 ${\mathrm{\β}}_1={\mathrm{\θ}}_0+\frac{1}{n}\mathrm{\α},{\mathrm{\β}}_2={\mathrm{\θ}}_0+\frac{2}{n}\mathrm{\α},$ 依次类推。

步骤54：根据各条直线上像素点灰度方差的大小，给各条直线倾斜角赋予权重 $Q_i=1-\frac{{\mathrm{\σ}}_i}{{\mathrm{\Σ}}_{k=1}^{2n}{\mathrm{\σ}}_k},$ 依照公式

${\mathrm{\θ}}_1=\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{\mathrm{\β}}_i\×Q_i=\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{\mathrm{\β}}_i\×(1-\frac{{\mathrm{\σ}}_i}{{\mathrm{\Σ}}_{k=1}^{2n}{\mathrm{\σ}}_k})$

计算得到的当前待识别摩尔条纹的新的角度；

步骤55：依据步骤52设定的迭代计算次数，重复步骤52到步骤54迭代计算过程，最终得到当前待识别摩尔条纹的精确位置，所述当前待识别摩尔条纹的精确位置包括当前待识别摩尔条纹的中心点坐标和角度，所述中心点坐标依然是步骤4中获取的初始位置信息中的中心点坐标，而角度为经过步骤51-55优化后的角度。

通过重复执行步骤5可以获得各摩尔条纹的精确位置信息，并且根据各摩尔条纹的中心点坐标能够得到摩尔条纹的周期。

步骤6：仪器校准达标判定。

对于计算得到的各条摩尔条纹，再以其中心点为基准，以相邻中心点的横坐标相减，计算它们的周期。依据实际需要，对角度和周期设定容许误差，当仪器调整至某一状态时，若计算得到的角度和周期都已经达标，就可以认为仪器对准已经完成。

本发明的上述摩尔条纹的识别方法可通过LabVIEW编程实现。

本发明利用LabVIEW编写软件，实现上述算法(注：LabVIEW上应该已安装有“视觉与运动”模块)。所有的图像处理过程都是基于图像灰度值进行的，所以首先应使用“IMAQ ImageToArray”函数将图像转换成灰度值二维矩阵。

图像光照修正部分。光源到传感器距离以及光照中心位置为输入量，像素点位置坐标以及像素点灰度值从图像读入，按之前所述带入公式 $a_{\mathrm{ij}}^{\′}=a_{\mathrm{ij}}\·{(I^2/({I^2+(i-p)}^2+{(j-q)}^2))}^{-\frac{3}{2}},$ 即可得到修正后的图像灰度矩阵，利用“IMAQ ArrayToImage”函数可得到修正后的图像。另外我们可以用“曲面”函数绘制图像灰度矩阵的三维图像，直观地观察修正效果。

图像滤波这部分用到了LabVIEW和MATLAB的混合编程技术。LabVIEW和MATLAB编程都有自己独到的优点，LabVIEW的图形化编程方式，可以让开发者更专注于算法本身，而MATLAB的图形处理工具箱则集成了诸多图像处理函数，给图像处理带来了极大的便利，因此将LabVIEW和MATLAB结合起来编程能大大提高我们的编程效率。LabVIEW和MATLAB的混合编程实现方法多种多样，本发明中采用的方法为在LabVIEW中使用“MATLAB script”节点调用MATLAB，其位于“数学＞脚本与公式＞脚本节点＞MATLAB脚本”，***该节点后在其左侧添加输入，右侧添加输出，节点内输入MATLAB处理代码，即可使用，简单方便。

因此滤波这一部分中，我们使用MATLAB script节点，修正后的图像灰度矩阵利用MATLAB中的“fft2”函数实现快速二维傅里叶变换，从前面板输入滤波中心坐标以及滤波半径大小，完成滤波函数构建。图像经过滤波后使用“ifft2”实现傅里叶反变换，得到滤波后的图像。

二值化部分。滤波后的图像在该部分实现二值化，该部分我们集成了“局部阈值法”“自动阈值法”和“手动阈值法”三种方法。其中“局部阈值法”使用“IMAQ Local Threshold”函数实现，“自动阈值法”使用“IMAQAutoBThreshold 2”函数实现，而“手动阈值法”使用“IMAQ Threshold”函数实现。具体输入参数根据各个函数需要以及图片状况由操作者输入。

图像细化其实是由细化和条纹初始位置确定两部分组成。图像的细化使用“IMAQ Skeleton”函数实现，条纹图像利用该函数处理后就能得到条纹的脉络。但是脉络图像可能会出现弯曲，毛刺多等问题，而且其具***置信息也无法确定，因此我们紧接着使用“IMAQ Find Straight Edges 2”函数来实现条纹位置的初步确定。该函数本是用于查找图像直线边缘，但是我们将其应用于细化后的图像，通过设置恰当的参数，它便可以实现查找图像中直线，并返回直线位置和角度的函数，因此可以用于条纹位置的初步确定。

条纹位置精确计算部分主要用到“IMAQ Line Profile”函数，该函数可以对图像上某一指定的直线返回该直线上像素点灰度值方差大小。获得方差后具体计算按照先前描述的去编程实现即可。

条纹精确位置确定后再前面板输入角度和方差的容许误差，将各个条纹位置与标准比照，即可判定角度和周期是否达标，并在前面板设置两盏布尔灯分别指示角度和周期判定结果，以灯亮表示达标。

图4(a)-(e)示出了本发明实施例中提出的一种X射线光栅相衬成像对准过程中莫尔条纹的角度和周期的自动识别方法与软件实现，本发明提出的上述方法也适用于任意含条纹图像的条纹识别。具体实施时候分为五个步骤完成。

如图4(a)步骤一所示，第一步为光照不均匀修正。先测量光源到探测器的距离并输入。读入原始摩尔条纹图像，将其转换为三维灰度图观察其灰度特征，依据观察到的结果，依次设置光照中心的横坐标和纵坐标，运行软件，图像的光照不均匀即可得到修正。修正后的图像以及其三维灰度图也显示于软件中，帮助我们观察修正效果。另外我们可以通过点击“存储修正参数”按钮来保存当前的修正设置。

如图4(b)步骤二所示，第二步为图像滤波。光照修正后的图像为输入，输入图像先进行二维快速傅里叶变换，输出变换后的频域图像便于我们观察需要的频率成分所处位置。依据观察到的图像结果，依次确定滤波器三个滤波点的横坐标纵坐标和滤波半径大小，设置完成之后运行软件输出滤波后的频域图像和其反变换，可以观图像的滤波效果，参数设置可以再微调。另外我们可以通过点击“存储滤波参数”按钮来保存当前的修正设置。

如图4(c)步骤三所示，第三步为图像二值化。输入滤波后的图像，选择一种二值化方法。这里我们选择的方法为局部阈值法中的Niblack方法，具体参数设置为，Niblack偏离系数为0.2，计算窗口大小为32像素乘以32像素。实际操作中我们也可以根据图片质量选择其他的二值化方法以取得最好的二值化效果。

如图4(d)步骤四所示，第四步为图像细化。二值化之后的图像作为输入，图像中可以划定需要细化和测量的范围。直线查找的选项较多，“细化方向”依据条纹方向而定，可以左右方向，也可以上下方向，“kernel size”输入最小值，直线选项中“number of lines”尽量大些以包括所有的直线，其他参数用默认值就可以。实际操作中也可以根据需要调整其他参数设置以得到更好的效果。检测到的条纹数目输出在软件界面的右下角。

如图4(e)步骤五所示，第五步为条纹计算和合格判定。步骤四查找到的直线信息如条纹中心点位置信息、角度和线段长度信息等输入到该步，依据直线信息在步骤二滤波后得到的图像上进行划线计算。划线的角度和数目由软件操作者输入，计算完成后，最终得到的条纹线将以红线形式显示在滤波后的图像中。我们输出每条条纹的角度信息以及各个条纹之间的周期信息，并将其以统计图的形式显示在程序中，方便观察。我们再设置角度的容许误差δ(角度±δ°以内为合格)和周期的容许误差Δ(若周期的平均值为周期像素以内为合格)，各条条纹的角度与周期若合格将点亮布尔指示灯。最后汇总各条条纹信息，再判定***总体的角度和周期是否达标，若达标点亮布尔灯，两盏灯都亮起即***对准完成，若不是，依据条纹角度和周期的分布情况选择合适的校准方法，继续调整光栅位置，重复各步工作。

以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种X射线光栅相位衬度成像中的摩尔条纹的识别方法，其包括：

步骤1：对待识别的摩尔条纹图像进行光照不均匀修正；

步骤2：对经过光照不均匀修正后的摩尔条纹图像进行滤波；

步骤3：对滤波后的摩尔条纹图像进行二值化，获得二值化图像；

步骤4：对所述二值化图像进行细化，以提取二值化图像中各摩尔条纹的初始位置信息；

步骤5：根据所提取的各摩尔条纹的初始位置信息识别各摩尔条纹的精确位置。

2.如权利要求1所述的方法，其中，所述步骤1中根据点光源照射度的距离平方反比定律和电光源照射度的距离平方反比余弦定律，对所述摩尔条纹图像进行光照不均匀修正。

3.如权利要求2所述的方法，其中，所述步骤1中修正后的摩尔条纹图像上的任意位置像素点灰度值如下计算：

$a_{\mathrm{ij}}^{\′}=a_{\mathrm{ij}}\·{(I^2/(I^2+{(i-p)}^2+{(j-q)}^2))}^{-\frac{3}{2}}$ 其中，a_ij和a′_ij分别为所述摩尔条纹图像上任意像素点(i，j)修正前和修正后的像素灰度值，I为点光源到摩尔条纹图像探测器照射位置的距离，p和q分别为光照中心位置像素点(p，q)的像素坐标。

4.如权利要求1-3任一项所述的方法，其中，步骤2具体包括：

步骤21：对所述经过光照不均匀修正后的摩尔条纹图像进行快速傅里叶变换；

步骤22：从快速傅里叶变换后的结果选取频域滤波函数进行滤波；

步骤23：将频域滤波后的图像进行傅里叶反变换，得到滤波后的摩尔条纹图像。

5.如权利要求4所述的方法，其中，步骤22中采用的频域滤波函数如下表示：

$H(u,v)=\left\{\begin{array}{c}1,u^2+v^2\≤r_1\\ 1,{(u-u_0)}^2+{(v-v_0)}^2\≤r_2\\ 0,\mathrm{others}\end{array}\right.$

其中r₁为低通滤波器滤波半径，r₂为带通滤波器滤波半径，u，v为滤波平面上点的坐标，u₀，v₀为带通滤波中心坐标。

6.如权利要求1-3、5任一项所述的方法，其中步骤5具体包括：

步骤51：获取当前待识别摩尔条纹的初始位置信息，设当前待识别摩尔条纹的中心点为(x₀，y₀)，长度为l，角度为θ₀；

步骤52：设定计算偏置角度α和划线数量n，以(x₀，y₀)为中心，l为长度，向左右侧各倾斜α划n条直线；

步骤53：设各直线上像素点的方差依次为σ₁，σ₂，σ₃，……，σ_2n，则计算各直线对应的直线偏角β_i；

步骤54：根据各条直线上像素点灰度方差的大小，给各条直线倾斜角赋予权重并依照公式

${\mathrm{\θ}}_1=\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{\mathrm{\β}}_i\×Q_i=\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{\mathrm{\β}}_i\×(1-\frac{{\mathrm{\σ}}_i}{{\mathrm{\Σ}}_{k=1}^{2n}{\mathrm{\σ}}_k})$

计算得到的当前待识别摩尔条纹的新的角度；

步骤55：依据步骤52设定的迭代计算次数，重复步骤52到步骤54迭代计算过程，最终得到当前待识别摩尔条纹的精确位置。

7.如权利要求1所述的方法，其中，步骤5中所述当前待识别摩尔条纹的精确位置包括当前待识别摩尔条纹的中心点和角度。

8.如权利要求7所述的方法，其中，该方法还包括：

步骤6：根据计算得到的各摩尔条纹的角度以及摩尔条纹周期判断是否达到仪器校准的标准。

9.一种X射线光栅相位衬度成像中的摩尔条纹的识别装置，其包括：

修正模块：对待识别的摩尔条纹图像进行光照不均匀修正；

滤波模块：对经过光照不均匀修正后的摩尔条纹图像进行滤波；

二值化模块：对滤波后的摩尔条纹图像进行二值化，获得二值化图像；

细化模块：对所述二值化图像进行细化，以提取二值化图像中各摩尔条纹的初始位置信息；

识别模块：根据所提取的各摩尔条纹的初始位置信息识别各摩尔条纹的精确位置。

10.如权利要求9所述的装置，其中，所述识别模块如下识别各摩尔条纹的精确位置：

首先获取当前待识别摩尔条纹的初始位置信息，设当前待识别摩尔条纹的中心点为(x₀，y₀)，长度为l，角度为θ₀；

设定迭代计算次数，计算偏置角度α和划线数量n，以(x₀，y₀)为中心，l为长度，向左右侧各倾斜α划n条直线；设各直线上像素点的方差依次为，σ₁，σ₂，σ₃，……，σ_2n，则计算各直线对应的直线偏角β_i；

根据各条直线上像素点灰度方差的大小，给各条直线倾斜角赋予权重 $Q_i=1-\frac{{\mathrm{\σ}}_i}{{\mathrm{\Σ}}_{k=1}^{2n}{\mathrm{\σ}}_k},$ 并依照公式

${\mathrm{\θ}}_1=\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{\mathrm{\β}}_i\×Q_i=\underset{i=1}{\overset{2n}{\mathrm{\Σ}}}{\mathrm{\β}}_i\×(1-\frac{{\mathrm{\σ}}_i}{{\mathrm{\Σ}}_{k=1}^{2n}{\mathrm{\σ}}_k})$

计算得到的当前待识别摩尔条纹的新的角度；

根据上述设定的迭代计算次数，重复上述迭代计算过程，最终得到当前待识别摩尔条纹的精确位置。