CN112233156B - 微纳米ct投影数据的中心切片对齐方法 - Google Patents

微纳米ct投影数据的中心切片对齐方法 Download PDF

Info

Publication number
CN112233156B
CN112233156B CN202011094025.4A CN202011094025A CN112233156B CN 112233156 B CN112233156 B CN 112233156B CN 202011094025 A CN202011094025 A CN 202011094025A CN 112233156 B CN112233156 B CN 112233156B
Authority
CN
China
Prior art keywords
projection
image
matching
projection image
corrected
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
CN202011094025.4A
Other languages
English (en)
Other versions
CN112233156A (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.)
Lima precision measurement technology (Suzhou) Co.,Ltd.
Original Assignee
Lima Precision Measurement Technology Suzhou Co ltd
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 Lima Precision Measurement Technology Suzhou Co ltd filed Critical Lima Precision Measurement Technology Suzhou Co ltd
Priority to CN202011094025.4A priority Critical patent/CN112233156B/zh
Publication of CN112233156A publication Critical patent/CN112233156A/zh
Application granted granted Critical
Publication of CN112233156B publication Critical patent/CN112233156B/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/30—Determination of transform parameters for the alignment of images, i.e. image registration
    • G06T7/38—Registration of image sequences
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/20—Analysis of motion
    • G06T7/262—Analysis of motion using transform domain methods, e.g. Fourier domain methods
    • 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/10056—Microscopic image
    • G06T2207/10061—Microscopic image from scanning electron microscope
    • 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/10—Image acquisition modality
    • G06T2207/10116—X-ray image

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Databases & Information Systems (AREA)
  • Algebra (AREA)
  • Software Systems (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Multimedia (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Operations Research (AREA)
  • Probability & Statistics with Applications (AREA)
  • Image Analysis (AREA)
  • Apparatus For Radiation Diagnosis (AREA)

Abstract

本发明公开了一种微纳米CT投影数据的中心切片对齐方法，该方法包括如下步骤：步骤1，获取样品在投影角
Figure DDA0002723049570000011
下的待校正投影图像序列为
Figure DDA0002723049570000012
取待校正投影图像序列
Figure DDA0002723049570000013
中的第一待校正投影图像PA和第二待校正投影图像PB，并获取第一待校正投影图像PA对应的第一参考投影图像PAstd以及第二待校正投影图像PB对应第二参考投影图像PBstd；步骤2，将待校正投影图像序列
Figure DDA0002723049570000014
中的第t张待校正投影图像Pt(t＝1，2，...，NA，t≠A，B)按照预设的搜索角序列
Figure DDA0002723049570000015
中的每一个角度进行旋转，得到一组旋转后的投影图像序列
Figure DDA0002723049570000016
再与第一参考投影图像PAstd进行二次投影的频域匹配校正。本发明能够在迭代重投影的基础上利用邻近匹配进行修正，从而加快了迭代速度，有利于提高计算效率。

Description

微纳米CT投影数据的中心切片对齐方法
技术领域
本发明涉及X射线纳米CT技术领域，特别是关于一种微纳米CT投影数据的中心切片对齐方法。
背景技术
X射线纳米CT技术是指由X射线透视成像技术与纳米显微镜技术结合的三维成像技术，其空间分辨率可达到10纳米左右。在生命科学、化学以及材料科学等领域有着广阔应用。
然而，正是因为纳米CT的高精度特性，纳米CT在运转中对各元器件的精度有着较高的要求，实际中往往不能满足这些要求，使得成像的图像质量降低，进而使得纳米CT无法大范围推广使用。其中比较关键的一个问题在于图像分辨率的降低，常用于纳米CT的转台在运转过程中，往往会发生100纳米以上的径跳和端跳，1角度以上的摇摆，这使得纳米CT所测量样品的空间分辨率难以维持在10纳米左右。为了提高纳米CT的空间分辨率，我们需要对数据作适当的校正。
目前国内常用的校正方法主要有：
(1)硬件上利用电容传感器等精密仪器进行检测校正。但是，电容传感器价格昂贵，而且需要恒温恒湿等较好的测试环境，成本较为高昂，CCD测量法需要稳定长效的无指向的散射光光源，同时无法得到转轴的轴心轨迹。
(2)软件上的图像处理算法进行计算校正。但是，迭代重投影的方法只能校正径跳的误差，而且在一些转台径跳变化剧烈时收敛速度变慢。中心切片的空间域配准方法计算量大，导致了计算效率低。提取投影数据的特征点要求投影数据具有良好的对比度，否则存在较大的误差。而纳米CT数据，比如细胞数据，往往不具有良好的对比度，因此该方法的抗噪性不足。
发明内容
本发明的目的在于提供一种微纳米CT投影数据的中心切片对齐方法来克服或至少减轻现有技术的上述缺陷中的至少一个。
为实现上述目的，本发明提供一种微纳米CT投影数据的中心切片对齐方法，该方法包括如下步骤：
步骤1，获取样品在投影角
Figure GDA0003441648180000021
下的待校正投影图像序列为
Figure GDA0003441648180000022
取所述待校正投影图像序列
Figure GDA0003441648180000023
中的第一待校正投影图像PA和第二待校正投影图像PB，并获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd；
步骤2，将所述待校正投影图像序列
Figure GDA0003441648180000024
中的第t张待校正投影图像Pt(t＝1,2,…,NA,t≠A,B)按照预设的搜索角序列
Figure GDA0003441648180000025
中的每一个角度进行旋转，得到一组旋转后的投影图像序列
Figure GDA0003441648180000026
再与所述第一参考投影图像PAstd按照如下子步骤进行中心切片频域匹配校正：
步骤21，计算第一参考投影图像PAstd和投影图像序列
Figure GDA0003441648180000027
的二维傅里叶变换，分别得到复数第一参考投影图像f(PAstd)、复数图像序列
Figure GDA0003441648180000028
分别截取f(PAstd)和
Figure GDA0003441648180000029
中每张复数图像的中心复数值f(PAstd)mid、值
Figure GDA00034416481800000210
步骤22，按照下式(1)计算投影图像Ptk到PAstd的匹配移动量mtk,将Ptk向y轴负向移动mtk个单位，得到P′tk，计算P′tk和PAstd的均方根误差errk,称之为投影图像Ptk和PAstd的匹配误差，PAstd与
Figure GDA00034416481800000211
之间的匹配误差则为
Figure GDA00034416481800000212
其中i为虚数单位：
Figure GDA00034416481800000213
步骤23，获取匹配误差
Figure GDA00034416481800000214
中的最小值索引kmin，得到摇摆校正角
Figure GDA00034416481800000215
和校正摇摆角后的第t张投影图像
Figure GDA00034416481800000216
步骤24，利用中心切片频域匹配获取
Figure GDA00034416481800000217
到PAstd的匹配移动量
Figure GDA00034416481800000218
记为
Figure GDA00034416481800000219
步骤25，将投影图像
Figure GDA00034416481800000220
向上平移
Figure GDA00034416481800000221
得到第t张投影图像的校正端跳误差后的投影图像P′t,std；
步骤26，取不同的待校正投影图像，返回步骤2，直至校正完毕所述待校正投影图像序列
Figure GDA0003441648180000031
中所有第t张待校正投影图像的摇摆和端跳，得到校正摇摆和端跳误差后的投影图像
Figure GDA0003441648180000032
进一步地，步骤24中的“频域匹配”为“中心切片频域匹配”，其具体包括：
Figure GDA0003441648180000033
和PAstd分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N,分别对
Figure GDA0003441648180000034
和PAstd进行二维傅里叶变换，得到两个二维的复值图像
Figure GDA0003441648180000035
和f(PAstd)，
Figure GDA0003441648180000036
和f(PAstd)的中心项即第(M/2，N/2)项为低频项，分别记为
Figure GDA0003441648180000037
和f(PAstd)mid，利用下式(2)计算y轴方向的匹配移动量
Figure GDA0003441648180000038
Figure GDA0003441648180000039
进一步地，该方法还包括：
步骤3，将步骤2获得的
Figure GDA00034416481800000310
记为校正摇摆和端跳误差后的迭代第1次后的投影图像
Figure GDA00034416481800000311
按照如下子步骤，利用迭代重投影和邻近匹配修正技术校正径跳误差：
步骤31，设定最大迭代次数Niter,令iter＝0；
步骤32，令iter＝iter+1，对校正摇摆和端跳误差后的迭代第iter次后的投影图像
Figure GDA00034416481800000312
进行邻近匹配修正，修正后的结果记为
Figure GDA00034416481800000313
步骤33，对投影图像
Figure GDA00034416481800000314
直接进行迭代重建，得到重建图像
Figure GDA00034416481800000315
步骤34，将
Figure GDA00034416481800000316
进行投影角
Figure GDA00034416481800000317
的投影，得到第iter次迭代重投影中重建图
Figure GDA00034416481800000318
的投影图像序列
Figure GDA00034416481800000319
步骤35，将不同投影角下的投影图像
Figure GDA00034416481800000320
和
Figure GDA00034416481800000321
的各行进行频域匹配，找出
Figure GDA0003441648180000041
到
Figure GDA0003441648180000042
的各行匹配移动量的平均值
Figure GDA0003441648180000043
并将投影图像
Figure GDA0003441648180000044
向右移动
Figure GDA0003441648180000045
得到投影图像
Figure GDA0003441648180000046
步骤36，如果iter<Niter，进入步骤32，否则进行步骤37；
步骤37，迭代终止，得到径跳校正完成后的投影图像
Figure GDA0003441648180000047
进一步地，步骤32中的“
Figure GDA0003441648180000048
进行邻近匹配修正”具体包括：
针对按照投影角度先后次序排列的投影图像
Figure GDA0003441648180000049
从第一投影图像
Figure GDA00034416481800000410
和第二投影图像
Figure GDA00034416481800000411
开始，每相邻两张投影图像对应的行数据进行x方向进行频域匹配，找出每行的匹配移动量，并计算其均值，将其中的前一张图像沿着x方向按平均的移动匹配量进行移动，直到遍历完所有图像。
进一步地，“每相邻两张投影图像对应的行数据进行x方向进行频域匹配”为第t张投影图像
Figure GDA00034416481800000412
和第t+1张投影图像
Figure GDA00034416481800000413
对应的行数据进行x方向进行频域匹配，其方法具体包括：
Figure GDA00034416481800000414
和
Figure GDA00034416481800000415
中每行对应的数组分别为x轴方向的两个一维N元实数数组，记为
Figure GDA00034416481800000416
和
Figure GDA00034416481800000417
其坐标均为{1,2,…,N},row＝1,2,…,Nr，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure GDA00034416481800000418
和
Figure GDA00034416481800000419
Figure GDA00034416481800000420
和
Figure GDA00034416481800000421
的中间项即第
Figure GDA00034416481800000422
项为低频项，分别记为
Figure GDA00034416481800000423
和
Figure GDA00034416481800000424
利用下式(3)和(4)计算二者的匹配移动量
Figure GDA00034416481800000425
Figure GDA0003441648180000051
Figure GDA0003441648180000052
进一步地，“修正后的结果记为
Figure GDA0003441648180000053
”的方法具体包括：
步骤321，t＝0；
步骤322，令t＝t+1，若t≥NA,则转到步骤326，否则进入步骤323；
步骤323，将图像
Figure GDA0003441648180000054
向x轴正方向移动
Figure GDA0003441648180000055
个单位，得到
Figure GDA0003441648180000056
步骤324，将图像
Figure GDA0003441648180000057
向x轴负方向移动
Figure GDA0003441648180000058
个单位，得到
Figure GDA0003441648180000059
步骤325，返回步骤322；
步骤326，修正完成，得到
Figure GDA00034416481800000510
进一步地，步骤35中的“将不同投影角下的投影图像
Figure GDA00034416481800000511
和
Figure GDA00034416481800000512
的各行进行频域匹配”为“将不同投影角下的投影图像
Figure GDA00034416481800000513
中的每一幅投影图像
Figure GDA00034416481800000514
和
Figure GDA00034416481800000515
中的每一幅投影图像
Figure GDA00034416481800000516
的相对应行进行频域匹配”的方法具体包括：
Figure GDA00034416481800000517
和
Figure GDA00034416481800000518
分别为投影图像
Figure GDA00034416481800000519
和
Figure GDA00034416481800000520
的第row行，row＝1,2,…,Nr，也就是x轴方向的两个一维N元实数数组，其坐标均为{1,2,…,N},分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure GDA00034416481800000521
和
Figure GDA00034416481800000522
Figure GDA00034416481800000523
和
Figure GDA00034416481800000524
的中间项即第N/2项为低频项，分别记为
Figure GDA0003441648180000061
和
Figure GDA0003441648180000062
利用下式(5)(6)计算二者的匹配移动量
Figure GDA0003441648180000063
Figure GDA0003441648180000064
Figure GDA0003441648180000065
进一步地，根据式(5)和(6)计算出来的
Figure GDA0003441648180000066
计算得到
Figure GDA0003441648180000067
到
Figure GDA0003441648180000068
的各行匹配移动量的平均值
Figure GDA0003441648180000069
的方法具体包括：
记第t组投影图像
Figure GDA00034416481800000610
到
Figure GDA00034416481800000611
的各行匹配移动量为
Figure GDA00034416481800000612
则各行匹配移动量的均值
Figure GDA00034416481800000613
依次计算t＝1,2,…,NA时对应的
Figure GDA00034416481800000614
得到
Figure GDA00034416481800000615
进一步地，步骤1中的“获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd”的方法具体包括：
步骤11，将PA,PB分别按照预设的搜索角序列
Figure GDA00034416481800000616
中的每一个角度进行旋转，得到两组旋转后的投影图像序列
Figure GDA00034416481800000617
步骤12，分别计算投影图像序列
Figure GDA00034416481800000618
和
Figure GDA00034416481800000619
的二维傅里叶变换，得到复数图像序列
Figure GDA00034416481800000620
分别截取复数图像序列
Figure GDA00034416481800000621
Figure GDA00034416481800000622
中每张复数图像中心的复数值
Figure GDA00034416481800000623
步骤13，按照下式计算
Figure GDA0003441648180000071
和
Figure GDA0003441648180000072
两两之间的匹配移动量mAi,Bj，将PAi向y轴正向移动mAi,Bj个单位，得到P′Ai，计算P′Ai和PBj的均方根误差erri,j,称之为投影图像PAi和PBj的匹配误差，
Figure GDA0003441648180000073
与
Figure GDA0003441648180000074
两两之间的匹配误差则为
Figure GDA0003441648180000075
其中i为虚数单位：
Figure GDA0003441648180000076
步骤14，找出匹配误差
Figure GDA0003441648180000077
中最小值的索引imin,jmin；
步骤15，利用频域匹配得到
Figure GDA0003441648180000078
到
Figure GDA0003441648180000079
的匹配移动量，记为
Figure GDA00034416481800000710
步骤16，将投影图像
Figure GDA00034416481800000711
向上平移
Figure GDA00034416481800000712
得到所述第一参考投影图像PAstd，投影图像
Figure GDA00034416481800000713
向下平移
Figure GDA00034416481800000714
得到第二参考投影图像PBstd。
进一步地，步骤13具体包括：
步骤131，i＝0
步骤132，i＝i+1，j＝0，若i>Ns，转到步骤135；
步骤133，(3)j＝j+1，若j>Ns，转到步骤132；
步骤134，对PAi和PBj进行中心切片频域匹配，并计算匹配误差errij，转到步骤133；步骤135，两两频域匹配过程完成，得到
Figure GDA00034416481800000715
进一步地，步骤15中的“利用中心切片频域匹配计算中心切片
Figure GDA00034416481800000716
到
Figure GDA00034416481800000717
的匹配移动量
Figure GDA00034416481800000718
”的方法具体包括：
Figure GDA00034416481800000719
和
Figure GDA00034416481800000720
分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N,分别对
Figure GDA00034416481800000721
和
Figure GDA00034416481800000722
进行二维傅里叶变换，得到两个二维的复值图像
Figure GDA00034416481800000723
和
Figure GDA00034416481800000724
Figure GDA00034416481800000725
和
Figure GDA00034416481800000726
的中心项即第(M/2，N/2)项为低频项，分别记为
Figure GDA00034416481800000727
和
Figure GDA00034416481800000728
利用下式(7)计算y轴方向的匹配移动量
Figure GDA00034416481800000729
Figure GDA0003441648180000081
进一步地，步骤134中的“对PAi和PBj进行中心切片频域匹配，并计算匹配误差errij”的方法具体包括：
PAi和PBj分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N，分别对PAi和PBj进行二维傅里叶变换，得到两个二维的复值图像f(PAi)和f(PBj)，f(PAi)和f(PBj)的中心项即第(M/2，N/2)项为低频项，分别记为f(PAi)mid和f(PBj)mid，利用下式(2)计算y轴方向的匹配移动量mAi,Bj：
Figure GDA0003441648180000082
将PAi向y轴正向移动mAi,Bj个单位，得到P′Ai，计算P′Ai和PBj的均方根误差，记为erri，j,称之为投影图像PAi和PBj的匹配误差。
本发明所涉及的校正方法属于软件层面上的方法，在估计端跳和摇摆信息方面，本发明采用中心切片的频域匹配方法，通过利用图像变换在频域的数学性质，为加快匹配的速度和精度，进而为提高计算效率提供了有利条件。在估计径跳信息方面，本发明提高采用加邻近匹配的迭代重投影方法在迭代重投影的基础上利用邻近匹配进行修正，从而加快了迭代速度，有利于提高计算效率。利用这些信息实现对投影数据的校正，重建出空间分辨率更高的图像。
附图说明
图1的最左侧图是四棱柱形状样品S自虚线区域示意的初始位置到在y轴方向上存在端跳
Figure GDA0003441648180000083
后的灰色填充区域示意的待校正位置示意图，中间图是两位置样品在扫描角度A下各自的投影图像PA。
图2为最左侧图是四棱柱形状样品S自虚线区域示意的初始位置到在y轴方向上存在端跳
Figure GDA0003441648180000084
后的灰色填充区域示意的待校正位置示意图，中间图是两位置样品在扫描角度A+90度下各自的投影图像PB。
图3为本发明频域匹配后所得的匹配移动量mAi，Bj，沿着y轴方向将投影PA和PB平移的端跳校正过程示意图。
具体实施方式
下面结合附图和实施例对本发明进行详细的描述。
如图1至图2所示，为了叙述简便，下面以长方体的样品S的中心切片原理为例，说明本发明实施例所在的XYZ坐标系，将样品S的一个表面置于纸面所在的平面内，视为XY平面，X轴与Y轴垂直。其中，X轴的正向如图中的中心切片到X轴上LB的延伸方向所示，即在纸面内向右的方向。Y轴在纸面内向上的方向为下文中的“向上”，反之为“向下”。Z轴垂直于XY平面，正方向为指向纸面内地方向，Z轴正方向为射线扫描方向。端跳方向定义为沿Y轴方向的跳动，即在Y轴方向的平移量。
本发明实施例提供的微纳米CT投影数据的中心切片对齐方法包括如下步骤：
步骤1，获取样品在投影角
Figure GDA0003441648180000091
下的待校正投影图像序列为
Figure GDA0003441648180000092
取所述待校正投影图像序列
Figure GDA0003441648180000093
中的第一待校正投影图像PA和第二待校正投影图像PB，并获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd。其中，任取待校正投影图像序列为
Figure GDA0003441648180000094
中，间隔90°投影角的一组投影图像，分别取PA,PB。
在一个实施例中，步骤1中，“获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd”的方法具体包括：
步骤11，将PA,PB分别按照预设的搜索角序列
Figure GDA0003441648180000095
中的每一个角度进行旋转，得到两组旋转后的投影图像序列
Figure GDA0003441648180000096
步骤12，分别计算投影图像序列
Figure GDA0003441648180000097
和
Figure GDA0003441648180000098
的二维傅里叶变换，得到复数图像序列
Figure GDA0003441648180000099
分别截取复数图像序列
Figure GDA00034416481800000910
Figure GDA00034416481800000911
中每张复数图像中心的复数值
Figure GDA00034416481800000912
步骤13，按照下式计算
Figure GDA0003441648180000101
和
Figure GDA0003441648180000102
两两之间的匹配移动量mAi,Bj，将PAi向y轴正向移动mAi,Bj个单位，得到P′Ai，计算P′Ai和PBj的均方根误差erri,j,称之为投影图像PAi和PBj的匹配误差，
Figure GDA0003441648180000103
与
Figure GDA0003441648180000104
两两之间的匹配误差则为
Figure GDA0003441648180000105
其中i为虚数单位：
Figure GDA0003441648180000106
步骤14，找出匹配误差
Figure GDA0003441648180000107
中最小值的索引imin,jmin。
步骤15，利用频域匹配得到
Figure GDA0003441648180000108
到
Figure GDA0003441648180000109
的匹配移动量,记为
Figure GDA00034416481800001010
步骤16，将投影图像
Figure GDA00034416481800001011
向上平移
Figure GDA00034416481800001012
得到所述第一参考投影图像PAstd，投影图像
Figure GDA00034416481800001013
向下平移
Figure GDA00034416481800001014
得到第二参考投影图像PBstd。
在一个实施例中，步骤13中的“两两”可以理解为：每固定第一个序列中的一个信号，第二个序列信号取遍，需要匹配Ns次，第一个序列要取Ns次信号固定，所以一共进行了
Figure GDA00034416481800001015
次匹配，其具体子步骤包括：
步骤131，i＝0；
步骤132，i＝i+1，j＝0，若i>Ns，转到步骤135；
步骤133，(3)j＝j+1，若j>Ns，转到步骤132；
步骤134，对PAi和PBj进行频域匹配，并计算匹配误差errij，转到步骤133；
步骤135，两两频域匹配过程完成，得到
Figure GDA00034416481800001016
进一步地，步骤134中的“PAi和PBj进行频域匹配”的方法包括：
PAi和PBj为y轴有两个一维N元实数数组，其坐标均为{1，2，…，N}，分别进行离散傅里叶变换，得到两个一维N元复数数组f(PAi)和f(PBj)。f(PAi)和f(PBj)的中间项即第N/2项，为低频项，分别记为f(PAi)mid和f(PBj)mid，计算二者的相角差
Figure GDA0003441648180000111
将PAi数组沿着y轴正向移动mAi，Bj个单位，仍取坐标{1,2,…,N}上的数所构成的数组，记为PAi′,PBj′与PBj的近似程度最高，PAi′和PBj达成频域匹配，PAi′和PBj之间的均方根误差RMSE(PAi′,PBj)，记为匹配误差ErrAi，Bj。
在一个实施例中，步骤15中的“利用中心切片频域匹配计算中心切片
Figure GDA0003441648180000112
到
Figure GDA0003441648180000113
的匹配移动量
Figure GDA0003441648180000114
”的方法具体包括；
Figure GDA0003441648180000115
和
Figure GDA0003441648180000116
分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N,分别对
Figure GDA0003441648180000117
和
Figure GDA0003441648180000118
进行二维傅里叶变换，得到两个二维的复值图像
Figure GDA0003441648180000119
和
Figure GDA00034416481800001110
Figure GDA00034416481800001111
和
Figure GDA00034416481800001112
的中心项即第(M/2，N/2)项为低频项，分别记为
Figure GDA00034416481800001113
和
Figure GDA00034416481800001114
利用下式(7)计算y轴方向的匹配移动量
Figure GDA00034416481800001115
Figure GDA00034416481800001116
步骤2，将所述待校正投影图像序列
Figure GDA00034416481800001117
中的第t张待校正投影图像为Pt(t＝1,2,…,NA,t≠A,B)与所述第一参考投影图像PAstd按照如下子步骤进行中心切片的频域匹配校正：
步骤21，计算第一参考投影图像PAstd和投影图像序列
Figure GDA00034416481800001118
的二维傅里叶变换，分别得到复数第一参考投影图像f(PAstd)、复数图像序列
Figure GDA00034416481800001119
分别截取f(PAstd)和
Figure GDA00034416481800001120
中每张复数图像的中心复数值f(PAstd)mid、值
Figure GDA00034416481800001121
步骤22，按照下式计算投影图像Ptk到PAstd的匹配移动量mtk,将Ptk向y轴负向移动mtk个单位，得到P′tk，计算P′Ai和PAstd的均方根误差errk,称之为投影图像Ptk和PAstd的匹配误差，PAstd与
Figure GDA00034416481800001122
之间的匹配误差则为
Figure GDA00034416481800001123
其中i为虚数单位：
Figure GDA00034416481800001124
步骤23，找出匹配误差
Figure GDA00034416481800001125
中的最小值索引kmin，得到摇摆校正角
Figure GDA0003441648180000121
和校正摇摆角后的第t张投影图像
Figure GDA0003441648180000122
步骤24，利用频域匹配计算
Figure GDA0003441648180000123
到PAstd的匹配移动量
Figure GDA0003441648180000124
步骤25，将投影图像
Figure GDA0003441648180000125
向上平移
Figure GDA0003441648180000126
得到第t张投影图像的校正端跳误差后的投影图像P′t,std。
步骤26，取不同的待校正投影图像，返回步骤2，直至校正完毕所述待校正投影图像序列
Figure GDA0003441648180000127
中所有第t张待校正投影图像的摇摆和端跳，得到校正摇摆和端跳误差后的投影图像
Figure GDA0003441648180000128
在一个实施例中，步骤22具体包括：
步骤221，
Figure GDA0003441648180000129
与PAstd分别为y轴方向的两个一维N元实数数组，其坐标均为{1,2,…,N},分别进行离散傅里叶变换，得到两个一维N元复数数组f(Ptk)和f(PAstd)，f(Ptk)和f(PAstd)的中间项即第N/2项为低频项，分别记为f(Ptk)mid和f(PAstd)mid，利用下式(1)计算二者的相角差mtk。
Figure GDA00034416481800001210
步骤222，将
Figure GDA00034416481800001211
数组沿y轴方向正向移动
Figure GDA00034416481800001212
个单位，仍取坐标{1,2,…,N}上的数所构成的数组P′tk，得到数组序列
Figure GDA00034416481800001213
步骤223，计算
Figure GDA00034416481800001214
与PAstd之间的均方根误差，得到匹配误差
Figure GDA00034416481800001215
在一个实施例中，步骤24中的“频域匹配”为“中心切片频域匹配”，其具体包括：
Figure GDA00034416481800001216
和PAstd分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N,分别对
Figure GDA00034416481800001217
和PAstd进行二维傅里叶变换，得到两个二维的复值图像
Figure GDA00034416481800001218
和f(PAstd)，
Figure GDA00034416481800001219
和f(PAstd)的中心项即第(M/2，N/2)项为低频项，分别记为
Figure GDA00034416481800001220
和f(PAstd)mid，利用下式(2)计算二者的匹配移动量
Figure GDA00034416481800001221
Figure GDA0003441648180000131
在一个实施例中，本发明实施例提供的微纳米CT投影数据的中心切片对齐方法还包括：
步骤3，将步骤2获得的
Figure GDA0003441648180000132
记为校正摇摆和端跳误差后的迭代第1次后的投影图像
Figure GDA0003441648180000133
按照如下子步骤，利用迭代重投影和邻近匹配修正技术校正径跳误差：
步骤31，设定最大迭代次数Niter,令iter＝0。
步骤32，令iter＝iter+1，对校正摇摆和端跳误差后的迭代第iter次后的投影图像
Figure GDA0003441648180000134
进行邻近匹配修正，修正后的结果记为
Figure GDA0003441648180000135
步骤33，对投影图像
Figure GDA0003441648180000136
直接进行迭代重建，得到重建图像
Figure GDA0003441648180000137
步骤34，将
Figure GDA0003441648180000138
进行投影角
Figure GDA0003441648180000139
的投影，得到第iter次迭代重投影中重建图
Figure GDA00034416481800001310
的投影图像序列
Figure GDA00034416481800001311
步骤35，将不同投影角下的投影图像
Figure GDA00034416481800001312
和
Figure GDA00034416481800001313
的各行进行频域匹配，找出
Figure GDA00034416481800001314
到
Figure GDA00034416481800001315
的各行匹配移动量的平均值
Figure GDA00034416481800001316
并将投影图像
Figure GDA00034416481800001317
向右移动
Figure GDA00034416481800001318
得到投影图像
Figure GDA00034416481800001319
步骤36，如果iter<Niter，进入步骤32，否则进行步骤37。
步骤37，迭代终止，得到径跳校正完成后的投影图像
Figure GDA00034416481800001320
在一个实施例中，步骤32中的“
Figure GDA0003441648180000141
进行邻近匹配修正”具体包括：针对按照投影角度先后次序排列的投影图像
Figure GDA0003441648180000142
从第一投影图像
Figure GDA0003441648180000143
和第二投影图像
Figure GDA0003441648180000144
开始，每相邻两张投影图像对应的行数据进行x方向进行频域匹配，找出每行的匹配移动量，并计算其均值，将其中的前一张图像沿着x方向按平均的移动匹配量进行移动，直到遍历完所有图像。
在一个实施例中，步骤32中的“每相邻两张投影图像对应的行数据进行x方向进行频域匹配移动量的均值计算方法包括：
记两张投影图像各行的匹配移动量为
Figure GDA0003441648180000145
则匹配移动量的均值
Figure GDA0003441648180000146
表示为下式：
Figure GDA0003441648180000147
在一个实施例中，“每相邻两张投影图像对应的行数据进行x方向进行频域匹配”为第t张投影图像
Figure GDA0003441648180000148
和第t+1张投影图像
Figure GDA0003441648180000149
对应的行数据进行x方向进行频域匹配，其方法具体包括：
Figure GDA00034416481800001410
和
Figure GDA00034416481800001411
中每行对应的数组分别为x轴方向的两个一维N元实数数组，记为
Figure GDA00034416481800001412
和
Figure GDA00034416481800001413
其坐标均为{1,2,…,N},row＝1,2,…,Nr，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure GDA00034416481800001414
和
Figure GDA00034416481800001415
Figure GDA00034416481800001416
和
Figure GDA00034416481800001417
的中间项即第N/2项为低频项，分别记为
Figure GDA00034416481800001418
和
Figure GDA00034416481800001419
利用下式(3)和(4)计算二者的匹配移动量
Figure GDA00034416481800001420
Figure GDA00034416481800001421
Figure GDA00034416481800001422
下面以第一投影图像
Figure GDA0003441648180000151
和第二投影图像
Figure GDA0003441648180000152
为例，说明每相邻两张投影图像对应的行数据进行x方向进行频域匹配的方法：
Figure GDA0003441648180000153
和
Figure GDA0003441648180000154
中每行对应的数组分别为x轴方向的两个一维N元实数数组，记为
Figure GDA0003441648180000155
和
Figure GDA0003441648180000156
其坐标均为{1,2,…,N},row＝1,2,…,Nr，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure GDA0003441648180000157
和
Figure GDA0003441648180000158
Figure GDA0003441648180000159
和
Figure GDA00034416481800001510
的中间项即第N/2项为低频项，分别记为
Figure GDA00034416481800001511
和
Figure GDA00034416481800001512
利用下式计算二者的匹配移动量
Figure GDA00034416481800001513
Figure GDA00034416481800001514
Figure GDA00034416481800001515
在一个实施例中，步骤32中的“修正后的结果记为
Figure GDA00034416481800001516
”的方法具体包括：
步骤321，t＝0；
步骤322，令t＝t+1，若t≥NA,则转到步骤326，否则进入步骤323；
步骤323，将图像
Figure GDA00034416481800001517
向x轴正方向移动
Figure GDA00034416481800001518
个单位，得到
Figure GDA00034416481800001519
步骤324，将图像
Figure GDA00034416481800001520
向x轴负方向移动
Figure GDA00034416481800001521
个单位，得到
Figure GDA00034416481800001522
步骤325，返回步骤322；
步骤326，修正完成，得到
Figure GDA00034416481800001523
在一个实施例中，步骤35中的“将不同投影角下的投影图像
Figure GDA0003441648180000161
和
Figure GDA0003441648180000162
的各行进行频域匹配”为“将不同投影角下的投影图像
Figure GDA0003441648180000163
中的每一幅投影图像
Figure GDA0003441648180000164
和
Figure GDA0003441648180000165
中的每一幅投影图像
Figure GDA0003441648180000166
的相对应行进行频域匹配”的方法具体包括：
Figure GDA0003441648180000167
和
Figure GDA0003441648180000168
分别为投影图像
Figure GDA0003441648180000169
和
Figure GDA00034416481800001610
的第row行，row＝1,2,…,Nrow，也就是x轴方向的两个一维N元实数数组，其坐标均为{1,2,…,N},分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure GDA00034416481800001611
和
Figure GDA00034416481800001612
Figure GDA00034416481800001613
和
Figure GDA00034416481800001614
的中间项即第N/2项为低频项，分别记为
Figure GDA00034416481800001615
和
Figure GDA00034416481800001616
利用下式(5)和(6)计算二者的匹配移动量
Figure GDA00034416481800001617
Figure GDA00034416481800001618
Figure GDA00034416481800001619
在一个实施例中，根据式(5)和(6)计算出来的
Figure GDA00034416481800001620
计算得到
Figure GDA00034416481800001621
到
Figure GDA00034416481800001622
的各行匹配移动量的平均值
Figure GDA00034416481800001623
的方法具体包括：
记第t组投影图像
Figure GDA00034416481800001624
到
Figure GDA00034416481800001625
的各行匹配移动量为
Figure GDA00034416481800001626
则各行匹配移动量的均值
Figure GDA00034416481800001627
依次计算t＝1,2,…,NA时对应的
Figure GDA00034416481800001628
得到
Figure GDA00034416481800001629
本发明提供的校正方法，相比硬件校正方法，易于操作推广，可以节省大量的硬件成本。相比现有的软件算法校正方法，校正精度更高，计算速度更快。
最后需要指出的是：以上实施例仅用以说明本发明的技术方案，而非对其限制。本领域的普通技术人员应当理解：可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (9)

1.一种微纳米CT投影数据的中心切片对齐方法，其特征在于，包括如下步骤：
步骤1，获取样品在投影角
Figure FDA0003423321550000011
下的待校正投影图像序列为
Figure FDA0003423321550000012
取所述待校正投影图像序列
Figure FDA0003423321550000013
中的第一待校正投影图像PA和第二待校正投影图像PB，并获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd；
步骤1中的“获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd”的方法具体包括：
步骤11，将PA,PB分别按照预设的搜索角序列
Figure FDA0003423321550000014
中的每一个角度进行旋转，得到两组旋转后的投影图像序列
Figure FDA0003423321550000015
步骤12，分别计算投影图像序列
Figure FDA0003423321550000016
和
Figure FDA0003423321550000017
的二维傅里叶变换，得到复数图像序列
Figure FDA0003423321550000018
分别截取复数图像序列
Figure FDA0003423321550000019
Figure FDA00034233215500000110
中每张复数图像中心的复数值
Figure FDA00034233215500000111
步骤13，按照下式计算
Figure FDA00034233215500000112
和
Figure FDA00034233215500000113
两两之间的匹配移动量mAi,Bj，将PAi向y轴正向移动mAi,Bj个单位，得到P′Ai，计算P′Ai和PBj的均方根误差erri,j,称之为投影图像PAi和PBj的匹配误差，
Figure FDA00034233215500000114
与
Figure FDA00034233215500000115
两两之间的匹配误差则为
Figure FDA00034233215500000116
其中i为虚数单位：
Figure FDA00034233215500000117
步骤14，找出匹配误差
Figure FDA0003423321550000021
中最小值的索引imin,jmin；
步骤15，利用频域匹配得到
Figure FDA0003423321550000022
到
Figure FDA0003423321550000023
的匹配移动量,记为
Figure FDA0003423321550000024
步骤16，将投影图像
Figure FDA0003423321550000025
向上平移
Figure FDA0003423321550000026
得到所述第一参考投影图像PAstd，投影图像
Figure FDA0003423321550000027
向下平移
Figure FDA0003423321550000028
得到第二参考投影图像PBstd；
步骤2，将所述待校正投影图像序列
Figure FDA0003423321550000029
中的第t张待校正投影图像Pt(t＝1,2,…,NA,t≠A,B)按照预设的搜索角序列
Figure FDA00034233215500000210
中的每一个角度进行旋转，得到一组旋转后的投影图像序列
Figure FDA00034233215500000211
再与所述第一参考投影图像PAstd按照如下子步骤进行中心切片频域匹配校正：
步骤21，计算第一参考投影图像PAstd和投影图像序列
Figure FDA00034233215500000212
的二维傅里叶变换，分别得到复数第一参考投影图像f(PAstd)、复数图像序列
Figure FDA00034233215500000213
分别截取f(PAstd)和
Figure FDA00034233215500000214
中每张复数图像的中心复数值f(PAstd)mid、值
Figure FDA00034233215500000215
步骤22，按照下式(1)计算投影图像Ptk到PAstd的匹配移动量mtk,将Ptk向y轴负向移动mtk个单位，得到P′tk，计算P′tk和PAstd的均方根误差errk,称之为投影图像Ptk和PAstd的匹配误差，PAstd与
Figure FDA00034233215500000216
之间的匹配误差则为
Figure FDA00034233215500000217
其中i为虚数单位：
Figure FDA00034233215500000218
步骤23，获取匹配误差
Figure FDA00034233215500000219
中的最小值索引kmin，得到摇摆校正角
Figure FDA00034233215500000220
和校正摇摆角后的第t张投影图像
Figure FDA00034233215500000221
步骤24，利用中心切片频域匹配获取
Figure FDA00034233215500000222
到PAstd的匹配移动量
Figure FDA00034233215500000223
记为
Figure FDA00034233215500000224
步骤25，将投影图像
Figure FDA00034233215500000225
向上平移
Figure FDA00034233215500000226
得到第t张投影图像的校正端跳误差后的投影图像P′t,std；
步骤26，取不同的待校正投影图像，返回步骤2，直至校正完毕所述待校正投影图像序列
Figure FDA00034233215500000227
中所有第t张待校正投影图像的摇摆和端跳，得到校正摇摆和端跳误差后的投影图像
Figure FDA0003423321550000031
2.如权利要求1所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤24中的“频域匹配”为“中心切片频域匹配”，其具体包括：
Figure FDA0003423321550000032
和PAstd分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N,分别对
Figure FDA0003423321550000033
和PAstd进行二维傅里叶变换，得到两个二维的复值图像
Figure FDA0003423321550000034
和f(PAstd)，
Figure FDA0003423321550000035
和f(PAstd)的中心项即第(M/2，N/2)项为低频项，分别记为
Figure FDA0003423321550000036
和f(PAstd)mid，利用下式(2)计算y轴方向的匹配移动量
Figure FDA0003423321550000037
Figure FDA0003423321550000038
3.如权利要求1或2所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，还包括：
步骤3，将步骤2获得的
Figure FDA0003423321550000039
记为校正摇摆和端跳误差后的迭代第1次后的投影图像
Figure FDA00034233215500000310
按照如下子步骤，利用迭代重投影和邻近匹配修正技术校正径跳误差：
步骤31，设定最大迭代次数Niter,令iter＝0；
步骤32，令iter＝iter+1，对校正摇摆和端跳误差后的迭代第iter次后的投影图像
Figure FDA00034233215500000311
进行邻近匹配修正，修正后的结果记为
Figure FDA00034233215500000312
步骤33，对投影图像
Figure FDA00034233215500000313
直接进行迭代重建，得到重建图像
Figure FDA00034233215500000314
步骤34，将
Figure FDA00034233215500000315
进行投影角
Figure FDA00034233215500000316
的投影，得到第iter次迭代重投影中重建图
Figure FDA00034233215500000317
的投影图像序列
Figure FDA00034233215500000318
步骤35，将不同投影角下的投影图像
Figure FDA00034233215500000319
和
Figure FDA00034233215500000320
的各行进行频域匹配，找出
Figure FDA0003423321550000041
到
Figure FDA0003423321550000042
的各行匹配移动量的平均值
Figure FDA0003423321550000043
并将投影图像
Figure FDA0003423321550000044
向右移动
Figure FDA0003423321550000045
得到投影图像
Figure FDA0003423321550000046
步骤36，如果iter<Niter，进入步骤32，否则进行步骤37；
步骤37，迭代终止，得到径跳校正完成后的投影图像
Figure FDA0003423321550000047
4.如权利要求3所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤32中的
Figure FDA0003423321550000048
进行邻近匹配修正”具体包括：
针对按照投影角度先后次序排列的投影图像
Figure FDA0003423321550000049
从第一投影图像
Figure FDA00034233215500000410
和第二投影图像
Figure FDA00034233215500000411
开始，每相邻两张投影图像对应的行数据进行x方向进行频域匹配，找出每行的匹配移动量，并计算其均值，将其中的前一张图像沿着x方向按平均的移动匹配量进行移动，直到遍历完所有图像。
5.如权利要求4所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，“每相邻两张投影图像对应的行数据进行x方向进行频域匹配”为第t张投影图像
Figure FDA00034233215500000412
和第t+1张投影图像
Figure FDA00034233215500000413
对应的行数据进行x方向进行频域匹配，其方法具体包括：
Figure FDA00034233215500000414
和
Figure FDA00034233215500000415
中每行对应的数组分别为x轴方向的两个一维N元实数数组，记为
Figure FDA00034233215500000416
和
Figure FDA00034233215500000417
其坐标均为{1,2,…,N},row＝1,2,…,Nr，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure FDA00034233215500000418
和
Figure FDA00034233215500000419
Figure FDA00034233215500000420
和
Figure FDA00034233215500000421
的中间项即第
Figure FDA00034233215500000422
项为低频项，分别记为
Figure FDA00034233215500000423
和
Figure FDA0003423321550000051
利用下式(3)和(4)计算二者的匹配移动量
Figure FDA0003423321550000052
Figure FDA0003423321550000053
Figure FDA0003423321550000054
6.如权利要求5所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，“修正后的结果记为
Figure FDA0003423321550000055
的方法具体包括：
步骤321，t＝0；
步骤322，令t＝t+1，若t≥NA,则转到步骤326，否则进入步骤323；
步骤323，将图像
Figure FDA0003423321550000056
向x轴正方向移动
Figure FDA0003423321550000057
个单位，得到
Figure FDA0003423321550000058
步骤324，将图像
Figure FDA0003423321550000059
向x轴负方向移动
Figure FDA00034233215500000510
个单位，得到
Figure FDA00034233215500000511
步骤325，返回步骤322；
步骤326，修正完成，得到
Figure FDA00034233215500000512
7.如权利要求6所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤35中的“将不同投影角下的投影图像
Figure FDA00034233215500000513
和
Figure FDA00034233215500000514
的各行进行频域匹配”为“将不同投影角下的投影图像
Figure FDA00034233215500000515
中的每一幅投影图像
Figure FDA00034233215500000516
和
Figure FDA00034233215500000517
中的每一幅投影图像
Figure FDA00034233215500000518
的相对应行进行频域匹配”的方法具体包括：
Figure FDA0003423321550000061
和
Figure FDA0003423321550000062
分别为投影图像
Figure FDA0003423321550000063
和
Figure FDA0003423321550000064
的第row行，row＝1,2,…,Nr，也就是x轴方向的两个一维N元实数数组，其坐标均为{1,2,…,N},分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure FDA0003423321550000065
和
Figure FDA0003423321550000066
Figure FDA0003423321550000067
和
Figure FDA0003423321550000068
的中间项即第N/2项为低频项，分别记为
Figure FDA0003423321550000069
和
Figure FDA00034233215500000610
利用下式(5)(6)计算二者的匹配移动量
Figure FDA00034233215500000611
Figure FDA00034233215500000612
Figure FDA00034233215500000613
8.如权利要求7所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，根据式(5)和(6)计算出来的
Figure FDA00034233215500000614
计算得到
Figure FDA00034233215500000615
到
Figure FDA00034233215500000616
的各行匹配移动量的平均值
Figure FDA00034233215500000617
的方法具体包括：
记第t组投影图像
Figure FDA00034233215500000618
到
Figure FDA00034233215500000619
的各行匹配移动量为
Figure FDA00034233215500000620
则各行匹配移动量的均值
Figure FDA00034233215500000621
依次计算t＝1,2,…,NA时对应的
Figure FDA00034233215500000622
得到
Figure FDA00034233215500000623
9.如权利要求1所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤13具体包括：
步骤131，i＝0
步骤132，i＝i+1，j＝0，若i>Ns，转到步骤135；
步骤133，j＝j+1，若j>Ns，转到步骤132；
步骤134，对PAi和PBj进行中心切片频域匹配，并计算匹配误差errij，转到步骤133；步骤135，两两频域匹配过程完成，得到
Figure FDA0003423321550000071
CN202011094025.4A 2020-10-14 2020-10-14 微纳米ct投影数据的中心切片对齐方法 Active CN112233156B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011094025.4A CN112233156B (zh) 2020-10-14 2020-10-14 微纳米ct投影数据的中心切片对齐方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011094025.4A CN112233156B (zh) 2020-10-14 2020-10-14 微纳米ct投影数据的中心切片对齐方法

Publications (2)

Publication Number Publication Date
CN112233156A CN112233156A (zh) 2021-01-15
CN112233156B true CN112233156B (zh) 2022-02-15

Family

ID=74113351

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011094025.4A Active CN112233156B (zh) 2020-10-14 2020-10-14 微纳米ct投影数据的中心切片对齐方法

Country Status (1)

Country Link
CN (1) CN112233156B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117557623B (zh) * 2023-11-14 2024-05-14 上海月新生科信息科技有限公司 一种冷冻电镜图像序列的精准快速对齐方法

Citations (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101011258A (zh) * 2005-11-15 2007-08-08 Ge医疗系统环球技术有限公司 X射线ct设备和x射线ct透视设备
WO2010127241A2 (en) * 2009-04-30 2010-11-04 The Regents Of The University Of California System and methods for fast implementation of equally-sloped tomography
CN102768759A (zh) * 2012-07-04 2012-11-07 深圳安科高技术股份有限公司 一种术中ct图像射束硬化伪影校正方法及装置
WO2015080911A1 (en) * 2013-11-26 2015-06-04 Saudi Arabian Oil Company Automated workflow for 3d core digital modeling from computerized tomography scanner (cts) images
CN105092617A (zh) * 2015-09-18 2015-11-25 重庆大学 一种基于x射线能谱ct和x射线荧光ct技术的双模态分子成像系统
CN106651929A (zh) * 2016-11-21 2017-05-10 中国科学院西安光学精密机械研究所 一种高精度亚像元位移产生方法
CN106875334A (zh) * 2017-03-20 2017-06-20 中国科学院地质与地球物理研究所 页岩ct成像方法及装置
CN107239708A (zh) * 2017-04-24 2017-10-10 广东工业大学 一种基于量子混沌映射和分数域变换的图像加密方法
CN108169257A (zh) * 2017-12-08 2018-06-15 中国科学院高能物理研究所 具有三种定量成像机制的纳米ct三维成像方法
CN109729029A (zh) * 2018-12-26 2019-05-07 西安烽火电子科技有限责任公司 一种时域干扰消除的频域判决反馈均衡方法与装置
CN109886170A (zh) * 2019-02-01 2019-06-14 长江水利委员会长江科学院 一种钉螺智能检测识别与统计系统
CN110782489A (zh) * 2019-10-21 2020-02-11 科大讯飞股份有限公司 影像数据的匹配方法、装置、设备及计算机可读存储介质
CN111553849A (zh) * 2020-03-28 2020-08-18 中国人民解放军战略支援部队信息工程大学 基于局部特征匹配的锥束ct几何伪影去除方法及装置

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8634626B2 (en) * 2010-06-29 2014-01-21 The Chinese University Of Hong Kong Registration of 3D tomography images
TWI669942B (zh) * 2016-10-14 2019-08-21 聯發科技股份有限公司 用以移除漣漪效應之平滑化濾波方法和裝置
CN107203002B (zh) * 2017-06-12 2019-05-24 中国科学院地质与地球物理研究所 反演速度模型的建立方法和地下结构的像的获得方法
CN111103585B (zh) * 2019-11-29 2023-05-26 西安电子科技大学 基于双通道联合处理的合成孔径宽带信号源侦察成像方法
CN111142105B (zh) * 2020-01-06 2023-06-13 西南电子技术研究所（中国电子科技集团公司第十研究所） 复杂运动目标isar成像方法
CN111220546A (zh) * 2020-03-30 2020-06-02 中北大学 一种数字剪切散斑的同轴笼式结构无损检测系统

Patent Citations (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101011258A (zh) * 2005-11-15 2007-08-08 Ge医疗系统环球技术有限公司 X射线ct设备和x射线ct透视设备
WO2010127241A2 (en) * 2009-04-30 2010-11-04 The Regents Of The University Of California System and methods for fast implementation of equally-sloped tomography
CN102768759A (zh) * 2012-07-04 2012-11-07 深圳安科高技术股份有限公司 一种术中ct图像射束硬化伪影校正方法及装置
WO2015080911A1 (en) * 2013-11-26 2015-06-04 Saudi Arabian Oil Company Automated workflow for 3d core digital modeling from computerized tomography scanner (cts) images
CN105092617A (zh) * 2015-09-18 2015-11-25 重庆大学 一种基于x射线能谱ct和x射线荧光ct技术的双模态分子成像系统
CN106651929A (zh) * 2016-11-21 2017-05-10 中国科学院西安光学精密机械研究所 一种高精度亚像元位移产生方法
CN106875334A (zh) * 2017-03-20 2017-06-20 中国科学院地质与地球物理研究所 页岩ct成像方法及装置
CN107239708A (zh) * 2017-04-24 2017-10-10 广东工业大学 一种基于量子混沌映射和分数域变换的图像加密方法
CN108169257A (zh) * 2017-12-08 2018-06-15 中国科学院高能物理研究所 具有三种定量成像机制的纳米ct三维成像方法
CN109729029A (zh) * 2018-12-26 2019-05-07 西安烽火电子科技有限责任公司 一种时域干扰消除的频域判决反馈均衡方法与装置
CN109886170A (zh) * 2019-02-01 2019-06-14 长江水利委员会长江科学院 一种钉螺智能检测识别与统计系统
CN110782489A (zh) * 2019-10-21 2020-02-11 科大讯飞股份有限公司 影像数据的匹配方法、装置、设备及计算机可读存储介质
CN111553849A (zh) * 2020-03-28 2020-08-18 中国人民解放军战略支援部队信息工程大学 基于局部特征匹配的锥束ct几何伪影去除方法及装置

Non-Patent Citations (7)

* Cited by examiner, † Cited by third party
Title
"A GPU-based multi-resolution approach to iterative reconstruction algorithms in x-ray 3D dual spectral computed tomography";Zhang Huitao等;《Neurocomputing》;20161231;第215卷(第26期);第71-81页 *
"A Phase Field Approach to Limited-angle Tomographic Reconstruction";Caty,Olivier等;《Fundamenta Informaticae》;20200208;第172卷(第2期);第203-219页 *
"Automated registration for augmenting micro-CT 3D images";Shane Jamie Latham等;《Proceedings Computation Techniques and Applications Conference》;20081231;第534-548页 *
"Reprojection Alignment for Trajectory Perturbation Estimation in Microtomography";Shane J. Latham等;《IEEE Transactions on Computation Imageing》;20180630;第4卷(第2期);第271-283页 *
"基于散射校正板的锥束微纳CT系统的散射校正";唐天旭等;《光学学报》;20190831;第39卷(第8期);第1-11页 *
"基于近邻主特征匹配的微纳米尺度位移测量";刘永俊等;《东北大学学报（自然科学版）》;20170228;第38卷(第2期);第195-199页 *
"解析-迭代重建算法在微纳CT中的应用";王成祥等;《中国体现学与图像分析》;20171231;第22卷(第3期);第299-306页 *

Also Published As

Publication number Publication date
CN112233156A (zh) 2021-01-15

Similar Documents

Publication Publication Date Title
US7950849B2 (en) Method and device for geometry analysis and calibration of volumetric imaging systems
WO2018218611A1 (zh) 锥形束计算机断层成像系统的几何参数确定方法
CN111383272B (zh) 一种双目视觉水果分拣并联机器人视觉盲区末端位姿检测方法
US9253449B2 (en) Mosaic picture generation
CN111612731B (zh) 基于双目显微视觉的测量方法、装置、系统及介质
CN106056645B (zh) 基于频域分析的ct图像平移运动伪影校正方法
CN111429494A (zh) 一种基于生物视觉的点云高精度自动配准方法
CN112233156B (zh) 微纳米ct投影数据的中心切片对齐方法
CN111899290A (zh) 一种联合偏振和双目视觉的三维重建方法
CN110057847A (zh) Tr层析扫描投影重排方法及装置
CN117557733B (zh) 基于超分辨率的自然保护区三维重建方法
CN115546289A (zh) 一种基于机器人的复杂结构件三维形貌测量方法
CN112233157B (zh) 微纳米ct系统转台误差的频域校正方法
CN112233158B (zh) 微纳米ct投影图像的二次投影配准方法
CN108364325B (zh) 规则样品x射线ct投影图位置平移偏差检测及校正方法
CN114638874B (zh) 基于因式分解和isea的空间目标三维重构方法
Guo et al. Iterative image reconstruction for limited‐angle CT using optimized initial image
CN109886872B (zh) 安检设备及其图像检测方法
CN114565720A (zh) 一种基于线结构光旋转扫描的主动式三维重建系统及方法
US7702180B2 (en) Imaging method and device for the computer-assisted evaluation of computer-tomographic measurements by means of direct iterative reconstruction
CN108663386B (zh) 基于特征纹理模板的锥束ct系统探测器角度偏差测量方法
CN113052929A (zh) 一种基于投影视角加权的直线扫描cl重建方法
US20070142723A1 (en) Correction of geometric distortions in magnetic resonance images
Luo et al. Geometric calibration based on a simple phantom for multi-lens microscopic CT
CN114092594B (zh) 一种锥束ct系统和轴对称外形样品的几何误差校正方法

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
TA01 Transfer of patent application right

Effective date of registration: 20210414

Address after: 215127 west of Yiye Road south of Huikai road Luzhi town Wuzhong District Suzhou City Jiangsu Province

Applicant after: Lima precision measurement technology (Suzhou) Co.,Ltd.

Address before: 100048 No. 105 West Third Ring Road North, Beijing, Haidian District

Applicant before: Capital Normal University

TA01 Transfer of patent application right
GR01 Patent grant
GR01 Patent grant