CN108714028B - 一种磁共振成像方法、装置及医学成像系统 - Google Patents

一种磁共振成像方法、装置及医学成像系统 Download PDF

Info

Publication number
CN108714028B
CN108714028B CN201810323015.XA CN201810323015A CN108714028B CN 108714028 B CN108714028 B CN 108714028B CN 201810323015 A CN201810323015 A CN 201810323015A CN 108714028 B CN108714028 B CN 108714028B
Authority
CN
China
Prior art keywords
information
susceptibility distribution
map
error
distribution map
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
CN201810323015.XA
Other languages
English (en)
Other versions
CN108714028A (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.)
Shanghai United Imaging Healthcare Co Ltd
Original Assignee
Shanghai United Imaging Healthcare 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 Shanghai United Imaging Healthcare Co Ltd filed Critical Shanghai United Imaging Healthcare Co Ltd
Priority to CN201810323015.XA priority Critical patent/CN108714028B/zh
Publication of CN108714028A publication Critical patent/CN108714028A/zh
Priority to US16/247,098 priority patent/US10950011B2/en
Priority to US17/188,918 priority patent/US11587269B2/en
Application granted granted Critical
Publication of CN108714028B publication Critical patent/CN108714028B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/05—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves 
    • A61B5/055—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves  involving electronic [EMR] or nuclear [NMR] magnetic resonance, e.g. magnetic resonance imaging
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • 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
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
    • G06V10/74—Image or video pattern matching; Proximity measures in feature spaces
    • G06V10/75—Organisation of the matching processes, e.g. simultaneous or sequential comparisons of image or video features; Coarse-fine approaches, e.g. multi-scale approaches; using context analysis; Selection of dictionaries
    • G06V10/751—Comparing pixel values or logical combinations thereof, or feature values having positional relevance, e.g. template matching
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/40—Detecting, measuring or recording for evaluating the nervous system
    • A61B5/4058—Detecting, measuring or recording for evaluating the nervous system for evaluating the central nervous system
    • A61B5/4064—Evaluating the brain
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/24—Arrangements or instruments for measuring magnetic variables involving magnetic resonance for measuring direction or magnitude of magnetic fields or magnetic flux
    • G01R33/243—Spatial mapping of the polarizing magnetic field
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/4818—MR characterised by data acquisition along a specific k-space trajectory or by the temporal order of k-space coverage, e.g. centric or segmented coverage of k-space
    • G01R33/4824—MR characterised by data acquisition along a specific k-space trajectory or by the temporal order of k-space coverage, e.g. centric or segmented coverage of k-space using a non-Cartesian trajectory
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00—Arrangements or instruments for measuring magnetic variables
    • G01R33/20—Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44—Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/48—NMR imaging systems
    • G01R33/54—Signal processing systems, e.g. using pulse sequences ; Generation or control of pulse sequences; Operator console
    • G01R33/56—Image enhancement or correction, e.g. subtraction or averaging techniques, e.g. improvement of signal-to-noise ratio and resolution
    • G01R33/565—Correction of image distortions, e.g. due to magnetic field inhomogeneities

Landscapes

  • Health & Medical Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • General Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Theoretical Computer Science (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Radiology & Medical Imaging (AREA)
  • General Physics & Mathematics (AREA)
  • Molecular Biology (AREA)
  • Biomedical Technology (AREA)
  • Pathology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • Biophysics (AREA)
  • Public Health (AREA)
  • Veterinary Medicine (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Computing Systems (AREA)
  • Artificial Intelligence (AREA)
  • Databases & Information Systems (AREA)
  • Evolutionary Computation (AREA)
  • Software Systems (AREA)
  • Multimedia (AREA)
  • Condensed Matter Physics & Semiconductors (AREA)
  • Signal Processing (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

本发明实施例提供了一种磁共振成像方法、装置及医学成像系统，涉及磁共振成像技术领域，在磁化率分布图的迭代过程中，动态更新误差限定信息，从而减小磁化率分布图的计算收敛误差。该方法包括根据初始场图和误差限定信息获取磁化率分布图；对误差限定信息与磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差；当更新后的磁化率分布图的计算收敛误差小于或等于误差阈值时，停止迭代过程，并输出该次更新后的磁化率分布图；其中，每次互相迭代更新包括：通过当前磁化率分布图更新当前误差限定信息；通过更新后的误差限定信息更新当前磁化率分布图。本发明实施例的技术方案适用磁化率分布定量计算过程。

Description

一种磁共振成像方法、装置及医学成像系统
【技术领域】
本发明涉及磁共振成像技术领域，尤其涉及一种磁共振成像方法、装置及医学成像系统。
【背景技术】
磁共振磁化率定量成像(Quantitative Susceptibility Mapping，QSM)技术是磁共振成像中一项新兴技术，用于计算组织或者成像物体的磁化率分布图。QSM技术的大致流程为：首先获得一组或多组磁共振相位信息或相位图像，然后利用磁共振成像的相位信息来获取组织局部的场图变化信息，最终通过场图与磁化率之间的物理关系反演出磁化率分布图。
QSM算法可以根据特征进行分类。其中，应用较为广泛的是基于贝叶斯统计的优化算法。目前，为了得到数值准确、伪影较小的磁化率分布图，基于贝叶斯统计的优化算法，通常是利用磁共振图像中的解剖结构先验信息作为磁化率分布的似然函数，为QSM计算所要解决的逆向偶极效应问题提供贝叶斯似然函数。这样做的假设前提是：生物组织磁化率分布的空间特性与解剖结构的空间特性高度相似，因此经过优化算法得到的解与先验信息一致的那部分将会得到强化，从而增强结果的准确度。目前此类解剖结构信息的获取通常从磁共振幅值图或者从初步估计的QSM图中获取，并且该先验信息在求解对应的逆向磁偶效应问题的计算中保持不变。但是实际应用中的一个问题是，如何获得最合适的先验信息。在QSM中，先验信息的作用是提供关于磁化率空间分布的特征，比如同一组织内部的分布连续性、不同组织间的分界等。在计算欠定逆向偶极问题时，可以以贝叶斯似然函数或者正则化项的形式，利用这些信息对某些解进行限制或强化，从而使得到的最优化解更接近真实解。理想情况下，如果真实的磁化率分布已知，如计算机模拟模型，所提取的先验信息将完全符合。然而现实情况下，真实的磁化率分布是未知的，因此通常使用其他的结构信息来对获得一个近似的先验信息，比如同一扫描下获得的幅值图中的结构信息，或者磁化率分布的初步估计结果。若是用幅值图中的结构信息，其空间特性与真实磁化率分布的空间特性存在差异，因此正则化项对解的约束可能会被弱化甚至是错误的，从而导致最优化解中出现误差。若是使用磁化率分布的初步估计结果来获取结构信息，则需要两次或以上的计算，增加计算时间，同时仍然可能存在空间特性估计上的差异。
另一方面，现有的QSM方法都需要解决计算相关的放射条纹伪影。该伪影是由于逆向偶极核函数中存在奇点而造成的(在数学上不可避免)，并且在物体外部和内部都会产生，其影响包括造成物体边界的模糊、与真实的组织结构混叠、以及使得磁化率数值产生误差。从数学表征上看，放射条纹伪影的出现是由于逆向偶极函数中存在无限大值奇点造成的，因此去除奇点、降低奇点及其附近的点的数值，可以降低伪影的强度，这就是k空间截断方法的基本原理。但是由于截断操作改变了逆向偶极函数的数学表达，因此，计算得到的结果数值将偏低。截断越多，伪影越少，数值越偏低估。而即使截断操作可以降低伪影强度，从原理上仍然不可能实现完全去除。
另一方面，从贝叶斯统计的角度看，解剖结构上的先验信息与磁化率分布的空间共位性可以认为是一种稀疏的信息，如组织结构之间的边界信息。因此将这种稀疏的信息作为正则化约束条件加入到最优化计算中，可以有效的降低放射条状伪影。然而，先验信息与真实磁化率分布之间的差异可能导致对某些解的过正则化或欠正则化，从而影响结果的准确度。另外，对于具有较强磁化率特性的组织或区域，如空气-组织交界处往往具有很强的局部场差异，因此相对其他组织来说，这些强磁化率区域的残余放射条状伪影仍然可能造成较强的影响。
鉴于此，有必要对现有QSM成像方法进行改进，以减小磁化率分布图的计算收敛误差。
【发明内容】
有鉴于此，本发明实施例提供了一种磁共振成像方法、装置及医学成像系统，在迭代计算过程中，通过动态更新先验信息和/或伪影信息，可以减小磁化率分布图的计算收敛误差。
第一方面，本发明实施例提供一种磁共振成像方法，所述方法包括：
根据初始场图和误差限定信息，获取磁化率分布图，所述误差限定信息包括先验信息、伪影信息中的至少一种；
对所述误差限定信息与所述磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差；
当更新后的磁化率分布图的计算收敛误差小于或等于误差阈值时，停止互相迭代更新过程，并输出该次更新后的磁化率分布图；
其中，每次互相迭代更新的过程包括：
通过当前磁化率分布图对当前误差限定信息进行更新；
通过更新后的误差限定信息对当前磁化率分布图进行更新。
如上所述的方面和任一可能的实现方式，进一步提供一种实现方式，当所述误差限定信息包括先验信息时，所述通过当前磁化率分布图对当前误差限定信息进行更新，包括：
提取当前磁化率分布图的特征，得到与当前磁化率分布图对应的特性图；
在所述特性图中提取更新后的先验信息。
如上所述的方面和任一可能的实现方式，进一步提供一种实现方式，当所述特性图为梯度特性图时，所述提取当前磁化率分布图的特征，得到与当前磁化率分布图对应的特性图，包括：
提取当前磁化率分布图沿参考方向的梯度特征，得到所述梯度特性图，所述参考方向包括X方向、Y方向或Z方向的一种或多种。
如上所述的方面和任一可能的实现方式，进一步提供一种实现方式，当所述误差限定信息包括先验信息时，所述通过当前磁化率分布图对当前误差限定信息进行更新，包括：
基于正向偶极效应模型方程，将当前磁化率分布图转换为先验场图；
根据所述先验场图，生成更新后的先验信息。
如上所述的方面和任一可能的实现方式，进一步提供一种实现方式，当所述误差限定信息包括伪影信息时，所述通过当前磁化率分布图对当前误差限定信息进行更新，包括：
基于正向偶极效应模型方程，将当前磁化率分布图转换为估计场图；
根据所述初始场图以及所述估计场图，生成差值场图；
基于正向偶极核函数，根据所述差值场图生成更新后的伪影信息。
如上所述的方面和任一可能的实现方式，进一步提供一种实现方式，当所述误差限定信息包括伪影信息时，所述通过更新后的误差限定信息对当前磁化率分布图进行更新，包括：
以正则化的形式将更新后的伪影信息带入所述互相迭代更新过程，得到更新后的磁化率分布图；或者，
在当前磁化率分布图中去除更新后的伪影信息，得到更新后的磁化率分布图。
如上所述的方面和任一可能的实现方式，进一步提供一种实现方式，所述先验信息包括空间平滑性、空间梯度、边缘信息、加权蒙版中的至少一种。
如上所述的方面和任一可能的实现方式，进一步提供一种实现方式，以以下多种迭代方法中的至少一种迭代方法，对所述误差限定信息与所述磁化率分布图进行至少一次互相迭代更新：
共轭梯度法、预处理共轭梯度法、二次最小法与分裂-布雷格曼法。
第二方面，本发明实施例提供一种磁共振成像装置，所述装置包括：
获取单元，用于根据初始场图和误差限定信息，获取磁化率分布图，所述误差限定信息包括先验信息、伪影信息中的至少一种；
更新单元，用于对所述误差限定信息与所述磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差；
处理单元，用于当更新后的磁化率分布图的计算收敛误差小于或等于误差阈值时，停止互相迭代更新过程，并输出该次更新后的磁化率分布图；
其中，每次互相迭代更新的过程包括：
通过当前磁化率分布图对当前误差限定信息进行更新；
通过更新后的误差限定信息对当前磁化率分布图进行更新。
第三方面，本发明实施例提供一种医学成像系统，所述医学成像系统包括处理器以及存储器；所述存储器用于存储指令，所述指令被所述处理器执行时，导致所述医学成像系统实现如上任一方面及任一可能的实现方式所述的方法。
本发明实施例提供了一种磁共振成像方法、装置及医学成像系统，首先根据初始场图以及误差限定信息，获取磁化率分布图，然后对误差限定信息与磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差。在互相迭代更新过程中，通过动态的更新先验信息和/或伪影信息，可以使先验信息或伪影信息逐步接近真实值，从而在迭代过程中逐步减小磁化率分布图的计算收敛误差，可以得到计算收敛误差小于或等于误差阈值、精确度较高的磁化率分布图。
【附图说明】
为了更清楚地说明本发明实施例的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其它的附图。
图1是本发明实施例提供的一种磁共振成像方法流程图；
图2(1)是本发明实施例提供的另一种磁共振成像方法流程图；
图2(2)是本发明实施例提供的另一种磁共振成像方法流程图；
图2(3)是本发明实施例提供的另一种磁共振成像方法流程图；
图3(1)是本发明实施例提供的另一种磁共振成像方法流程图；
图3(2)是本发明实施例提供的另一种磁共振成像方法流程图；
图4是本发明实施例提供的另一种磁共振成像方法流程图；
图5是本发明实施例提供的另一种磁共振成像方法流程图；
图6(1)是本发明实施例提供的另一种磁共振成像方法流程图；
图6(2)是本发明实施例提供的另一种磁共振成像方法流程图；
图7(a)是本发明实施例提供的一种基于先验信息动态更新的磁化率分布图；
图7(b)是本发明实施例提供的另一种基于先验信息动态更新的磁化率分布图；
图7(c)是本发明实施例提供的另一种基于先验信息动态更新的磁化率分布图；
图8(1)是本发明实施例提供的一种基于伪影信息动态更新的示例图；
图8(2)是本发明实施例提供的另一种基于伪影信息动态更新的示例图；
图9(1)是本发明实施例提供的磁共振成像系统的结构框图；
图9(2)是本发明实施例的图像生成器的示例性架构；
图10是本发明实施例的处理器结构框图。
【具体实施方式】
为了更好的理解本发明的技术方案，下面结合附图对本发明实施例进行详细描述。
应当明确，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其它实施例，都属于本发明保护的范围。
在本发明实施例中使用的术语是仅仅出于描述特定实施例的目的，而非旨在限制本发明。在本发明实施例和所附权利要求书中所使用的单数形式的“一种”、“所述”和“该”也旨在包括多数形式，除非上下文清楚地表示其他含义。
取决于语境，如在此所使用的词语“如果”可以被解释成为“在……时”或“当……时”或“响应于确定”或“响应于检测”。类似地，取决于语境，短语“如果确定”或“如果检测(陈述的条件或事件)”可以被解释成为“当确定时”或“响应于确定”或“当检测(陈述的条件或事件)时”或“响应于检测(陈述的条件或事件)”。
应当理解，本文中使用的术语“和/或”仅仅是一种描述关联对象的关联关系，表示可以存在三种关系，例如，A和/或B，可以表示：单独存在A，同时存在A和B，单独存在B这三种情况。另外，本文中字符“/”，一般表示前后关联对象是一种“或”的关系。
本发明实施例提供了一种磁共振成像方法，适用于QSM过程中，如图1所示，为本发明实施例提供的磁共振成像方法流程图，包括：
101、根据初始场图和误差限定信息，获取磁化率分布图，所述误差限定信息包括先验信息、伪影信息中的至少一种。
需要说明的是，可以采集受检者目标部位的磁共振信号，根据受检者目标部位的磁共振信号的相位信息确定初始场图。
在QSM中，先验信息的作用是提供关于磁化率空间分布的特征，比如同一组织内部的分布连续性、不同组织间的分界等。在计算欠定逆向偶极问题时，可以以贝叶斯似然函数或者正则化项的形式，利用这些信息对某些解进行限制或强化，从而使得到的最优化解更接近真实解。在具体的实现方式中，先验信息可以是与被检者目标部位对应的空间平滑性、空间梯度、边缘信息、加权蒙版中的一种或多种形式的组合，且上述先验信息可以从受检者目标部位的幅值图、相位图磁化率分布图、R2*(T2*)图或者磁敏感加权成像(Susceptibility weighted imaging，SWI)图上获取。
具体的，伪影信息可以是放射条状伪影。由于逆向偶极问题的欠定特性，QSM求解过程中会产生计算相关的放射条状伪影。从数学表征上看，放射条状伪影的出现是由于逆向偶极函数中存在无限大值奇点造成的。去除磁化率分布图中的残余放射条状伪影，可以有效降低伪影的强度，提高磁化率分布图的准确度。
具体的，可以根据初始场图和初始的先验信息，或者根据初始场图和初始的伪影信息，或者根据初始场图和初始的先验信息、伪影信息，来获取初始的磁化率分布图。
102、对所述误差限定信息与所述磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差。
其中，每次互相迭代更新的过程包括：
1021、通过当前磁化率分布图对当前误差限定信息进行更新。
在上述迭代过程中，在一种可能的实现方式中，可以通过当前磁化率分布图的特性图，来动态更新先验信息；在另一种可能的实现方式中，可以基于当前磁化率分布图得到先验场图，通过先验场图来动态更新先验信息。但需要说明的是，本发明实施例提供的技术方案中，先验信息的更新并不只限于通过磁化率分布图进行动态更新，在其他可选的实现方式中，先验信息还可以来源于幅值图、相位图、SWI图或者R2*(T2*)图等，因此先验信息还可以通过这些类型的医学图像进行逐步逼近真实值的动态更新。或者，先验信息可以结合这些类型的医学图像和中间迭代计算过程的磁化率分布图来进行动态更新。
在上述迭代过程中，可以将当前磁化率分布图与初始场图结合，来动态更新伪影信息。
1022、通过更新后的误差限定信息对当前磁化率分布图进行更新。
上述迭代过程，将更新后的先验信息和/或伪影信息带入计算过程，对当前磁化率分布图进行更新。
在具体的实现过程中，对所述误差限定信息与所述磁化率分布图进行互相迭代更新时，可以通过以下任意基于迭代和逆向偶极模型的最优化算法，比如共轭梯度法、预处理共轭梯度法、二次最小法、分裂-布雷格曼法等。
需要说明的是，在每次误差限定信息与磁化率分布图互相迭代更新之后，均对更新后的磁化率分布图的计算收敛误差进行确定，若更新后的磁化率分布图的计算收敛误差小于或等于误差阈值时，则停止互相迭代更新过程，并输出该次更新后的磁化率分布图(即执行步骤103)；若更新后的磁化率分布图的误差大于误差阈值，则继续互相迭代更新过程(即循环步骤1021和步骤1022)。
103、当更新后的磁化率分布图的计算收敛误差小于或等于误差阈值时，停止互相迭代更新过程，并输出该次更新后的磁化率分布图。
其中，误差阈值可以根据理想磁化率分布图的误差信息进行预设，当磁化率分布图的计算收敛误差小于或等于误差阈值时，认为磁化率分布图的计算收敛误差符合要求，磁化率分布图满足精确度标准，可以进行输出。
本发明实施例提供的技术方案通过互相迭代方式优化计算磁化率分布图。在中间计算过程中，动态更新先验信息和/或伪影信息，使先验信息和/或伪影信息逐步接近真实值，从而在迭代计算磁化率分布图的过程中逐步提升磁化率分布图的精确度，直至磁化率分布图的计算收敛误差小于或等于误差阈值时，输出最终结果。
上述发明实施例提供了一种磁共振成像方法，首先根据初始场图以及误差限定信息，获取磁化率分布图，然后对误差限定信息与磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差。在互相迭代更新过程中，通过动态的更新先验信息和/或伪影信息，可以使先验信息或伪影信息逐步接近真实值，从而在迭代过程中逐步减小磁化率分布图的计算收敛误差，可以得到计算收敛误差小于或等于误差阈值、精确度较高的磁化率分布图。
进一步来说，结合前述方法流程，当误差限定信息包括先验信息时，本发明实施例的另一种可能的实现方式针对步骤1021中先验信息的更新，还提供了以下两种可能的实现方法。
第一种实现方法，如图2(1)所示，步骤1021具体包括：
201、提取当前磁化率分布图的特征，得到与当前磁化率分布图对应的特性图。
202、在所述特性图中提取更新后的先验信息。
在一种具体的实现方式中，当前磁化率分布图的特征可以是磁化率分布图沿参考方向的梯度特征，则上述实现方法可以得到当前磁化率分布图对应的梯度特性图，基于梯度特性图进行先验信息的更新。其中，参考方向可以是三维空间的X、Y、Z中的任意一个或多个方向，X、Y、Z三个方向可以互相垂直或有一定夹角(不相互平行即可)。
第二种实现方法，如图2(2)所示，步骤1021具体包括：
203、基于正向偶极效应模型方程，将当前磁化率分布图转换为先验场图。
204、根据所述先验场图，生成更新后的先验信息。
第二种实现方法中，在上次迭代更新磁化率分布图之后，可以基于正向偶极效应模型方程，通过当前磁化率分布图转换得到场图，将场图作为先验场图来提取新的先验信息，以达成先验信息的更新。
上述两种实现方法在迭代计算过程中，基于中间迭代过程的磁化率分布图来对先验信息进行动态提取和更新，将更新后的先验信息应用到下一次磁化率分布图的迭代计算中。这样不仅可以使得先验信息逐渐接近真实的磁化率分布的空间特性，提高结果的准确性，同时可以进一步加快最优化算法的迭代收敛速度。
进一步来说，结合前述方法流程，当误差限定信息包括伪影信息时，本发明实施例的另一种可能的实现方式针对步骤1021中伪影信息的更新，还提供了以下具体的方法流程，如图2(3)所示，步骤1021具体包括：
205、基于正向偶极效应模型方程，将当前磁化率分布图转换为估计场图。
需要说明的是，这里的估计场图与步骤203中的先验场图均是通过当前磁化率分布图转换得到的场图，仅因为所处环境不同，所以用估计场图和先验场图两种不同的方式进行命名。
206、根据所述初始场图以及所述估计场图，生成差值场图。
207、基于正向偶极核函数，根据所述差值场图生成更新后的伪影信息。
具体的，可以根据部分正向偶极核函数函数和差值场图来计算得到更新后的伪影信息。
其中，部分正向偶极核函数只包含了原正向偶极核函数中含有和接近零值点的部分,例如可以通过对原正向偶极核函数进行阈值计算获得部分正向偶极核函数。
上述实现方式在迭代计算过程中，基于中间迭代过程的磁化率分布图动态地更新伪影信息，并可以将其以正则化项的形式应用到迭代计算过程中，或直接在磁化率分布图中去除更新后的伪影信息。该方式可以使磁化率分布图的伪影强度降到可忽略的水平，同时不影响其他计算过程，从而在降低条状伪影的同时，还保证了计算结果的准确度。
进一步来说，结合前述方法流程，当误差限定信息包括伪影信息时，本发明实施例的另一种可能的实现方式针对步骤1022中磁化率分布图的更新，还提供了以下两种可能的实现方式。
第一种实现方式，如图3(1)所示，步骤1022包括：
301、以正则化的形式将更新后的伪影信息带入所述互相迭代更新过程，得到更新后的磁化率分布图。
第二种实现方式，如图3(2)所示，步骤1022包括：
302、在当前磁化率分布图中去除更新后的伪影信息，得到更新后的磁化率分布图。
上述两种基于伪影信息对磁化率分布图进行更新的方法，为磁化率分布图的更新提供了多样化的实现方式。
根据以上技术方案，本发明实施例给出了基于以下公式迭代计算磁化率分布图的具体示例：
Figure BDA0001625800320000121
第一部分
Figure BDA0001625800320000122
是逆向偶极效应模型方程的二次最小化方程。其中，逆向偶极效应模型方程是场图与对应磁化率分布图进行互相转换的计算方程。其中，F和F-1分别是傅里叶与反傅里叶算子，D是原偶极核函数，χ是磁化率分布，φ是磁场分布。如果只利用这一部分进行逆向偶极方程的数值求解，计算得到的磁化率分布图将包含严重的放射条状伪影和数值误差。
第二部分
Figure BDA0001625800320000123
是先验信息在三个方向上梯度的规范式。其中，α是该部分的正则化因子，其数值选择可以调节计算结果在数据连续性与空间平滑性之间的平衡。P是一个加权蒙板(先验信息)，可以动态地从最优化迭代计算过程中的初始或中间χ中得到。
Figure BDA0001625800320000124
是梯度算子。第二部分是对先验信息的动态更新，并将其应用于最优化迭代计算磁化率分布图中的具体实现。
第三部分
Figure BDA0001625800320000125
是残余放射条状伪影(伪影信息的一种具体形式)的规范式。其中，β是该部分的正则化因子。D′是部分正向偶极核函数，只包含了原正向偶极核函数D中含有和接近零值点的部分，例如可以通过对D进行阈值计算获得：D′＝D·(|D|<thresh)。Δφ是由估计场图(由初始磁化率分布图或中间过程迭代更新的磁化率分布图计算得到)与初始场图两者得到的差值场图。第三部分是对残余放射条状伪影的动态更新，并将其应用于最优化迭代计算磁化率分布图中的具体实现。
上述示例的具体方程式是实现本技术方案的方式之一，本技术方案还可以应用于不同的正则化项，或者不同的最优化迭代形式。
基于本发明实施例提供的技术方案，本发明实施例还提供了以下几种的磁共振成像方法以及误差限定参数更新方法的具体流程。
如图4为本发明另一实施例提供的磁共振成像方法流程图，该方法包括如下步骤：
401、获取/采集对象目标部位的磁共振信号。
在一些实施例中，对象可包括物质、组织、样本、身体等，或者其任意组合。对象可包括头、胸、肺、胸膜、纵膈、腹、大肠、小肠、膀胱、胆囊、三焦、盆腔、骨干、四肢、骨架、血管等，或者其任意组合。在一些实施例中，对象目标部位可受激产生自旋信号，该自旋信号可经过相位编码处理而形成磁共振信号。
402、根据对象目标部位的磁共振信号确定初始场图。
在一些实施例中，可提取磁共振信号的相位信息，该相位信息可对应初始场图。
403、确定初始先验信息。
该初始先验信息可以是与受检者目标部位对应的空间平滑性、空间梯度、边缘信息、加权蒙版中的至少一种，这些信息在一定程度上能够反映目标部位的轮廓。可选地，上述信息可以通过如下过程获得：根据目标部位对应的磁共振信号重建幅值图、相位图、磁化率分布图、R2*(T2*)图或者SWI图，从重建的图像中提取先验信息。
404、根据初始先验信息和初始场图获取初始磁化率分布图。
可选地，初始磁化率分布图的确定过程可前述参考公式1。
405、对初始磁化率分布图进行一次或多次迭代优化，获取目标部位对应的磁化率分布图，且在所述磁化率分布图迭代优化的过程中，更新先验信息。
可选地，可根据第一次迭代过程中产生的第一磁化率分布图确定第一先验信息；根据第一先验信息更新第一磁化率分布图而产生第二磁化率分布图。进一步地，可根据第二磁化率分布图确定第二先验信息；根据第二先验信息更新第二磁化率分布图而产生第三磁化率分布图。当然需要说明的是，本申请中对于迭代优化次数并没有具体限制，如可经过三次、四次或者更多次迭代获得满足要求的磁化率分布图。
如图5为本发明另一实施例提供的先验信息的迭代更新确定过程，包括如下步骤：
501、输入先验图/映射图(map)。
该先验图/映射图可以是幅值图、相位图、磁化率分布图、R2*(T2*)图或者SWI图。在此实施例中先验图为前一次迭代更新后的磁化率分布图。
502、确定先验图/映射图的特征图。
在此实施例中，特征图为前一次迭代更新后的磁化率分布图分别沿X、Y和Z方向的梯度图。
503、根据特征图产生加权蒙版(property mask)，该加权蒙版为当前次迭代的先验信息。
在上述过程中，当在磁化率分布图迭代优化的过程中基于迭代计算过程中的中间磁化率分布图来进行先验信息的动态提取和更新，并将更新后的该先验信息应用到下一次迭代中，这样不仅可以使得先验信息逐渐接近真实的磁化率分布的空间特性，同时可以进一步加快最优化算法的迭代收敛速度和提高结果的准确性。
如图6(1)为本发明另一实施例提供的磁共振成像方法流程图，该方法包括如下步骤：
601、获取/采集对象目标部位的磁共振信号。
602、根据对象目标部位的磁共振信号确定初始场图。
在一些实施例中，可提取磁共振信号的相位信息，该相位信息可对应初始场图。
603、确定初始先验信息。
该初始先验信息可以是与受检者目标部位对应的空间平滑性、空间梯度、边缘信息、加权蒙版中的至少一种，这些信息在一定程度上能够反映目标部位的轮廓。可选地，上述信息可以通过如下过程获得：根据目标部位对应的磁共振信号重建幅值图、相位图、磁化率分布图、R2*(T2*)图或者SWI图，从重建的图像中提取先验信息。
604、根据初始先验信息和初始场图获取初始磁化率分布图。
可选地，初始磁化率分布图的确定过程可前述参考公式1。
605、对初始磁化率分布图进行一次或多次迭代优化，获取目标部位对应的磁化率分布图，且在所述磁化率分布图迭代优化的过程中，更新场图对应的残余放射条状伪影。
如图6(2)为本发明另一实施例提供的场图对应的残余放射条状伪影的迭代更新确定过程，包括如下步骤：
611、获取初始场图和估计场图。
其中，估计场图可基于前一次迭代更新后的磁化率分布图得到。
需要说明的是，磁场场分布图(场图)与其对应的磁化率分布图可使用磁偶效应模型或者逆向磁偶效应模型进行互相转换计算，可选地分别使用磁偶核函数或者逆向磁偶核函数进行计算。偶极模型方程为：
ΔB(r)＝B0·FT-1{g(k)·χ(k)} (公式2)
其中，磁偶核函数为：
Figure BDA0001625800320000151
逆向偶极模型方程为：
χ(r)＝FT-1{g-1(k)·FT[ΔB(r)]}/B0 (公式4)
其中，逆向磁偶核函数为：
Figure BDA0001625800320000161
其中，ΔBz(r)表示场图；Δχ(r)表示磁化率分布图；FT和FT-1分别表示傅里叶与反傅里叶变换矩阵；χ是磁化率分布；k表示傅里叶域中的坐标向量；kx表示沿x方向的坐标，ky表示沿y方向的坐标，kz表示沿z方向的坐标；r表示空间域中的坐标向量；B0表示外磁场；g和g-1分别表示正向磁偶核变换矩阵和逆向磁偶核变换矩阵。
612、根据初始场图和估计场图，产生差值场图。
在此实施例中，差值场图中的每个像素点为初始场图的像素点和估计场图的对应像素点的像素差值。
613、根据差值场图确定场图对应的残余放射条状伪影。
在此实施例中，可获取偶极方程D，联合利用偶极方程D和差值场图确定场图对应的残余放射条状伪影。
本实施例中，在磁化率分布图迭代优化的过程中动态地估计残余放射条状伪影，并将该残余放射条状伪影以正则化项的形式应用到最优化计算过程中。该方式可使残余放射条状伪影强度降到可忽略的水平，同时不影响其他正则化项对计算的效果，从而在降低伪影的同时保证计算结果的准确度。
如图7所示，为本技术方案中对先验信息动态更新的具体实验结果示例。本示例中先验信息为根据图像空间梯度特征而获得的加权蒙版。图像空间梯度在x、y、z三个方向上独立计算，得到三幅对应的梯度图Gx，Gy和Gz，对每一副梯度图均进行蒙版提取，分别以mask(Gx)、mask(Gy)、mask(Gz)标识。本实例中通过阈值的方法进行像素设置，当梯度图中某像素的值大于阈值，则对应的加权蒙版的像素值设为0，反之则设为1。
图7(a)中，初始加权蒙版(即第一次迭代之前)由幅值图计算获得；然后在后面磁化率分布图的迭代计算过程中，该加权蒙版则仅根据上次迭代所得到的中间磁化率分布图计算获得，并在更新后参与此次迭代计算。
图7(b)中，初始加权蒙版由大脑蒙版计算获得，即初始加权蒙版不包含有用的先验知识；然后在后面磁化率分布图的迭代计算过程中，该加权蒙版则根据上次迭代所得到的中间磁化率分布图计算获得，并在更新后参与此次迭代计算。
图7(c)中，初始加权蒙版由幅值图计算获得；然后在后面磁化率分布图的迭代计算过程中，根据上次迭代所得到的中间磁化率分布图计算获得更新的加权蒙版，并同时叠加初始加权蒙版的信息获得最终加权蒙版，最终加权蒙版参与此次迭代计算。与图7(a)的不同之处在于，图7(c)的最终加权蒙版既包含了中间磁化率分布图的空间先验知识，又包含了幅值图的空间先验知识。
在首次迭代(迭代英文全称iteration，缩写iter，首次迭代即iter＝1)时，输入初始加权蒙版，并最终迭代计算9次(last，iter＝9)之后，对比图7(a)和图7(b)、图7(c)可以看出，通过结合幅值图和中间迭代计算过程中的磁化率分布图来动态更新先验信息，可以得到最为清晰(结果最为精确)的磁化率分布图。
如图8(1)所示，为基于残余放射条状伪影，本技术方案与传统方法的实验结果对比示例。图中，左列为计算机模拟生成的理想磁化率分布图；中列为使用动态更新残余放射条状伪影并进行伪影去除方法(本技术方案)得到的磁化率分布图与理想磁化率分布图之间的差别图；右列为使用传统方法得到的磁化率分布图与理想磁化率分布图之间的差别图。其中，ppm(part per million)，表示百万分之一。-0.3-0.3ppm是指左列图像灰度对应的数值，最亮的灰度对应0.3ppm或更高的值，反之亦然。-0.05-0.05ppm是指中列和右列图像灰度对应的数值，最亮的灰度对应0.05ppm或更高的值，反之亦然。对比图8(1)中的中列和右列可以看出，使用本技术方案对残余放射条状伪影进行动态更新和伪影去除方法得到的磁化率分布图更接近于真实的磁化率分布，两者间差别更小、更均匀。
如图8(2)所示，为基于本技术方案，对残余放射条状伪影动态更新和伪影去除的真实实验结果示例，其中，本示例为大脑数据对降低残余放射条状伪影效果的演示。具体的，(a)为使用了动态更新残余放射条状伪影与伪影去除得到的磁化率分布图；(b)为传统方法得到的磁化率分布图；(c)为(a)与(b)的差值图；(d)为估计得到的残余放射条状伪影的示例。由图8(2)中可以看出，使用了本发明的动态更新残余放射条状伪影与伪影去除方法，可以有效降低磁化率分布图中的伪影程度。
本发明实施例还提供一种计算可读存储介质，包括计算机可读指令，当计算机读取并执行所述计算机可读指令时，使得计算机执行如上任一实施例所述的方法。
本发明实施例提供了一种医学成像系统，该医学成像系统可以执行上述图1至图6中的一种或多种流程，进行磁化率分布定量计算。可选地，该医学成像系统可以是磁共振成像系统(magnetic resonance，MR)，也可以是正电子放射断层照相磁共振成像(PET-MRI)多模态系统。在一个实施例中，以单模态磁共振成像系统为例说明，图9(1)是根据本申请的一些实施例的磁共振成像(MRI)系统的框图。如所解说，MRI系统900可包括MRI扫描器910、控制模块920、处理模块930，以及显示模块940。MRI扫描器910可包括磁体模块911以及射频(RF)模块912。在一些实施例中，MRI扫描器910可对一对象执行扫描，以获取/采集受检者目标部位的磁共振信号。在一些实施例中，该扫描可以是用于生成磁共振(MR)图像的成像扫描，或是用于校准MRI系统900的预扫描。磁体模块911可包括主磁场生成器和/或梯度磁场生成器。主磁场生成器在扫描期间可产生静态磁场B0。主磁场生成器可以是各种类型的，包括例如永磁体、超导电磁体、电阻性电磁体等。梯度磁场生成器可包括生成分别在“X”、“Y”、“Z”方向上的磁场梯度Gx、Gy、Gz。如本文中所使用的，X、Y和Z方向可代表坐标系中的X、Y和Z轴。仅作为示例，X轴和Z轴可以在水平面中，X轴和Y轴可以在垂直面中，Z轴可以沿机架的旋转轴。在一些实施例中，X轴、Y轴和Z轴可以由梯度磁场生成器(即，在梯度磁场生成器中的梯度线圈)来指定。梯度磁场可编码和/或读出位于MRI扫描器910内的对象的空间信息。在一些实施例中，磁体模块911在扫描期间可生成沿特定方向上的磁场梯度。仅作为示例，磁体模块111可生成第一方向上的第一磁场梯度、第二方向上的第二磁场梯度，以及第三方向上的第三磁场梯度。在一些实施例中，第一、第二和第三方向可分别沿X轴、Y轴和Z轴。在一些实施例中，沿X轴、Y轴和/或Z轴的磁场梯度可对应于k空间中的不同编码/读出方向(例如，Kx轴的方向、Ky轴的方向、Kz轴的方向或任何其他方向)。RF模块912可包括RF发射线圈和/或接收线圈。这些RF线圈可向/从感兴趣的对象发射或接收RF信号。在一些实施例中，磁体模块911和/或RF模块912的功能、大小、类型、几何、位置、量、和/或幅值可根据一个或多个具体条件来确定或改变。例如，根据功能和/或大小上的不同，RF线圈可被分类为容积线圈和局部线圈。在一些实施例中，容积线圈可包括鸟笼线圈、横电磁线圈、表面线圈、鞍形线圈等。在本申请的一些实施例中，局部线圈可包括鸟笼线圈、螺线管线圈、鞍形线圈、柔性线圈等。在一些实施例中，磁体模块911和射频(RF)模块912可被设计成围绕对象以形成隧道类型MRI扫描器(即，闭孔MRI扫描器，参见图9(2))或开放MRI扫描器。
控制模块920可控制MRI扫描器910的磁体模块911和/或RF模块912、处理模块930、和/或显示模块940。仅作为示例，控制模块920可控制X方向、Y方向和Z方向上的磁场梯度。在一些实施例中，控制模块920可从/向MRI扫描器910、处理模块930、和/或显示模块940接收或发送信息。根据一些实施例，控制模块920可从显示模块940接收由例如用户提供的命令，并根据收到的命令来调节磁体模块911和/或RF模块912以拍摄感兴趣对象的图像。仅作为示例，该命令可涉及用户选择的成像应用，该成像应用可以是磁敏感加权成像或磁化率定量成像等。
处理模块930可处理接收自不同模块的不同种类的信息。为进一步理解本公开，以下给出了数个示例，但这些示例并不限定本申请的范围。例如，在一些实施例中，处理模块930可处理接收自RF模块912的MR信号，并且基于这些信号生成一个或多个MR图像并将这些图像递送给显示模块940。在一些实施例中，处理模块930可处理用户或操作者经由显示模块940所作的数据输入并将该数据变换成具体命令，以及将这些命令供应给控制模块920。在一个实施例中，处理模块930可执行如图1至图6中所列举的流程方法，对接收的磁共振信号进行处理，以生成幅值图、相位图、磁化率分布图、R2*(T2*)图或者SWI图等。
显示模块940可接收输入和/或显示输出信息。输入和/或输出信息可包括程序、软件、算法、数据、文本、数字、图像、声音或类似物等或其任何组合。例如，用户或操作者可以输入初始MR参数或条件来发起扫描。作为另一示例，一些信息可以从外部源导入，诸如从软盘、硬盘、无线终端或类似物等或其任何组合导入。在一些实施例中，控制模块920、处理模块930、和/或显示模块940可被集成到图像生成器中。用户可设置MR扫描中的参数，控制成像规程，观看通过图像生成器生产出的图像。在一个实施例中，显示模块940可以设置为CRT显示屏或LCD显示屏，用以显示幅值图、相位图、磁化率分布图、R2*(T2*)图或者SWI图等。
图9(2)为根据本申请的一些实施例的图像生成器的示例性架构。在一些实施例中，控制模块920、处理模块930、和/或显示模块940或其部分或其组合可经由其硬件、软件程序、固件或其组合来被实现在图像生成器上。
图像生成器可包括内部通信总线、中央处理单元(CPU)、I/O接口、COM端口，以及一个或多个存储器器件。内部通信总线可在图像生成器的诸组件之间传送数据。中央处理单元(CPU)可执行计算机指令。计算机指令可涉及例程、程序、对象、组件、数据结构、规程、模块等。在一些实施例中，CPU可处理接收来自MRI扫描器910、控制模块920或MRI系统900的任何其他组件的数据或信息。在一些实施例中，CPU可包括一个或多个处理器。处理器可包括微控制器、微处理器、精简指令集计算机(RISC)、应用专用集成电路(ASIC)、应用专用指令集处理器(ASIP)、中央处理单元(CPU)、图形处理单元(GPU)、物理处理单元(PPU)、微控制器单元、数字信号处理器(DSP)、现场可编程门阵列(FPGA)、高级RISC机(ARM)、可编程逻辑器件(PLD)、能够执行一个或多个功能的任何电路或处理器或类似物等或其任何组合。例如，处理器可包括用以处理接收自MRI扫描器110的MRI数据以供图像重构的微控制器。
一个或多个存储器器件可存储接收来自MRI扫描器910的数据或信息。在一些实施例中，存储器器件可包括盘、随机存取存储器(RAM)、只读存储器(ROM)或类似物等或其任何组合。盘可由例如磁盘、光盘、软盘、光碟或zip盘等来实现。RAM可由例如动态RAM(DRAM)、双倍数据率同步动态RAM(DDR SDRAM)、静态RAM(SRAM)、晶闸管RAM(T-RAM)、以及零电容器RAM(Z-RAM)等来实现。ROM可由例如掩模ROM(MROM)、可编程ROM(PROM)、可擦式可编程ROM(PEROM)、电可擦式可编程ROM(EEPROM)、压缩盘ROM(CD-ROM)、以及数字多用盘ROM等来实现。在一些实施例中，存储器器件可存储一个或多个程序和/或指令以执行本公开中描述的示例性方法。例如，ROM可存储用于基于MR数据来重构MR图像的程序或算法。
图像生成器可包括连接到网络的一个或多个COM端口以供应数据通信。这些通信端口(COM端口)可经由网络向/从MRI扫描器910传送或接收信息。在一些实施例中，通信端口可包括有线端口(例如，通用串行总线(USB)、高清多媒体接口(HDMI)、无线端口(诸如蓝牙端口、红外接口、以及WiFi端口)或类似物等或其任何组合。
I/O接口可支持图像生成器与一个或多个外围设备之间的信息输入或输出。在一些实施例中，外围设备可包括终端、键盘、触摸屏、光标控制设备、遥控器或类似物等或其任何组合。终端可包括例如移动设备(例如，智能电话、智能手表、微型计算机或类似物等)、个人计算机或类似物等或其任何组合。例如，终端可由计算机来实现，计算机可为通用计算机或专门设计的计算机。光标控制设备可包括鼠标、轨迹球或光标方向键以将方向信息和命令选择传达给例如处理模块930或控制显示器设备上的光标移动。
经由I/O接口的信息输入和/或输出可包括程序、软件、算法、数据、文本、数字、图像、声音或类似物等或其任何组合。例如，用户可输入某些初始参数或条件以发起MRI数据处理。在一些实施例中，经由I/O接口的信息输入可经由键盘、触摸屏、语音传感器、运动传感器、脑监视系统或任何其他设备来输入。经由I/O接口输出的信息可被传送到显示模块940、扬声器、打印机、计算设备或类似物等或其组合。
本领域技术人员将认识到，本文的教导顺应于各种修改和/或增强。例如，尽管本文中描述的各种组件的实现可被实施在硬件设备中，但是它也可被实现为仅软件方案，例如安装在现存服务器上。另外，如本文中所公开的图像处理设备可被实现为固件、固件/软件组合、固件/硬件组合或硬件/固件/软件组合。
如图10所示为处理模块930的结构框图，该处理模块930可以是中央处理单元(CPU)、应用专用集成电路(ASIC)、应用专用指令集处理器(ASIP)、图形处理器(GPU)、物理处理单元(PPU)、微控制器、数字信号处理器(DSP)、现场可编程门阵列(FPGA)、ARM或类似物等或其任何组合，其包括：
获取单元1001，用于根据初始场图和误差限定信息，获取磁化率分布图，所述误差限定信息包括先验信息、伪影信息中的至少一种。在一个实施例中，获取单元1001可接收MRI扫描器910采集的磁共振信号，并根据受检者目标部位的磁共振信号确定初始场图。可选地，获取单元1001可根据磁共振信号幅值图、相位图、磁化率分布图、R2*(T2*)图或者SWI图中的一种或多种，并从上述图像中确定空间平滑性、空间梯度、边缘信息、加权蒙版以作为初始先验信息。在另一实施例中，获取单元1001可根据初始先验信息和初始场图获取初始磁化率分布图。
更新单元1002，用于对所述误差限定信息与所述磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差。更新单元1002可从获取单元1001确定初始磁化率分布图，对初始磁化率分布图进行一次或多次迭代优化获取目标部位对应的磁化率分布图，且在磁化率分布图迭代优化的过程中，更新先验信息和/或场图对应的残余放射条状伪影。在一个实施例中，更新单元1002可首先更新先验信息，并基于更新后的先验信息重新确定磁化率分布图。在一个实施例中，更新单元1002可首先更新场图对应的残余放射条状伪影，并基于更新后的场图对应的残余放射条状伪影重新确定磁化率分布图。
处理单元1003，从更新单元1002获取更新后的磁化率分布图，用于当更新后的磁化率分布图的计算收敛误差小于或等于误差阈值时，停止互相迭代更新过程，并输出该次更新后的磁化率分布图。
其中，每次互相迭代更新的过程包括：
通过当前磁化率分布图对当前误差限定信息进行更新；
通过更新后的误差限定信息对当前磁化率分布图进行更新。
上述发明实施例提供了一种磁共振成像装置，首先根据初始场图以及误差限定信息，获取磁化率分布图，然后对误差限定信息与磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差。在互相迭代更新过程中，通过动态的更新先验信息和/或伪影信息，可以使先验信息或伪影信息逐步接近真实值，从而在迭代过程中逐步减小磁化率分布图的计算收敛误差，可以得到计算收敛误差小于或等于误差阈值、精确度较高的磁化率分布图。
所属领域的技术人员可以清楚地了解到，为描述的方便和简洁，上述描述的系统，装置和单元的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。
在本发明所提供的几个实施例中，应该理解到，所揭露的系统，装置和方法，可以通过其它的方式实现。例如，以上所描述的装置实施例仅仅是示意性的，例如，所述单元的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如，多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。另一点，所显示或讨论的相互之间的耦合或直接耦合或通信连接可以是通过一些接口，装置或单元的间接耦合或通信连接，可以是电性，机械或其它的形式。
所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本实施例方案的目的。
另外，在本发明各个实施例中的各功能单元可以集成在一个处理单元中，也可以是各个单元单独物理存在，也可以两个或两个以上单元集成在一个单元中。上述集成的单元既可以采用硬件的形式实现，也可以采用硬件加软件功能单元的形式实现。
上述以软件功能单元的形式实现的集成的单元，可以存储在一个计算机可读取存储介质中。上述软件功能单元存储在一个存储介质中，包括若干指令用以使得一台计算机装置(可以是个人计算机，服务器，或者网络装置等)或处理器(Processor)执行本发明各个实施例所述方法的部分步骤。而前述的存储介质包括：U盘、移动硬盘、只读存储器(Read-Only Memory，ROM)、随机存取存储器(Random Access Memory，RAM)、磁碟或者光盘等各种可以存储程序代码的介质。
以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明保护的范围之内。

Claims (10)

1.一种磁共振成像方法，其特征在于，所述方法包括：
根据初始场图和误差限定信息，获取磁化率分布图，所述误差限定信息包括先验信息或伪影信息，所述先验信息提供关于磁化率空间分布的特征，所述伪影信息为QSM求解过程中产生的放射条状伪影；
对所述误差限定信息与所述磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差；
当更新后的磁化率分布图的计算收敛误差小于或等于误差阈值时，停止互相迭代更新过程，并输出该次更新后的磁化率分布图；
其中，每次互相迭代更新的过程包括：
通过当前磁化率分布图对当前误差限定信息进行更新；
通过更新后的误差限定信息对当前磁化率分布图进行更新。
2.根据权利要求1所述的方法，其特征在于，当所述误差限定信息包括先验信息时，所述通过当前磁化率分布图对当前误差限定信息进行更新，包括：
提取当前磁化率分布图的特征，得到与当前磁化率分布图对应的特性图；
在所述特性图中提取更新后的先验信息。
3.根据权利要求2所述的方法，其特征在于，当所述特性图为梯度特性图时，所述提取当前磁化率分布图的特征，得到与当前磁化率分布图对应的特性图，包括：
提取当前磁化率分布图沿参考方向的梯度特征，得到所述梯度特性图，所述参考方向包括X方向、Y方向或Z方向的一种或多种。
4.根据权利要求1所述的方法，其特征在于，当所述误差限定信息包括先验信息时，所述通过当前磁化率分布图对当前误差限定信息进行更新，包括：
基于正向偶极效应模型方程，将当前磁化率分布图转换为先验场图；
根据所述先验场图，生成更新后的先验信息。
5.根据权利要求1所述的方法，其特征在于，当所述误差限定信息包括伪影信息时，所述通过当前磁化率分布图对当前误差限定信息进行更新，包括：
基于正向偶极效应模型方程，将当前磁化率分布图转换为估计场图；
根据所述初始场图以及所述估计场图，生成差值场图；
基于正向偶极核函数，根据所述差值场图生成更新后的伪影信息。
6.根据权利要求1所述的方法，其特征在于，当所述误差限定信息包括伪影信息时，所述通过更新后的误差限定信息对当前磁化率分布图进行更新，包括：
以正则化的形式将更新后的伪影信息带入所述互相迭代更新过程，得到更新后的磁化率分布图；或者，
在当前磁化率分布图中去除更新后的伪影信息，得到更新后的磁化率分布图。
7.根据权利要求1至6中任意一项所述的方法，其特征在于，所述先验信息包括空间平滑性、空间梯度、边缘信息、加权蒙版中的至少一种。
8.根据权利要求1所述的方法，其特征在于，以以下多种迭代方法中的至少一种迭代方法，对所述误差限定信息与所述磁化率分布图进行至少一次互相迭代更新：
共轭梯度法、二次最小法与分裂-布雷格曼法。
9.一种磁共振成像装置，其特征在于，所述装置包括：
获取单元，用于根据初始场图和误差限定信息，获取磁化率分布图，所述误差限定信息包括先验信息或伪影信息，所述先验信息提供关于磁化率空间分布的特征，所述伪影信息为QSM求解过程中产生的放射条状伪影；
更新单元，用于对所述误差限定信息与所述磁化率分布图进行至少一次互相迭代更新，并在每次互相迭代更新后确定磁化率分布图的计算收敛误差；
处理单元，用于当更新后的磁化率分布图的计算收敛误差小于或等于误差阈值时，停止互相迭代更新过程，并输出该次更新后的磁化率分布图；
其中，每次互相迭代更新的过程包括：
通过当前磁化率分布图对当前误差限定信息进行更新；
通过更新后的误差限定信息对当前磁化率分布图进行更新。
10.一种医学成像系统，其特征在于，所述医学成像系统包括处理器以及存储器；所述存储器用于存储指令，所述指令被所述处理器执行时，导致所述医学成像系统实现如权利要求1至8任一项所述的方法。
CN201810323015.XA 2018-04-11 2018-04-11 一种磁共振成像方法、装置及医学成像系统 Active CN108714028B (zh)

Priority Applications (3)

Application Number Priority Date Filing Date Title
CN201810323015.XA CN108714028B (zh) 2018-04-11 2018-04-11 一种磁共振成像方法、装置及医学成像系统
US16/247,098 US10950011B2 (en) 2018-04-11 2019-01-14 Method and system for determining magnetic susceptibility distribution
US17/188,918 US11587269B2 (en) 2018-04-11 2021-03-01 Method and system for determining magnetic susceptibility distribution

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810323015.XA CN108714028B (zh) 2018-04-11 2018-04-11 一种磁共振成像方法、装置及医学成像系统

Publications (2)

Publication Number Publication Date
CN108714028A CN108714028A (zh) 2018-10-30
CN108714028B true CN108714028B (zh) 2022-02-25

Family

ID=63898927

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810323015.XA Active CN108714028B (zh) 2018-04-11 2018-04-11 一种磁共振成像方法、装置及医学成像系统

Country Status (2)

Country Link
US (2) US10950011B2 (zh)
CN (1) CN108714028B (zh)

Families Citing this family (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108829639B (zh) * 2018-07-06 2023-10-27 上海联影医疗科技股份有限公司 一种磁共振成像方法和设备
US10890640B2 (en) 2019-03-18 2021-01-12 Uih America, Inc. Systems and methods for signal representation determination in magnetic resonance imaging
US11069063B2 (en) 2019-03-18 2021-07-20 Uih America, Inc. Systems and methods for noise analysis
CN110174632B (zh) * 2019-06-10 2021-06-01 上海东软医疗科技有限公司 磁共振成像方法、装置、成像设备及磁共振成像系统
US11796618B2 (en) 2019-07-12 2023-10-24 Shanghai United Imaging Healthcare Co., Ltd. Systems and methods for magnetic resonance imaging
US11423593B2 (en) * 2019-12-19 2022-08-23 Shanghai United Imaging Intelligence Co., Ltd. Systems and methods for reconstructing a medical image using meta learning
CN113050010B (zh) * 2019-12-26 2023-03-21 上海联影医疗科技股份有限公司 噪音分析的系统、方法
WO2021092600A2 (en) * 2020-12-14 2021-05-14 Innopeak Technology, Inc. Pose-over-parts network for multi-person pose estimation
CN112860267B (zh) * 2021-04-23 2021-07-30 武汉深之度科技有限公司 一种内核裁剪方法及计算设备
CN113655415B (zh) * 2021-08-16 2023-01-17 东北大学 一种磁场分布的增强现实在线可视化方法
CN118067158A (zh) * 2022-11-23 2024-05-24 中国科学院深圳先进技术研究院 基于转换区域和局部场图迭代的化学位移编码成像方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103293498A (zh) * 2012-02-27 2013-09-11 医影量有限责任公司 提供磁化率定量成像的系统和方法
CN103764025A (zh) * 2011-10-12 2014-04-30 株式会社日立制作所 磁共振成像装置以及磁化率强调图像生成方法
CN107072592A (zh) * 2014-11-11 2017-08-18 株式会社日立制作所 磁共振成像装置以及定量性磁化率匹配方法

Family Cites Families (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8781197B2 (en) * 2008-04-28 2014-07-15 Cornell University Tool for accurate quantification in molecular MRI
US8422756B2 (en) 2010-04-27 2013-04-16 Magnetic Resonance Innovations, Inc. Method of generating nuclear magnetic resonance images using susceptibility weighted imaging and susceptibility mapping (SWIM)
WO2012037067A1 (en) 2010-09-14 2012-03-22 Massachusetts Institute Of Technology Multi-contrast image reconstruction with joint bayesian compressed sensing
US9285449B2 (en) 2011-06-15 2016-03-15 Chunlei Liu Systems and methods for imaging and quantifying tissue magnetism with magnetic resonance imaging
US8886283B1 (en) * 2011-06-21 2014-11-11 Stc.Unm 3D and 4D magnetic susceptibility tomography based on complex MR images
US9664761B2 (en) * 2013-06-26 2017-05-30 Medimagemetric LLC Joint estimation of chemical shift and quantitative susceptibility map using MRI signal
US9612300B2 (en) * 2013-11-25 2017-04-04 Wisconsin Alumni Research Foundation System and method for object-based initialization of magnetic field inhomogeneity in magnetic resonance imaging
US9542763B2 (en) * 2014-04-25 2017-01-10 The General Hospital Corporation Systems and methods for fast reconstruction for quantitative susceptibility mapping using magnetic resonance imaging
US10126397B2 (en) 2014-05-09 2018-11-13 The General Hospital Corporation Systems and methods for fast magnetic resonance image reconstruction using a heirarchically semiseparable solver

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103764025A (zh) * 2011-10-12 2014-04-30 株式会社日立制作所 磁共振成像装置以及磁化率强调图像生成方法
CN103293498A (zh) * 2012-02-27 2013-09-11 医影量有限责任公司 提供磁化率定量成像的系统和方法
CN107072592A (zh) * 2014-11-11 2017-08-18 株式会社日立制作所 磁共振成像装置以及定量性磁化率匹配方法

Also Published As

Publication number Publication date
US20190318511A1 (en) 2019-10-17
CN108714028A (zh) 2018-10-30
US11587269B2 (en) 2023-02-21
US10950011B2 (en) 2021-03-16
US20210183118A1 (en) 2021-06-17

Similar Documents

Publication Publication Date Title
CN108714028B (zh) 一种磁共振成像方法、装置及医学成像系统
CN111513716B (zh) 使用扩展灵敏度模型和深度神经网络进行磁共振图像重建的方法和系统
Cruz et al. Accelerated motion corrected three‐dimensional abdominal MRI using total variation regularized SENSE reconstruction
JP5843876B2 (ja) 磁気共鳴イメージング装置および磁化率強調画像生成方法
Zhao et al. Image reconstruction from highly undersampled (k, t)-space data with joint partial separability and sparsity constraints
JP5902317B2 (ja) 磁気共鳴イメージング装置および定量的磁化率マッピング法
Sumpf et al. Fast T2 mapping with improved accuracy using undersampled spin-echo MRI and model-based reconstructions with a generating function
US8217652B2 (en) Spatial intensity correction for RF shading non-uniformities in MRI
US11348230B2 (en) Systems and methods for generating and tracking shapes of a target
US8942445B2 (en) Method and system for correction of lung density variation in positron emission tomography using magnetic resonance imaging
US9430854B2 (en) System and method for model consistency constrained medical image reconstruction
US9709641B2 (en) Magnetic resonance imaging apparatus, image processing apparatus, and susceptibility map calculation method
US9928596B2 (en) Motion corrected imaging system
US10895622B2 (en) Noise suppression for wave-CAIPI
EP2526437A1 (en) Electric properties tomography imaging method and system
Roeloffs et al. Model‐based reconstruction for T1 mapping using single‐shot inversion‐recovery radial FLASH
CN107290700B (zh) 一种相位校正方法、装置及磁共振系统
CN113924503B (zh) 时域磁共振的参数图确定
van der Heide et al. High‐resolution in vivo MR‐STAT using a matrix‐free and parallelized reconstruction algorithm
Tian et al. Evaluation of pre‐reconstruction interpolation methods for iterative reconstruction of radial k‐space data
Zou et al. Quantitative MR relaxation using MR fingerprinting with fractional-order signal evolution
Zhang et al. Improved dynamic contrast-enhanced mri using low rank with joint sparsity
US11022667B2 (en) System and method for image reconstruction
Mickevicius et al. On the use of low-dimensional temporal subspace constraints to reduce reconstruction time and improve image quality of accelerated 4D-MRI
Duan et al. Adaptive transform learning and joint sparsity based PLORAKS parallel magnetic resonance image reconstruction

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
CB02 Change of applicant information

Address after: 201807 Shanghai city Jiading District Industrial Zone Jiading Road No. 2258

Applicant after: Shanghai Lianying Medical Technology Co., Ltd

Address before: 201807 Shanghai city Jiading District Industrial Zone Jiading Road No. 2258

Applicant before: SHANGHAI UNITED IMAGING HEALTHCARE Co.,Ltd.

CB02 Change of applicant information
GR01 Patent grant
GR01 Patent grant