CN108710092B - 基于数值模拟技术的仿真核磁共振波谱分析方法及分析系统 - Google Patents

Abstract

本发明公开了一种基于数值模拟技术的仿真核磁共振波谱分析方法及系统，所述方法包括构建模拟样品分子信息模型及模拟样品分子信息库；构建波谱仪硬件参数；构建波谱序列物理模型；基于前述构建的所述模拟样品分子信息模型及模拟样品分子信息库、波谱仪仪相关硬件参数和波谱序列物理模型，构建核磁共振波谱原始数据采集物理数学模型，进行模拟原始数据采集算法和程序实现，完成原始数据的模拟采集；将前述模拟采集的原始数据经傅立叶变换得到波谱图。本发明方法及系统可应用于核磁共振波谱应用技术的配套虚拟实验实训平台以及相关技术领域应用人员自学等多方面，具有广泛的应用前景。

Description

基于数值模拟技术的仿真核磁共振波谱分析方法及分析系统

技术领域

本发明涉及核磁共振技术领域，具体涉及一种基于数值模拟技术的仿真核磁共振波谱分析方法及分析系统。

背景技术

核磁共振波谱是化学系研究分子和物质结构的重要工具，是国产化程度最低的分析仪器；设备昂贵，台套数有限，科研机时安排紧凑，教学实验机时非常有限；担心学生实验操作会退化设备性能，调坏相关参数；新手易出现操作失误导致设备损坏；实验流程复杂，样品配备，匀场、锁场和调谐匹配需要缓慢仔细调节，占用时间多等原因，导致波谱实验开展薄弱。

核磁共振波谱分析技术在有机化学相关领域得到越来越广泛的应用，比如制药、石油化工、有机材料、橡胶材料、生物医药等；越来越多的科研人员需要使用核磁共振波谱分析仪开展相关科研和测试工作。但是核磁共振波谱分析技术涉及到众多的核磁共振技术基础，尤其是序列技术，加上待测样品的复杂多变性，使得仪器的使用和测试数据结果的解释使用出现了很多问题。这些都是源于各领域的从业人员对于核磁共振波谱分析技术的基础原理不了解，对波谱仪的操作流程以及相关参数设置对谱图的影响规律不清楚。为了使相关从业人员更好的快速理解掌握核磁共振波谱分析技术的原理，更好的设置采集参数和硬件调节，从而快速高效准确的操作实验，采集数据并处理谱图，用以开展真机操作前的模拟培训，本发明可以满足这个要求。

发明内容

本发明提出了一种基于数值模拟技术的仿真核磁共振波谱分析方法及分析系统，克服了现有技术存在的上述技术问题。

本发明提出的基于数值模拟技术的仿真核磁共振波谱分析方法，其包括：

构建模拟样品分子信息模型及模拟样品分子信息库；

构建波谱仪硬件参数；

构建波谱序列物理模型；

基于前述构建的所述模拟样品分子信息模型及模拟样品分子信息库、波谱仪仪相关硬件参数、和波谱序列物理模型，结合起来，构建核磁共振波谱原始数据采集物理数学模型，进行模拟原始数据采集算法和程序实现，完成原始数据的模拟采集；

将前述模拟采集的原始数据，进行FT变换，得到波谱图。即获得信号频率的一维/二维或三维分布图，也就是一维/二维或三维波谱图。

进一步地，本发明方法还包括：将上一步得到的波谱图进行谱图处理。所述谱图处理包括相位校正、基线校正、寻峰、积分、测距、化学位移校正等谱图处理方法。

进一步地，本发明方法还包括：构建所述波谱仪硬件参数的校正模型及模拟调节，用以行硬件系统参数校正，为保证与真机系统操作一致。所述模拟调节包括wobb(自动调谐匹配)、shim(匀场)、lock(锁场)、P1调节等。所述校正模型的构建涉及到物理数学模型的构建。现有技术中在测试样品之前需要将仪器系统校正到最优状态，在每次实验之前都需要进行校正，才能保证最终谱图是有效的，而这些步骤是波谱实验中最花费时间、难度很大。本发明则通过基于数值模拟技术对硬件参数进行校正模型的构建及模拟调节，解决前述提到的这些技术问题。

其中，步骤1中，构建模拟样品分子信息模型及模拟样品分子信息库，构建常用待分析样品分子的信息模型以及信息库。所述样品分子信息按照自旋体系来确定，可以是两个自旋体系：一个为AmBnCd，A、B、C为相互之间存在偶合的三种自旋核，核子数目分别为m，n，d。另一个为XpYqZw，X、Y、Z为相互之间存在偶合的另外三种自旋核，核子数目分别为p，q，w。本发明中，样品分子信息模型不局限于前述两个自旋体系，还可以是三个、四个或更多自旋体系。每个自旋体系内的核子种类不局限于三种，还可以是四种、五种或更多核子种类。优选地，只是自旋体系之间的核是没有J间接偶合的。

具体地，对于AmBnCd，构建其化学位移信息为Wa，Wb，Wc，分别表示A，B，C三种核的化学位移值，量纲为ppm；构建两两之间的J偶合信息为Jab，Jbc，Jac，分别表示AB核，BC核，AC核之间的J偶合值，单位为Hz。构建两两之间的D偶合信息为Dab，Dbc，Dac，分别表示AB核，BC核，AC核之间的直接相互偶合值，单位为Hz。XYZ或其他自旋体系依次类推。

本发明中，为便于扩展分子信息库，也可以由用户自行根据上述规则进行编辑新的分子信息，增加进信息库后，即可成为一种新的样品分子。自行构建分子信息可以通过EXCEl数据库添加，或者通过软件界面参数添加或其他形式添加。主要包括自旋体系数(1，2，3，4或更多)，自旋种类(1，2，3或4)，每种自旋的数目、化学位移、J偶合和D值；需要注意的是，两个不同自旋体系之间的自旋核之间是独立的，即不存在间接偶合J或直接偶合D。

其中，分子样品可以为液体样品分子，固体样品分子。固体样品分子主要在于其D值差异。

其中，所述弛豫分析相关硬件参数包括主磁场B0、主磁场不均匀性deltaB0、射频场B1、梯度场G，PL1W，P1(或称为P90)，P2(或称为P180)等。

其中，所述波谱序列包括单脉冲序列、去偶序列(包括dept45，dept90，dept35等)、同核化学位移相关二维谱COSY序列，同核化学位移相关overhause增强二维谱NOESY序列，异核化学位移相关谱(1H-13C HSQC序列)、异核化学位移相关谱(1H-13C HMBC序列)等。

本发明中序列构建则是通过给出不同的数据采集模式的控制参数来实现。除了采用上述常用的波谱序列，还可以自行构建新的波谱序列方法和模型。上述序列为本领域基本理论。而现有技术中实际仪器采用序列发生器(硬件)产生一系列控制信号驱动相关硬件完成相关工作。

其中，所述波谱图为一维谱图、二维谱图或三维谱图。

本发明还提出了一种基于数值模拟技术的仿真核磁共振波谱分析系统，用于实现本发明所述的基于数值模拟技术的仿真核磁共振波谱分析方法。所述系统包括：模拟样品分子信息模型及模拟样品分子信息库的构建模块，用于构建模拟样品分子信息模型及模拟样品分子信息库；波谱仪相关硬件参数的构建模块，用于构建波谱仪相关硬件参数；波谱序列物理模型的构建模块，用于构建不同波谱序列物理模型；核磁共振FID数据采集物理数学模型的构建模块，用于构建核磁共振FID数据采集物理数学模型，用于原始数据的模拟采集；以及波谱图模块，用于得到信号频率的一维、二维或三维波谱图。

进一步，还包括谱图处理模块，用于对前述所得的波谱图进行谱图处理，包括进行相位校正、基线校正、寻峰、积分、测距、化学位移校正等谱图处理。

进一步，还包括波谱硬件参数的模拟调节校正模型的构建模块，用于构建波谱硬件参数的模拟调节校正模型。所述模拟调节校正模块包括自动调谐匹配wobb、匀场shim、锁场lock、P1调节等模块进行调节校正。

本发明还提出了基于数值模拟技术的仿真核磁共振波谱分析系统的用途，所述分析系统应用于核磁共振波谱应用技术的配套虚拟实验实训平台中，还可以应用于相关波谱分析技术领域应用人员自学技术原理中。本发明还提出了基于本发明内容可以开发出核磁共振弛豫谱分析技术仿真分析仪，该仿真分析仪可作为一般核磁共振弛豫谱应用技术的配套虚拟实验实训平台或者作为相关弛豫谱分析技术领域应用人员用于自学技术原理。本发明分析方法及系统为仿真分析仪即为一款软件(除了软件运行的计算机载体之外，不含任何硬件)。本发明可以实现与一台波谱分析仪(含硬件和软件)完全相同的数据采集和处理功能。

本发明模拟Bruker的核磁共振波谱仪的TopSpin2.0的基本功能,可以实现常见小有机分子的液体核磁共振波谱实验的基本过程。其步骤及功能包括：1，采用常见的有机小分子的化学位移和J耦合信息作为样品信息；2，可模拟样品的进入(inject)和弹出(eject)磁体腔功能；3，可模拟实现探头的手动/自动调谐匹配过程；4，可模拟锁场/匀场调节过程；5，可模拟设置采样参数进行样品信号采集的过程；6，可模拟进行基本谱图处理和分析功能。7，可完成一维谱/二维谱的数据采集与谱图处理。

本发明优势及有益效果包括：无需庞大的硬件配套下，可获得与核磁共振波谱仪相同的实验效果，避免硬件成像教学仪对硬件的高要求。可用作实际样品的测量效果的参考比对以及参数设置了解。本发明有益效果还包括，其技术特点可用于直观快速的学习了解核磁共振波谱技术的原理；可用于波谱相关专业人员的早期自学原理和操作，熟悉参数的设置规则；也可以用于实际测试样品的信号和波谱结果比对参考。

本发明采用数值模拟仿真技术，与一般仿真仪器操作的表象不同，仿真的是本质同时也包括操作过程。虚拟仿真技术是当前的发展热点。虚拟技术应用于大型昂贵设备和高新技术的教学培训是必然趋势。与通常意义上的基于结构和过程的形象化仿真虚拟技术不同，本发明项目将极其抽象晦涩的核磁共振波谱理论全部用数学公式进行表达，通过程序进行实现。最终的软件可以完全摆脱硬件，而实现与硬件仪器一模一样的实验效果。换言之，通常的虚拟技术主要用基于动画设计等技术。本发明是核磁共振理论本质和计算机程序开发技术的完美结合。本发明数字模拟技术仿真，仿真效果上注重本质相同，而非结构形象化。

本发明应用之一，可用于相关领域人才的实验教学和实践培训；不需要硬件(硬件单台价格在200万到1000万之间)，不受制于硬件，实现低成本批量化的实验操作与培训，充分掌握数据采集参数和波谱处理过程对不同样品分子波谱的影响规律，从而真正掌握波谱分析技术理论，提高分析技术分析和操作能力，为后续工作奠定基础。

本发明又一应用，可用于仪器开发和测试校准；实际波谱仪器受到众多的因素影响，比如电子学噪声、涡流、系统死时间等，导致无法评价最终的波谱是否受到设备硬件的影响，还是真实反应了样品的特性。应用本虚拟系统，可将标准样品得出理想硬件条件下的波谱作为测试标准，从而评价或制定硬件指标。

本发明又一应用，可用做同种样品分子在不同设备下谱图结果参照；由于波谱仪设备昂贵，所以难以得到同种样品在不同设备(场强不同，均匀性也可能不同)下的谱图表现。采用本虚拟谱仪，可以快速给出同一个分子在不同场强(60MHz到1.2GHz)设备下不同均匀性(1ppm到1ppb)的谱图表现。用以产生理想谱图和不同设备下的谱图参照物。

附图说明

图1A为本发明方法流程示意图,图1B表示模拟样品分子信息，包括J偶合，化学位移等信息，通过构建不同的样品分子信息作为输入。图1C表示按图1A流程得到的输出谱图，即图1B表示的分子信息通过图1A方法虚拟数据采集和FT后可以得到图1C表示的谱图输出。图1D为本发明实施例的流程示意图。

图2为本发明系统结构示意图。

图3为本发明方法中原始数据模拟采集流程示意图。

图4为本发明实例的VMRS主界面，通过计算机编程软件(本实例为matlab软件)进行界面编程后的仿真系统界面。

图5为本发明实施例中wobb界面示意图。

图6为本发明实施例中shim界面示意图。

图7为本发明实施例中PL1W调节界面示意图。

图8为本发明实施例中章动谱效果示意图。

图9为本发明实施例中相位校正成功后的谱图示意图。

图10为本发明实施例中谱图积分效果图示意图。

图11为本发明实施例中谱图自动寻峰效果图示意图。

图12为应用本发明虚拟采集得到的二维谱图contour格式显示的示意图。

图13为应用本发明虚拟采集得到的二维谱图mesh格式显示的示意图。

图14为应用本发明虚拟采集得到的二维谱数据经一维FT后的数据(部分)，可看到不同谱峰的相位变化。本图只显示了25条一维FT后的数据(一维数组)，本实例总数据为100条。数据总条数为num_t1。

具体实施方式

结合以下具体实施例和附图，对本发明作进一步的详细说明。实施本发明的过程、条件、实验方法等，除以下专门提及的内容之外，均为本领域的普遍知识和公知常识，本发明没有特别限制内容。

本发明基于数值模拟技术的仿真核磁共振波谱分析方法，参考图1A、图1D所示，包括步骤：

步骤1，根据待分析样品分子的分子信息模型和分子信息库；分子信息模型构建完全由本发明完成。具体分子信息库，可由使用者按照信息模型要求自行编辑后添加。

步骤2，构建波谱仪相关的硬件参数，包括主磁场B0，射频场B1，梯度场G，主磁场不均匀性deltaB0，射频场功率PL1W等；

步骤3，构建不同波谱脉冲序列的物理数学模型；

步骤4，将上述步骤1、2、3所构建的样品分子信息模型和分子信息库，波谱仪设备硬件参数、不同波谱脉冲序列的物理数学模型序列参数结合在一起，构建核磁共振波谱原始数据采集物理数学模型，进行模拟数据采集算法和程序实现，完成原始FID数据的模拟采集。

步骤5，对步骤4模拟采集的原始数据进行经傅立叶变换(FT变换)，得到一维/二维/三维波谱图。

进一步，步骤6，对步骤5得到的波谱图进行相位校正、基线校正、寻峰、积分、测距、化学位移校正等谱图处理；

进一步，构建所述波谱仪硬件参数的校正及模拟调节模型；所述校正及模拟调节包括自动调谐匹配wobb、匀场shim、锁场lock、P1调节等。

模拟样品分子的一些信息，包括J偶合，化学位移等信息，通过构建不同的样品分子信息作为输入，如图1B所示，如图1B表示的分子信息通过如图1A、图1D所示的本发明方法虚拟数据采集和FT后可得到如图1C表示的谱图输出。

如图2所示，本发明还提出了一种基于数值模拟技术的仿真核磁共振波谱分析系统，用于实现本发明所述的基于数值模拟技术的仿真核磁共振波谱分析方法。所述系统包括：模拟样品分子信息模型及模拟样品分子信息库的构建模块，用于构建模拟样品分子信息模型及模拟样品分子信息库；波谱仪相关硬件参数的构建模块，用于构建波谱仪相关硬件参数；波谱序列物理模型的构建模块，用于构建不同波谱序列物理模型；核磁共振FID数据采集物理数学模型的构建模块，用于构建核磁共振FID数据采集物理数学模型，用于原始数据的模拟采集；以及波谱图模块，用于得到信号频率的一维、二维或三维波谱图。

进一步，还包括谱图处理模块，用于对前述所得的波谱图进行谱图处理，包括进行相位校正、基线校正、寻峰、积分、测距、化学位移校正等谱图处理。

进一步，还包括波谱硬件参数的模拟调节校正模型的构建模块，用于构建波谱硬件参数的模拟调节校正模型。所述模拟调节校正模块包括自动调谐匹配wobb、匀场shim、锁场lock、P1调节等模块进行调节校正。

具体地，对于步骤1，本发明给出分子信息构建模型，但不限于这种模型，其余以实现本发明目的的样品模型构建都基于本发明构思，属于本发明保护范围。

构建常用待分析样品分子信息模型以及分子信息库。所述分子信息按照自旋体系来确定，可以是两个自旋体系，其中一个为AmBnCd，另一个为XpYqZw；ABC为相互之间存在偶合的三种自旋核，核子数目分别为下标m，n，d。XYZ为相互之间存在偶合的另外三种自旋核，核子数目分别为下标p，q，w。具体的分子信息模型不局限于两个自旋体系，还可以是三个，四个或更多。每个自旋体系内的核子种类也可以不局限于三种，可以是四种、五种或更多。只是自旋体系之间的核可以认为是没有J偶合的。

对于AmBnCd，构建其化学位移信息为Wa，Wb，Wc，分别表示A，B，C三种核的化学位移值，量纲为ppm；构建两两之间的J偶合信息为Jab，Jbc，Jac，分别表示AB核，BC核，AC核之间的J偶合值，单位为Hz。构建两两之间的D偶合信息为Dab，Dbc，Dac，分别表示AB核，BC核，AC核之间的直接相互偶合值，单位为Hz。同样地，XYZ或其他自旋体系依次类推。

进一步，为便于扩展分子信息库，还可以由用户自行根据上述规则进行编辑新的分子信息，增加进信息库后，即可成为一种新的样品分子信息。自行构建分子信息可以通过EXCEl数据库添加，或者通过软件界面参数添加或其他形式添加。样品分子信息库主要包括自旋体系数(1，2，3，4或更多)，自旋种类(1，2，3或4)，每种自旋的数目、化学位移、J偶合和D值。需要注意的是，两个不同自旋体系之间的自旋核之间是独立的，即不存在间接偶合J或直接偶合D。

步骤2中，构建弛豫分析仪相关的硬件参数，包括主磁场B0，主磁场不均匀性deltaB0，射频场B1，梯度场G，PL1W(射频功率)等。

步骤3中，构建不同波谱序列的物理模型，包括单脉冲序列、去偶序列(dept45，dept90，dept35等)、同核化学位移相关二维谱COSY序列，同核化学位移相关overhause增强二维谱NOESY序列，异核化学位移相关谱(1H-13C HSQC序列)、异核化学位移相关谱(1H-13CHMBC序列)等。不限于这几种序列；具体波谱序列模型包括脉冲种类、脉冲个数、脉冲时间时间、脉冲幅度、间隔时间、梯度施加时间和间隔时间、梯度幅度、梯度步进值、步进次数、数据采集时间、采样间隔、采样点数，反转时间间隔和次数以及各事件的时间次序等。

上述波谱序列为基本理论，实际仪器采用序列发生器(硬件)产生一系列控制信号驱动相关硬件完成相关工作。本发明的序列构建则是通过给出不同的数据采集模式的控制参数来实现。本发明不保护基本的序列理论，但是基于本发明目的的构建用以虚拟采集的序列(包括参数)模式是本发明保护内容。

比如单脉冲序列，分别构建以下参数来实现：D0(重复时间)、SW(采样带宽)、TD(采样点数)、TDeff(有效采样点数)、DS(空采次数)、P90(90度脉冲宽度)、P180(180度脉冲宽度)，NS(累加次数)，DE，D1(延时采样时间)，NUC1(H通道)；这些参数都与单脉冲序列数据采集相关，后续综合的数据采集中都会用到这些参数。

比如二维COSY序列，需要构建的参数为：t1-step(第一维时间步进值)、t1-num(第一维个数)、D0(重复时间)、SW(采样带宽)、TD(采样点数)、P90(90度脉冲宽度)、P180(180度脉冲宽度)，t2(t2维时间数组)。

比如HSQC序列，与上述参数相同，增加NUC2(C通道)。

比如DEPT去偶序列，增加去偶频率，去偶脉宽，去偶脉冲波形等。

步骤4中，将上述1，2，3所构建的模拟样品分子信息及分子信息库，硬件参数和波谱序列参数结合在一起，构建核磁共振波谱FID原始数据采集物理数学模型，进行模拟数据采集算法和程序实现，完成原始数据的模拟采集。本发明创新在于提出了以上步骤，与现有技术实践中需要使用核磁共振波谱仪硬件设备相比较，具有实质性的显著区别且具有预料不到的有益效果。

现有技术中的实际波谱仪所用的真实的硬件部件，按照序列发生器产生的指令对真实样品进行信号采集，但具体采集的信号是什么样的，信号遵循什么样的规律，对于现有技术中的硬件而言都是不需要关心的。

而，本发明的信号是通过物理数学模型虚拟产生的。本发明中，不同的波谱序列，虚拟产生的原始数据也是不同。且，本发明既可以模拟一维谱(单脉冲，去偶dept90等)，也可以模拟二维谱(同核化学位移相关谱，异核化学位移相关谱)，或者三维谱。

如图3所示，步骤4中，构建综合原始数据模拟采集物理数学模型，具体地：

一维谱采集模式，对于每一个自旋体系内的每一种自旋核，分别计算每一个自旋体系内的所有的自旋核J偶合因素导致的谱线裂分的幅值项和频率项(即J偶合裂分峰幅值项、J偶合裂分峰频率项)，自身化学位移导致的频率项(即化学位移项)，以及旋转坐标系数据采集项、T2弛豫项和T1弛豫项，磁化矢量翻转项(磁化矢量反转)，将前述各项连乘，分别获得每一种自旋核信号，然后，将每一个自旋核体系下的多个自旋核信号进行连加，得到累加的每一个自旋核体系信号。然后，将多个自旋核体系信号进行连加，再叠加上模拟随机噪声，即可得到模拟采集的核磁共振波谱原始数据。获得的原始数据得以显示并存储。

例如，在第一个自旋核体系(ABC)下，将各项连乘，分别得到自旋核信号A、自旋核信号B、自旋核信号C，将其连加得到累加的第一个自旋核体系信号。同样，在第二个自旋核体系(XYZ)下，将各项连乘，分别得到自旋核信号X、自旋核信号Y、自旋核信号Z，进行连加，得到累加的第二个自旋核体系信号。然后，再将第一个自旋核体系信号、第二个自旋核体系信号连加，再叠加上模拟随机噪声，即可得到模拟采集的核磁共振波谱原始数据。

其中，自旋核体系、自旋核的数量均可以是1-n的正整数，n大于1。

二维谱采集模式下，在一维谱的基础上，考虑每一个自旋体系内的每一种自旋核灾t1维演化期导致的对角峰频率项(对应量子算符Ix)和交叉峰频率项(对应量子算符IzIy)；将对角锋频率项仍然作为该自旋核的乘积因子；将交叉峰频率项则作为J偶合的另一方自旋核的乘积因子，即极化转移的体现。

对于去偶，可通过去偶频率及频率范围，将对应去偶频率或在去偶频率范围内的J偶合置0。

如图3中，J偶合导致的裂分峰的幅度项Na可由下列双重循环得到(以AmBnCd自旋体系的A核为例，A，B，C自旋核的数目分别为m，n，d，相邻J偶合分别为Jab，Jac，Jbc)：

For i＝1:n+1

For j＝1:d+1

Na＝m/n/d*A0(n+1，i)*A0(d+1，j)；

End

End

如图3中，J偶合导致的裂分峰的频率项Ja可由下列双重循环得到(以AmBnCd自旋体系的A核为例，A，B，C自旋核的数目分别为m，n，d，相邻J偶合分别为Jab，Jac，Jbc)

For i＝1:n+1

For j＝1:d+1

exp(-i*Ja*t)＝exp(-i*(I0(n+1，I)*Jab)+I0(d+1，j)*Jac)*t)；

End

End

I0，A0分别为前述构建的7*7的二维矩阵；

上述计算中，如果自旋核种类超过3种，比如4种，则需要进行3重(4-1重)循环；依次类推。

具体地，如一维谱图测量的单脉冲序列的模拟数据采集：如图3所示，

第一步：首先构建数据采样模型，即一个一维离散的时间数组，包含采样点数，空扫次数，延迟采样时间，采样带宽，回波时间，采样时间t等多个参数在内，如t＝(DE:1/SW:(DE+TD/SW)；为序列确定的采样规律，其中包括采样带宽SW，延迟时间DE，采样点数TD；模拟采样数据算法实例为：

第二步，分别构建不同偶合核数目的谱线裂分的峰数数组I0，I0为一个二维数组(矩阵)，行数对应的是不同的J偶合数目，每行的数值表示谱峰裂分的数目和频率偏移相对量。一个具体的I0实例矩阵为7*7的二维矩阵：

I0＝【0，0，0，0，0，0，0；

-1/2，1/2，0，0，0，0，0；

-1，0，1，0，0，0，0，0；

-3/2，-1/2，1/2，3/2，0，0，0；

-2，-1，0，1，2，0，0；

-5/2，-3/2，-1/2，1/2，3/2，5/2，0；

-3，-2，-1，0，1，2，3；】

第一行表示周围没有J偶合核，频率裂分数为0，频率偏移值为0；第二行表示周围有1个J偶合核，裂分数为2，频率偏移相对量分别为-1/2和1/2；第三行表示周围有2个J偶合核，裂分数为3，频率偏移相对量分别为-1，0，1；依次类推。该矩阵最多可以模拟周围有6个J偶合核的情况，基本涵盖了常用的分子情况；

第三步，分别构建不同偶合核数目的谱线裂分的峰的幅值数组A0，A0为一个二维数组(矩阵)，行数对应的是不同的J偶合数目，每行的数值表示谱峰裂分的数目和幅值。一个具体的A0实例矩阵为7*7的二维矩阵：

I0＝【1，0，0，0，0，0，0；

1，1，0，0，0，0，0；

1，2，1，0，0，0，0，0；

1，3，3，1，0，0，0；

1，4，6，4，1，0，0；

1，5，10，10，5，1，0；

1，6，15，20，15，6，1；】

第一行表示周围没有J偶合核，频率裂分数为0，原始峰的幅值为1；第二行表示周围有1个J偶合核，裂分数为2，裂分峰的幅值分别为1和1；第三行表示周围有2个J偶合核，裂分数为3，裂分峰的幅值分别为1，2，1；依次类推。该矩阵最多可以模拟周围有6个J偶合核的情况，基本涵盖了常用的分子情况；

第四步，对每一个自旋体系内的每一种自旋核(如前述之A，B，C等)，计算其基础幅值信息，考虑射频激励效应、样品T1弛豫效应、样品本征T2弛豫效应、样品实际T2*弛豫效应等因素，最终幅值为各项因素之积；T1弛豫效应为:1-exp(-TR/T1).

其中涉及到自旋核数目(如前述之m，n，d等)，TR，TE，T1，T2，B1，tau，deltaB0等前述构建的模型参数；

第五步，对每一个自旋体系内的每一种自旋核(如前述之A，B，C等)的频率裂分信息，具体为化学位移的频率效应，周围具有J偶合效应自旋核对该核的频率裂分效应Jab，即I0*A0*Jab，逐个按照周围的自旋核的数目对应的I0和A0矩阵内的元素进行点积再与Jab相乘；

其中涉及到I0，A0，J，W等模型和分子信息；

第六步，对4-5步的效应进行相乘后再考虑旋转坐标系下的数据采集过程，得到某一种自旋核的信号；

第七步，重复4-6步，对全部自旋核信号都进行重复计算后累加，得到一个自旋体系内的信号；

第八步，重复4-7步，得到另一个自旋体系的采集信号后进行叠加。同时增加进随机噪声；即可得到单脉冲序列的最终模拟采集数据；

第九步，如果有超过2个自旋体系的分子，则重复第八步。

具体地，如COSY二维谱的模拟数据采集模型：

第一步：首先构建t2维数据采样模型，即一个一维离散的时间数组，包含采样点数，空扫次数，延迟采样时间，采样带宽，回波时间，采样时间t等多个参数在内，如t＝(DE:1/SW:(DE+TD/SW)；为序列确定的采样规律，其中包括采样带宽SW，延迟时间DE，采样点数TD；

第二步，构建t1维(演化期)的数据采集模型，也是一个离散的时间数组；包涵t1步进值，t1步进数；一个具体实例为，t1＝【0，step_t1，num_t1*step_t1】；其中，step_t1为演化期步进值；num_t1为演化期步进数；一个具体的实例为step_t1＝1ms；num_t1＝128；其他数值都可以。

第三步，分别构建不同偶合核数目的谱线裂分的峰数数组I0，I0为一个二维数组(矩阵)，行数对应的是不同的J偶合数目，每行的数值表示谱峰裂分的数目和频率偏移相对量。一个具体的I0实例矩阵为7*7的二维矩阵：

I0＝【0，0，0，0，0，0，0；

-1/2，1/2，0，0，0，0，0；

-1，0，1，0，0，0，0，0；

-3/2，-1/2，1/2，3/2，0，0，0；

-2，-1，0，1，2，0，0；

-5/2，-3/2，-1/2，1/2，3/2，5/2，0；

-3，-2，-1，0，1，2，3；】

第一行表示周围没有J偶合核，频率裂分数为0，频率偏移值为0；第二行表示周围有1个J偶合核，裂分数为2，频率偏移相对量分别为-1/2和1/2；第三行表示周围有2个J偶合核，裂分数为3，频率偏移相对量分别为-1，0，1；依次类推。该矩阵最多可以模拟周围有6个J偶合核的情况，基本涵盖了常用的分子情况；

第四步，分别构建不同偶合核数目的谱线裂分的峰的幅值数组A0，A0为一个二维数组(矩阵)，行数对应的是不同的J偶合数目，每行的数值表示谱峰裂分的数目和幅值。一个具体的A0实例矩阵为7*7的二维矩阵：

I0＝【1，0，0，0，0，0，0；

1，1，0，0，0，0，0；

1，2，1，0，0，0，0，0；

1，3，3，1，0，0，0；

1，4，6，4，1，0，0；

1，5，10，10，5，1，0；

1，6，15，20，15，6，1；】

第一行表示周围没有J偶合核，频率裂分数为0，原始峰的幅值为1；第二行表示周围有1个J偶合核，裂分数为2，裂分峰的幅值分别为1和1；第三行表示周围有2个J偶合核，裂分数为3，裂分峰的幅值分别为1，2，1；依次类推。该矩阵最多可以模拟周围有6个J偶合核的情况，基本涵盖了常用的分子情况；

第五步，计算每个自旋体系内的每种自旋核自演化(对应着对角峰)模型以及对每另一种核的极化转移(对应着交叉峰)模型；

第六步，对每一个自旋体系内的每一种自旋核(如前述之A，B，C等)，计算其基础幅值信息，考虑射频激励效应、样品T1弛豫效应、样品本征T2弛豫效应、样品实际T2*弛豫效应等因素，最终幅值为各项因素之积；T1弛豫效应为:1-exp(-TR/T1).

其中涉及到自旋核数目(如前述之m，n，d等)，TR，TE，T1，T2，B1，tau，deltaB0，回波个数等前述构建的模型参数；

第七步，对每一个自旋体系内的每一种自旋核(如前述之A，B，C等)的频率裂分信息，具体为化学位移的频率效应，周围具有J偶合效应自旋核对该核的频率裂分效应J，即I0*A0*J，逐个按照周围的自旋核的数目对应的I0和A0矩阵内的元素进行点积再与J相乘；

其中涉及到I0，A0，J，W等模型和分子信息；

第八步，对六-七步的效应进行相乘后再考虑旋转坐标系下的数据采集过程，得到某一种自旋核的t2维采集信号；

第九步，将第五步得到的t1维自演化信息和其他种自旋核对该核的极化转移信息，进行复数表示后，与第八步的结果进行相乘，作为某一种自旋核的最终信号。

第十步，重复第六-九步，对全部自旋核信号都进行重复计算后累加，得到一个自旋体系内的信号；

第十一步，重复第六-十步，得到另一个自旋体系的采集信号后进行叠加。同时增加进随机噪声；即可得到单脉冲序列的最终模拟采集数据；

第十二步，如果有超过2个自旋体系的分子，则重复第十步。

步骤5，对前述步骤4中模拟采集的原始数据，进行FT变换，可得到频谱实部、虚部或模谱，即波谱输出。由于实部的底部基线更窄，分辨更好，故一般都输出实部谱图。本发明可以得出一维谱图/二维谱图或三维谱图。

进一步包括步骤6，实部谱图一般信息不够，需要进一步处理才能给出更多有用信息提示。比如包括进行相位校正、基线校正、寻峰、积分、测距、化学位移校正等谱图处理。谱图处理可采用本领域现有技术来实现。

进一步包括，硬件参数调节校正步骤。核磁共振波谱的有效结果依赖于高度均匀和稳定的主磁场，探头线圈的良好匹配调谐状态，以及射频角度的准确性。任何一个环节的不足都会对结果造成巨大的影响。因此实际操作实验中，都会分别对上述各环节进行硬件调节和校正。本发明同样可以采用数值模拟仿真技术实现调节校正功能，包括wobb，shim，lock和PL1W调节等模块及调节校正。

wobb模块：

第一步，构建探头线圈谐振电路模型；可以是串联谐振模型，也可以是并联谐振模型；本发明的一个实例为串联谐振模型；串联谐振电路模型里包括线圈等效电感L(6nH)，等效电阻r(0.2欧姆)，匹配电容cm(13pf可调)，调谐电容ct(120pf可调)；

第二步，在调谐频率范围(580MHz～620MHz之间)计算电路模型的功率反射系数；

第三步，显示反射系数与拟调谐工作频率之间的差距；

第四步，分别调节匹配电容cm和调谐电容ct，使得反射系数的峰值频率与工作频率一致(调谐成功)；且反射系数的值最小(匹配最佳)；这个步骤与实际谱仪的手动调谐匹配过程一致，需要花费一定时间调节；

第五步，也可以跳过第四步，通过自动调谐匹配理论计算模型，直接计算得到最佳调谐匹配所需的调谐电容值和匹配电容值，系统一下子就达到最佳调谐匹配状态；这个步骤与实际谱仪的自动调谐匹配过程完全一致，调节时间短；

本发明实施例中wobb实例效果，如图5所示，可设置目标调谐频率后，目标频率处会显示一条红色竖线；分别调节(粗调/细调)调谐电容和匹配电容，可观察到谐振垫偏移和反射曲线锐度调整。仔细缓慢调节2个电容，可以使得反射曲线尖峰点在红线上，同时曲线尖锐度最高，显示的分贝值最小(如-140dB等。)也可以点击自动调谐按钮，即可自动计算出最佳的调谐电容和匹配电容，并且反射曲线显示为最佳状态；

shim模块：如图6所示，

第一步，构建主磁场非均匀场模型；构建模型有两种，一种是按照磁场空间位置，构建不同阶次的非均匀项；但这种模型在实际匀场操作时需要进行逐点计算，时间很长，与实际仪器的操作过程不吻合。故本发明不采用这种建模方法。本发明采用的构建模型是，通过随机函数构建20个[-0.5，0.5]之间的随机数，同时按照实际仪器中的不均匀权重乘上权重因子，采用每次随机构建不均匀项的方法与实际设备工作较为吻合，因为不同样品放置入磁体后，会形成不同的磁场不均匀性。

第二步，分别构建20组匀场项的权重因子。比如Z项乘以200，表示最大偏差为+-100Hz；具体为：

[x,y,z,z2,z3,z4,z5,z6,xz,yz,xy,z2y2,xyz,x3,y3,xz2,yz2,x2y2z,yz3,xz3]＝[50,50,200,100,60,20,10,2,16,16,20,24,10,6,6,2,2,2,2,2]；等式左边数组表示20组非均匀匀项，分别对应着20组匀场线圈调节电流；等式右边为各非均匀项的非均匀权重。

第三步，调节20组匀场线圈的其中一组电流，对应给出其产生的磁场分量，逐次与原始随机差生的非均匀项进行对应累加；对剩余的非均匀性进行求模后，输出显示FID信号，同时输出FID信号的积分值。当调节电流产生的磁场与原设定的非均匀性值相反时，叠加后会使磁场变均匀，FID信号拖尾变长，积分面积变大。

第四步，依次反复调节20组线圈电流，每组线圈调节原则是使得积分面积变大；具体调节设置了粗调/细调/微调三个档次，可先粗调，再细调，最后微调。

第五步，第四步调节最佳后，对FID信号进行FT变换，测量并评价其半高宽，即完成主磁场均匀性调节和评价。

具体实例效果，如图6所示的shim界面示意图：匀场前，可单击Acquire按钮，采集FID信号，再FT后，观察频谱宽度较宽，表示均匀性较差；如果sample-in按钮显示不为绿色，则表示没有放置标样，FID信号为0；此时可单击sample-in，以置入样品。

单击匀场线圈z，z按钮变绿色，表示为正在调节线圈；同时界面下面的调节滑动条显示为z，实际调节值也显示为z线圈的实际调节电流。此时调节滑动条(可以粗调/细调/微调三档切换)，随着电流调节，FID信号拖尾会轻微变化，同时FID信号的积分面积数值会出现相应增减变化。数值增加表示磁场均匀在优化方向，则继续该方向电流调节；否则则表示磁场均匀性退化，则反向调节。调节到z电流无法使积分值增加为止，表示z已达到暂时最优。

再单击x，x按钮变绿色，表示调节x线圈电流有效。与上面调节类似，使FID积分面积最大；此后再依次调节其余18组线圈。再此回头调节x，y，z等；直至积分值最大后，再次FT观察频谱展宽，分辨率得以提高。

如果单击spin on/off使其变绿后，x，y方向的不均匀性值直接为0了。任何x，y调节都会退化磁场。

lock模块：

第一步，构建锁场样品信息，包括旋磁比gama，T1，T2等；

第二步，构建磁场时间漂移曲线；

第三步，构建锁场电流信息，包括电流调节精度；

第四步，设定锁场频率值，输出中心频率漂移；同时计算二者之间的频率偏差，并通过电流磁场比计算出所需锁场电流，进行频率追踪；

第五步，锁场状态下，按照0.5s每次不断重复步骤4，使得实际频率稳定在所需频率值上。

P1调节模块：

第一步，构建PL1W射频功率值调节范围(0-200W)以及射频脉宽P1调节范围(0-100us)；

第二步，构建射频功率到线圈电流到射频场B1之间的物理数学模型，建立射频功率和射频场强度B1之间的转换关系，本发明采用的最终关系是B1＝2.6*sqrt(PL1W)/1000(Gs)。系数2.6可以根据不同的线圈直径，线圈导线直径以及线圈材料差异而不同；

第三步，构建射频翻转角度计算数学模型，theta＝gama*B1*p1；gama为旋磁比；

第四步，构建信号采集数学模型，S＝sin(theta)*exp(-1i*w*t)*exp(-t/T2)；

第五步，分别调节PL1W或p1，观察信号幅值变化规律；依照相关理论确定不同PL1W下的90度，180度脉宽。

第六步，为更真实的模拟探头线圈中的射频场不均匀性，可以构建B1不均匀性衰减项，本实例为B1(p1)＝exp(-p1/0.0005)；表示B1的非均匀衰减指数时间为500us，该值可根据探头线圈具体情况更改。

第七步，设定自动步进p1，依次采集信号，FT后显示峰值，可得到章动谱图(如图8所示)；可一次性得到不同角度的p1值。也可以观察到B1的不均匀程度；

实例如图7所示为PL1W和P1调节效果：调节射频场功率PL1W滑动条，其功率值出现变化，右侧信号频谱幅值也出现变化，按照射频翻转角度与信号幅度之间的变换规律变化；

调节射频脉冲宽度P1滑动条，其脉宽值出现变化，右侧信号频谱幅值也出现变化，按照射频翻转角度与信号幅度之间的变换规律变化；

本发明实施例中章动谱效果，如图8所示：设定min-p1，max-p1以及步进值后，单击go按钮。右侧会分别出现不同p1值下的谱图峰值。峰值变化与射频翻转角度与信号幅度之间的变换规律变化规律一致，遵循正弦变化规律，即为章动谱；

整体调节过程与真机调节过程和参数设置一摸一样。

本发明实施例中相位校正效果，如图9所示：对于乙醇分子进行数据采集和FT后，其实部谱图出现色散形成分，不全为吸收线形。在主界面上单击相位校正按钮后，出现一条十字线，用途是确定待调节的pivot(0阶校正主峰)。在某一组峰处单击左键，即完成pivot设置；

此后，增/减调节零阶相位PL0，观察pivot的峰形变化，直至pivot峰形体现出标准的吸收线型；此后再增/减调节一阶相位PL1，观察pivot之外的峰形变化，直至所有峰形体现出标准的吸收线型。PL1条结过程中，pivot峰形保持不变；

pivot设置原则，一般习惯是设置最高峰附近。但实际上可以设置不同的峰作为pivot，都能实现最终的完美相位校正。

本系统上实现的调节过程和调节效果与真实波谱仪一摸一样。

对经过相位校正后的谱图，单击自动积分按钮，则出现谱图的积分谱效果，如图10所示的本发明实施例中谱图自动积分效果。

对经过相位校正后的谱图，单击自动寻峰按钮，则出现谱图的积自动寻峰和峰位标示效果，如图11所示本发明实施例中谱图自动寻峰效果。

实施例中，对乙醇(还可以是针对分子库里的任何其他分子)进行二维数据采集，并进行2DFT后得到乙醇分子的二维谱图显示效果，如图12所示为采集得到的二维谱图contour显示效果。实施例中应用本发明虚拟采集得到的二维谱图mesh显示效果，如图13所示。

图14为图12二维谱图的中间效果，只对t2维数据进行FT后的效果显示。图14为应用本发明虚拟采集得到的二维谱数据经一维FT后的数据(部分)，可看到不同谱峰的相位变化。图14只显示了25条一维FT后的数据(一维数组)，本实例总数据为100条。数据总条数为num_t1。

本发明的保护内容不局限于以上实施例。在不背离本发明构思的精神和范围下，本领域技术人员能够想到的变化和优点都被包括在本发明中，并且以所附权利要求书为保护范围。

Claims (5)

1.一种基于数值模拟技术的仿真核磁共振波谱分析方法，其特征在于，所述方法包括：

构建模拟样品分子信息模型及模拟样品分子信息库；所述样品分子信息按两个及两个以上自旋体系来确定，每个自旋体系内的核子种类为三种、四种、五种或五种以上；所述自旋体系包括：自旋体系AmBnCd以及XpYqZw，其中，ABC为相互之间存在偶合的三种自旋核，核子数目分别为m，n，d；XYZ为相互之间存在偶合的另三种自旋核，核子数目分别为p，q，w；

构建波谱仪硬件参数；所述波谱仪硬件参数包括主磁场B0、offset、主磁场不均匀性deltaB0、射频场B1、梯度场G、射频功率PL1W、P1；

硬件参数调节校正步骤；采用数值模拟仿真技术实现调节校正，所述调节校正方法包括自动调谐匹配wobb、匀场shim、锁场lock、P1调节；

构建波谱序列物理模型；

基于前述构建的所述样品分子信息模型及信息库、所述波谱仪硬件参数、所述波谱序列物理模型，构建核磁共振波谱数据采集物理数学模型，进行模拟数据采集算法和程序实现，完成原始数据的模拟采集；所述核磁共振波谱数据采集物理数学模型构建原始数据的步骤中，一维谱采集模式下，对每一个自旋体系内的每一种自旋核，分别计算每一个自旋体系内的所有自旋核J偶合因素导致的谱线裂分的幅值项和频率项、自身化学位移导致的频率项、旋转坐标系下的数据采集项、T2弛豫项、T1弛豫项，磁化矢量翻转项，将上述各项连乘，获得每一个自旋核信号，然后对每一种自旋核以及然后将每一个自旋核体系内的所有自旋核信号进行连加，得到累加的每一个自旋核体系信号；然后将所有自旋体系的自旋体系信号连加，再叠加上模拟随机噪声，即得到模拟采集的核磁共振波谱原始数据；所述构建核磁共振波谱数据采集物理数学模型的步骤中，二维谱采集模式下，在一维谱的基础上，将每一个自旋体系内的每一种自旋核在t1维演化期导致的对角锋的频率项作为该自旋核的乘积因子；将每一个自旋体系内的每一种自旋核在t1维演化期导致的交叉峰频率项作为J偶合的另一方自旋核的乘积因子；

将上一步骤所模拟采集的原始数据，经傅立叶变换得到波谱图。

2.如权利要求1所述的仿真核磁共振波谱分析方法，其特征在于，进一步地，将如所述权利要求1方法获得的波谱图进行谱图处理，所述谱图处理方法包括相位校正、基线校正、寻峰、积分、测距、化学位移校正。

3.如权利要求1所述的仿真核磁共振波谱分析方法，其特征在于，所述波谱序列包括单脉冲序列、去偶序列、同核化学位移相关二维谱COSY序列，同核化学位移相关overhause增强二维谱NOESY序列、异核化学位移相关谱1H-13C HSQC序列、异核化学位移相关谱1H-13CHMBC序列。

4.一种基于数值模拟技术的仿真核磁共振波谱分析系统，其特征在于，采用如权利要求1-3之任一项所述的基于数值模拟技术的仿真核磁共振波谱分析方法，所述系统包括：

模拟样品分子信息模型及分子信息库的构建模块，用于构建模拟样品分子信息模型及分子信息库；

波谱仪硬件参数的构建模块，用于构建波谱仪相关硬件参数；

波谱仪硬件参数调节校正模块，用以构建波谱仪硬件参数的模拟调节及校正模型；所述调节校正模块包括自动调谐匹配wobb模块、匀场shim模块、锁场lock模块、P1调节模块；

波谱序列物理模型的构建模块，用于构建不同波谱序列的物理模型；

核磁共振波谱数据采集物理数学模型的构建模块，用于构建核磁共振波谱数据采集物理数学模型，进行模拟数据采集算法和程序实现，实现原始数据的模拟采集；

波谱反演模块，用于将前述模拟采集的原始数据用于FT变换得到波谱图。

5.如权利要求4所述的仿真核磁共振波谱分析系统，其特征在于，进一步包括谱图处理模块，用于对获得的所述波谱图进行谱图处理，包括进行相位校正、基线校正、寻峰、积分、测距、化学位移校正等谱图处理。

Images