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

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

Info

Publication number
CN112233156A
CN112233156A CN202011094025.4A CN202011094025A CN112233156A CN 112233156 A CN112233156 A CN 112233156A CN 202011094025 A CN202011094025 A CN 202011094025A CN 112233156 A CN112233156 A CN 112233156A
Authority
CN
China
Prior art keywords
projection
image
matching
projection image
astd
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.)
Granted
Application number
CN202011094025.4A
Other languages
English (en)
Other versions
CN112233156B (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
Capital Normal University
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 Capital Normal University filed Critical Capital Normal University
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 BDA0002723049550000021
下的待校正投影图像序列为
Figure BDA0002723049550000022
取所述待校正投影图像序列
Figure BDA0002723049550000023
中的第一待校正投影图像PA和第二待校正投影图像PB，并获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd；
步骤2，将所述待校正投影图像序列
Figure BDA0002723049550000024
中的第t张待校正投影图像Pt(t＝1，2，...，NA，t≠A，B)按照预设的搜索角序列
Figure BDA0002723049550000025
中的每一个角度进行旋转，得到一组旋转后的投影图像序列
Figure BDA0002723049550000026
再与所述第一参考投影图像PAstd按照如下子步骤进行二次投影的频域匹配校正：
步骤21，计算第一参考投影图像PAstd和投影图像序列
Figure BDA0002723049550000027
的二维傅里叶变换，分别得到复数第一参考投影图像f(PAstd)、复数图像序列
Figure BDA0002723049550000028
分别截取f(PAstd)和
Figure BDA0002723049550000029
中每张复数图像的中心复数值f(PAstd)mid、值
Figure BDA00027230495500000210
步骤22，按照下式(1)计算投影图像Ptk到PAstd的匹配移动量mtk，将Ptk向y轴负向移动mtk个单位，得到P′tk，计算P′Ai和PAstd的均方根误差errk，称之为投影图像Ptk和PAstd的匹配误差，PAstd与
Figure BDA00027230495500000211
之间的匹配误差则为
Figure BDA00027230495500000212
其中i为虚数单位：
Figure BDA00027230495500000213
步骤23，获取匹配误差
Figure BDA00027230495500000214
中的最小值索引kmin，得到摇摆校正角
Figure BDA00027230495500000215
和校正摇摆角后的第t张投影图像
Figure BDA00027230495500000216
步骤24，利用中心切片频域匹配获取
Figure BDA00027230495500000217
到PAstd的匹配移动量
Figure BDA00027230495500000218
记为
Figure BDA00027230495500000219
步骤25，将投影图像
Figure BDA00027230495500000220
向上平移
Figure BDA00027230495500000221
得到第t张投影图像的校正端跳误差后的投影图像P′t，std；
步骤26，取不同的待校正投影图像，返回步骤2，直至校正完毕所述待校正投影图像序列
Figure BDA0002723049550000031
中所有第t张待校正投影图像的摇摆和端跳，得到校正摇摆和端跳误差后的投影图像
Figure BDA0002723049550000032
进一步地，步骤24中的“频域匹配”为“中心切片频域匹配”，其具体包括：
Figure BDA0002723049550000033
和PAstd分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N，分别对
Figure BDA0002723049550000034
和PAstd进行二维傅里叶变换，得到两个二维的复值图像
Figure BDA0002723049550000035
和f(PAstd)，
Figure BDA0002723049550000036
和f(PAstd)的中心项即第(M/2，N/2)项为低频项，分别记为
Figure BDA0002723049550000037
和f(PAstd)mid，利用下式(2)计算y轴方向的匹配移动量
Figure BDA0002723049550000038
Figure BDA0002723049550000039
进一步地，该方法还包括：
步骤3，将步骤2获得的
Figure BDA00027230495500000310
记为校正摇摆和端跳误差后的迭代第1次后的投影图像
Figure BDA00027230495500000311
按照如下子步骤，利用迭代重投影和邻近匹配修正技术校正径跳误差：
步骤31，设定最大迭代次数Niter，令iter＝0；
步骤32，令iter＝iter+1，对校正摇摆和端跳误差后的迭代第iter次后的投影图像
Figure BDA00027230495500000312
进行邻近匹配修正，修正后的结果记为
Figure BDA00027230495500000313
步骤33，对投影图像
Figure BDA00027230495500000314
直接进行迭代重建，得到重建图像
Figure BDA00027230495500000315
步骤34，将
Figure BDA00027230495500000316
进行投影角
Figure BDA00027230495500000317
的投影，得到第iter次迭代重投影中重建图
Figure BDA00027230495500000318
的投影图像序列
Figure BDA00027230495500000319
步骤35，将不同投影角下的投影图像
Figure BDA00027230495500000320
和
Figure BDA00027230495500000321
的各行进行频域匹配，找出
Figure BDA0002723049550000041
到
Figure BDA0002723049550000042
的各行匹配移动量的平均值
Figure BDA0002723049550000043
并将投影图像
Figure BDA0002723049550000044
向右移动
Figure BDA0002723049550000045
得到投影图像
Figure BDA0002723049550000046
步骤36，如果iter＜Niter，进入步骤32，否则进行步骤37；
步骤37，迭代终止，得到径跳校正完成后的投影图像
Figure BDA0002723049550000047
进一步地，步骤32中的“
Figure BDA0002723049550000048
进行邻近匹配修正”具体包括：
针对按照投影角度先后次序排列的投影图像
Figure BDA0002723049550000049
从第一投影图像
Figure BDA00027230495500000410
和第二投影图像
Figure BDA00027230495500000411
开始，每相邻两张投影图像对应的行数据进行x方向进行频域匹配，找出每行的匹配移动量，并计算其均值，将其中的前一张图像沿着x方向按平均的移动匹配量进行移动，直到遍历完所有图像。
进一步地，“每相邻两张投影图像对应的行数据进行x方向进行频域匹配”为第t张投影图像
Figure BDA00027230495500000412
和第t+1张投影图像
Figure BDA00027230495500000413
对应的行数据进行x方向进行频域匹配，其方法具体包括：
Figure BDA00027230495500000414
和
Figure BDA00027230495500000415
中每行对应的数组分别为x轴方向的两个一维N元实数数组，记为
Figure BDA00027230495500000416
和
Figure BDA00027230495500000417
其坐标均为{1，2，...，N}，row＝1，2，...，Nr，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure BDA00027230495500000418
和
Figure BDA00027230495500000419
和
Figure BDA00027230495500000420
的中间项即第
Figure BDA00027230495500000421
项为低频项，分别记为
Figure BDA00027230495500000422
和
Figure BDA00027230495500000423
利用下式(3)和(4)计算二者的匹配移动量
Figure BDA00027230495500000424
Figure BDA0002723049550000051
Figure BDA0002723049550000052
进一步地，“修正后的结果记为
Figure BDA0002723049550000053
”的方法具体包括：
步骤321，t＝0；
步骤322，令t＝t+1，若t≥NA，则转到步骤326，否则进入步骤323；
步骤323，将图像
Figure BDA0002723049550000054
向x轴正方向移动
Figure BDA0002723049550000055
个单位，得到
Figure BDA0002723049550000056
步骤324，将图像
Figure BDA0002723049550000057
向x轴负方向移动
Figure BDA0002723049550000058
个单位，得到
Figure BDA0002723049550000059
步骤325，返回步骤322；
步骤326，修正完成，得到
Figure BDA00027230495500000510
进一步地，步骤35中的“将不同投影角下的投影图像
Figure BDA00027230495500000511
和
Figure BDA00027230495500000512
的各行进行频域匹配”为“将不同投影角下的投影图像
Figure BDA00027230495500000513
中的每一幅投影图像
Figure BDA00027230495500000514
和
Figure BDA00027230495500000515
中的每一幅投影图像Pt (iter)的相对应行进行频域匹配”的方法具体包括：
Figure BDA00027230495500000516
和Pt (iter)(row)分别为投影图像
Figure BDA00027230495500000517
和Pt (iter)的第row行，row＝1，2，...，Nr，也就是x轴方向的两个一维N元实数数组，其坐标均为{1，2，...，N}，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure BDA00027230495500000518
和f(Pt (iter)(row))，
Figure BDA00027230495500000519
和f(Pt (iter)(row))的中间项即第N/2项为低频项，分别记为
Figure BDA0002723049550000061
和f(Pt (iter)(row))mid，利用下式(5)(6)计算二者的匹配移动量
Figure BDA0002723049550000062
Figure BDA0002723049550000063
Figure BDA0002723049550000064
进一步地，根据式(5)和(6)计算出来的
Figure BDA0002723049550000065
计算得到
Figure BDA0002723049550000066
到
Figure BDA0002723049550000067
的各行匹配移动量的平均值
Figure BDA0002723049550000068
的方法具体包括：
记第t组投影图像
Figure BDA0002723049550000069
到Pt (iter)的各行匹配移动量为
Figure BDA00027230495500000610
则各行匹配移动量的均值
Figure BDA00027230495500000611
依次计算t＝1，2，...，NA时对应的
Figure BDA00027230495500000612
得到
Figure BDA00027230495500000613
进一步地，步骤1中的“获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd”的方法具体包括：
步骤11，将PA，PB分别按照预设的搜索角序列
Figure BDA00027230495500000614
中的每一个角度进行旋转，得到两组旋转后的投影图像序列
Figure BDA00027230495500000615
步骤12，分别计算投影图像序列
Figure BDA00027230495500000616
和
Figure BDA00027230495500000617
的二维傅里叶变换，得到复数图像序列
Figure BDA00027230495500000618
分别截取复数图像序列
Figure BDA00027230495500000619
中每张复数图像中心的复数值
Figure BDA00027230495500000620
步骤13，按照下式计算
Figure BDA0002723049550000071
和
Figure BDA0002723049550000072
两两之间的匹配移动量mAi，Bj，将PAi向y轴正向移动mAi，Bj个单位，得到P′Ai，计算P′Ai和PBj的均方根误差erri，j，称之为投影图像PAi和PBj的匹配误差，
Figure BDA0002723049550000073
与
Figure BDA0002723049550000074
两两之间的匹配误差则为
Figure BDA0002723049550000075
其中i为虚数单位：
Figure BDA0002723049550000076
步骤14，找出匹配误差
Figure BDA0002723049550000077
中最小值的索引imin，jmin；
步骤15，利用频域匹配得到
Figure BDA0002723049550000078
到
Figure BDA0002723049550000079
的匹配移动量，记为
Figure BDA00027230495500000710
步骤16，将投影图像
Figure BDA00027230495500000711
向上平移
Figure BDA00027230495500000712
得到所述第一参考投影图像PAstd，投影图像
Figure BDA00027230495500000713
向下平移
Figure BDA00027230495500000714
得到第二参考投影图像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 BDA00027230495500000715
进一步地，步骤15中的“利用中心切片频域匹配计算二次投影
Figure BDA00027230495500000716
到
Figure BDA00027230495500000717
的匹配移动量
Figure BDA00027230495500000718
”的方法具体包括：
Figure BDA00027230495500000719
和
Figure BDA00027230495500000720
分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N，分别对
Figure BDA00027230495500000721
和
Figure BDA00027230495500000722
进行二维傅里叶变换，得到两个二维的复值图像
Figure BDA00027230495500000723
和
Figure BDA00027230495500000724
和
Figure BDA00027230495500000725
的中心项即第(M/2，N/2)项为低频项，分别记为
Figure BDA00027230495500000726
和
Figure BDA00027230495500000727
利用下式(7)计算y轴方向的匹配移动量
Figure BDA00027230495500000728
Figure BDA0002723049550000081
进一步地，步骤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轴方向的匹配移动量
Figure BDA0002723049550000082
Figure BDA0002723049550000083
将PAi向y轴正向移动
Figure BDA0002723049550000084
个单位，得到P′Ai，计算P′Ai和PBj的均方根误差，记为erri，j，称之为投影图像PAi和PBj的匹配误差。
本发明所涉及的校正方法属于软件层面上的方法，在估计端跳和摇摆信息方面，本发明采用二次投影的频域匹配方法，通过利用图像变换在频域的数学性质，为加快匹配的速度和精度，进而为提高计算效率提供了有利条件。在估计径跳信息方面，本发明提高采用加邻近匹配的迭代重投影方法在迭代重投影的基础上利用邻近匹配进行修正，从而加快了迭代速度，有利于提高计算效率。利用这些信息实现对投影数据的校正，重建出空间分辨率更高的图像。
附图说明
图1为本发明实施例提供的基于二次投影频域匹配的微纳米CT***转台误差校正方法中的两个一维N元数组Ltk和
Figure BDA0002723049550000085
的坐标与对应数值示意图。
图2为本发明实施例提供的基于二次投影频域匹配的微纳米CT***转台误差校正方法中的利用傅立叶变换找出数组Ltk和
Figure BDA0002723049550000086
的频域中间的相角差m的示意图。
图3为本发明实施例提供的基于二次投影频域匹配的微纳米CT***转台误差校正方法中的一维N元数组Ltk沿y轴正方向移动m个单位的示意图。
图4为图3中的数组Ltk和
Figure BDA0002723049550000087
达成匹配的状态示意图。
图5的最左侧图是四棱柱形状样品S自虚线区域示意的初始位置到在y轴方向上存在端跳
Figure BDA0002723049550000088
后的灰色填充区域示意的待校正位置示意图，中间图是两位置样品在扫描角度A下各自的投影图像PA，最右侧图是扫描角度A下的两投影沿着x方向的进行二次投影的各自得到的一维数组结果。
图6为最左侧图是四棱柱形状样品S自虚线区域示意的初始位置到在y轴方向上存在端跳
Figure BDA0002723049550000091
后的灰色填充区域示意的待校正位置示意图，中间图是两位置样品在扫描角度A+90度下各自的投影图像PB，最右侧图是扫描角度B下的两投影沿着x方向的进行二次投影的各自得到的一维数组结果。
图7为本发明根据二次投影LA和LB频域匹配后所得的匹配移动量myAB，沿着y轴方向将投影PA和PB平移的端跳校正过程示意图。
具体实施方式
下面结合附图和实施例对本发明进行详细的描述。
如图4至图6所示，为了叙述简便，下面以长方体的样品S的二次投影原理为例，说明本发明实施例所在的XYZ坐标系，将样品S的一个表面置于纸面所在的平面内，视为XY平面，X轴与Y轴垂直。其中，X轴的正向如图中的二次投影到X轴上LB的延伸方向所示，即在纸面内向右的方向。Y轴在纸面内向上的方向为下文中的“向上”，反之为“向下”。Z轴垂直于XY平面，正方向为指向纸面内地方向，Z轴正方向为射线扫描方向。端跳方向定义为沿Y轴方向的跳动，即在Y轴方向的平移量。
本发明实施例提供的微纳米CT投影数据的中心切片对齐方法包括如下步骤：
步骤1，获取样品在投影角
Figure BDA0002723049550000092
下的待校正投影图像序列为
Figure BDA0002723049550000093
取所述待校正投影图像序列
Figure BDA0002723049550000094
中的第一待校正投影图像PA和第二待校正投影图像PB，并获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd。其中，任取待校正投影图像序列为
Figure BDA0002723049550000095
中，间隔90°投影角的一组投影图像，分别取PA，PB。
在一个实施例中，步骤1中，“获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd”的方法具体包括：
步骤11，将PA，PB分别按照预设的搜索角序列
Figure BDA0002723049550000096
中的每一个角度进行旋转，得到两组旋转后的投影图像序列
Figure BDA0002723049550000101
步骤12，计算
Figure BDA0002723049550000102
的二次投影，得到两组二次投影信号序列
Figure BDA0002723049550000103
其中，二次投影技术为现有技术，在此不展开描述。
步骤13，对
Figure BDA0002723049550000104
两两进行频域匹配，并计算匹配误差
Figure BDA0002723049550000105
步骤14，找出匹配误差
Figure BDA0002723049550000106
中最小值的索引imin，jmin。
步骤15，利用频域匹配计算二次投影
Figure BDA0002723049550000107
到
Figure BDA0002723049550000108
的匹配移动量
Figure BDA0002723049550000109
步骤16，将投影图像
Figure BDA00027230495500001010
向上平移
Figure BDA00027230495500001011
得到所述第一参考投影图像PAstd，投影图像
Figure BDA00027230495500001012
向下平移
Figure BDA00027230495500001013
得到第二参考投影图像PBstd。
在一个实施例中，步骤13中的“两两”可以理解为：每固定第一个序列中的一个信号，第二个序列信号取遍，需要匹配Ns次，第一个序列要取Ns次信号固定，所以一共进行了
Figure BDA00027230495500001014
次匹配，其具体子步骤包括：
步骤131，i＝0；
步骤132，i＝i+1，j＝0，若i＞Ns，转到步骤135；
步骤133，(3)j＝j+1，若j＞Ns，转到步骤132；
步骤134，对LAi和LBj进行频域匹配，并计算匹配误差errij，转到步骤133；
步骤135，两两频域匹配过程完成，得到
Figure BDA00027230495500001015
进一步地，步骤134中的“LAi和LBj进行频域匹配”的方法包括：
LAi和LBj为y轴有两个一维N元实数数组，其坐标均为{1，2，...，N}，分别进行离散傅里叶变换，得到两个一维N元复数数组f(LAi)和f(LBj)。f(LAi)和f(LBj)的中间项即第N/2项，为低频项，分别记为f(LAi)mid和f(LBj)mid，计算二者的相角差mAi，
Figure BDA0002723049550000111
将LAi数组沿着y轴正向移动mAi，Bj个单位，仍取坐标{1，2，...，N}上的数所构成的数组，记为LAi′，LBj′与LBj的近似程度最高，LAi′和LBj达成频域匹配，LAi′和LBj之间的均方根误差RMSE(LAi′，LBj)，记为匹配误差ErrAi，Bj。
在一个实施例中，步骤15中的“利用中心切片频域匹配计算二次投影
Figure BDA0002723049550000112
到
Figure BDA0002723049550000113
的匹配移动量
Figure BDA0002723049550000114
”的方法具体包括；
Figure BDA0002723049550000115
和
Figure BDA0002723049550000116
分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N，分别对
Figure BDA0002723049550000117
和
Figure BDA0002723049550000118
进行二维傅里叶变换，得到两个二维的复值图像
Figure BDA0002723049550000119
和
Figure BDA00027230495500001110
Figure BDA00027230495500001111
和
Figure BDA00027230495500001112
的中心项即第(M/2，N/2)项为低频项，分别记为
Figure BDA00027230495500001113
和
Figure BDA00027230495500001114
利用下式(7)计算y轴方向的匹配移动量
Figure BDA00027230495500001115
Figure BDA00027230495500001116
步骤2，将所述待校正投影图像序列
Figure BDA00027230495500001117
中的第t张待校正投影图像为Pt(t＝1，2，...，NA，t≠A，B)与所述第一参考投影图像PAstd按照如下子步骤进行二次投影的频域匹配校正：
步骤21，计算第一参考投影图像PAstd和投影图像序列
Figure BDA00027230495500001118
的二维傅里叶变换，分别得到复数第一参考投影图像f(PAstd)、复数图像序列
Figure BDA00027230495500001119
分别截取f(PAstd)和
Figure BDA00027230495500001120
中每张复数图像的中心复数值f(PAstd)mid、值
Figure BDA00027230495500001121
步骤22，按照下式计算投影图像Ptk到PAstd的匹配移动量mtk，将Ptk向y轴负向移动mtk个单位，得到P′tk，计算P′Ai和PAstd的均方根误差errk，称之为投影图像Ptk和PAstd的匹配误差，PAstd与
Figure BDA00027230495500001122
之间的匹配误差则为
Figure BDA00027230495500001123
其中i为虚数单位：
Figure BDA00027230495500001124
步骤23，找出匹配误差
Figure BDA00027230495500001125
中的最小值索引kmin，得到摇摆校正角
Figure BDA0002723049550000121
和校正摇摆角后的第t张投影图像
Figure BDA0002723049550000122
步骤24，计算
Figure BDA0002723049550000123
的二次投影，得到校正摇摆角后的第t张投影图的二次投影信号
Figure BDA0002723049550000124
并利用频域匹配计算
Figure BDA0002723049550000125
到
Figure BDA0002723049550000126
的匹配移动量
Figure BDA0002723049550000127
步骤25，将投影图像
Figure BDA0002723049550000128
向上平移
Figure BDA0002723049550000129
得到第t张投影图像的校正端跳误差后的投影图像P′t，std。
步骤26，取不同的待校正投影图像，返回步骤2，直至校正完毕所述待校正投影图像序列
Figure BDA00027230495500001210
中所有第t张待校正投影图像的摇摆和端跳，得到校正摇摆和端跳误差后的投影图像
Figure BDA00027230495500001211
在一个实施例中，步骤22具体包括：
步骤221，如图1和图2所示，
Figure BDA00027230495500001212
与
Figure BDA00027230495500001213
分别为y轴方向的两个一维N元实数数组，其坐标均为{1，2，...，N}，分别进行离散傅里叶变换，得到两个一维N元复数数组f(Ltk)和
Figure BDA00027230495500001214
f(Ltk)和
Figure BDA00027230495500001215
的中间项即第N/2项为低频项，分别记为f(Ltk)mid和
Figure BDA00027230495500001216
利用下式(1)计算二者的相角差
Figure BDA00027230495500001217
Figure BDA00027230495500001218
步骤222，如图3所示，将
Figure BDA00027230495500001219
数组沿y轴方向正向移动
Figure BDA00027230495500001220
个单位，仍取坐标{1，2，...，N}上的数所构成的数组L′tk，得到数组序列
Figure BDA00027230495500001221
步骤223，计算
Figure BDA00027230495500001222
与
Figure BDA00027230495500001223
之间的均方根误差，得到匹配误差
Figure BDA00027230495500001224
在一个实施例中，步骤24中的“频域匹配”为“中心切片频域匹配”，其具体包括：
Figure BDA00027230495500001225
和PAstd分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N，分别对
Figure BDA00027230495500001226
和PAstd进行二维傅里叶变换，得到两个二维的复值图像
Figure BDA00027230495500001227
和f(PAstd)，
Figure BDA00027230495500001228
和f(PAstd)的中心项即第(M/2，N/2)项为低频项，分别记为
Figure BDA00027230495500001229
和f(PAstd)mid，利用下式(2)计算二者的匹配移动量
Figure BDA0002723049550000131
Figure BDA0002723049550000132
在一个实施例中，本发明实施例提供的微纳米CT投影数据的中心切片对齐方法还包括：
步骤3，将步骤2获得的
Figure BDA0002723049550000133
记为校正摇摆和端跳误差后的迭代第1次后的投影图像
Figure BDA0002723049550000134
按照如下子步骤，利用迭代重投影和邻近匹配修正技术校正径跳误差：
步骤31，设定最大迭代次数Niter，令iter＝0。
步骤32，令iter＝iter+1，对校正摇摆和端跳误差后的迭代第iter次后的投影图像
Figure BDA0002723049550000135
进行邻近匹配修正，修正后的结果记为
Figure BDA0002723049550000136
步骤33，对投影图像
Figure BDA0002723049550000137
直接进行迭代重建，得到重建图像
Figure BDA0002723049550000138
步骤34，将
Figure BDA0002723049550000139
进行投影角
Figure BDA00027230495500001310
的投影，得到第iter次迭代重投影中重建图
Figure BDA00027230495500001311
的投影图像序列
Figure BDA00027230495500001312
步骤35，将不同投影角下的投影图像
Figure BDA00027230495500001313
和
Figure BDA00027230495500001314
的各行进行频域匹配，找出
Figure BDA00027230495500001315
到
Figure BDA00027230495500001316
的各行匹配移动量的平均值
Figure BDA00027230495500001317
并将投影图像
Figure BDA00027230495500001318
向右移动
Figure BDA00027230495500001319
得到投影图像
Figure BDA00027230495500001320
步骤36，如果iter＜Niter，进入步骤32，否则进行步骤37。
步骤37，迭代终止，得到径跳校正完成后的投影图像
Figure BDA00027230495500001321
在一个实施例中，步骤32中的“
Figure BDA0002723049550000141
进行邻近匹配修正”具体包括：针对按照投影角度先后次序排列的投影图像
Figure BDA0002723049550000142
从第一投影图像
Figure BDA0002723049550000143
和第二投影图像
Figure BDA0002723049550000144
开始，每相邻两张投影图像对应的行数据进行x方向进行频域匹配，找出每行的匹配移动量，并计算其均值，将其中的前一张图像沿着x方向按平均的移动匹配量进行移动，直到遍历完所有图像。
在一个实施例中，步骤32中的“每相邻两张投影图像对应的行数据进行x方向进行频域匹配移动量的均值计算方法包括：
记两张投影图像各行的匹配移动量为
Figure BDA0002723049550000145
则匹配移动量的均值
Figure BDA0002723049550000146
表示为下式：
Figure BDA0002723049550000147
在一个实施例中，“每相邻两张投影图像对应的行数据进行x方向进行频域匹配”为第t张投影图像
Figure BDA0002723049550000148
和第t+1张投影图像
Figure BDA0002723049550000149
对应的行数据进行x方向进行频域匹配，其方法具体包括：
Figure BDA00027230495500001410
和
Figure BDA00027230495500001411
中每行对应的数组分别为x轴方向的两个一维N元实数数组，记为
Figure BDA00027230495500001412
和
Figure BDA00027230495500001413
其坐标均为{1，2，...，N}，row＝1，2，...，Nr，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure BDA00027230495500001414
和
Figure BDA00027230495500001415
和
Figure BDA00027230495500001416
的中间项即第N/2项为低频项，分别记为
Figure BDA00027230495500001417
和
Figure BDA00027230495500001418
利用下式(3)和(4)计算二者的匹配移动量
Figure BDA00027230495500001419
Figure BDA00027230495500001420
Figure BDA00027230495500001421
下面以第一投影图像
Figure BDA0002723049550000151
和第二投影图像
Figure BDA0002723049550000152
为例，说明每相邻两张投影图像对应的行数据进行x方向进行频域匹配的方法：
Figure BDA0002723049550000153
和
Figure BDA0002723049550000154
中每行对应的数组分别为x轴方向的两个一维N元实数数组，记为
Figure BDA0002723049550000155
和
Figure BDA0002723049550000156
其坐标均为{1，2，...，N}，row＝1，2，...，Nr，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure BDA0002723049550000157
和
Figure BDA0002723049550000158
和
Figure BDA0002723049550000159
的中间项即第N/2项为低频项，分别记为
Figure BDA00027230495500001510
和
Figure BDA00027230495500001511
利用下式计算二者的匹配移动量
Figure BDA00027230495500001512
Figure BDA00027230495500001513
Figure BDA00027230495500001514
在一个实施例中，步骤32中的“修正后的结果记为
Figure BDA00027230495500001515
”的方法具体包括：
步骤321，t＝0；
步骤322，令t＝t+1，若t≥NA，则转到步骤326，否则进入步骤323；
步骤323，将图像
Figure BDA00027230495500001516
向x轴正方向移动
Figure BDA00027230495500001517
个单位，得到
Figure BDA00027230495500001518
步骤324，将图像
Figure BDA00027230495500001519
向x轴负方向移动
Figure BDA00027230495500001520
个单位，得到
Figure BDA00027230495500001521
步骤325，返回步骤322；
步骤326，修正完成，得到
Figure BDA00027230495500001522
在一个实施例中，步骤35中的“将不同投影角下的投影图像
Figure BDA0002723049550000161
和
Figure BDA0002723049550000162
的各行进行频域匹配”为“将不同投影角下的投影图像
Figure BDA0002723049550000163
中的每一幅投影图像
Figure BDA0002723049550000164
和
Figure BDA0002723049550000165
中的每一幅投影图像Pt (iter)的相对应行进行频域匹配”的方法具体包括：
Figure BDA0002723049550000166
和Pt (iter)(row)分别为投影图像
Figure BDA0002723049550000167
和Pt (iter)的第row行，row＝1，2，...，Nrow，也就是x轴方向的两个一维N元实数数组，其坐标均为{1，2，...，N}，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure BDA0002723049550000168
和f(Pt (iter)(row))，
Figure BDA0002723049550000169
和f(Pt (iter)(row))的中间项即第N/2项为低频项，分别记为
Figure BDA00027230495500001610
和f(Pt (iter)(row))mid，利用下式(5)和(6)计算二者的匹配移动量
Figure BDA00027230495500001611
Figure BDA00027230495500001612
Figure BDA00027230495500001613
在一个实施例中，根据式(5)和(6)计算出来的
Figure BDA00027230495500001614
计算得到
Figure BDA00027230495500001615
到
Figure BDA00027230495500001616
的各行匹配移动量的平均值
Figure BDA00027230495500001617
的方法具体包括：
记第t组投影图像
Figure BDA00027230495500001618
到Pt (iter)的各行匹配移动量为
Figure BDA00027230495500001619
则各行匹配移动量的均值
Figure BDA00027230495500001620
依次计算t＝1，2，...，NA时对应的
Figure BDA00027230495500001621
得到
Figure BDA00027230495500001622
本发明提供的校正方法，相比硬件校正方法，易于操作推广，可以节省大量的硬件成本。相比现有的软件算法校正方法，校正精度更高，计算速度更快。
最后需要指出的是：以上实施例仅用以说明本发明的技术方案，而非对其限制。本领域的普通技术人员应当理解：可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (10)

1.一种微纳米CT投影数据的中心切片对齐方法，其特征在于，包括如下步骤：
步骤1，获取样品在投影角
Figure FDA0002723049540000011
下的待校正投影图像序列为
Figure FDA0002723049540000012
取所述待校正投影图像序列
Figure FDA0002723049540000013
中的第一待校正投影图像PA和第二待校正投影图像PB，并获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd；
步骤2，将所述待校正投影图像序列
Figure FDA0002723049540000014
中的第t张待校正投影图像Pt(t＝1，2，...，NA，t≠A，B)按照预设的搜索角序列
Figure FDA0002723049540000015
中的每一个角度进行旋转，得到一组旋转后的投影图像序列
Figure FDA0002723049540000016
再与所述第一参考投影图像PAstd按照如下子步骤进行二次投影的频域匹配校正：
步骤21，计算第一参考投影图像PAstd和投影图像序列
Figure FDA0002723049540000017
的二维傅里叶变换，分别得到复数第一参考投影图像f(PAstd)、复数图像序列
Figure FDA0002723049540000018
分别截取f(PAstd)和
Figure FDA0002723049540000019
中每张复数图像的中心复数值f(PAstd)mid、值
Figure FDA00027230495400000110
步骤22，按照下式(1)计算投影图像Ptk到PAstd的匹配移动量mtk，将Ptk向y轴负向移动mtk个单位，得到P′tk，计算P′Ai和PAstd的均方根误差errk，称之为投影图像Ptk和PAstd的匹配误差，PAstd与
Figure FDA00027230495400000111
之间的匹配误差则为
Figure FDA00027230495400000112
其中i为虚数单位：
Figure FDA00027230495400000113
步骤23，获取匹配误差
Figure FDA00027230495400000114
中的最小值索引kmin，得到摇摆校正角
Figure FDA00027230495400000115
和校正摇摆角后的第t张投影图像
Figure FDA00027230495400000116
步骤24，利用中心切片频域匹配获取
Figure FDA00027230495400000117
到PAstd的匹配移动量
Figure FDA00027230495400000118
记为
Figure FDA00027230495400000119
步骤25，将投影图像
Figure FDA00027230495400000120
向上平移
Figure FDA00027230495400000121
得到第t张投影图像的校正端跳误差后的投影图像P′t，std；
步骤26，取不同的待校正投影图像，返回步骤2，直至校正完毕所述待校正投影图像序列
Figure FDA0002723049540000021
中所有第t张待校正投影图像的摇摆和端跳，得到校正摇摆和端跳误差后的投影图像
Figure FDA0002723049540000022
2.如权利要求1所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤24中的“频域匹配”为“中心切片频域匹配”，其具体包括：
Figure FDA0002723049540000023
和PAstd分别为两个相同大小的二维投影图像，在x轴和y轴对应的像素单元数分别为M和N，分别对
Figure FDA0002723049540000024
和PAstd进行二维傅里叶变换，得到两个二维的复值图像
Figure FDA0002723049540000025
和f(PAstd)，
Figure FDA0002723049540000026
和f(PAstd)的中心项即第(M/2，N/2)项为低频项，分别记为
Figure FDA0002723049540000027
和f(PAstd)mid，利用下式(2)计算y轴方向的匹配移动量
Figure FDA0002723049540000028
Figure FDA0002723049540000029
3.如权利要求1或2所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，还包括：
步骤3，将步骤2获得的
Figure FDA00027230495400000210
记为校正摇摆和端跳误差后的迭代第1次后的投影图像
Figure FDA00027230495400000211
按照如下子步骤，利用迭代重投影和邻近匹配修正技术校正径跳误差：
步骤31，设定最大迭代次数Niter，令iter＝0；
步骤32，令iter＝iter+1，对校正摇摆和端跳误差后的迭代第iter次后的投影图像
Figure FDA00027230495400000212
进行邻近匹配修正，修正后的结果记为
Figure FDA00027230495400000213
步骤33，对投影图像
Figure FDA00027230495400000214
直接进行迭代重建，得到重建图像
Figure FDA00027230495400000215
步骤34，将
Figure FDA00027230495400000216
进行投影角
Figure FDA00027230495400000217
的投影，得到第iter次迭代重投影中重建图
Figure FDA00027230495400000218
的投影图像序列
Figure FDA00027230495400000219
步骤35，将不同投影角下的投影图像
Figure FDA0002723049540000031
和
Figure FDA0002723049540000032
的各行进行频域匹配，找出
Figure FDA0002723049540000033
到
Figure FDA0002723049540000034
的各行匹配移动量的平均值
Figure FDA0002723049540000035
并将投影图像
Figure FDA0002723049540000036
向右移动
Figure FDA0002723049540000037
得到投影图像
Figure FDA0002723049540000038
步骤36，如果iter＜Niter，进入步骤32，否则进行步骤37；
步骤37，迭代终止，得到径跳校正完成后的投影图像
Figure FDA0002723049540000039
4.如权利要求3所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤32中的“
Figure FDA00027230495400000310
进行邻近匹配修正”具体包括：
针对按照投影角度先后次序排列的投影图像
Figure FDA00027230495400000311
从第一投影图像
Figure FDA00027230495400000312
和第二投影图像
Figure FDA00027230495400000313
开始，每相邻两张投影图像对应的行数据进行x方向进行频域匹配，找出每行的匹配移动量，并计算其均值，将其中的前一张图像沿着x方向按平均的移动匹配量进行移动，直到遍历完所有图像。
5.如权利要求4所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，“每相邻两张投影图像对应的行数据进行x方向进行频域匹配”为第t张投影图像
Figure FDA00027230495400000314
和第t+1张投影图像
Figure FDA00027230495400000315
对应的行数据进行x方向进行频域匹配，其方法具体包括：
Figure FDA00027230495400000316
和
Figure FDA00027230495400000317
中每行对应的数组分别为x轴方向的两个一维N元实数数组，记为
Figure FDA00027230495400000318
和
Figure FDA00027230495400000319
其坐标均为{1，2，...，N}，row＝1，2，...，Nr，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure FDA00027230495400000320
和
Figure FDA00027230495400000321
Figure FDA00027230495400000322
和
Figure FDA00027230495400000323
的中间项即第
Figure FDA00027230495400000324
项为低频项，分别记为
Figure FDA0002723049540000041
和
Figure FDA0002723049540000042
利用下式(3)和(4)计算二者的匹配移动量
Figure FDA0002723049540000043
Figure FDA0002723049540000044
Figure FDA0002723049540000045
6.如权利要求5所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，“修正后的结果记为
Figure FDA0002723049540000046
”的方法具体包括：
步骤321，t＝0；
步骤322，令t＝t+1，若t≥NA，则转到步骤326，否则进入步骤323；
步骤323，将图像
Figure FDA0002723049540000047
向x轴正方向移动
Figure FDA0002723049540000048
个单位，得到
Figure FDA0002723049540000049
步骤324，将图像
Figure FDA00027230495400000410
向x轴负方向移动
Figure FDA00027230495400000411
个单位，得到
Figure FDA00027230495400000412
步骤325，返回步骤322；
步骤326，修正完成，得到
Figure FDA00027230495400000413
7.如权利要求6所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤35中的“将不同投影角下的投影图像
Figure FDA00027230495400000414
和
Figure FDA00027230495400000415
的各行进行频域匹配”为“将不同投影角下的投影图像
Figure FDA00027230495400000416
中的每一幅投影图像
Figure FDA00027230495400000417
和
Figure FDA00027230495400000418
中的每一幅投影图像
Figure FDA00027230495400000419
的相对应行进行频域匹配”的方法具体包括：
Figure FDA0002723049540000051
和
Figure FDA0002723049540000052
分别为投影图像
Figure FDA0002723049540000053
和
Figure FDA0002723049540000054
的第row行，row＝1，2，...，Nr，也就是x轴方向的两个一维N元实数数组，其坐标均为{1，2，...，N}，分别进行离散傅里叶变换，得到两个一维N元复数数组
Figure FDA0002723049540000055
和
Figure FDA0002723049540000056
Figure FDA0002723049540000057
和
Figure FDA0002723049540000058
的中间项即第N/2项为低频项，分别记为
Figure FDA0002723049540000059
和
Figure FDA00027230495400000510
利用下式(5)(6)计算二者的匹配移动量
Figure FDA00027230495400000511
Figure FDA00027230495400000512
Figure FDA00027230495400000513
8.如权利要求7所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，根据式(5)和(6)计算出来的
Figure FDA00027230495400000514
计算得到
Figure FDA00027230495400000515
到
Figure FDA00027230495400000516
的各行匹配移动量的平均值
Figure FDA00027230495400000517
的方法具体包括：
记第t组投影图像
Figure FDA00027230495400000518
到
Figure FDA00027230495400000519
的各行匹配移动量为
Figure FDA00027230495400000520
则各行匹配移动量的均值
Figure FDA00027230495400000521
依次计算t＝1，2，...，NA时对应的
Figure FDA00027230495400000522
得到
Figure FDA00027230495400000523
9.如权利要求1所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤1中的“获取所述第一待校正投影图像PA对应的第一参考投影图像PAstd以及所述第二待校正投影图像PB对应第二参考投影图像PBstd”的方法具体包括：
步骤11，将PA，PB分别按照预设的搜索角序列
Figure FDA0002723049540000061
中的每一个角度进行旋转，得到两组旋转后的投影图像序列
Figure FDA0002723049540000062
步骤12，分别计算投影图像序列
Figure FDA0002723049540000063
和
Figure FDA0002723049540000064
的二维傅里叶变换，得到复数图像序列
Figure FDA0002723049540000065
分别截取复数图像序列
Figure FDA0002723049540000066
中每张复数图像中心的复数值
Figure FDA0002723049540000067
步骤13，按照下式计算
Figure FDA0002723049540000068
和
Figure FDA0002723049540000069
两两之间的匹配移动量mAi，Bj，将PAi向y轴正向移动mAi，Bj个单位，得到P′Ai，计算P′Ai和PBj的均方根误差erri，j，称之为投影图像PAi和PBj的匹配误差，
Figure FDA00027230495400000610
与
Figure FDA00027230495400000611
两两之间的匹配误差则为
Figure FDA00027230495400000612
其中i为虚数单位：
Figure FDA00027230495400000613
步骤14，找出匹配误差
Figure FDA00027230495400000614
中最小值的索引imin，jmin；
步骤15，利用频域匹配得到
Figure FDA00027230495400000615
到
Figure FDA00027230495400000616
的匹配移动量，记为
Figure FDA00027230495400000617
步骤16，将投影图像
Figure FDA00027230495400000618
向上平移
Figure FDA00027230495400000619
得到所述第一参考投影图像PAstd，投影图像
Figure FDA00027230495400000620
向下平移
Figure FDA00027230495400000621
得到第二参考投影图像PBstd。
10.如权利要求9所述的微纳米CT投影数据的中心切片对齐方法，其特征在于，步骤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 FDA0002723049540000071
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 true CN112233156A (zh) 2021-01-15
CN112233156B 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)

Cited By (1)

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

Citations (19)

* 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
US20110317898A1 (en) * 2010-06-29 2011-12-29 Lin Shi Registration of 3D tomography images
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成像方法及装置
CN107203002A (zh) * 2017-06-12 2017-09-26 中国科学院地质与地球物理研究所 反演速度模型及其建立方法和地下结构的像的获得方法
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 西安烽火电子科技有限责任公司 一种时域干扰消除的频域判决反馈均衡方法与装置
CN109845266A (zh) * 2016-10-14 2019-06-04 联发科技股份有限公司 用以移除涟漪效应的平滑化滤波方法和装置
CN109886170A (zh) * 2019-02-01 2019-06-14 长江水利委员会长江科学院 一种钉螺智能检测识别与统计***
CN110782489A (zh) * 2019-10-21 2020-02-11 科大讯飞股份有限公司 影像数据的匹配方法、装置、设备及计算机可读存储介质
CN111103585A (zh) * 2019-11-29 2020-05-05 西安电子科技大学 基于双通道联合处理的合成孔径宽带信号源侦察成像方法
CN111142105A (zh) * 2020-01-06 2020-05-12 西南电子技术研究所（中国电子科技集团公司第十研究所） 复杂运动目标isar成像方法
CN111220546A (zh) * 2020-03-30 2020-06-02 中北大学 一种数字剪切散斑的同轴笼式结构无损检测***
CN111553849A (zh) * 2020-03-28 2020-08-18 中国人民解放军战略支援部队信息工程大学 基于局部特征匹配的锥束ct几何伪影去除方法及装置

Patent Citations (19)

* 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
US20110317898A1 (en) * 2010-06-29 2011-12-29 Lin Shi Registration of 3D tomography images
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技术的双模态分子成像***
CN109845266A (zh) * 2016-10-14 2019-06-04 联发科技股份有限公司 用以移除涟漪效应的平滑化滤波方法和装置
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 广东工业大学 一种基于量子混沌映射和分数域变换的图像加密方法
CN107203002A (zh) * 2017-06-12 2017-09-26 中国科学院地质与地球物理研究所 反演速度模型及其建立方法和地下结构的像的获得方法
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 科大讯飞股份有限公司 影像数据的匹配方法、装置、设备及计算机可读存储介质
CN111103585A (zh) * 2019-11-29 2020-05-05 西安电子科技大学 基于双通道联合处理的合成孔径宽带信号源侦察成像方法
CN111142105A (zh) * 2020-01-06 2020-05-12 西南电子技术研究所（中国电子科技集团公司第十研究所） 复杂运动目标isar成像方法
CN111553849A (zh) * 2020-03-28 2020-08-18 中国人民解放军战略支援部队信息工程大学 基于局部特征匹配的锥束ct几何伪影去除方法及装置
CN111220546A (zh) * 2020-03-30 2020-06-02 中北大学 一种数字剪切散斑的同轴笼式结构无损检测***

Non-Patent Citations (7)

* Cited by examiner, † Cited by third party
Title
CATY,OLIVIER等: ""A Phase Field Approach to Limited-angle Tomographic Reconstruction"", 《FUNDAMENTA INFORMATICAE》 *
SHANE J. LATHAM等: ""Reprojection Alignment for Trajectory Perturbation Estimation in Microtomography"", 《IEEE TRANSACTIONS ON COMPUTATION IMAGEING》 *
SHANE JAMIE LATHAM等: ""Automated registration for augmenting micro-CT 3D images"", 《PROCEEDINGS COMPUTATION TECHNIQUES AND APPLICATIONS CONFERENCE》 *
ZHANG HUITAO等: ""A GPU-based multi-resolution approach to iterative reconstruction algorithms in x-ray 3D dual spectral computed tomography"", 《NEUROCOMPUTING》 *
刘永俊等: ""基于近邻主特征匹配的微纳米尺度位移测量"", 《东北大学学报（自然科学版）》 *
唐天旭等: ""基于散射校正板的锥束微纳CT***的散射校正"", 《光学学报》 *
王成祥等: ""解析-迭代重建算法在微纳CT中的应用"", 《中国体现学与图像分析》 *

Cited By (2)

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

Also Published As

Publication number Publication date
CN112233156B (zh) 2022-02-15

Similar Documents

Publication Publication Date Title
US9253449B2 (en) Mosaic picture generation
CN106056645B (zh) 基于频域分析的ct图像平移运动伪影校正方法
CN111612731B (zh) 基于双目显微视觉的测量方法、装置、***及介质
CN110133014B (zh) 一种芯片内部缺陷检测方法及***
CN111899290A (zh) 一种联合偏振和双目视觉的三维重建方法
CN112233156B (zh) 微纳米ct投影数据的中心切片对齐方法
CN117557733B (zh) 基于超分辨率的自然保护区三维重建方法
CN112233157B (zh) 微纳米ct***转台误差的频域校正方法
CN112233158B (zh) 微纳米ct投影图像的二次投影配准方法
CN110308165B (zh) 使用Pi线优化的计算机断层扫描投影的几何对齐、样本移动校正以及强度归一化
CN108364325B (zh) 规则样品x射线ct投影图位置平移偏差检测及校正方法
Huang et al. Calibrating sensing drift in tomographic inversion
CN116205993A (zh) 一种用于3d aoi的双远心镜头高精度标定方法
CN109886872B (zh) 安检设备及其图像检测方法
US7702180B2 (en) Imaging method and device for the computer-assisted evaluation of computer-tomographic measurements by means of direct iterative reconstruction
CN114565720A (zh) 一种基于线结构光旋转扫描的主动式三维重建***及方法
Donné et al. Robust plane-based calibration for linear cameras
US20070142723A1 (en) Correction of geometric distortions in magnetic resonance images
CN113052929A (zh) 一种基于投影视角加权的直线扫描cl重建方法
Zemek et al. Automatic marker-free estimation methods for the axis of rotation in sub-micron X-ray computed tomography
Luo et al. Geometric calibration based on a simple phantom for multi-lens microscopic CT
Kishore et al. Implicit Reconstructions from Deformed Projections for CryoET
CN108663386A (zh) 基于特征纹理模板的锥束ct***探测器角度偏差测量方法
CN114092594B (zh) 一种锥束ct***和轴对称外形样品的几何误差校正方法
Nakini et al. Distortion correction in 3d-modeling of roots for plant phenotyping

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