CN106462984A - 用于差分相位对比成像中的谱相位展开的无偏置正则化 - Google Patents

用于差分相位对比成像中的谱相位展开的无偏置正则化 Download PDF

Info

Publication number
CN106462984A
CN106462984A CN201580029369.8A CN201580029369A CN106462984A CN 106462984 A CN106462984 A CN 106462984A CN 201580029369 A CN201580029369 A CN 201580029369A CN 106462984 A CN106462984 A CN 106462984A
Authority
CN
China
Prior art keywords
regularization
critical point
object function
phase shift
phase
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.)
Pending
Application number
CN201580029369.8A
Other languages
English (en)
Inventor
T·克勒
R·普罗克绍
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.)
Koninklijke Philips NV
Original Assignee
Koninklijke Philips Electronics NV
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 Koninklijke Philips Electronics NV filed Critical Koninklijke Philips Electronics NV
Publication of CN106462984A publication Critical patent/CN106462984A/zh
Pending legal-status Critical Current

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/002D [Two Dimensional] image generation
    • G06T11/003Reconstruction from projections, e.g. tomography
    • G06T11/005Specific pre-processing for tomographic reconstruction, e.g. calibration, source positioning, rebinning, scatter correction, retrospective gating
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/0002Inspection of images, e.g. flaw detection
    • G06T7/0012Biomedical image inspection
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/48Diagnostic techniques
    • A61B6/484Diagnostic techniques involving phase contrast X-ray imaging
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10116X-ray image
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2211/00Image generation
    • G06T2211/40Computed tomography
    • G06T2211/408Dual energy

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Radiology & Medical Imaging (AREA)
  • Quality & Reliability (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Apparatus For Radiation Diagnosis (AREA)

Abstract

一种方法和相关的装置(SP),用于针对缠绕伪影来校正相位偏移图像数据。在包括干涉测量装备(G0、G1、G2)的成像系统(IM)的探测器(D)处探测所述数据。在涉及利用正则化优化目标函数的相位展开方法中,提出了两阶段方案。在具有正则化的一个阶段中和没有正则化的另一阶段中处理所测量的数据。这允许改进经校正的(相位展开的)相位偏移数据的准确度,因为能够避免由所述正则化所引起的不期望的偏置。

Description

用于差分相位对比成像中的谱相位展开的无偏置正则化
技术领域
本发明涉及用于处理测量的相位偏移数据的方法,涉及信号处理装置,涉及计算机程序单元并且涉及计算机可读介质。
背景技术
在包括医学领域的技术应用的范围中,相位对比成像因为其极好的软组织对比度已经被认为优于仅基于吸收的成像技术。在一些应用中,相位对比成像需要在成像装置中安装的干涉仪。所述干涉测量装置允许提取干涉图样的特定相位偏移,其已经被发现与期望的折射相关,并且因此与由要被成像的对象所引起的相位偏移相关。遗憾的是,经由干涉仪获得的相位偏移数据可能是模糊的。这是因为干涉仪仅具有有限的角范围。超出该范围的任何测量结果将被“缠”绕到该范围中,使得事实上不是相位偏移自身被观测到,而仅仅是测量的相位偏移值“对所述范围取模”。相位缠绕是不期望的,因为其可以在相位对比图像中引起可视的伪影。为了防止这种不明确,在过去已经提出了若干种“相位展开”算法。然而,已经发现,通过现有算法可获得的展开相位数据的准确性偶尔是不准确的。
发明内容
因此,存在改善在处理测量的相位偏移数据时的准确度的需求。
本发明的目的由独立权利要求的主题来解决,其中,在从属权利要求中并入了另外的实施例。
应当注意,下文所描述的本发明的各方面同样适用于信号、信号处理装置、计算机程序单元以及计算机可读介质。
根据本发明的第一方面,提供了一种用于展开相位偏移数据的方法,包括:
-接收由探测器测量的相位偏移数据;
-将所述相位偏移数据作为参数包括到没有正则化项的非凸目标函数或非凹目标函数;
-找到所述目标函数的临界点;
-仅在所述临界点之中找到正则化的目标函数的目标临界点,所述正则化的目标函数是根据所述目标函数和正则化项形成的;并且
-输出所述目标临界点作为针对展开的相位偏移数据的估计。
根据第二方面,提供了一种用于展开相位偏移数据的备选方法,包括:
-接收由探测器测量的相位偏移数据;
-将所述相位偏移数据作为参数包括到根据非凸目标函数或非凹目标函数和正则化项形成的正则化的目标函数中;
-找到所述正则化的目标函数的临界点,所述临界点包括引导临界点;
-找到没有所述正则化项的所述目标函数的非正则化的临界点;
-从这样找到的非正则化的临界点之中选取目标点,所述目标点比至少一个其他非正则化临界点中的至少一个更接近所述引导临界点;并且
-输出所选取的目标点作为针对展开的相位偏移数据的估计。
换言之,所述方法针对相位缠绕伪影来校正输入的相位偏移数据,并且输出如此校正的数据,作为展开的相位偏移数据。
在本文中所使用的目标函数是如从最优化理论的领域获知的。换言之,所述目标函数是实值函数,其将所述函数的域的元素映射到实数上。所述目标函数取决于背景有时也被称为效用函数或代价函数。那些函数的临界点是函数的域中的其中目标函数假设局部(或全局)最大值或最小值的参数或位置。因此,所述临界点根据背景还被称为最大化者或最小化者。在本文中设想了目标函数的最大化和最小化中的任一个。在一个实施例中,在所述方法中的任一种方法中搜索临界值等于通过其他已知的数值算法来求解针对所述临界点的目标函数。
所提出的方法允许避免或者至少缓解由正则化项的存在而引入的偏置效应。正则化子或正则化类似地是实值函数,其通常仅仅是目标函数的域中的元素的函数。其允许加强期望的性质,例如对要搜索的临界点的平滑或“缩小”。已经观测到,如果目标函数具有多个局部最小值(亦即,目标函数既不是凸的,也不是凹的),则所述正则化项施加偏置,从而将所述目标函数的全局最小值偏移到域的正确范围中。这是期望的,但目标函数还朝向所述范围在所述目标函数的局部极值上引入不期望的偏置(例如,小的绝对值)。由于所搜索的被展开的相位数据在那些局部极值之中,因此该偏置贯穿于基于正则化的优化的许多相位展开算法的输出。
为了克服该偏置效应,在本文中提出了两种备选方法。根据一种方法,首先搜索没有正则化的目标函数,并且然后,结果被应用到正则化的函数。在备选实施例中,首先搜索正则化的函数,并且然后,结果被应用到未正则化的目标函数。换言之,反转搜索步骤的顺序。在这两种方法中,输出(亦即,展开的相位偏移数据)最终是未正则化的目标函数的点。换言之,正则化的目标函数的临界点不是这样返回作为展开的相位偏移数据。相反,正则化的目标函数的临界点被用作引导手段,以从未正则化的目标函数的临界点中进行选择。以这种方式,由正则化项的存在引起的在搜索步骤中的偏置能够被避免或者至少被缓解,并且再者是能够考虑正则化要求。又换言之,所述方法允许在观察对具有正则化和避免或缓解由所述正则化项这样引起的偏置的需求之间寻求平衡。
根据一个实施例,目标函数是关于测量的数据的似然函数。在所述似然函数潜在的概率密度是以下中的任意一个:缠绕的高斯密度、冯米塞斯密度、或者适于对诸如角数据的“方向性”数据的概率建模的其他合适的概率密度。
但是,并不是说所述目标函数必须具有概率来源。所提出的方法可应用于使用对目标函数的优化的任何相位展开算法,而不管目标函数是如何获得的,其通过统计/概率或其他信号建模方法。
在根据第二方面的方法中,根据一个实施例,所述目标点与引导点最接近(当由合适的距离量度来测量时)。然而,这不是强制的,因为设想了备选实施例,其中,当考虑所有未正则化的临界点时,目标点仅仅处在引导点的预定义附近内是足够的,从而可能不是与引导点“最接近”。
在一个实施例中,所述正则化项是吉洪诺夫正则化子。这允许解释如下事实:由于探测器的有限的空间分辨率,所测量的数据可能是欠采样的。
根据一个实施例,在所述方法中的任一种方法中,在给定像素处的测量的数据独立于在(任意或所有)邻近像素处测量的数量来处理。换言之,其中,当所述方法被应用于在所述给定像素处测量的数据时,所述方法不依赖于在探测器D的给定像素的邻近像素处测量的数据。
“引导”临界点是这样临界点:在所述临界点处,当由(正则化的)目标函数来评价时,所述目标函数返回比至少一个其他临界点中的至少一个“更好”的值。限定词“更好”取决于手头上的优化任务。例如,在最小化问题中,值越低越好。相反地,在最大化问题中,值越大越好。具体而言并且根据一个实施例,所述引导临界点可以是这样的一个临界点:在所述临界点处,所述目标函数假设所有其他临界点的最大值或最小值。但这不是强制性的,因为设想了备选实施例,在所述备选实施例中,引导点仅仅返回比预定义阈值更好的值或者返回比预定义数量的临界点完成的更好的值是足够的。
如在本文中所使用的“临界点”是这样的点:其中,所述目标函数假设局部(全局)最小值或最大值。如在本文中所使用的,针对任意给定的优化问题,所述“临界点”遍及局部最小值或局部最大值。在本文中未设想混合的情形,其中,针对给定的优化,临界点一次是局部最小值,并且然后是局部最大值。
附图说明
图1示出了用于谱差分相位对比成像的成像布置;
图2示出了由相位缠绕引起的图像噪声或伪影的范例;
图3示出了根据不同实施例的信号处理器的框图;
图4是正则化的和未正则化的目标函数的曲线图;
图5示出了根据处理相位偏移数据的不同实施例的方法的流程图。
具体实施方式
图1示出了具有相位对比成像能力,尤其是差分相位对比成像(DPCI)的成像系统IM的基本部件。存在X射线源XR,其用于生成X射线辐射波XB,穿过检查区域之后的X射线辐射波能由探测器D的探测器像素px探测。所述相位对比成像能力是由在X射线源XR与辐射敏感探测器D之间布置的干涉仪来实现的。
所述干涉仪(其在一个非限制性实施例中为Talbot类型或者为Talbot-Lau类型的)包括两个G1、G2(Talbot类型)或者更多个,优选三个光栅G0、G1和G2(Talbot-Lau类型)。在X射线源侧的第一衰减光栅G0具有周期p0,以匹配和引起在X射线源XR处发射的X射线辐射波前的空间相干。
吸收光栅G1(具有周期p1)被放置在距X射线源的距离D处,并且引起更下游的具有周期p2的干涉图样。所述干涉图样能够由探测器D探测。现在,当(要被成像的)样本被引入到在X射线源与探测器之间的检查区域中时,那么干涉图样的相位被偏移。由于沿着通过样本的各自路径的累积的反射(因此,名为DCPI),该干涉图样偏移(如在他处已经报告的,例如在F M Epple等人,Unwrapping differential X-ray phase contrast imagesthrough phase estimation from multiple energy data,OPTICS EXPRESS,2013年12月2日,第21卷,第24期中)与由于沿通过样本的各自的路径的积累的折射的相位偏移的梯度ΦΔ成比例。换言之,如果然后是测量干涉的相位改变,则这将允许提取由样本中的折射所引起的相位偏移的偏移(或梯度)。
遗憾的是,干涉图样的相位偏移通常过小而不能直接空间地分辨。大多数X射线探测器的分辨率能力不会允许这种情况。因此,为了对该干涉图样相位偏移进行“采样”,与干涉图样具有相同周期p2的第二衰减光栅G2被放置在距光栅G1距离l处。根据在本文中所有预想的不同的实施例,对干涉图样相位偏移(并且因此,由样本引起的相位梯度的)的实际提取能够以多种不同的方式来实现。
一般而言,针对差分相位提取所要求的是在探测器D与光栅的至少一个光栅之间的相对运动。在一个实施例中,这能够通过使用致动器来跨不同的、离散的光栅位置横向地(亦即,沿着与光栅垂直的x方向)移动例如分析器光栅G2,并且然后在每个光栅位置处测量在每个像素px处的强度。在每个像素处的强度将被发现为以正弦方式震荡。换言之,在分析器光栅G2的运动期间,每个像素记录作为时间的函数(或者更好地,作为不同光栅位置的函数)的一时间序列的不同强度(在各自的像素处)。该方法(“相位步进”)已经由F.Pfeiffer等人在“Phase retrieval and differential phase-contrast imaging with low-brilliance X-ray sources”Nature Phys.Lett.2,258–261(2006)中进行了描述。
如例如在早先引用的参考Epple中的在29104页上的等式(1a)、(1b)所描述的,在每个像素px处的震荡强度信号在其他量之中对干涉图样的期望的相位偏移(相位偏移,下文简单地称为)进行“编码”。
干涉图样偏移能够通过由数据采集电路(未示出)将每像素px的探测器信号处理成数值形式来获得,所述数据采集电路具体而言包括A/D转换电路。所述数值然后由傅里叶分析器模块(未示出)来处理,或者由被配置为执行合适的曲线拟合操作的模块来处理,以获得每像素的测量的相位偏移这些参数的集合然后能够由处理单元来处理。
在图1的成像系统IM中的探测器D是光子计数类型的。探测器D包括电路ASIC,其赋予探测器布置的光子计数和能量分辨能力。该能力允许谱成像,因为每个能量响应能够与被认为构成被成像对象的特定基材料(骨骼、水、脂肪、对比剂等)相关联。将在任意给定像素px处拾取的脉冲与特定的能量阈值(也被称为“分箱”)进行比较,所述能量阈值表示不同能量水平Ej=E1、E2、E3。如果脉冲超过特定阈值,那么设置针对所述能量阈值的计数器,并且例如由傅里叶分析器将所述信号处理为针对特定能量E1、E2或E3的相位偏移
换言之,在探测器的光子能量分辨实施例中,不同能量水平Ej处测量的相位偏移数据通过不同的能量通道E1到E3(在附图中仅示出了三个这样的能量水平E1、E2和E3,但这仅仅是出于说明的目的,并且在本文中设想了大于二的任意其他数量的能量水平)来提供。
然后,由信号处理器SP每像素地处理针对每个能量通道的干涉图样的相位偏移角如将在下文进一步详细描述的,以展开相位模糊。单个处理器SP能够运行在诸如工作DPU的数字处理单元上。然后,可以通过合适的成像软件将展开的相位角绘制成相位对比图像以用于在屏幕MT上显示,或者其可以以其他方式处理或存储。
相位缠绕,如简要提到的,是这样的一种现象:其尤其可以在差分相位对比成像内部发生,因为干涉仪的动态范围通常被限于特定角范围,例如,在-π与π之间。换言之,干涉仪(亦即,三个(或者至少两个)光栅G0到G2的系统)实际上不允许测量这样的干涉图样的相位偏移,而是仅真实相位偏移“以2π求模”:干涉仪对范围外的相位偏移的响应对应于将真实的相位偏移(表达为以弧度为单位的长度)围绕单位圆许多次,只要长度允许,并且然后读出在单位圆的周长上的角位置,其中,长度的缠绕的端部处于所述角位置处。关于缠绕的数量的信息丢失。仅仅端部的所述角位置可用作为干涉仪的响应。相位展开算法旨在根据缠绕的相位和额外的信息或模型假设来恢复丢失的信息。信息的这种丢失是不期望的,因为其在相位对比图像数据中引入“伪影”不连续,这继而导致图像噪声,所述图像噪声能够潜在地掩盖影像中的关键细节。在图2的窗格A)和B)中的示范性相位对比影像示出了会出现什么样的错误。在左侧上的图2A)示出了以25keV成像的体模体(在这种情况下,为人脚),其结构引起折射梯度,所述折射梯度将超出干涉仪的动态范围。窗格A)示出了常规吸收X射线图像。右侧窗格B)示出了对应的相位对比图像。箭头对指向具有由相位缠绕不连续所引起的噪声累积或伪影的图像部分。在图像B)中的其他图像部分看上去足够平滑,因此,此处不存在相位缠绕。相位缠绕在诸如重建的后续处理期间或者在相位缠绕的数据与吸收对比图像数据组合时导致甚至更多的扰动伪影。
因此,期望关于相位偏移来分辨模糊,以展开所测量的(角)相位偏移数据。已经发现,在许多基于探测器的系统中,不能够保证差分相位对比图像是根据香农采样条件实际采样的。探测器的空间分辨率(亦即,在邻近像素之间的内部空间)通常不足够小而不能够满足(空间)采样定理。这一事实已经被申请人发现,如果使用标准的相位展开算法,则严重地破坏结果的保真度。
参考图3,示出了用于分辨相位模糊的信号处理器SP的两个实施例A、B,亦即,这两者被配置为操作为“相位展开器”。更为具体地,如在本文中所提出的信号处理器旨在用于与现有的相位展开算法(其涉及求解最优问题)相组合。
在先前引用的Epple参考文献的第29105页上的等式(5)中描述了这样的相位展开算法的范例。类似于其他相位展开算法,Epple描述了发现特定目标函数f的最大化者或最小化者。所述目标函数总体而言是代价函数或效用函数,其将预定义可允许范围的元素(所述函数的“域”,在这种情况下为以弧度为单位测量的角度值)映射为数值。所述目标函数是标量值函数,其定义2D中的曲线或者3维或更高维度中的表面,其通常包括多个局部最小值或最大值。搜寻的最优的点是其中发生局部最小值或最大值(统称为局部“极值”)的那些点。那些点将被称为最大化者(针对最大化问题)或最小化者(针对最小化问题),或者统称为“临界点”,但是,针对给定优化问题,所述临界点或者是局部最小值,或者是最大值,但不是这两者。点是否是“临界”的点不仅取决于函数做什么,而且取决于特定的域,亦即,所述函数被应用于或者在其处被评价的元素的范围上。所述域可以是连续的(例如,实值),但是当将所述函数限制到离散的子域(所有的整数表示)时,那么可以获得关于哪个点是否是临界的点的不同的结果。针对临界点的充分但不必要的条件是,其中,所述目标函数是可微分的,并且所述函数的其梯度或导数消失。对所述目标函数的邻近的简单采样允许决定所找的临界点是否是最小化者或最大化者。
在下文中,当考虑最大化问题(亦即,找到、搜索或求解最大化者)时,应当理解,所有所述的内容同样适用于相反的问题,亦即,最小化问题,如每个最大化问题能够被重写为最小化问题,并且反之亦然。换言之,如本文中所提出的方法和信号处理器被设想用于最大化或最小化中的任一个。
返回参考相位展开算法的具体情况,所测量的和潜在缠绕的相位角被应用于考虑中的特定类型的目标函数,并且然后,数值算法搜索临界点。所述临界点然后可以被输出(可能在一些重标度之后),作为展开的相位,其然后能够被传递到重建器以用于CT图像重建或者一些其他的图像处理后端。
又更为具体地并且如在Epple中示范性描述的,所提出的信号处理器SP在其上操作的目标函数可以基于信号建模方法,例如,统计学方法,诸如最大似然(ML)或其他考虑。
所述目标函数根据所采用的信号建模方法来公式化,并且不同能量水平E1、E2、E3处的所测量的(并且可能缠绕的)相位偏移被作为参数应用到所述目标函数。对应的展开的相位偏移数据被包括在目标函数中,作为希望通过合适的数值算法或技术来求解或“搜索”的自由变量。换言之,所述目标函数将所测量的、可能缠绕的相位偏移数据与(还未知的)展开的相位偏移数据相关。信号处理装置SP操作于通过以要在下文在图3和图5处更为详细解释的方式应用合适的数值算法来计算展开的相位偏移数据。
作为对上文的扩展,多个已知相位展开算法(诸如在Epple中一个)不仅仅求解目标函数的临界点,而是求解正则化的目标函数的临界点。该正则化通过使用正则化子项R产生。正则化子项或者有时被称为惩罚函数是这样的项,其与在考虑中的目标函数算数地组合,以形成正则化的目标函数。惩罚或正则化项被用于加强针对优化问题的解(亦即,临界点)的特定期望的性质。例如,可以希望满足特定的平滑度条件,其在使用无正则化项的目标函数时将排除其他可能的解决方案。在其他实施例中,可能希望超过大的一些地支持针对优化问题的“小的”解,其中,通过合适的范数(诸如在惩罚项中并入的L1或L2)来测量小的程度。
在图4中图解地显示了正则化子项R对目标函数f的效应,其中,考虑了最小化问题。图4示出了具有不同局部最小值和最大值的未正则化的目标函数f的曲线图,换言之,函数f是非凸的,并且其结果是,如在本文中所提出的方法和信号处理器是针对这样的非凸/非凹函数的具体应用。如能够看到的,在函数的域(其在这种情况下被认为是相位角的范围)的特定元素处出现函数f的局部最小值。如果现在将该函数与正则化项R进行组合(例如,通过形成两个的函数和)以形成正则化的目标函数f+R,其曲线图/表面被朝向原点有效地挤压。使用正则化项R的意图(在最小化的背景下)是与其距原点的距离成比例地“惩罚”局部最小值。换言之,支持“小的”最小值,亦即,那些更接近原点的。但是,遗憾的是,还存在另一不想要的效应:正则化子R的应用还引起偏置b,其中,所述临界点自身被偏移到原点。该偏置已经被观察到有时引起不准确。如在本文中所提出的根据图3的信号处理器操作于移除该偏置,并且因此保留具有正则化函数的益处,但具有返回的解针对原始问题更准确的增加的益处。再次地,已经在图4中关于最小化所述的所有内容经过必要的修正同样适用于最大化,并且在本文中设想了非凹函数和这些实施例。
简要地,并且如在图3中能够看到的,信号处理器SR的实施例A)、B)包括合适的接口IN和OUT,以分别接收每像素的(可能相位缠绕损坏的)所测量的相位角并且输出展开的相位角估计(针对所述像素),其中,针对一些或所有测量的值,相位缠绕损坏已经被移除或者至少被减轻。
信号处理器SP包括至少两个优化模块R和NR。一个模块R操作为优化器以搜索正则化的目标函数的临界点,而另一模块NR求解或操作为搜索未正则化的目标函数的临界点。
宽泛地,如在图3A)、B)中示意性示出的,能够反转在处理流水线中的两个优化模块NR和R的顺序。换言之,在图A的实施例中,相位角被首先应用到未被正则化的目标函数,并且进一步在下游,该第一阶段优化的结果在第二阶段处被进一步处理,在第二阶段中,所述目标函数被正则化。在图3的实施例B)中,该处理流水线被反转,其中,所接收的相位角被首先应用于被正则化的目标函数,并且来自于其的输出然后被馈送到第二模块,在第二模块中,所述目标函数未被正则化。简言之,在本文中所提出的是两阶段信号处理模块,其包括两个优化器,一个用于关于无正则化项的目标函数的优化,并且一个用于关于具有正则化项的目标函数的优化。两个优化器模块能够以任一顺序来应用。现在将首先参考图5中的流程图A)来解释图3的信号处理器SP的操作。
在步骤S305处,(潜在相位缠绕的)能量分辨的相位角被接收作为通过成像系统IM的干涉成像装备G0-G2在探测器D处“看到”的数据。
在步骤S310处,这样接收的相位偏移数据被应用到无正则化项的目标函数F,亦即,利用所测量的数据来填充针对参数的占位符。如在之前已经提示的,关于由正则化项的存在所引起的朝向原点的临界点的偏移的复杂性能够仅在目标函数的表面中存在两个或更多个最小值时应用。因此,下文本发明尤其应用于非凸的(或者非凹的,取决于优化器是否操作用于执行最大化或最小化)目标函数。
返回参考图5的流程图A),在步骤S315处,信号处理器SP的第一模块NR搜索所述(非凸或非凹)目标函数的临界点。注意,在处理流水线中的该阶段处不涉及正则化项R。在S315处的搜索或求解操作能够通过针对那些临界点进行有效求解的现有数值最小化或最大化算法中的任意算法来实现。在本文中设想的范例是下山单纯形方法、(非线性)共轭梯度方法、牛顿-拉夫逊方法等,合适的方法能够被选取以适用于关于手头的特定目标函数(其结构一般由所选取的信号模型来确定)的不适定性。示范性实施例是ML方法(给定合适类型的概率密度)的统计学方法,其中,出于数值原因,通常优化似然函数的负对数。在本发明的探测器的背景下,根据一个实施例,(缠绕的)高斯或冯·米塞斯(von Mises)类型的密度被建模为信号模型。
这样找到的临界点然后被转发到阶段S320至第二优化器R,其中,现在关于组合的目标函数,亦即,正则化的目标函数,来处理所述临界点,所述正则化的目标函数是根据在模块NR的阶段S315处的先前目标函数和正则化项来形成的。然而,与步骤S310相对,该时间,先前连续的域现在是离散的,亦即,被限制到在步骤S315处找到的临界点的集合的谨慎空间。换言之,在步骤S320中,仅在临界点之中进行搜索,以揭示正则化的目标函数的目标临界点。实现在步骤S320中的现在正则化的优化的一种方式是在临界点的离散并且有限的集合处评价正则化的目标函数,并且找到这样的一个:其中所述目标函数假设或“返回”最小值,或者至少其中所述目标函数假设小于可预定义的可接受阈值的值。如果存在满足该条件的更多的点,亦即,如果超过一个的临界点返回正则化的目标函数的落在所述阈值之下的值,那么可以实行随机实验以返回所述点中的单个点,其作为量化点在用户接口中被呈现给人类用户,以供最终的选择。例如,能够通过利用量化临界点中的相应的一个替换先前的、测量的相位偏移值,来针对量化临界点中的每个生成样本相位对比图像,并且然后,用户能够针对视觉平滑度来检查样本图像,以如此支持最终的选择。
在步骤S325处,目标点针对相关像素被输出,作为展开的相位偏移数据,或者直接地,或者在一些会话或缩放之后,这取决于所使用的具体优化算法和/或目标函数。
根据一个实施例,针对每个单像素数据,独立地处理(但可能并行地或者以序列化方式)以上步骤S305-S325。换言之,在本方法中不涉及对针对任何给定像素或过程的邻近数据的评价。这不同于诸如在申请人的WO 2012/038857中公开的其他方法,在所述其他方法中,为了达到针对像素的展开的相位值,需要针对所述像素的邻近像素值的评价。
现在参考图5中的流程图B),其对应于由根据图3的实施例b)的信号处理单元可实施的方法。
该方法类似于在实施例A)中的先前的方法,但是现在,步骤S315和S320的顺序在处理流水线中实质上被反转。
更为详细的,在步骤S405处,接收针对任意给定像素的相位偏移数据。
所测量的数据然后在步骤S410处被应用到由模块R实施的正则化的目标函数。正则化的目标函数由非凸目标函数和正则化项来形成,类似于上文在步骤S310处已经解释的。
不同于步骤S315,在步骤S415处进行搜索以首先在处理流水线中使用模块R来首先找到正则化的目标函数的临界点。在临界点之中找到的是“引导”临界点,亦即,“最佳地”满足正则化的目标函数或者至少好于如此找到的其他临界点中的一些或所有的一个临界点。例如,当所述正则化的目标函数在所述引导临界点处被评价时,所述引导临界点可以返回最小值,或者在这样评价时返回最高值。换言之,所述引导临界点在两个意义上“引导”。其关于吸引最低的代价或者具有最高的利用率方面是领先的,取决于优化任务如何被公式化。但是,所述引导点“引导”或指引对期望的展开的相位偏移数据的搜索,如现在将参考接下来的两个步骤S420、S425解释的。在一些实施例中,所述引导点可能不是全局最佳的,但可以是仅仅要求是较好的,亦即,对正则化的目标函数返回比大多数其余临界点更有利的响应,或者可以返回比预定义阈值更好的响应。
处理流程现在前进到步骤S420处的模块NR,以找到“非正则化的”临界点,亦即,未正则化的目标函数的临界点,亦即,现在,在步骤S415中的目标函数在没有正则化项的情况下被孤立地优化。
然后,在步骤S425中,未正则化的优化问题是通过从在步骤S420中找到的非正则化的临界点之中选取这样的一个来解决的:其最接近当求解正则化的优化问题时在步骤S415中找到返回的引导临界点。接近度是由诸如L1或L2的一些合适的范数来测量的。在一些实施例中,可以存在超过一个“最接近的”点,因为仅仅要求非正则化的临界点比可预定义的距离阈值更接近所述引导临界点。在后一种情况下,随机算法可以运行到绘制算法决定性。备选地,呈现了用户接口,用于如上文针对实施例A)所描述的最终选择。
然后,在步骤S430处,可能地在应用合适的转换或缩放操作之后,在步骤S425中选取的点然后可以被输出作为展开的相位偏移数据。
如将从两个实施例意识到的,在这两种情况下,输出的点总是未正则化的优化问题的那些问题。然而,正则化的优化问题仍被实行以确保满足正则化要求。但是,通过优化目标函数的正则化的版本找到的临界点仅被用作指引手段或“线索”,以从关于未正则化的目标函数搜索的未正则化的临界点中进行选择。这确保了在步骤S325或S430处输出的最终结果均没有由于正则化项的干扰的偏置。
如上文简要提到的,正则化项和目标函数的特定形式将取决于所使用的特定的相位展开算法和/或潜在信号模型。本方法因此可以被认为如添加到现有的相位展开算法,其涉及关于非凸目标函数或非凹目标函数的优化问题。
同样地,在两个模块R和NR中的被用于找到临界点的数值技术可以是相同的或者可以不是相同的。
例如,根据一个实施例,并且类似于已经由Epple所描述的内容,具有正则化项R的目标函数可以具有如下结构:
其中:
w指示不同的能量分箱,
是所测量的潜在缠绕的相位偏移,
M是要估计的展开的相位偏移
Ew是能量分箱w的有效能量,
是所测量的相位偏移的标准偏差,
R是正则化函数,并且
λ是正则化参数以控制正则化子R的响应的强度。
模型参数Ew、λ可以通过校准测量来获得,并且在测试测量中针对最佳的结果调节λ。
在(1)中的潜在信号模型利用针对展开问题的相位梯度的公知的能量相关性,尽管以上模型仍然能够与在Epple中所描述的单色辐射一起使用。
依据(1)的以上目标函数是根据最大似然方法导出的。展开的相位偏移数据是使用最大似然函数的估计的相位t。更为具体地,通过对所述最大似然函数的对数的最小化来估计所述相位梯度数据。
已经假设,通过冯·米塞斯密度来支配相位数据测量探测器。然而,在本文中也可以包括其他密度,例如,替换地,可以假设缠绕的高斯密度。T Weber等人在“Measurementsand simulations analyzing the noise behavior of grating-based X-ray phase-contrast imaging”,Nuclear Instruments and Methods in Physics Research A648(2011),第S273–S275页中证实了(参见Weber的图4),冯·米塞斯以足够的准确度描述相位偏移测量结果,并且T Weber等人在“Noise in x-ray grating-based phase-contrastimaging”,Med.Phys.38(7),2011年7月,第4133-4140页中讨论了,冯·米塞斯分布的各自的变化,其能够被用于估计例如上文(1)中的H Gudbjartsson等人在“The RicianDistribution of Noisy MRI Data”,Magn.Reson.Med.,1995年12月,34(6),第910-914页中在等式(5)中描述了针对合适的密度的其他选项。
应当意识到,取决于被用于对所测量的数据的统计进行建模的密度或分布类型,所述似然函数将具有(1)不同的结构,并且因此可以调用其他数值技术来找到临界值。然而,本领域技术人员将意识到,以上两阶段的正则化和非正则化方法能够被应用到任意的目标函数,不管其是否是从最大似然方法导出的。例如,还在本文中设想了基于其他数值方法的目标函数。
如能够看到的,正则化项仅仅是在考虑中的特定样本点M的函数。根据一个实施例,所述正则化项是吉洪诺夫(Tikhonov)正则化子。这允许解释如下事实,即,在图2中的相位梯度(参见窗格B)的整个相位对比图像可能经历混叠,其使得通过在图像域中的平滑化约束的正则化有缺点。相反地,提出了吉洪诺夫正则化,其惩罚M的平方,因此支持针对相位梯度M的小值。
作为范例,将图3的方法A)的步骤S315应用到目标函数(1)要求针对未正则化的对数似然函数的物理自由变量找到在可允许的范围/域中的(在该情况下)局部最小值:
注意,所述正则化其在步骤S315中不存在。
数值地找到的局部最小值的位置由M1、M2、…、Mn来指代。因为在(2)中没有应用正则化,因此针对M的真实值的最佳估计(其在这些值之中)不经历偏置。在随后步骤S320中应用正则化,其中,从离散最小化问题获得适当的“目标”局部最小值:
该优化问题非常类似于原始正则化问题(2),除了现在仅考虑未正则化的对数似然函数的无偏置的局部最小值M1、M2、…、Mn。换言之,在对(2)的优化中考虑的连续域不被限制到由子集{M1、M2、…、Mn}定义的离散域。因此,相位展开的最终的值Mk将不经历或者将经历较少的正则化诱导的偏置。
在备选中,如果要应用如依据图3中的B)的备选方法,首先在针对临界点的连续域中求解等式(1),并且通过在那些临界点处进行评价来建立引导临界点。然后,再次在连续域中求解(2),以求解未正则化的临界点。然后,从所述未正则化的临界点选取这样的点作为输出,即所述点最接近或者至少足够接近(相对于距离阈值)所述引导临界点。
在一个实施例中,图像数据处理装置SP被编程在合适的科学计算平台(诸如)中,并且可以被翻译为适于运行在计算系统(诸如成像器的工作站DPU)上的C++或C进程。
在本发明的另一示范性实施例中,提供了一种计算机程序或一种计算机程序单元,其特征在于适于在适当的系统上运行根据前面的实施例之一所述的方法的方法步骤。
因此,所述计算机程序单元可以被存储在计算机单元上,所述计算机单元也可以是本发明的实施例的部分。该计算单元可以适于执行以上描述的方法的步骤或诱发以上描述的方法的步骤的执行。此外,其可以适于操作以上描述的装置的部件。所述计算单元能够适于自动地操作和/或运行用户的命令。计算机程序可以被加载到数据处理器的工作存储器中。所述数据处理器由此可以被装备为执行本发明的方法。
本发明的该示范性实施例涵盖从一开始就使用本发明的计算机程序或借助于更新将现有程序转变为使用本发明的程序的计算机程序两者。
更进一步地,所述计算机程序单元能够提供实现如以上所描述的方法的示范性实施例的流程的所有必需步骤。
根据本发明的另一示范性实施例,提出了一种计算机可读介质,例如CD-ROM,其中,所述计算机可读介质具有存储在所述计算机可读介质上的计算机程序单元,其中,所述计算机程序单元由前面部分描述。
计算机程序可以存储和/或分布在与其他硬件一起提供或作为其他硬件的部分提供的诸如光学存储介质或固态介质的适当的介质上,但是计算机程序也可以以其他的形式分布,例如经由因特网或其他有线或无线的远程电信系统分布。
然而,所述计算机程序也可以存在于诸如万维网的网络上并能够从这样的网络中下载到数据处理器的工作存储器中。根据本发明的另一示范性实施例,提供了一种用于使得计算机程序单元能够被下载的介质,其中,所述计算机程序单元被布置为执行根据本发明的之前描述的实施例之一所述的方法。
必须指出,本发明的实施例参考不同主题加以描述。具体而言,一些实施例参考方法类型的权利要求加以描述,而其他实施例参考设备类型的权利要求加以描述。然而,本领域技术人员将从以上和下面的描述中了解到,除非另行指出,除了属于一种类型的主题的特征的任何组合之外,涉及不同主题的特征之间的任何组合也被认为由本申请公开。然而,所有特征能够被组合以提供超过特征的简单加和的协同效应。
尽管已经在附图和前面的描述中详细说明和描述了本发明,但这样的说明和描述被认为是说明性或示范性的而非限制性的。本发明不限于所公开的实施例。通过研究附图、说明书和从属权利要求,本领域的技术人员在实践所主张的本发明时能够理解和实现所公开的实施例的其他变型。
在权利要求中,词语“包括”不排除其他单元或步骤,并且,词语“一”或“一个”并不排除多个。单个处理器或其他单元可以履行权利要求书中记载的若干项目的功能。在互不相同的从属权利要求中记载了特定措施并不指示不能有利地使用这些措施的组合。权利要求中的任何附图标记不应被解释为对范围的限制。

Claims (15)

1.一种用于展开相位偏移数据的方法,包括:
接收(S305)由探测器(D)测量的相位偏移数据;
将所述相位偏移数据作为参数包括(S310)到没有正则化项的非凸目标函数或非凹目标函数中;
找到(S315)所述目标函数的临界点;
仅在所述临界点之中找到(S320)正则化的目标函数的目标临界点,所述正则化的目标函数是根据所述目标函数和正则化项而形成的;并且
输出(S325)所述目标临界点作为针对展开的相位偏移数据的估计。
2.一种用于展开相位偏移数据的方法,包括:
接收(S405)由探测器(D)测量的相位偏移数据;
将所述相位偏移数据作为参数包括(S410)到根据非凸目标函数或非凹目标函数和正则化项而形成的正则化的目标函数中;
找到(S415)所述正则化的目标函数的临界点,所述临界点包括引导临界点;
找到(S420)没有所述正则化项的所述目标函数的非正则化的临界点;
从如此找到的所述非正则化的临界点之中选取(S425)目标点,所述目标点比至少一个其他非正则化的临界点中的至少一个更接近所述引导临界点;并且
输出(S430)所选取的目标点作为针对展开的相位偏移数据的估计。
3.根据权利要求1或2所述的方法,其中,所述目标函数是能够根据关于所测量的数据的似然函数导出的。
4.根据权利要求3所述的方法,其中,针对所述似然函数的潜在概率密度是基于缠绕的高斯密度的。
5.根据权利要求3所述的方法,其中,针对所述似然函数的潜在概率密度是基于冯·米塞斯密度的。
6.根据前述权利要求中的任一项所述的方法,其中,所述正则化项是吉洪诺夫正则化子。
7.根据前述权利要求中的任一项所述的方法,其中,在所述方法中,在给定像素处测量的所述相位偏移数据是独立于在邻近像素处测量的相位偏移数据来处理的。
8.一种被配置为执行前述权利要求中的任一项所述的方法的信号处理装置(SP)。
9.根据权利要求8所述的信号处理装置,被布置用于根据关于所测量的数据的似然函数导出所述目标函数。
10.根据权利要求9所述的信号处理装置,被布置用于将针对所述似然函数的潜在概率密度基于缠绕的高斯密度。
11.根据权利要求8、9或10所述的信号处理装置,被布置用于使用吉洪诺夫正则化子作为正则化项。
12.根据权利要求8至11中的任一项所述的信号处理装置,被布置用于独立于在邻近像素处测量的相位偏移数据来处理在给定像素处测量的相位偏移数据。
13.根据权利要求8至12中的任一项所述的信号处理装置,其中,所述探测器(D)是能量分辨的。
14.一种用于控制根据权利要求8至中的任一项所述的信号处理装置(SP)的计算机程序单元,所述计算机程序单元在由数据处理单元(PU)运行时适于执行根据权利要求1至7中的任一项所述的方法的步骤。
15.一种在其上存储有根据权利要求14所述的程序单元的计算机可读介质。
CN201580029369.8A 2014-06-02 2015-04-21 用于差分相位对比成像中的谱相位展开的无偏置正则化 Pending CN106462984A (zh)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
EP14170740 2014-06-02
EP14170740.6 2014-06-02
PCT/EP2015/058566 WO2015185259A1 (en) 2014-06-02 2015-04-21 Biais-free regularization for spectral phase-unwrapping in differential phase contrast imaging.

Publications (1)

Publication Number Publication Date
CN106462984A true CN106462984A (zh) 2017-02-22

Family

ID=50943073

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201580029369.8A Pending CN106462984A (zh) 2014-06-02 2015-04-21 用于差分相位对比成像中的谱相位展开的无偏置正则化

Country Status (5)

Country Link
US (1) US10037600B2 (zh)
EP (1) EP2953097B1 (zh)
JP (1) JP6495945B2 (zh)
CN (1) CN106462984A (zh)
WO (1) WO2015185259A1 (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108537842A (zh) * 2017-12-29 2018-09-14 南京理工大学 差分相衬显微成像中背景非均匀性的校正与补偿方法
CN109470729A (zh) * 2017-09-06 2019-03-15 株式会社岛津制作所 放射线相位差摄影装置

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR20140111818A (ko) * 2013-03-12 2014-09-22 삼성전자주식회사 엑스선 영상 장치 및 그 제어 방법
WO2019056309A1 (en) * 2017-09-22 2019-03-28 Shenzhen United Imaging Healthcare Co., Ltd. METHOD AND SYSTEM FOR GENERATING A PHASE CONTRAST IMAGE

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102217934A (zh) * 2011-04-08 2011-10-19 中国科学院深圳先进技术研究院 磁共振成像方法及系统
WO2012038857A1 (en) * 2010-09-20 2012-03-29 Koninklijke Philips Electronics N.V. Phase gradient unwrapping in differential phase contrast imaging
US20120185801A1 (en) * 2011-01-18 2012-07-19 Savant Systems, Llc Remote control interface providing head-up operation and visual feedback when interacting with an on screen display
WO2013091078A1 (en) * 2011-12-23 2013-06-27 Liu Junmin Method for phase unwrapping
WO2013104966A1 (en) * 2012-01-12 2013-07-18 Koninklijke Philips N.V. Generating attenuation image data and phase image data in an x-ray system
WO2013083893A4 (fr) * 2011-12-05 2013-08-22 Solystic Dispositif d'empilage pour objets plats empilés sur chant et machine de tri postal
CN103632345A (zh) * 2013-11-27 2014-03-12 中国科学技术大学 一种基于正则化的mri图像非均匀性校正方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6741357B2 (en) * 2001-08-14 2004-05-25 Seagate Technology Llc Quadrature phase shift interferometer with unwrapping of phase
US7372393B2 (en) * 2006-07-07 2008-05-13 Mitsubishi Electric Research Laboratories, Inc. Method and system for determining unwrapped phases from noisy two-dimensional wrapped-phase images
US8433739B2 (en) * 2007-06-27 2013-04-30 L-3 Communications Integrated Systems, L.P. Methods and systems for detecting repetitive synchronized signal events
US8855265B2 (en) 2009-06-16 2014-10-07 Koninklijke Philips N.V. Correction method for differential phase contrast imaging
JP6105586B2 (ja) * 2011-08-31 2017-03-29 コーニンクレッカ フィリップス エヌ ヴェKoninklijke Philips N.V. エネルギー高感度検出の微分位相コントラストイメージング

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2012038857A1 (en) * 2010-09-20 2012-03-29 Koninklijke Philips Electronics N.V. Phase gradient unwrapping in differential phase contrast imaging
US20120185801A1 (en) * 2011-01-18 2012-07-19 Savant Systems, Llc Remote control interface providing head-up operation and visual feedback when interacting with an on screen display
CN102217934A (zh) * 2011-04-08 2011-10-19 中国科学院深圳先进技术研究院 磁共振成像方法及系统
WO2013083893A4 (fr) * 2011-12-05 2013-08-22 Solystic Dispositif d'empilage pour objets plats empilés sur chant et machine de tri postal
WO2013091078A1 (en) * 2011-12-23 2013-06-27 Liu Junmin Method for phase unwrapping
WO2013104966A1 (en) * 2012-01-12 2013-07-18 Koninklijke Philips N.V. Generating attenuation image data and phase image data in an x-ray system
CN103632345A (zh) * 2013-11-27 2014-03-12 中国科学技术大学 一种基于正则化的mri图像非均匀性校正方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
王彦飞: "地震波干涉偏移及预条件正则化最小二乘偏移成像方法对比", 《地球物理学报》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109470729A (zh) * 2017-09-06 2019-03-15 株式会社岛津制作所 放射线相位差摄影装置
CN108537842A (zh) * 2017-12-29 2018-09-14 南京理工大学 差分相衬显微成像中背景非均匀性的校正与补偿方法

Also Published As

Publication number Publication date
US10037600B2 (en) 2018-07-31
JP2017520294A (ja) 2017-07-27
EP2953097B1 (en) 2016-10-26
US20170091933A1 (en) 2017-03-30
JP6495945B2 (ja) 2019-04-03
EP2953097A1 (en) 2015-12-09
WO2015185259A1 (en) 2015-12-10

Similar Documents

Publication Publication Date Title
CN110050281A (zh) 学习图像中的对象的注释
CN107407714A (zh) 用于根据b0图和b1图来计算导出值的mri方法
CN102955142B (zh) 用于迭代的mr重建方法的采样模式
CN106462984A (zh) 用于差分相位对比成像中的谱相位展开的无偏置正则化
CN109060849A (zh) 一种确定辐射剂量调制线的方法、系统和装置
Liu et al. Method for B0 off‐resonance mapping by non‐iterative correction of phase‐errors (B0‐NICE)
CN104599301B (zh) 一种pet图像的重建方法和装置
Aja-Fernández et al. Effective noise estimation and filtering from correlated multiple-coil MR data
WO2017092973A1 (en) Removal of image artifacts in sense-mri
Tao et al. Robust motion correction for myocardial T1 and extracellular volume mapping by principle component analysis‐based groupwise image registration
JP6700316B2 (ja) B0不均一性マップ及び被検体磁気感受性マップを用いる骨mri
CN108802648B (zh) 一种基于梯度回波的磁共振定量成像方法和装置
US10297049B2 (en) Statistically weighted regularization in multi-contrast imaging
US12067652B2 (en) Correction of magnetic resonance images using multiple magnetic resonance imaging system configurations
Soliman et al. Max‐IDEAL: a max‐flow based approach for IDEAL water/fat separation
CN109124666A (zh) 一种确定辐射剂量调制线的方法、系统和装置
Koppers et al. Reconstruction of diffusion anisotropies using 3D deep convolutional neural networks in diffusion imaging
Diao et al. Parameter estimation for WMTI‐Watson model of white matter using encoder–decoder recurrent neural network
CN112204411A (zh) Cest磁共振成像中的运动检测
CN116563096A (zh) 用于图像配准的形变场的确定方法、装置以及电子设备
CN113966471A (zh) 螺旋k空间采样磁共振图像的重建
CN109791187A (zh) B0校正灵敏度编码磁共振成像
CN110073234B (zh) 医学仪器及其方法
CN107810518B (zh) 图像处理系统和方法
CN104173049A (zh) 用于借助磁共振断层成像系统记录图像数据组的方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
WD01 Invention patent application deemed withdrawn after publication

Application publication date: 20170222

WD01 Invention patent application deemed withdrawn after publication