CN111932573B - 一种光学系统空间分辨率的自动测试方法 - Google Patents

一种光学系统空间分辨率的自动测试方法 Download PDF

Info

Publication number
CN111932573B
CN111932573B CN202010635729.1A CN202010635729A CN111932573B CN 111932573 B CN111932573 B CN 111932573B CN 202010635729 A CN202010635729 A CN 202010635729A CN 111932573 B CN111932573 B CN 111932573B
Authority
CN
China
Prior art keywords
point
image
gray
edge
fitting
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010635729.1A
Other languages
English (en)
Other versions
CN111932573A (zh
Inventor
齐子诚
倪培君
郑颖
郭智敏
付康
张荣繁
左欣
张维国
乔日东
徐向群
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Weapon Science Academy Ningbo Branch
Original Assignee
China Weapon Science Academy Ningbo Branch
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China Weapon Science Academy Ningbo Branch filed Critical China Weapon Science Academy Ningbo Branch
Priority to CN202010635729.1A priority Critical patent/CN111932573B/zh
Publication of CN111932573A publication Critical patent/CN111932573A/zh
Application granted granted Critical
Publication of CN111932573B publication Critical patent/CN111932573B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/136—Segmentation; Edge detection involving thresholding
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B15/00—Measuring arrangements characterised by the use of electromagnetic waves or particle radiation, e.g. by the use of microwaves, X-rays, gamma rays or electrons
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N23/00—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00
    • G01N23/02—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material
    • G01N23/04—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material
    • G01N23/046—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material using tomography, e.g. computed tomography [CT]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/20—Image enhancement or restoration using local operators
    • G06T5/30—Erosion or dilatation, e.g. thinning
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/40—Image enhancement or restoration using histogram techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/0002—Inspection of images, e.g. flaw detection
    • G06T7/0004—Industrial image inspection
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/11—Region-based segmentation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/13—Edge detection
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/60—Analysis of geometric attributes
    • G06T7/66—Analysis of geometric attributes of image moments or centre of gravity
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2223/00—Investigating materials by wave or particle radiation
    • G01N2223/03—Investigating materials by wave or particle radiation by transmission
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2223/00—Investigating materials by wave or particle radiation
    • G01N2223/10—Different kinds of radiation or particles
    • G01N2223/101—Different kinds of radiation or particles electromagnetic radiation
    • G01N2223/1016—X-ray
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2223/00—Investigating materials by wave or particle radiation
    • G01N2223/60—Specific applications or type of materials
    • G01N2223/646—Specific applications or type of materials flaws, defects
    • G01N2223/6462—Specific applications or type of materials flaws, defects microdefects
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10072—Tomographic images
    • G06T2207/10081—Computed x-ray tomography [CT]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30108—Industrial image inspection

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Biochemistry (AREA)
  • Pulmonology (AREA)
  • Radiology & Medical Imaging (AREA)
  • Electromagnetism (AREA)
  • Chemical & Material Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • General Health & Medical Sciences (AREA)
  • Immunology (AREA)
  • Pathology (AREA)
  • Geometry (AREA)
  • Quality & Reliability (AREA)
  • Analysing Materials By The Use Of Radiation (AREA)

Abstract

一种光学系统空间分辨率的自动测试方法，包括：对贴合在一起的两测试块进行CT扫描，并利用自动阈值分割方法将其转换为二值化图像，且计算出二值化图像中质心的位置；对二值化图像进行形态学处理，并计算边缘图像中距离质心最近的点N；将质心与点N连成线段，将边缘图像中每个像素值为1的边缘点按照位于线段CN的左、右侧进行分类；之后，对分类后得到的第一集合和第二集合中所有边缘点进行拟合得到双圆的圆心坐标；利用两圆心坐标和计算出两圆心对应的中轴线；接着，提取中轴线上的灰度分布，对每个灰度值进行归一化处理后则进行线性拟合；最后，则计算MTF曲线。该方法的测试模体结构简单，加工难度小且成本低，测量精度更高，速度快。

Description

一种光学系统空间分辨率的自动测试方法
技术领域
本发明涉及光学系统性能测试技术领域，特别涉及一种光学系统空间分辨率的自动测试方法。
背景技术
工业CT检测技术是一种在X射线检测技术上发展起来的实用无损检测手段，是一种特殊的光学系统，具有成像直观，定量、定位、定性准确，广泛应用于工业无损探伤、医疗卫生等领域。在CT设备研发、生产、验收、使用和调试维修等过程中均需要涉及对工业CT系统的性能评价，因此对生产厂家和用户迫切需要一种测试便捷、精度高、实用性强的工业CT系统的性能评价方法。
空间分辨率是工业CT系统的性能评价重要指标之一，它表示CT系统分辨微小细节的量度，一般已能够分辨的两个细节之间的最小间距来定量表示。目前，空间分辨率测试方法主要分为条纹(圆孔)模体测量法和调制传递函数法。其中条纹(圆孔)模体测量法是采用一系列不同周期性的结构(线对、圆孔、方孔等)，对模体进行扫描重建，观察CT图像以能分辨条纹或圆孔图案的最小周期作为极限分辨率，可以方便直观地得到CT系统的分辨率，但是测量结果带有主观性，测量结果不准确，测量结果的精度根据加工模体的周期结构间隔而定，加工难度大，成本极高；调制传递函数法是调制度与输入空间频率的关系曲线，在数值上等于点扩散函数的傅里叶变换，一般情况下把MTF曲线上调制度10％对应的线对数作为CT系统的极限分辨率。调制传递函数法可分为点扩散函数法(PSF)、线扩散函数法(LSF)和边缘扩散函数法(ESF)，该方法对模体要求不高，噪声影响较大，可用于比较设备之间及工艺条件下的空间分辨率，给出归一化的量化结果是不直观，需要等效转换才能与实际线对卡结果进行比较。
因此需要进一步改进。
发明内容
本发明所要解决的技术问题是针对现有技术的现状，提供一种测试模体加工难度低，测试结果直观且测试精度高、速度快的光学系统空间分辨率的自动测试方法。
本发明解决上述技术问题所采用的技术方案为：一种光学系统空间分辨率的自动测试方法，其特征在于：用于自动测量线阵工业CT系统或面阵工业CT系统中的空间分辨率，包括以下步骤：
步骤1、通过机械化加工手段制造待测试模体；其中，待测试模体为两个材质与尺寸均相同的圆柱形或球形测试块；
步骤2、将步骤1中的两测试块贴合在一起，并对贴合在一起的两测试块进行CT扫描，获得两测试块截面对应的CT图像；
步骤3、利用自动阈值分割方法将步骤2中的CT图像转换为二值化图像b(x,y)，并计算出二值化图像b(x,y)中质心C的位置(xc,yc)，其中，x为该二值化图像区域的横坐标，y为该二值化图像区域的纵坐标，该二值化图像区域的尺寸为m×n像素，xc和yc分别对应为质心C对应在二值化图像中的横坐标和纵坐标；
步骤4、对二值化图像b(x,y)进行形态学处理，得到边缘图像e(x,y)；
步骤5、计算出边缘图像e(x,y)中距离质心C最近的点N(xn,yn)，xn和yn分别为点N的横坐标和纵坐标；
步骤6、将质心C与步骤5中的点N连成线段CN，将边缘图像e(x,y)中每个像素值为1的边缘点按照位于线段CN的左、右侧进行分类，并将位于线段CN左侧的边缘点组成第一集合SL，将位于线段CN右侧的边缘点组成第二集合SR；
步骤7、分别对第一集合SL和第二集合SR中的所有边缘点进行最小二乘法拟合得到双圆的圆心坐标，分别对应为L(xL,yL)和R(xR,yR)；
步骤8、利用两圆心坐标L(xL,yL)和R(xR,yR)计算出两圆心对应的中轴线l；
步骤9、设两圆心L(xL,yL)和R(xR,yR)的连线与中轴线l的交点为O(xo,yo)，设中轴线l上任意一点n，设该点n距离交点O(xo,yo)的距离为h，设任一圆上的点m，该点m和点n之间的连线平行于两圆心线段，则该点n与点m之间的间距a计算公式为：
Figure BDA0002568975780000021
其中，D为测试块的直径；
步骤10、提取中轴线l上的灰度分布，并对中轴线l上的每个灰度值进行归一化处理，且对中轴线l上每个归一化后的灰度值进行分段线性拟合，获取中轴线l上归一化且线性拟合后的灰度分布曲线，将其设为f(a)；
其中，归一化处理的方法为：将中轴线l上的每个灰度值分别除以材料灰度值均值T；
步骤11、计算
Figure BDA0002568975780000031
设线对卡频率为t，t与a的关系式可表示为/>
Figure BDA0002568975780000032
将/>
Figure BDA0002568975780000033
变换为关于t的函数，即/>
Figure BDA0002568975780000034
步骤12、计算MTF曲线，其表达式为：
Figure BDA0002568975780000035
具体的，所述步骤3中质心C的位置(xc,yc)计算公式为：
Figure BDA0002568975780000036
其中，num为计算满足条件的x或y的数量。
在本方案中，所述步骤4中边缘图像e(x,y)的计算公式为：
e(x,y)＝b(x,y)-b(x,y)ΘB；
其中，Θ表示腐蚀运算，b(x,y)ΘB表示二值化图像b(x,y)被结构元素B所腐蚀，
Figure BDA0002568975780000037
作为改进，所述步骤6中任一个边缘点的分类具体步骤为：
假设边缘图像e(x,y)中某个像素值为1的边缘点记为M(x1,y1)，当点M(x1,y1)使以下公式成立时，则该边缘点M(x1,y1)属于第一集合SL；反之，当点M(x1,y1)使以下公式不成立时，则该边缘点M(x1,y1)属于第二集合SR；
Figure BDA0002568975780000038
进一步的，所述步骤10中材料灰度值均值T获得的具体步骤为：
步骤10-1、在CT图像中，选取两测试块内部区域进行图像灰度直方图统计，设该区域内图像的灰度直方图序列为hf：其中，hf＝M(j)，j为选取的区域图像中灰度级对应的像素值，即灰度值，j∈[0,L-1]，M(j)为选取的区域图像中灰度值为j的像素数量，L-1为最大灰度值；
步骤10-2、对输出序列直方图hf进行高斯拟合，得到拟合后的材料灰度高斯曲线函数Hf(j)，计算公式为：
Figure BDA0002568975780000041
其中，a为材料灰度拟合后高斯曲线的幅值，T为材料灰度拟合后高斯曲线的灰度值均值，σ为材料灰度拟合后高斯曲线的灰度值标准偏差；
步骤10-3、根据上述输出序列直方图hf与拟合后的材料灰度高斯曲线函数Hf(j)之间的对应关系，得到材料灰度拟合后高斯曲线的材料灰度值均值T。
与现有技术相比，本发明的优点在于：该测试方法的测试模体结构简单，加工难度小且成本低，并通过计算获得两测试块之间缝隙的灰度分布，通过自动测量变狭缝工业CT可识别的长度来定量判断检测设备的空间分辨率，其测量精度更高，速度快，且测试结果为尺寸物理量，更为直观。
附图说明
图1为本发明实施例中两测试块截面对应的CT图像；
图2为在图1的CT图像中画出质心位置的示意图；
图3为对图1的CT图像提取边缘图像的示意图；
图4为对图1的CT图像中画出中轴线的示意图；
图5为本发明实施例中测试结果示意图；
图6为本发明实施例中双圆模体示意图。
具体实施方式
以下结合附图实施例对本发明作进一步详细描述。
如图1～4所示，一种光学系统空间分辨率的自动测试方法，用于自动测量线阵工业CT系统或面阵工业CT系统中的空间分辨率，包括以下步骤：
步骤1、通过机械化加工手段制造待测试模体；其中，待测试模体为两个材质与尺寸均相同的圆柱形或球形测试块；
步骤2、将步骤1中的两测试块贴合在一起，并对贴合在一起的两测试块进行CT扫描，获得两测试块截面对应的CT图像；
步骤3、利用自动阈值分割方法将步骤2中的CT图像转换为二值化图像b(x,y)，并计算出二值化图像b(x,y)中质心C的位置(xc,yc)，其中，x为该二值化图像区域的横坐标，y为该二值化图像区域的纵坐标，该二值化图像区域的尺寸为m×n像素，xc和yc分别对应为质心C对应在二值化图像中的横坐标和纵坐标；
质心C的位置(xc,yc)计算公式为：
Figure BDA0002568975780000051
Figure BDA0002568975780000052
其中，num为计算满足条件的x或y的数量；
步骤4、对二值化图像b(x,y)进行形态学处理，得到边缘图像e(x,y)；
通过先对二值化图像b(x,y)进行形态学腐蚀计算，并利用二值化图像b(x,y)减去腐蚀后的图像，则得到该边缘图像e(x,y)；
其计算公式为：
e(x,y)＝b(x,y)-b(x,y)ΘB；
其中，Θ表示腐蚀运算，b(x,y)ΘB表示二值化图像b(x,y)被结构元素B所腐蚀，
Figure BDA0002568975780000053
步骤5、计算出边缘图像e(x,y)中距离质心C最近的点N(xn,yn)，xn和yn分别为点N的横坐标和纵坐标；
计算流程为：计算边缘图像e(x,y)中每个边缘点与质心C的距离，并提取出距离最小的点，则该点为点N；
步骤6、将质心C与步骤5中的点N连成线段CN，将边缘图像e(x,y)中每个像素值为1的边缘点按照位于线段CN的左、右侧进行分类，并将位于线段CN左侧的边缘点组成第一集合SL，将位于线段CN右侧的边缘点组成第二集合SR；
其中，任一个边缘点的分类具体步骤为：
假设边缘图像e(x,y)中某个像素值为1的边缘点记为M(x1,y1)，当点M(x1,y1)使以下公式成立时，则该边缘点M(x1,y1)属于第一集合SL；反之，当点M(x1,y1)使以下公式不成立时，则该边缘点M(x1,y1)属于第二集合SR；
Figure BDA0002568975780000061
并通过使用上述相同的方法依次对边缘图像e(x,y)中的每个边缘点均进行分类，从而能得到第一集合SL和第二集合SR；
步骤7、分别对第一集合SL和第二集合SR中的所有边缘点进行最小二乘法拟合得到双圆的圆心坐标，分别对应为L(xL,yL)和R(xR,yR)；
步骤8、利用两圆心坐标L(xL,yL)和R(xR,yR)计算出两圆心对应的中轴线l；
中轴线l上的点坐标(xl,yl)可表示为：
Figure BDA0002568975780000062
步骤9、设两圆心L(xL,yL)和R(xR,yR)的连线与中轴线l的交点为O(xo,yo)，设中轴线l上任意一点n，设该点n距离交点O(xo,yo)的距离为h(单位为毫米)，设任一圆上的点m，该点m和点n之间的连线平行于两圆心线段，则该点n与点m之间的间距a(单位为毫米)计算公式为：
Figure BDA0002568975780000063
其中，D为测试块的直径(单位为毫米)；
步骤10、提取中轴线l上的灰度分布，并对中轴线l上的每个灰度值进行归一化处理，且对中轴线l上每个归一化后的灰度值进行分段线性拟合，获取中轴线l上归一化且线性拟合后的灰度分布曲线，将其设为f(a)；
其中，归一化处理的方法为：将中轴线l上的每个灰度值分别除以材料灰度值均值T；
步骤11、计算
Figure BDA0002568975780000071
设线对卡频率为t，t与a的关系式可表示为/>
Figure BDA0002568975780000072
将/>
Figure BDA0002568975780000073
变换为关于t的函数，即/>
Figure BDA0002568975780000074
步骤12、计算MTF曲线，其表达式为：
Figure BDA0002568975780000075
本实施例中，材料灰度值均值T获得的具体步骤为：
步骤10-1、在CT图像中，选取两测试块内部区域进行图像灰度直方图统计，设该区域内图像的灰度直方图序列为hf：其中，hf＝M(j)，j为选取的区域图像中灰度级对应的像素值，即灰度值，j∈[0,L-1]，M(j)为选取的区域图像中灰度值为j的像素数量，L-1为最大灰度值；
步骤10-2、对输出序列直方图hf进行高斯拟合，得到拟合后的材料灰度高斯曲线函数Hf(j)，计算公式为：
Figure BDA0002568975780000076
其中，a为材料灰度拟合后高斯曲线的幅值，T为材料灰度拟合后高斯曲线的灰度值均值，σ为材料灰度拟合后高斯曲线的灰度值标准偏差；
步骤10-3、根据上述输出序列直方图hf与拟合后的材料灰度高斯曲线函数Hf(j)之间的对应关系，得到材料灰度拟合后高斯曲线的材料灰度值均值T；
通过灰度直方图统计及该灰度直方图对应的高斯拟合曲线能确定出最优的灰度值均值T，因此该方法能有效减少噪声干扰，使该灰度值均值T的结果更加准确，提高该自动测试结果的准确性。
自动测试方法能得到空间分辨率结果的证明过程如下：
调制度MT可表示为
Figure BDA0002568975780000077
周期性条形模体可表示为矩形脉冲函数，则线对卡模体中波峰p(x)、波谷q(x)可以分别表示为：
Figure BDA0002568975780000081
Figure BDA0002568975780000082
其中：n∈Z，2ai为线对卡模体中波峰(波谷)矩形脉冲宽度，ε(x)为单位阶跃信号，定义为
Figure BDA0002568975780000083
工业CT成像系统是线性非时变系统，成像过程导致对比度下降模型可近似为高斯退化函数模型g(x)，则成像过程简化为输入函数与高斯退化函数卷积：
Figure BDA0002568975780000084
同理，波谷成像后表示为
Figure BDA0002568975780000085
不考虑噪声的情况下，则最大波峰Imax必然在x＝0处：
Figure BDA0002568975780000086
将(6)代入(1)式，即可获得线对卡(5线对)CT成像后调制度MT表达式如下：
Figure BDA0002568975780000087
使用双圆变间距模型进行分析：使用两个具有相同公称直径(D)的校准柱体(球体)，使其在一条线(一个点)上相互点接触，由于CT重建造成接触区域图像畸变，该畸变区域的范围随着CT系统空间分辨率的提高而减小。因此，可将该模体作为一种狭缝宽度连续变化的矩形波测试卡。
如图6所示，两个相邻圆形模体平行于圆心线段的间距可写为：
Figure BDA0002568975780000091
与圆心线段距离为h的灰度分布函数可以表示为
Figure BDA0002568975780000092
式中，ah为与圆心线段距离为h的双圆间波谷宽度。需要注意，f(x,ah)表示归一化的灰度值，且0≤f(x,ah)≤1。
将x＝0代入式(9)可得不同间距下中轴线上的灰度分布，有
Figure BDA0002568975780000093
由于正态分布函数g(x)的积分无解析解，引入g(x)积分近似计算公式
Figure BDA0002568975780000094
将(11)代入(10)式并进行简化，可得
Figure BDA0002568975780000095
其中W31；
根据上述双圆变间距与周期性结构(线对)等效分析：设工业CT成像过程的高斯退化函数模型g(x)为：
Figure BDA0002568975780000096
根据标准正态分布可知，在3倍标准偏差的区间内，g(x)的积分函数约占整体的99.73％。则当线对宽度ai等于标准差
Figure BDA0002568975780000097
时，在3倍标准偏差的区间内包含3组线对，此时的MT值约为1.1％。当1.1％≤MT≤100％时，最大波峰Imax可近似为：
Figure BDA0002568975780000101
将其代入公式(7)，引入g(x)积分近似计算公式并简化，则MT可表示为：
Figure BDA0002568975780000102
因此通过以上证明过程中的公式(12)和公式(13)，可以说明步骤11中通过计算
Figure BDA0002568975780000103
之后再根据行业惯例，将计算空间分辨率的频率用1mm几个线对来表示，其中，一个线对的宽度为4a，则几个线对就是1mm除以4a，因此MIF曲线的线对卡的频率t与a的关系式可表示为：/>
Figure BDA0002568975780000104
可以说明该自动测试方法的测试结果准确，且更加方便测试。
在高能工业CT系统上进行验证该测试方法的准确性，前期研究表明CT系统扫描空间分辨率与微动次数密切相关，微动次数越大空间分辨率越高。本发明中分别采用1、2、5、10微动扫描测试模体，并采用本发明专利方法进行测试空间分辨率，测试结果如图5所示，测试结果及趋势与线对卡基本一致。
以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (5)

1.一种光学系统空间分辨率的自动测试方法，其特征在于：用于自动测量线阵工业CT系统或面阵工业CT系统中的空间分辨率，包括以下步骤：
步骤1、通过机械化加工手段制造待测试模体；其中，待测试模体为两个材质与尺寸均相同的圆柱形或球形测试块；
步骤2、将步骤1中的两测试块贴合在一起，并对贴合在一起的两测试块进行CT扫描，获得两测试块截面对应的CT图像；
步骤3、利用自动阈值分割方法将步骤2中的CT图像转换为二值化图像b(x,y)，并计算出二值化图像b(x,y)中质心C的位置(xc,yc)，其中，x为该二值化图像区域的横坐标，y为该二值化图像区域的纵坐标，该二值化图像区域的尺寸为m×n像素，xc和yc分别对应为质心C对应在二值化图像中的横坐标和纵坐标；
步骤4、对二值化图像b(x,y)进行形态学处理，得到边缘图像e(x,y)；
步骤5、计算出边缘图像e(x,y)中距离质心C最近的点N(xn,yn)，xn和yn分别为点N的横坐标和纵坐标；
步骤6、将质心C与步骤5中的点N连成线段CN，将边缘图像e(x,y)中每个像素值为1的边缘点按照位于线段CN的左、右侧进行分类，并将位于线段CN左侧的边缘点组成第一集合SL，将位于线段CN右侧的边缘点组成第二集合SR；
步骤7、分别对第一集合SL和第二集合SR中的所有边缘点进行最小二乘法拟合得到双圆的圆心坐标，分别对应为L(xL,yL)和R(xR,yR)；
步骤8、利用两圆心坐标L(xL,yL)和R(xR,yR)计算出两圆心对应的中轴线l；
步骤9、设两圆心L(xL,yL)和R(xR,yR)的连线与中轴线l的交点为O(xo,yo)，设中轴线l上任意一点n，设该点n距离交点O(xo,yo)的距离为h，设任一圆上的点m，该点m和点n之间的连线平行于两圆心线段，则该点n与点m之间的间距a计算公式为：
Figure FDA0002568975770000011
其中，D为测试块的直径；
步骤10、提取中轴线l上的灰度分布，并对中轴线l上的每个灰度值进行归一化处理，且对中轴线l上每个归一化后的灰度值进行分段线性拟合，获取中轴线l上归一化且线性拟合后的灰度分布曲线，将其设为f(a)；
其中，归一化处理的方法为：将中轴线l上的每个灰度值分别除以材料灰度值均值T；
步骤11、计算
Figure FDA0002568975770000021
设线对卡频率为t，t与a的关系式可表示为/>
Figure FDA0002568975770000022
将
Figure FDA0002568975770000023
变换为关于t的函数，即/>
Figure FDA0002568975770000024
步骤12、计算MTF曲线，其表达式为：
Figure FDA0002568975770000025
2.根据权利要求1所述的光学系统空间分辨率的自动测试方法，其特征在于：所述步骤3中质心C的位置(xc,yc)计算公式为：
Figure FDA0002568975770000026
其中，num为计算满足条件的x或y的数量。
3.根据权利要求1所述的光学系统空间分辨率的自动测试方法，其特征在于：所述步骤4中边缘图像e(x,y)的计算公式为：
e(x,y)＝b(x,y)-b(x,y)ΘB；
其中，Θ表示腐蚀运算，b(x,y)ΘB表示二值化图像b(x,y)被结构元素B所腐蚀，
Figure FDA0002568975770000027
4.根据权利要求1所述的光学系统空间分辨率的自动测试方法，其特征在于：所述步骤6中任一个边缘点的分类具体步骤为：
假设边缘图像e(x,y)中某个像素值为1的边缘点记为M(x1,y1)，当点M(x1,y1)使以下公式成立时，则该边缘点M(x1,y1)属于第一集合SL；反之，当点M(x1,y1)使以下公式不成立时，则该边缘点M(x1,y1)属于第二集合SR；
Figure FDA0002568975770000031
5.根据权利要求1所述的光学系统空间分辨率的自动测试方法，其特征在于：所述步骤10中材料灰度值均值T获得的具体步骤为：
步骤10-1、在CT图像中，选取两测试块内部区域进行图像灰度直方图统计，设该区域内图像的灰度直方图序列为hf：其中，hf＝M(j)，j为选取的区域图像中灰度级对应的像素值，即灰度值，j∈[0,L-1]，M(j)为选取的区域图像中灰度值为j的像素数量，L-1为最大灰度值；
步骤10-2、对输出序列直方图hf进行高斯拟合，得到拟合后的材料灰度高斯曲线函数Hf(j)，计算公式为：
Figure FDA0002568975770000032
其中，a为材料灰度拟合后高斯曲线的幅值，T为材料灰度拟合后高斯曲线的灰度值均值，σ为材料灰度拟合后高斯曲线的灰度值标准偏差；
步骤10-3、根据上述输出序列直方图hf与拟合后的材料灰度高斯曲线函数Hf(j)之间的对应关系，得到材料灰度拟合后高斯曲线的材料灰度值均值T。
CN202010635729.1A 2020-07-03 2020-07-03 一种光学系统空间分辨率的自动测试方法 Active CN111932573B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010635729.1A CN111932573B (zh) 2020-07-03 2020-07-03 一种光学系统空间分辨率的自动测试方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010635729.1A CN111932573B (zh) 2020-07-03 2020-07-03 一种光学系统空间分辨率的自动测试方法

Publications (2)

Publication Number Publication Date
CN111932573A CN111932573A (zh) 2020-11-13
CN111932573B true CN111932573B (zh) 2023-05-23

Family

ID=73312191

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010635729.1A Active CN111932573B (zh) 2020-07-03 2020-07-03 一种光学系统空间分辨率的自动测试方法

Country Status (1)

Country Link
CN (1) CN111932573B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116413290B (zh) * 2023-02-22 2023-11-28 奥影检测科技(上海)有限公司 一种工业ct散射校正方法
CN116879873B (zh) * 2023-09-08 2023-11-24 立臻精密智造(昆山)有限公司 成像质量检测方法、装置、电子设备及存储介质

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110060293A (zh) * 2019-04-24 2019-07-26 中国兵器科学研究院宁波分院 一种ct检测系统的缺陷检出性能极限评估方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DE10126546A1 (de) * 2001-05-30 2002-12-05 Wilfried Donner Vorrichtung und Verfahren zur Ermittlung eines ortsabhängigen Intensitäts- und Farbprofils und/oder Schärfeprofils optischer Linsensysteme

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110060293A (zh) * 2019-04-24 2019-07-26 中国兵器科学研究院宁波分院 一种ct检测系统的缺陷检出性能极限评估方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
陈津平 ; 崔芳芳 ; 邹晶 ; 胡晓东 ; 吕寒玉 ; .显微CT三维分辨率测试方法研究进展.纳米技术与精密工程.2017,(第02期),全文. *

Also Published As

Publication number Publication date
CN111932573A (zh) 2020-11-13

Similar Documents

Publication Publication Date Title
Kirk et al. Fractal parameters and computer image analysis applied to wear particles isolated by ferrography
CN111932573B (zh) 一种光学系统空间分辨率的自动测试方法
CN107796826B (zh) 基于齿心距离曲线分析的微型双联齿轮断齿缺陷检测方法
CN110231394B (zh) 基于交流电磁场的非铁磁性材料不规则裂纹成像方法
CN115018844B (zh) 一种基于人工智能的塑料薄膜质量评估方法
CN112215810B (zh) 一种疲劳测试裂纹监测方法及装置
CN110060293B (zh) 一种ct检测系统的缺陷检出性能极限评估方法
CN109596530B (zh) 光学表面麻点和灰尘缺陷分类的暗场偏振成像方法
CN109655234B (zh) 一种针对于相机的自动化测试方法
CN116109635B (zh) 复合悬式绝缘子表面质量检测方法、装置、设备及介质
CN110243923A (zh) 基于交流电磁场的腐蚀缺陷可视化成像及评估方法
CN105631876A (zh) 一种基于全局二值化的ct图像分辨率自动测试方法
CN110544243B (zh) 一种ct图像小缺陷自动检出、定量及可靠性评价方法
CN112767396B (zh) 缺陷的检测方法、装置和计算机可读存储介质
CN107300562A (zh) 一种测量继电器成品触点间距的x射线无损检测方法
CN116205923B (zh) 基于x-ray的汽车轮毂内部缺陷无损检测方法
CN110081832B (zh) 逻辑投影成像测量非透明物体结构参数的方法及装置
CN103282938A (zh) 用于通过x射线计算机断层扫描装置确定位于待检查物体中的结构位置的方法和评估设备
US7257248B2 (en) Non-contact measurement system and method
CN114549403B (zh) 一种机械零部件侧剖面多单体智能精密几何圆心检测方法
CN111879799B (zh) 一种光学系统空间分辨率的手动测试方法
CN108844977B (zh) 一种剔除角度倾斜影响的工业ct系统空间分辨率测试方法及评价方法
CN113554688B (zh) 一种基于单目视觉的o型密封圈尺寸测量方法
KR101385592B1 (ko) 영상인식 방법 및 그 시스템
CN113222880A (zh) 分析检测方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant