1、第 41 卷第 12 期大学物理Vol41 No122022 年 12 月COLLEGEPHYSICSDec 2022收稿日期:20220308;修回日期:20220618基金项目:国家自然科学基金(81971607);广州医科大学学科建设项目(024102206183)资助作者简介:傅洪波(1974)，男，湖南长沙人，广州医科大学生物医学工程学院数理系副教授，主要从事大学物理教学和数学建模和医学成像研究工作通信作者:谢国喜，Email:guoxixie 163com物理自然技术社会磁共振成像 K 空间信号的数理模型 以自旋回波序列为例傅洪波，陈小伟，谢国喜(广州医科大学 生物医学工程学院，广
2、东 广州511436)摘要:K 空间是磁共振成像信号空间编码的基础，准确理解 K 空间信号的表达式是掌握磁共振成像原理的关键和难点本研究以常见的自旋回波脉冲序列为例，结合梯度磁场的空间编码过程，详细地推导了 K 空间信号的表达式 在信号采集环节，我们利用一个简单的物理磁化模型，简化了法拉第电磁感应定律下的信号表达式 在数学处理上，我们从微元分析出发，逐步推导出 K 空间数据的二维积分表达式，详细地展示了推导过程中的相关细节和依据 本模型推导过程中有效地将一系列磁共振概念和理论联系起来，可以帮助对磁共振成像理论感兴趣的研究者深入系统地认知 K 空间的概念和相关理论关键词:磁共振成像;K 空间;自
3、旋回波中图分类号:4452文献标识码:A文章编号:1000-0712(2022)12-0039-04【DOI】1016854/jcnki1000-0712220120K 空间也称为傅里叶空间，是由带有空间编码信息的磁共振成像(Magnetic esonance Imaging，MI)信号填充而成的数据矩阵，每一组 MI 图像都有其相应的 K 空间数据矩阵，它们是傅里叶变换对的关系，对 K 空间数据进行反傅里叶变换，就能重建出 MI 图像 因此 K 空间是 MI 信号空间编码的基础，也是掌握磁共振成像理论的关键环节但是，K 空间信号表达式的数学物理推导过程鲜见于各类论著及论文中 前期调研发现，现
4、有国内外论著，包括大量教材，直接给出 K 空间的信号表达式1-3，再结合相应图示，以几何化的直观方式讲述 K 空间是什么和怎样形成，很少有论文(著)对K 空间信号表达式的来源进行讨论，而是更多是关注磁共振技术新发展、数值仿真和实际应用4-6 其中论文 7 虽然能针对成像数据与 K 空间之间的关系进行数理推导，但推导过程偏于简单，缺乏一些关键细节，比如针对法拉第感应定律的相关推导，仍采用磁化强度矢量，而非定律要求的磁感应强度矢量，对两者关系也未详细说明，导致整个推导过程中的数学物理过程难以理解 另外，我们在互联网上也了解到，比如，Albert Einstein College of Medici
5、ne 和Massachusetts Institute of Technology 开设相应 MI课程，涉及 K 空间的部分，也未曾给出一个逻辑完整的介绍，因此根据上述调研的结果，我们认为有必要对此问题展开研究，形成相应的教学支撑值得一提的是，Haack E Mark 和 Brown obertW 等人合作的 2 部著作8，9 从严密的物理模型和数学处理给出了磁共振三维信号一般表达的推导，但理解该推导过程难度甚高，既要掌握电动力学的相应物理模型又要熟悉矢量微分等数学处理，另外一些关键细节也被省略，没有和序列匹配起来，教与学两方面都殊为不易 但是深刻理解磁共振成像 K 空间信号表达式既是磁共振成
6、像理论学习的难点和关键，也是相关技术性研发的理论基石，特别是国内MI 这种关键医疗设备仍大量依赖国外产品的背景下，加强我国大学相应的教学，具有深远意义针对这一问题，本文在相关文献的基础上，以大学普通物理理论为基础，并从自旋回波序列这一最基本的成像序列出发，结合梯度磁场的空间编码过程、以实例的形式完成 K 空间信号表达式的数学物理推导，并给出推导过程中的相关细节和依据，试图给学习者提供一个既简单明了，又便于系统地掌握K 空间信号模型40大学物理第 41 卷1磁共振成像与 K 空间概述磁共振成像是利用核磁共振原理进行成像的技术，其成像过程简述如下:自旋原子核在静磁场 B0的作用下磁化产生宏观磁化强
7、度矢量(以下简称磁化矢量)，磁化矢量在射频场的作用下产生横向分量，横向分量在梯度磁场的作用下实现空间编码，对经过空间编码的磁共振信号进行采集、解调、滤波和放大等处理后，成为填充 K 空间的数据，对 K 空间数据进行逆傅里叶变换就可以重建出扫描图像 K空间信号 S(kx，ky)与重建图像(x，y)之间是傅里叶变换对的关系:S(kx，ky)=(x，y)ei2(kxx+kyy)dxdy(1)接下来，我们对式(1)给出相应模型和推导2微元模型的建立根据量子力学理论8，9，处于静磁场中，磁化的氢核产生纵向磁化矢量最大值为M0=14223kTB0，(0kT)(2)其中 B0为静磁场的磁场强度，0为质子的共
8、振频率，为磁旋比，k 为玻耳兹曼常量，T 为热力学温度，为约化普朗克常量，表示该体积内的自旋密度(质子密度)数为方便起见，假设成像层面垂直于静磁场，成像层面为 xy 平面，层厚为 z，选取该层面内一个体素微元 zdxdy 为研究对象，根据式(2)可构造微元:M0zdxdy=22B012kT0(x，y，z)zdxdy=D0(x，y)dxdy(3)其中，D=22B0/12kT，0(x，y，z)为体素微元的质子体密度，这样 0(x，y，z)zdxdy 表示该体素微元内的质子数，由于断层厚度 z 为常量，可以重新定义 0(x，y)为二维像素微元的质子面密度，这样乘以微元面积 dxdy 后仍等于对应质子
9、数如图 1 所示，在自旋回波序列 90脉冲激励下，纵向磁化矢量将翻转至 xy 平面，翻转到 xy 平面的横向磁化矢量 Mxy等于纵向磁化矢量，Mxy在静磁场下进动，并按 T*2进行指数衰减 在 180重聚脉冲作用下，衰减的横向磁化强度矢量在 TE时刻重聚，重聚信号最大值与初始横向磁化矢量之间是 T2指数衰减的关系，并随之在剩下的 TS/2 时间内按函数exp(t/T*2)衰减图 1自旋回波序列和空间编码梯度磁场为处理方便，我们仅讨论 TE时刻后采集的数据 结合式(3)和 Bloch 方程可知 Mxy在回波时间 TE后的信号大小变化的表达式为Mxy=K0(x，y)(1eT/T1)eTE/T2et
10、/T*2=C(x，y)et/T*2(4)式(4)引入变量 C(x，y)进行了简化，其中的0(x，y)、纵向弛豫时间 T1、横向弛豫时间 T2都是空间位置的变量 该表达式也可以说明 Mxy与 0(x，y)、T1和 T2有关，可以通过调节 T和 TE参数实现质子密度、T1和 T2的加权成像 特别指出的是，区别于文献 8，9，我们单独把 exp(t/T*2)作为一项而不是和 exp(TE/T2)合并，一方面可以使得最终K 空间信号表达式可以体现 3 种常见加权成像方式，同时可以简化后续数学处理当 90射频脉冲激发后，纵向磁化矢量翻转为横向磁化矢量 Mxy，受静磁场和梯度磁场作用，将以角速度 进动，其
11、大小与 Mxy所处位置磁场强度的大小成正比 设某时刻它与 x 轴正向的夹角为 t+，则断层微元内的横向磁化矢量 Mxy在 x 和 y 轴的投影大小分别为:dMx=C(x，y)et/T*2cos(t+)dxdy(5)dMy=C(x，y)et/T*2sin(t+)dxdy(6)式(5)、式(6)中的共有项 t、为:t=(B0+BGx)tx=B0+xBx()tx(7)=(B0+BGy)ty=B0+yBy()ty(8)其中 BGx、BGy分别为频率编码和相位编码的磁场梯度大小，tx为频率编码梯度磁场作用时间，ty为相位第 12 期傅洪波，等:磁共振成像 K 空间信号的数理模型 以自旋回波序列为例41编
12、码梯度磁场作用时间 根据图 1 可知，先进行相位编码再进行频率编码图 2xy 平面的横向磁化矢量和正交采集线圈示意图3磁共振信号的采集磁共振成像信号是根据法拉第电磁感应定律进行采集获得的 如图 2 所示，如果在频率编码梯度时间的同时，在 xy 轴分别放置一个有效接收面积为 S 的闭合线圈，旋转变化的横向磁化矢量 Mxy将使得穿过线圈的磁通量发生变换，变化的磁通量将在闭合线圈上产生感应电动势，由此检测到磁共振成像信号，采集得到的数据就是 K 空间数据受论著 8 启发，我们做出以下分析:磁共振成像样品置于静磁场 B0中，发生磁化后，产生纵向磁化矢量 M0，假设成像组织是由各向同性的顺磁物质构成的，
13、总的磁感应强度 B 为静磁场强度 B0和介质磁化电流产生的附加磁场 B的矢量和10:B=B0+B(9)根据磁化理论可知，真空磁导率为 0，在静磁场中的磁化矢量为 M0的介质，产生的附加磁场 B为B=0M0(10)因此在射频 90脉冲作用下，横向磁化矢量 Mxy在静磁场作用下产生的附加磁场最大值为B=0M0(11)该附加磁场在接收线圈处产生相应变化磁通量，即=BS(12)其中 S 为线圈的有效接收面积由于接收线圈中总的磁通量是所有微元 dxdy的附加磁场 B的累加贡献，结合定积分定义，根据式(5)、式(6)和式(12)，可得图 2 中 x 和 y 轴处的感应线圈中产生的磁通量分别为x(t)=S0
14、dMx=0SC(x，y)et/T*2cos(t+)dxdy(13)y(t)=S0dMy=0SC(x，y)et/T*2sin(t+)dxdy(14)其中 S 为线圈总的有效接收面积，为简单计，这里设两组线圈的 S 相等根据法拉第电磁感应定律，感应电动势为=ddt(15)将式(13)、式(14)对时间求导可得x(t)=0SC(x，y)et/T*2sin(t+)dxdy+1T*20SC(x，y)et/T*2cos(t+)dxdy(16)y(t)=0SC(x，y)et/T*2cos(t+)dxdy+1T*20SC(x，y)et/T*2sin(t+)dxdy(17)由于 频率相对于 1/T*2(也相当于
15、频率)有 4个数量级以上差别9，10，因此可以忽略式(16)和式(17)右边表达式第 2 项的影响 另外，根据图 2，exp(t/T*2)表达式中 t 的取值范围为 0，Ts/2，在核磁共振成像中 TS取值 1 ms 左右，TS/2T*2，因此exp(t/T*2)取值接近 1，可视为常数8，9 而且 是考虑了频率编码的梯度磁场叠加静磁场的结果，实际上静磁场的 0Gx，这是由于静磁场的磁感应强度是 T 级别，而梯度为 mT 级别，因此 往往用0替代8，9，最终式(16)、式(17)可进一步简化为x(t)=0S0C(x，y)sin(t+)dxdy(18)y(t)=0S0C(x，y)cos(t+)d
16、xdy(19)4K 空间的积分表达K 空间是复空间，K 空间数据是由 x、y 方向测得的信号强度分别作为复数信号的实部和虚部，即(t)=y+ix(20)根据欧拉公式得(t)=00SC(x，y)ei(t+)dxdy(21)成像断层的所有微元集成后的磁共振信号为(t)=00SC(x，y)ei(t+)dxdy(22)实际应用中在得到式(22)之前，还需借助信号“解调”技术消除上面式(18)、式(19)感应信号中类似载波的高频成分 0 解调是用频率接近或等于0的正弦或余弦函数与信号相乘 比如用参考信号sin(0t+0ty)乘以式(21)，借助三角恒等式，可得解调后信号为42大学物理第 41 卷sig=sin(0t+0ty)sin(t+)=12cos (0)t+0ty 12cos (+0)t+0ty=12cos(t+)12cos (+0)t+0ty(23)对于式(23)中对应高频信号部分，利用低通滤波去掉，只剩下 对应的低频部分，对于虚部通道也做类似处理，最终式(21)变为(t)=00C(x，y)ei(t+)dxdy(24)其中t=xt0Bxdt=xt0Gx(t)dt(25)=yt0Bydt=y
