CN108010099A - 一种x射线多能谱ct有限角扫描和图像迭代重建方法 - Google Patents
一种x射线多能谱ct有限角扫描和图像迭代重建方法 Download PDFInfo
- Publication number
- CN108010099A CN108010099A CN201711261688.9A CN201711261688A CN108010099A CN 108010099 A CN108010099 A CN 108010099A CN 201711261688 A CN201711261688 A CN 201711261688A CN 108010099 A CN108010099 A CN 108010099A
- Authority
- CN
- China
- Prior art keywords
- image
- energy spectrum
- energy
- scanning
- spectrum
- 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
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T11/00—2D [Two Dimensional] image generation
- G06T11/003—Reconstruction from projections, e.g. tomography
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T11/00—2D [Two Dimensional] image generation
- G06T11/003—Reconstruction from projections, e.g. tomography
- G06T11/005—Specific pre-processing for tomographic reconstruction, e.g. calibration, source positioning, rebinning, scatter correction, retrospective gating
-
- 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
- G06T2211/00—Image generation
- G06T2211/40—Computed tomography
- G06T2211/416—Exact reconstruction
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Apparatus For Radiation Diagnosis (AREA)
Abstract
本发明公开一种X射线多能谱CT有限角扫描模式和图像迭代重建方法，包括：应用多个能量谱的X射线分别在各自的有限角度范围内扫描物体，获得被测物体的在各个能量谱下的投影数据；基于这些投影数据分别重建各个能量谱下的被测物体图像。利用被测物体的一致性，以及各能量谱下重建图像的互补性，设计联合正则化方法，有效地去除各个能量谱下重建图像中由于有限角扫描而造成的伪影，计算得到高质量的重建图像；进而应用图像域分解方法或投影域分解方法计算得到被测物体的基材料图像或双效应图像。本发明适用于多能谱CT在有限角度扫描条件下的图像重建，不仅降低了X射线剂量，缩短扫描时间，且具有重建图像质量高，收敛速度快等优点。
Description
技术领域
本发明涉及X射线CT成像技术领域，具体而言，涉及一种X射线多能谱CT有限角扫描和图像迭代重建方法。
背景技术
X射线计算机断层成像技术(X ray Computed Tomography,简称X射线CT)已广泛应用于医疗诊断、工业无损检测、生物样品检测、安全检查等领域。目前常用的X射线管发出的是服从一定谱分布的多色X射线束，而不是理想的单能量X射线束，采集的投影数据为多色投影数据，因此传统CT也称为单能谱CT。单能谱CT仅使用一个X射线能谱对被测物体在不同角度下进行扫描，利用X射线与物质的相互作用原理，获得被测物体在这些角度下的投影数据，利用这些投影数据反求出物体内部结构信息。在单能谱CT成像中，有可能出现不同物质的射线衰减属性非常接近的情况，这些物质的CT值相近，难以区分。此外传统CT算法中忽略了X射线的多色性，导致在重建图像中产生硬化伪影。多能谱CT使用多个X射线能谱扫描被测物体，或者采用光子计数探测器，获得不同能谱段的光子穿过物体的衰减信息，能够测得比传统单能谱CT更多的被测物体信息。利用这些投影信息可以重建出被测物体的特定基材料的密度图像或双效应图像，进而得到等效原子序数和电子密度图像，基于这些信息对物质进行区分。由于多能谱CT具有更好的物质区分能力，因此具有广泛的应用前景。
虽然由完全投影数据可以唯一重建被测物体的内部信息图像，但在实际CT扫描中，由于受到检测剂量、检测环境、检测时间和成本，以及被测物体材质和结构等的限制，常常难以获得完整投影数据。在不完整投影情况下的进行CT图像重建，是一个病态反问题，利用传统方法重建的图像中会出现严重的伪影。怎样基于不完整投影实现高质量CT图像重建，是CT领域的一个研究热点，具有重要的实用价值。如果扫描角度受限，采集的投影数据的角度数小于180度，基于这些投影数据的CT成像研究被称为有限角CT图像重建，是不完全数据CT成像的一个主要研究内容。多能谱CT扫描中同样也会出现扫描角度受限的情况，基于有限角度的多能谱投影数据的CT图像重建，是本发明着重论述的内容。
从数学角度来分析，多能谱CT成像问题是一个高维数的非线性方程组求解的问题。传统的多能谱CT求解方法大致可以分为两大类：投影域分解法和图像域分解法，这两种分解方法各自都有直接法和迭代法的求解方法。在有限角的情况下，多能谱CT成像问题是一个欠定的非线性方程组求解问题。直接法虽然步骤简单，但难以解决有限角多能谱CT成像问题。迭代法利用数值方法或者优化方法构造迭代结构，通过对图像重建结果逐步修正求出多基材料密度图像等信息。该方法不仅计算量大，而且针对有限角度多能谱CT成像这种不适定问题求解时难以收敛。鉴于正则化方法是求解不适定反问题的一种有效方法，因此可以基于多能谱CT成像的特点，设计一种适用于有限角多能谱CT成像的正则化函数，应用合理的正则化方法，加速迭代收敛速度，得到高质量的稳定的多能谱CT图像。
发明内容
本发明提供一种X射线多能谱CT有限角扫描模式和图像迭代重建方法，用以解决现有技术中当扫描角度受限时，传统的多能谱CT重建算法难以重建被测物体的基材料密度图像及双效应图像等的问题。
本发明的方法包括以下步骤：
S1：用多个X射线能谱在各自有限角范围内扫描被测物体，获得被测物体的多能谱投影数据；
S2：利用所测得的多个能谱的投影数据，用迭代算法重建各个能谱的被测物体图像；
S3：基于各个能谱的被测物体图像，应用图像域分解方法或投影域分解方法计算得到被测物体的基材料图像或双效应图像。
进一步地，其中所述步骤S1中的多能谱CT系统，其射线源可以是能够发出一定能量谱的X射线的普通射线源，也可以是能够发出单能X射线的射线源。
进一步地，其中所述步骤S1中的多能谱CT有限角扫描模式，每个能谱均可以用不足180度的范围扫描被测物体，以适应特定CT系统的扫描限制。在扫描中，不同能谱的X射线的扫描角度范围尽量不同，但可以彼此重叠一部分。
进一步地，其中所述步骤S2中用于重建各个能谱下的被测物体图像μs的正则化约束模型函数为：
其中参数As是第s个能谱扫描时的系统投影矩阵；qs是测得的能谱多色投影；ws是方程的权值；U为“联合能谱图像”；αs和βs分别是后两项的权重因子。通过最小化该函数，实现了各个能谱下的被测物体图像的高质量重建。
进一步地，上述联合能谱图像定义为各个能谱图像的加权和：
这里参数λs为能谱图像μs的权重因子。
进一步地，所述步骤S2中的正则化约束的图像重建模型，可分为如下两个子问题交替迭代求解：
这里是第一个子问题在第n次迭代时得到的解。
进一步地，其中所述步骤S2包括：
S21：为待重建的被测物体的各个能谱图像及联合能谱图像赋初始值；
S22：求解上述正则化约束的图像重建模型中的两个子问题，更新当前重建的能谱被测物体图像；
S23：对当前重建的能谱被测物体图像加入像素值非负的约束条件；
S24：利用第n次迭代中生成的所有能谱图像，更新联合能谱图像。
S25：重复步骤S22-S24，直到满足迭代结束条件，重建出被测物体在各个能谱扫描下的被测物体图像。
进一步地，在所述步骤S3中，在计算被测物体的基材料图像或双效应图像时，无论用图像域分解方法或投影域分解方法，均可采用各自的直接标定分解法和迭代分解法。
进一步地，在所述步骤S3中，当应用投影域分解方法进行基材料图像或双效应图像重建时，需要首先将各个重建的能谱图像在完全角度下作线积分计算，估计出所需的完全角度的投影数据。
进一步地，无论是利用所测得的多能谱投影数据迭代重建各个能谱的被测物体图像，还是基于各个能谱的被测物体图像，应用图像域分解方法或投影域分解方法计算基材料图像或双效应图像，均可采用并行方法计算，并基于硬件并行计算平台加速实现。
与现有技术相比，本发明方法适用于多能谱CT在有限角度扫描条件下的图像重建，不仅降低了X射线剂量，缩短扫描时间，且具有重建图像质量高，抗噪性能强等优点。
附图说明
为了更清楚地说明本发明实施例的技术方案，下面将对实施例所需要使用的附图作简单地介绍。显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明一个实施例的X射线多能谱CT有限角扫描和图像迭代重建方法流程图；
图2为用作测试模体的FORBILD胸腔模型图像；
图3为实验所用的X射线源分别在140kV和80kV管电压下发出的高能谱和低能谱示意图；
图4为对被测模体进行高、低能谱有限角扫描的示意图；
图5a为本发明一个实施例的高能谱下采集到的110度范围的测试模体的多色投影数据图像；
图5b为本发明一个实施例的低能谱下采集到的110度范围的测试模体的多色投影数据图像；
图6a为基于图5a中110度范围的高能谱多色投影数据，采用传统ART算法重建的高能谱图像；
图6b为基于图5b中110度范围的低能谱多色投影数据，采用传统ART算法重建的低能谱图像；
图7a为基于图5a中110度范围的高能谱多色投影数据，采用本发明方法重建的高能谱图像；
图7b为基于图5b中110度范围的低能谱多色投影数据，采用本发明方法重建的低能谱图像；
图8a为基于图7a中重建的高能谱图像，采用本发明方法估算出的360°全角度的高能谱多色投影数据图像；
图8b为基于图7b中重建的低能谱图像，采用本发明方法估算出的360°全角度的低能谱多色投影数据图像；
图9a为本发明一个实施例的骨基材料密度图像重建的结果图像；
图9b为本发明一个实施例的水基材料密度图像重建的结果图像；
图9c为本发明一个实施例所重建的关于能量为70keV的光子的线性衰减系数图像。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行完整地描述。显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有付出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
在忽略X光子散射影响的情况下，X射线多能谱CT图像重建问题可表示为如下的数学模型：
其中，参数s表示能谱索引；ws(E)表示第s个归一化的X射线能谱在能量E处的值；μ(x,E)为被测物体在点x处关于能量E的线性衰减系数；qs,L表示沿着给定X射线路径L所探测到的多色投影值；由于是含有多种能量的X射线能谱的投影，而不是单一能量的X射线的投影，因此这里称为“多色投影”。Ωs是对应于第s个能谱的X射线扫描路径集合。在有限角扫描情况下，不同能谱的X射线的扫描路径尽量覆盖不同的角度范围，但彼此的角度范围可以重叠一部分，即也可以完全不重合，即不重合时，测得的多能谱多色投影在几何上不匹配。
多能谱CT图像重建问题可归结为利用测得的多色投影qs,L的集合求解μ(x,E)。为了减少未知数的个数，在求解方程(1)时，将μ(x,E)分解为与物质特性相关量的线性组合：
其中am(E)是一个能量依赖函数，而fm(x)是一个空间位置函数。在基材料分解模型中，am(E)表示第m种基材料的在能量E处的质量衰减系数，fm(x)表示被测物体中第m种基材料的密度空间分布，即第m种基材料的密度图像。将式(2)代入式(1)得到：
可进一步写为：
其中：
pm(L)＝R(fm(x))＝∫Lfm(x)dl,L∈Ωs,m＝1,2,…,M. (5)
这里pm(L)表示沿着射线L穿过fm(x)的线积分；R表示Radon变换算子。
在数值实现时，被重建的基材料密度图像被离散化为像素阵列。设fm＝(fm,0,fm,1,…,fm,J-1)τ表示fm(x)的离散形式，这里fm,j是在fm(x)的第j个像素上的采样值，参数J是每幅图像的像素总数，τ表示向量的转置。设RL＝(rL,0,rL,1,…,rL,J-1)是一个CT系统的投影向量，这里rL,j表示基材料密度图像的第j个像素对沿着射线L的投影的贡献权值。公式(5)能重写为如下向量形式：
由测得的各个能谱的多色投影实现基材料图像或者双效应图像分解的方法很多，大致可分为两大类：图像域分解方法和投影域分解方法。这两大类方法均包含各自的直接标定分解法和迭代分解法。
图像域分解方法首先由多色投影直接重建的被测物体图像，即：
μs＝R-1(qs) s＝1,2,…,Ns. (7)
这里R-1表示Radon逆变换算子，如可采用传统的波反投影(FBP)类型解析算法，或者代数迭代法(ART)等迭代型算法实现；μs表示由多色投影直接重建的被测物体图像，简称为能谱图像，以区别于基材料图像和双效应图像。得到能谱图像后，图像域分解方法将能谱图像分解为基材料图像或双效应图像:
m＝1,2,…,M. (8)
这里Dimg表示图像域基材料或双效应分解算子，可以用传统的图像域直接分解方法，如Brooks、Coleman等方法实现，也可以是图像域迭代分解方法实现。
投影域分解方法首先将多色投影分解为基材料图像或者双效应图像的线积分，即：
m＝1,2,…,M. (9)
这里Dproj表示投影域分解算子，可采用EART、MDIR等方法实现。得到各个角度的基材料图像或者双效应图像的线积分后，即可用传统单能图像重建方法实现基材料图像或者双效应图像的逐个重建，即：
fm＝R-1(pm)，m＝1,2,…,M. (10)
在扫描角度受限的情况下，如果直接利用所采集的有限角多色投影重建基材料图像，一方面由于投影数据的不足，另一方面由于多能分解过程不平衡性所导致的基材料图像噪声的放大，无论是利用传统ART或SART算法通过重建能谱图像再图像域分解的方法，还是利用EART等方法直接进行投影域分解的方法，所重建出来的图像都会失真，出现严重的条纹伪影和渐变伪影等。
为此，本发明的方法首先利用正则化约束的重建方法重建各个能谱的被测物体的图像，然后再基于这些图像进行基材料分解或者双效应分解。为了能够重建出在各个能谱有限角扫描下的高质量的被测物体图像，本发明利用了一个已证明的性质：如果采集的投影数据中包含与待重建图像的特征边界相切的射线的投影数据，这部分图像边缘容易被恢复出来；反之该特征边界就难以被重建出来。根据这个性质，在扫描时，本发明用不同的X射线能谱在不同有限角范围内扫描被测物体，获得被测物体的多能谱投影数据。由于扫描范围不同，直接用多色投影数据重建的各个能谱图像具有互补性，即重建出的被测物体的有效边缘不同，且图像中伪影的朝向和位置也不相同。进而利用被测物体及扫描几何参数的一致性，设计了如下的不同能谱重建图像的联合正则化约束下的图像重建方法，有效地去除各个能谱重建图像中由于有限角扫描而造成的伪影，获得高质量的重建图像。利用正则化约束的重建方法重建每一个能谱图像μs的具体方法如下：
上式中包含三项，第一项为数据保真项，第二项和第三项为图像正则化约束项，αs和βs是用于调节各项平衡的权重因子，可以在迭代中根据当前图像的有效像素值、噪声值以及伪影值之间的大小关系，动态调整这两个参数的值。
数据保真项通过最小化加权投影残差的L2范数的平方，使重建的图像尽可能满足Radon变换。这里As是第s个能谱扫描时的系统投影矩阵；qs是测得的能谱多色投影；ws是方程的权值，由其投影数据的方差决定。
第二项为图像μs的总变差最小约束项，用于在保持边界的情况下，去除图像μs中的噪声。
第三项利用被测物体的一致性，以及各能谱重建图像的互补性，通过最小化“联合能谱图像U”在Nd个方向上差分的L0范数之和，对能谱图像μs进行正则化，去除图像中的伪影。这里联合能谱图像定义为各个能谱图像的加权和：
其中参数λs为能谱图像μs的权值，用于调节各个能谱图像的幅度，使之处于同一水平。
第三项可进一步表示为：
这里每个能谱图像μs包含Nj个像素，每个像素处有Nd个差分方向。注意这里是逐个差分方向依次求最小化。
作为一个实施例，重建所有能谱图像的主要步骤如下：
1)设定每个能谱图像的迭代初值为0或者某个常数，然后计算联合能谱图像U(0)的迭代初值；
2)通过求解上面式(11)，迭代更新每个能谱图像
公式(11)分为两个子问题交替迭代求解：
子问题(14)和(15)可用惩罚加权最小二乘法(PWLS)、惩罚似然法(PL)、分裂Bregman算法等进行求解。这里子问题(15)中的是子问题(14)在第n次迭代时得到的解。注意对重建的图像加像素值非负的约束。
3)利用步骤(2)中第n次迭代中生成的所有能谱图像更新联合能谱图像:
4)重复步骤2)～3)，直到满足迭代结束条件，重建出能谱图像。
得到重建的所有能谱图像后，应用图像域分解方法或投影域分解方法，重建各个基材料图像或者双效应图像。在应用投影域分解方法进行图像重建时，只须将各个能谱重建图像在完全角度下作线积分计算，便可估计出所需的完全的投影数据：
q's＝R(μs) (17)
这里R表示Radon变换算子。然后再利用公式(9)的投影域分解方法分解出基材料图像或者双效应图像的线积分。最后利用公式(10)实现基材料图像或者双效应图像的逐个重建。在使用迭代重建算法进行多基材料或双效应图像重建时，也可以对待重建图像进行正则化约束，如图像的总变差最小化等，用于更好地去除图像中的噪声和伪影。
以下通过一个具体实施例来说明本发明的多能谱CT有限角扫描和图像迭代重建方法的实现过程。图1为本发明一个实施例的X射线多能谱CT有限角扫描和图像迭代重建方法流程图。首先对被测物体进行X射线多能谱有限角扫描，得到各个能谱的多色投影数据；然后利用所测得的多个能谱的投影数据，用迭代算法重建各个能谱的被测物体图像；最后基于各个能谱的被测物体图像，应用图像域分解方法或投影域分解方法计算得到被测物体的基材料图像或双效应图像。为提高重建的速度，图1实施例的被测物体的图像重建过程还可以采用并行方法计算，并基于硬件并行计算平台如图形处理器(GPU)等加速实现。
为了描述简单，本实施例选用目前最常用的双能谱扫描模式。图2所示为本实施例中双能谱CT系统发出X射线低、高能谱，这两个能谱是X光管的管电压分别为80kV和140kV(1mm铜滤波片)下发出的。在数值实现时，能谱的采样间隔设为1keV。图3为本实施例中所用的测试模体，该测试模体是一个标准的FORBILD胸腔模型，尺寸为400mm×200mm。分别选水和骨组织作为两种基材料，水的密度为1.0g/cm3，骨组织的密度为1.92g/cm3。在双能谱CT扫描中，采用扇束扫描。系统的扫描参数为：射线源到探测器的距离为1200mm，射线源到转台中心的距离为981mm，探测器由512个1.345mm的探测单元组成。如图4所示，射线源在140kV和80kV管电压下发出高、低能谱，每隔1°采集投影数据一次，分别对被测模体做110°范围的扫描，其中重合的扫描范围为20°。采集到的测试模体的多色投影数据如图5所示，其中图5a为高能谱下采集的多色投影数据；图5b为低能谱下采集的多色投影数据。基于这些投影数据，用联合代数迭代算法SART逐角度修正迭代一次后直接重建的被测模体的高能谱图像和低能谱图像分别如图6a和图6b所示，重建图像的尺寸为512×512，即Nj＝512×512个像素。每个像素处包含垂直和水平两个差分方向，即Nd＝2。由于只有有限角度的投影数据，数据量不足，采用传统迭代算法得到的图像中不仅有明显的条带状伪影和渐变伪影，而且部分边缘没有得到良好重建。
图7为采用本发明方法，图像初值置0迭代120次后重建的结果图像,其中图7a是高能谱图像，图7b是低能谱图像。比较图6和图7，可以看出本发明实施例的方法可以高质量地重建高、低能谱图像。不仅有效地消除了条带状伪影和渐变伪影，而且所有边缘得到了良好的重建。
图8a为基于图7a中重建的高能谱图像，采用公式(17)估算出的360°全角度的高能谱多色投影数据图像。图8b为基于图7b中重建的低能谱图像，采用公式(17)估算出的360°全角度的低能谱多色投影数据图像。
图9为在已知能谱的情况下，基于图8中计算得到的全角度多色投影，采用基于投影域分解的扩展的联合代数迭代算法(E-SART)，迭代3次后重建的结果图像，其中图9a是骨基材料密度图像，图9b是水基材料密度图像，图9c是关于能量为70keV的光子的线性衰减系数图像。这里线性衰减系数图像是通过将重建的水基图像与水对指定能量的X射线光子质量衰减系数的乘积以及重建的骨基图像与骨对指定能量的X射线光子质量衰减系数的乘积求和得到的。为了清楚显示重建的细节，图9a的显示灰度窗设为[0.4 1.0]，图9b的显示灰度窗设为[0.9 1.1]，图9c的显示灰度窗设为[-100HU 100HU]。由图9可以看出本发明实施例的方法可以很好地分解出双基材料密度图像，基材料密度图像中仅有一些轻微的采样伪影，在线性衰减系数图像中基本消除了射线硬化伪影。
进而用归一化均方距离(NMSB)和归一化平均绝对距离(NMAD)作量化测试。本发明方法重建图像的相应NMSB和NMAD值分别为0.090112和0.063286，由此可知虽然扫描角度受限，但本发明方法的重建的图像误差很小，进一步验证了本发明方法的优点。
本实施例说明本发明方法适用于多能谱CT在有限角度扫描条件下的图像重建，不仅降低了X射线剂量，缩短扫描时间，且具有重建图像质量高，抗噪性能强等优点。
本实施例为了绘图简单，仅给出了在X射线扇束有限角扫描情况下的图示，本发明方法同样适用于X射线锥形束有限角扫描情况。
本领域普通技术人员可以理解：附图只是一个实施例的示意图，附图中的模块或流程并不一定是实施本发明所必须的。
本领域普通技术人员可以理解：实施例中的装置中的模块可以按照实施例描述分布于实施例的装置中，也可以进行相应变化位于不同于本实施例的一个或多个装置中。上述实施例的模块可以合并为一个模块，也可以进一步拆分成多个子模块。
最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围。
Claims (10)
1.一种X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，包括以下步骤：
S1：用多能谱CT扫描模式通过多个能量谱的X射线在各自有限角范围内扫描物体，获得被测物体的多能谱投影数据；
S2：利用所测得的多能谱投影数据，用迭代算法重建各个能量谱下的被测物体图像；
S3：基于各个能量谱下的被测物体图像，应用图像域分解方法或投影域分解方法计算得到被测物体的基材料图像或双效应图像。
2.根据权利要求1所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，在步骤S1中所述的多个能量谱的X射线通过多能谱CT系统获得，其射线源包括能够发出宽能谱X射线的普通射线源，以及能够发出窄能谱X射线的射线源；所述多能谱CT扫描模式包括单射线源单探测器系统的不同X射线能谱的多次扫描模式、多射线源多探测器的同时扫描模式以及基于光子探测器的不同能谱段扫描模式中的至少一种。
3.根据权利要求1所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，在步骤S1中，多个能量谱的X射线在各自有限角范围内扫描物体是指每个能谱的圆周扫描范围均不足180度；在扫描中，不同能谱的X射线的扫描角度范围彼此不同或彼此重叠一部分。
4.根据权利要求1所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，在步骤S2中，设共采用Ns个能谱的X射线扫描物体，Ns为正整数，所获得的Ns组有限角投影数据，通过最小化下式中的正则化约束模型函数实现各个能谱下的被测物体图像μs的重建：
<mrow>
<msub>
<mi>&mu;</mi>
<mi>s</mi>
</msub>
<mo>=</mo>
<munder>
<mrow>
<mi>arg</mi>
<mi> </mi>
<mi>m</mi>
<mi>i</mi>
<mi>n</mi>
</mrow>
<mrow>
<msub>
<mi>&mu;</mi>
<mi>s</mi>
</msub>
<mo>&GreaterEqual;</mo>
<mn>0</mn>
</mrow>
</munder>
<mrow>
<mo>(</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mo>|</mo>
<mo>|</mo>
<msub>
<mi>w</mi>
<mi>s</mi>
</msub>
<mo>(</mo>
<mrow>
<msub>
<mi>A</mi>
<mi>s</mi>
</msub>
<msub>
<mi>&mu;</mi>
<mi>s</mi>
</msub>
<mo>-</mo>
<msub>
<mi>q</mi>
<mi>s</mi>
</msub>
</mrow>
<mo>)</mo>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
<mo>+</mo>
<msub>
<mi>&alpha;</mi>
<mi>s</mi>
</msub>
<mo>|</mo>
<mo>|</mo>
<msub>
<mi>&mu;</mi>
<mi>s</mi>
</msub>
<mo>|</mo>
<msub>
<mo>|</mo>
<mrow>
<mi>T</mi>
<mi>V</mi>
</mrow>
</msub>
<mo>+</mo>
<msub>
<mi>&beta;</mi>
<mi>s</mi>
</msub>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>d</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>N</mi>
<mi>d</mi>
</msub>
</munderover>
<mo>|</mo>
<mo>|</mo>
<msub>
<mo>&part;</mo>
<mi>d</mi>
</msub>
<mi>U</mi>
<mo>|</mo>
<msub>
<mo>|</mo>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
<mo>,</mo>
<mi>s</mi>
<mo>=</mo>
<mn>1</mn>
<mo>,</mo>
<mn>2</mn>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<msub>
<mi>N</mi>
<mi>s</mi>
</msub>
<mo>.</mo>
</mrow>
上式中的函数共包含三项，其中第一项为数据保真项，其参数As是第s个能谱扫描时的系统投影矩阵，qs是测得的能谱多色投影矩阵，ws是方程的权值，设为其投影数据的方差的倒数，第一项通过最小化加权投影残差的L2范数的平方，使重建的图像尽可能满足投影变换；第二项为图像μs的总变差最小约束项，用于在保持边界的情况下，去除图像μs中的噪声；第三项通过最小化“联合能谱图像U”，在Nd个方向上差分的L0范数之和，对能谱图像μs进行正则化，去除图像中的伪影；公式中αs和βs是用于平衡这三项的调节参数，0<αs<1,0<βs<1。
5.根据权利要求4所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，正则化约束的图像重建模型中，联合能谱图像U定义为各个能谱图像的加权和：
<mrow>
<mi>U</mi>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>s</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>N</mi>
<mi>s</mi>
</msub>
</munderover>
<msub>
<mi>&lambda;</mi>
<mi>s</mi>
</msub>
<msub>
<mi>&mu;</mi>
<mi>s</mi>
</msub>
</mrow>
其中参数λs为能谱图像μs的权值，用于调节各个能谱图像的幅度，使之处于同一水平。
6.根据权利要求4所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，正则化约束的图像重建模型中，其第三项进一步表示为：
<mrow>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>d</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>N</mi>
<mi>d</mi>
</msub>
</munderover>
<mo>|</mo>
<mo>|</mo>
<msub>
<mo>&part;</mo>
<mi>d</mi>
</msub>
<mi>U</mi>
<mo>|</mo>
<msub>
<mo>|</mo>
<mn>0</mn>
</msub>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>d</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>N</mi>
<mi>d</mi>
</msub>
</munderover>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>N</mi>
<mi>j</mi>
</msub>
</munderover>
<mo>#</mo>
<mo>{</mo>
<mi>j</mi>
<mo>|</mo>
<mo>|</mo>
<msub>
<mo>&part;</mo>
<mi>d</mi>
</msub>
<msub>
<mi>U</mi>
<mi>j</mi>
</msub>
<mo>|</mo>
<mo>&NotEqual;</mo>
<mn>0</mn>
<mo>}</mo>
</mrow>
这里Nj是每个能谱图像μs包含的像素数目，Nd如权利要求4中所述是差分方向个数，#表示计数，该项通过逐个差分方向最小化联合能谱图像在Nj个像素处的差分不为0的个数，实现对能谱图像μs中伪影的去除。
7.根据权利要求4所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，正则化约束的图像重建模型通过如下两个子问题交替迭代求解：
<mrow>
<msubsup>
<mi>&mu;</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>)</mo>
<mo>&prime;</mo>
</mrow>
</msubsup>
<mo>=</mo>
<munder>
<mrow>
<mi>arg</mi>
<mi> </mi>
<mi>m</mi>
<mi>i</mi>
<mi>n</mi>
</mrow>
<mrow>
<msubsup>
<mi>&mu;</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>)</mo>
<mo>&prime;</mo>
</mrow>
</msubsup>
<mo>&GreaterEqual;</mo>
<mn>0</mn>
</mrow>
</munder>
<mrow>
<mo>(</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mo>|</mo>
<mo>|</mo>
<msub>
<mi>w</mi>
<mi>s</mi>
</msub>
<mo>(</mo>
<mrow>
<msub>
<mi>A</mi>
<mi>s</mi>
</msub>
<msubsup>
<mi>&mu;</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</msubsup>
<mo>-</mo>
<msub>
<mi>q</mi>
<mi>s</mi>
</msub>
</mrow>
<mo>)</mo>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
<mo>+</mo>
<msub>
<mi>&alpha;</mi>
<mi>s</mi>
</msub>
<mo>|</mo>
<mo>|</mo>
<msubsup>
<mi>&mu;</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</msubsup>
<mo>|</mo>
<msub>
<mo>|</mo>
<mrow>
<mi>T</mi>
<mi>V</mi>
</mrow>
</msub>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<msubsup>
<mi>&mu;</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>)</mo>
</mrow>
</msubsup>
<mo>=</mo>
<munder>
<mrow>
<mi>arg</mi>
<mi> </mi>
<mi>m</mi>
<mi>i</mi>
<mi>n</mi>
</mrow>
<mrow>
<msubsup>
<mi>&mu;</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>)</mo>
</mrow>
</msubsup>
<mo>&GreaterEqual;</mo>
<mn>0</mn>
</mrow>
</munder>
<mrow>
<mo>(</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mo>|</mo>
<mo>|</mo>
<msubsup>
<mi>&mu;</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>)</mo>
</mrow>
</msubsup>
<mo>-</mo>
<msubsup>
<mi>&mu;</mi>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>)</mo>
<mo>&prime;</mo>
</mrow>
</msubsup>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
<mo>+</mo>
<msub>
<mi>&beta;</mi>
<mi>s</mi>
</msub>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>d</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<msub>
<mi>N</mi>
<mi>d</mi>
</msub>
</munderover>
<mo>|</mo>
<mo>|</mo>
<msub>
<mo>&part;</mo>
<mi>d</mi>
</msub>
<msup>
<mi>U</mi>
<mrow>
<mo>(</mo>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</msup>
<mo>|</mo>
<msub>
<mo>|</mo>
<mn>0</mn>
</msub>
<mo>)</mo>
</mrow>
</mrow>
这里是第一个子问题在第n次迭代时得到的解，n为正整数。
8.根据权利要求7所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，其中所述步骤S2包括：
S21：为待重建的被测物体的各个能谱图像及联合能谱图像赋初始值；
S22：分别交替求解权利要求7中的有限角多能谱图像重建模型的两个子问题，更新当前重建的能谱被测物体图像；
S23：对当前重建的能谱被测物体图像加入像素值非负的约束条件；
S24：利用第n次迭代中生成的所有能谱图像，更新联合能谱图像。
S25：重复步骤S22-S24，直到满足迭代结束条件，重建出被测物体在各个能谱扫描下的被测物体图像。
9.根据权利要求1所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，在所述步骤S3中，在计算被测物体的基材料图像或双效应图像时，图像域分解方法或投影域分解方法均可采用各自的直接标定分解法或迭代分解法；在采用直接标定分解法时，首先用标定模体标定出分解公式的参数；在采用迭代分解法时，预先标定出能谱或者投影公式。
10.根据权利要求1所述的X射线多能谱CT有限角扫描和图像迭代重建方法，其特征在于，在所述步骤S3中，当应用投影域分解方法进行基材料图像或双效应图像重建时，需要首先将各个重建的能谱图像在完全角度下作线积分计算，估计出所需的完全角度的投影数据；和/或
所述图像迭代重建方法适用于X射线平行束、扇束有限角扫描和锥形束有限角扫描中的任意一种；和/或
利用所测得的多能谱投影数据迭代重建各个能谱的被测物体图像，或基于各个能谱的被测物体图像，应用图像域分解方法或投影域分解方法计算基材料图像或双效应图像，均可采用并行方法计算，并基于硬件并行计算平台加速实现。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201711261688.9A CN108010099B (zh) | 2017-12-04 | 2017-12-04 | 一种x射线多能谱ct有限角扫描和图像迭代重建方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201711261688.9A CN108010099B (zh) | 2017-12-04 | 2017-12-04 | 一种x射线多能谱ct有限角扫描和图像迭代重建方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN108010099A true CN108010099A (zh) | 2018-05-08 |
| CN108010099B CN108010099B (zh) | 2020-12-25 |
Family
ID=62056401
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201711261688.9A Active CN108010099B (zh) | 2017-12-04 | 2017-12-04 | 一种x射线多能谱ct有限角扫描和图像迭代重建方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN108010099B (zh) |
Cited By (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109009181A (zh) * | 2018-06-07 | 2018-12-18 | 西安交通大学 | 双能量ct下同时估计x射线球管光谱和重建图像的方法 |
| CN109146994A (zh) * | 2018-09-17 | 2019-01-04 | 南京航空航天大学 | 一种面向多能谱x射线ct成像的金属伪影校正方法 |
| CN109697691A (zh) * | 2018-12-27 | 2019-04-30 | 重庆大学 | 一种基于l0范数和奇异值阈值分解的双正则项优化的有限角投影重建方法 |
| CN109875593A (zh) * | 2019-02-02 | 2019-06-14 | 中北大学 | Ct成像方法、存储介质和装置 |
| CN110942494A (zh) * | 2019-02-08 | 2020-03-31 | 苏州波影医疗技术有限公司 | 利用三个数据域迭代重建计算机断层扫描图像的系统 |
| CN111127477A (zh) * | 2019-12-12 | 2020-05-08 | 明峰医疗系统股份有限公司 | 一种ct多色谱体素仿真方法 |
| CN111968060A (zh) * | 2020-08-28 | 2020-11-20 | 首都师范大学 | 一种基于倾斜投影修正技术的多能谱ct快速迭代重建方法 |
| CN112200882A (zh) * | 2020-10-10 | 2021-01-08 | 首都师范大学 | 基于傅里叶变换微分性质的板状物ct图像快速重建方法 |
| CN112288762A (zh) * | 2020-10-15 | 2021-01-29 | 西北工业大学 | 一种有限角度ct扫描的离散迭代重建方法 |
| CN112581556A (zh) * | 2020-12-25 | 2021-03-30 | 上海联影医疗科技股份有限公司 | 多能ct图像硬化校正方法、装置、计算机设备和存储介质 |
| CN114916950A (zh) * | 2022-07-21 | 2022-08-19 | 中国科学院深圳先进技术研究院 | 基于多层平板探测器的高空间分辨能谱ct图像重建方法 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103559699A (zh) * | 2013-11-18 | 2014-02-05 | 首都师范大学 | 一种基于投影估计的多能谱ct图像重建方法 |
| CN103559729A (zh) * | 2013-11-18 | 2014-02-05 | 首都师范大学 | 一种双能谱ct图像迭代重建方法 |
| US20150363947A1 (en) * | 2014-06-16 | 2015-12-17 | The University Of Chicago | Spectral x-ray computed tomography reconstruction using a vectorial total variation |
| CN105488826A (zh) * | 2015-12-17 | 2016-04-13 | 首都师范大学 | 一种基于fbp的能谱ct迭代成像方法和成像系统 |
-
2017
- 2017-12-04 CN CN201711261688.9A patent/CN108010099B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103559699A (zh) * | 2013-11-18 | 2014-02-05 | 首都师范大学 | 一种基于投影估计的多能谱ct图像重建方法 |
| CN103559729A (zh) * | 2013-11-18 | 2014-02-05 | 首都师范大学 | 一种双能谱ct图像迭代重建方法 |
| US20150363947A1 (en) * | 2014-06-16 | 2015-12-17 | The University Of Chicago | Spectral x-ray computed tomography reconstruction using a vectorial total variation |
| CN105488826A (zh) * | 2015-12-17 | 2016-04-13 | 首都师范大学 | 一种基于fbp的能谱ct迭代成像方法和成像系统 |
Non-Patent Citations (1)
| Title |
|---|
| 赵云松 等: "双能谱CT的迭代重建模型及重建方法", 《电子学报》 * |
Cited By (20)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109009181B (zh) * | 2018-06-07 | 2024-04-05 | 西安交通大学 | 双能量ct下同时估计x射线球管光谱和重建图像的方法 |
| CN109009181A (zh) * | 2018-06-07 | 2018-12-18 | 西安交通大学 | 双能量ct下同时估计x射线球管光谱和重建图像的方法 |
| CN109146994A (zh) * | 2018-09-17 | 2019-01-04 | 南京航空航天大学 | 一种面向多能谱x射线ct成像的金属伪影校正方法 |
| CN109146994B (zh) * | 2018-09-17 | 2023-06-20 | 南京航空航天大学 | 一种面向多能谱x射线ct成像的金属伪影校正方法 |
| CN109697691A (zh) * | 2018-12-27 | 2019-04-30 | 重庆大学 | 一种基于l0范数和奇异值阈值分解的双正则项优化的有限角投影重建方法 |
| CN109875593B (zh) * | 2019-02-02 | 2023-03-21 | 中北大学 | X射线多谱分离成像方法、存储介质和装置 |
| CN109875593A (zh) * | 2019-02-02 | 2019-06-14 | 中北大学 | Ct成像方法、存储介质和装置 |
| CN110942494A (zh) * | 2019-02-08 | 2020-03-31 | 苏州波影医疗技术有限公司 | 利用三个数据域迭代重建计算机断层扫描图像的系统 |
| CN110942494B (zh) * | 2019-02-08 | 2024-02-20 | 苏州波影医疗技术有限公司 | 利用三个数据域迭代重建计算机断层扫描图像的系统 |
| CN111127477A (zh) * | 2019-12-12 | 2020-05-08 | 明峰医疗系统股份有限公司 | 一种ct多色谱体素仿真方法 |
| CN111127477B (zh) * | 2019-12-12 | 2023-10-20 | 明峰医疗系统股份有限公司 | 一种ct多色谱体素仿真方法 |
| CN111968060B (zh) * | 2020-08-28 | 2022-07-08 | 首都师范大学 | 一种基于倾斜投影修正技术的多能谱ct快速迭代重建方法 |
| CN111968060A (zh) * | 2020-08-28 | 2020-11-20 | 首都师范大学 | 一种基于倾斜投影修正技术的多能谱ct快速迭代重建方法 |
| CN112200882B (zh) * | 2020-10-10 | 2022-11-25 | 首都师范大学 | 基于傅里叶变换微分性质的板状物ct图像快速重建方法 |
| CN112200882A (zh) * | 2020-10-10 | 2021-01-08 | 首都师范大学 | 基于傅里叶变换微分性质的板状物ct图像快速重建方法 |
| CN112288762B (zh) * | 2020-10-15 | 2023-05-09 | 西北工业大学 | 一种有限角度ct扫描的离散迭代重建方法 |
| CN112288762A (zh) * | 2020-10-15 | 2021-01-29 | 西北工业大学 | 一种有限角度ct扫描的离散迭代重建方法 |
| CN112581556A (zh) * | 2020-12-25 | 2021-03-30 | 上海联影医疗科技股份有限公司 | 多能ct图像硬化校正方法、装置、计算机设备和存储介质 |
| CN114916950B (zh) * | 2022-07-21 | 2022-11-01 | 中国科学院深圳先进技术研究院 | 基于多层平板探测器的高空间分辨能谱ct图像重建方法 |
| CN114916950A (zh) * | 2022-07-21 | 2022-08-19 | 中国科学院深圳先进技术研究院 | 基于多层平板探测器的高空间分辨能谱ct图像重建方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN108010099B (zh) | 2020-12-25 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108010099B (zh) | 一种x射线多能谱ct有限角扫描和图像迭代重建方法 | |
| US10716527B2 (en) | Combined sinogram- and image-domain material decomposition for spectral computed tomography (CT) | |
| US10945695B2 (en) | Apparatus and method for dual-energy computed tomography (CT) image reconstruction using sparse kVp-switching and deep learning | |
| US11864939B2 (en) | Apparatus and method that uses deep learning to correct computed tomography (CT) with sinogram completion of projection data | |
| GB2551029B (en) | Method for performing material decoposition using a dual-energy X-Ray CT and corresponding Dual-energy X-Ray CT Apparatus | |
| CN108010098B (zh) | 一种双能谱ct基材料图像迭代重建方法 | |
| US10621756B2 (en) | Apparatus and method for correcting bias in low-count computed tomography projection data | |
| US10111638B2 (en) | Apparatus and method for registration and reprojection-based material decomposition for spectrally resolved computed tomography | |
| CN110189389B (zh) | 基于深度学习的双能谱ct投影域基材料分解方法及装置 | |
| US10360677B2 (en) | Apparatus and method for joint-edge-preserving regularization to reduce noise in four-dimensional computed tomography images | |
| Zhao et al. | Patient-specific scatter correction for flat-panel detector-based cone-beam CT imaging | |
| WO2019067524A1 (en) | TD MONOCHROMATIC IMAGE RECONSTRUCTION FROM CURRENT INTEGRATION DATA VIA AUTOMATIC LEARNING | |
| Zhang et al. | Noise correlation in CBCT projection data and its application for noise reduction in low‐dose CBCT | |
| US9875527B2 (en) | Apparatus and method for noise reduction of spectral computed tomography images and sinograms using a whitening transform | |
| Lin et al. | A fast poly-energetic iterative FBP algorithm | |
| Yao et al. | Multi-energy computed tomography reconstruction using a nonlocal spectral similarity model | |
| WO2019060688A1 (en) | SYSTEM AND METHOD FOR LOW DOSE MULTISPECTRAL X-RAY TOMOGRAPHY | |
| Ding et al. | A high-resolution photon-counting breast CT system with tensor-framelet based iterative image reconstruction for radiation dose reduction | |
| Sheng et al. | A sequential regularization based image reconstruction method for limited-angle spectral CT | |
| Qiu et al. | Evaluating iterative algebraic algorithms in terms of convergence and image quality for cone beam CT | |
| Hu et al. | A practical material decomposition method for x-ray dual spectral computed tomography | |
| Tang et al. | Analytical fan-beam and cone-beam reconstruction algorithms with uniform attenuation correction for SPECT | |
| Zhao et al. | Iterative dual energy material decomposition from spatial mismatched raw data sets | |
| Wang et al. | A framelet-based iterative maximum-likelihood reconstruction algorithm for spectral CT | |
| US10799192B2 (en) | Method and apparatus for partial volume identification from photon-counting macro-pixel measurements |
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 |