CN114577100B - 磁场目标定位计算方法 - Google Patents

Abstract

本发明公开了一种磁场目标定位计算方法。根据本发明的计算方法，用于产生磁场的励磁线圈被分割为励磁线圈子块。每个励磁线圈子块作为电流元，从而计算出每个励磁线圈子块在空间中任意一点P的磁感应强度。各个励磁线圈子块在空间中任意一点P的磁感应强度叠加起来，得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系。基于所采集的目标定位装置处的磁感应强度信号以及励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，得到所述目标定位装置的位置和方向。本发明的磁场目标定位技术可以用于医疗手术中，特别是介入类手术中。而且，本发明的定位计算方法更具普适性，不用将励磁线圈横截面的形状限定为圆形。

Description

磁场目标定位计算方法

技术领域

本发明涉及电磁场，更具体涉及磁场目标定位计算方法。

背景技术

现代医疗技术中，可以通过导管、鞘管等耗材介入生物体，对生物体组织进行治疗。但是在术中，需要对导管、导丝、导引器(鞘管)、探针等目标物进行精准的定位和跟踪。在不同的生物体组织进行介入治疗时，对定位的精度需求不同，通常是精度越高越好，定位越准确越好。

由于导管等目标物通常通过血管、消化道等介入生物体内，本身尺寸会被设计得较小，如果再额外添加一定尺寸的定位设备，会使得尺寸上不满足介入体内的要求。另外，尽管也可以通过例如X光、磁共振成像等图像的方式来观察到目标物的位置，但这样的位置通常不能满足手术级别的定位精度要求。

因此，在医疗应用中，特别是介入类手术中，为了能够在不过多占用目标物尺寸的情况下，尽量保证定位的精度，可能会采用磁场目标定位技术。

在运用磁场目标定位技术时，涉及到根据毕奥-萨伐尔定律(Biot-Savart Law)列方程组求解目标物的位置和方向。一般来说，当用于产生磁场的励磁线圈的截面形状为圆形时，可以将其等效为磁偶极子。而对于目标物而言，用于磁感应强度采集的目标定位装置通常是传感器线圈，也可以等效为磁偶极子。由于励磁线圈和传感器线圈都被等效为磁偶极子，会使得求解过程得到简化。然而，上述等效的前提是励磁线圈的截面形状必须为圆形，如果是其他形状则不能适用这样的等效，也就不能使用简化求解方式来进行磁场目标定位计算。

因此，希望能够提供一种更普遍适用的磁场目标定位计算方法，能够突破励磁线圈截面形状是圆形的限制，从而能够更广泛地应用于磁场目标定位，且仍然保证定位的精度。

发明内容

本发明针对励磁线圈的各种横截面形状提出了一种更普遍适用的定位计算方法。根据本发明的计算方法，任何横截面形状的励磁线圈都可以被分割为更小的子块，用每一个子块作为电流元，然后再将所有子块在空间中的磁感应强度进行叠加，得到空间中任意一点的磁感应强度，与采集得到的目标传感器线圈上的磁感应强度相对照，就可以列方程求解得到目标定位装置的位置和方向，例如三维坐标、俯仰角和旋转角。

根据本发明的第一方面，提供一种用于磁场目标定位计算的方法。该方法可以包括：将用于产生磁场的励磁线圈分割为励磁线圈子块；将每个励磁线圈子块作为电流元，计算出每个励磁线圈子块在空间中任意一点P的磁感应强度；将各个励磁线圈子块在空间中任意一点P的磁感应强度进行叠加，得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系；以及基于所采集的目标定位装置处的磁感应强度信号以及励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，得到所述目标定位装置的位置和方向。

在根据本发明第一方面的方法中，优选地，根据作为电流元的每个励磁线圈子块的位置和摆放方向，以及所述励磁线圈子块的电流强度，基于毕奥-萨伐尔定律，得到励磁线圈子块在空间中任意一点P的磁感应强度。

在根据本发明第一方面的方法中，优选地，将励磁线圈沿轴向分成M段，得到M个子线圈片，将子线圈片等效为子线圈片的轮廓后，再将所述轮廓分段。采用毕奥-萨伐尔定律计算出轮廓各段在磁场中任意一点P的磁感应强度分量。将轮廓各段在P的磁感应强度分量进行叠加，得到所述轮廓在P的磁感应强度；将M个轮廓在P的磁感应强度在轴向上进行叠加，得到励磁线圈在P的磁感应强度，从而得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系。根据电磁感应定律列出目标定位装置的传感器线圈在P法向向量方向上的磁感应强度，其中，法向向量是指在传感器线圈截面的法向单位向量；根据励磁线圈在P的磁感应强度和传感器线圈在P法向向量上的磁感应强度相等的原理列出磁感应电动势方程，从而求解得到所述目标定位装置的位置和方向。

在根据本发明第一方面的方法中，优选地，所述轮廓可以是封闭轮廓。

在根据本发明第一方面的方法中，优选地，所述轮廓可以由圆弧段和直线段构成。

在根据本发明第一方面的方法中，优选地，所述轮廓的形状可以是圆角矩形。

在根据本发明第一方面的方法中，优选地，所述轮廓的形状可以是圆形。

在根据本发明第一方面的方法中，优选地，所述目标定位装置可以位于医疗介入到生物体内的医学装置上。

在根据本发明第一方面的方法中，优选地，所述目标定位装置的位置和方向可以包括所述目标定位装置的三维坐标、俯仰角和旋转角。

在根据本发明第一方面的方法中，优选地，所述磁场是采用固定在固定装置上的多个磁场生成器来产生的。每个磁场生成器设置在不同位置或不同摆放方向以产生相应的磁场。每个磁场生成器包括所述励磁线圈。

在根据本发明第一方面的方法中，优选地，所述多个磁场生成器是至少6个磁场生成器。

在根据本发明第一方面的方法中，优选地，所采集的目标定位装置处的磁感应强度信号包括每个磁场生成器作用于所述目标定位装置的相应磁感应强度信号分量。

在根据本发明第一方面的方法中，优选地，基于每个磁场生成器作用于所述目标定位装置的相应磁感应强度信号分量以及每个磁场生成器的励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，列出磁感应电动势方程组，从而求解得到所述目标定位装置的位置和方向。

在根据本发明第一方面的方法中，优选地，计算所述信号分量的模量值；得到模量值最大和/或最小的信号分量，并将该信号分量对应的方程去除；将剩下的方程组成优化超定方程组，求解所述目标定位装置的位置和方向。

在根据本发明第一方面的方法中，优选地，将信号分量均分为多个组(由于位置采用XYZ三维坐标来表述的原因，而且每三个线圈ⅩYZ组合成一组，因此分为四组，以每组三个模量之和判断取舍为最优)，并计算各组的信号模量之和；比较各组的信号模量之和，得到模量之和值最大和/或最小的信号分量组，并将该信号分量组对应的方程去除；将剩下的方程组成优化超定方程组，求解所述目标定位装置的位置和方向。

在根据本发明第一方面的方法中，优选地，根据Levenberg-Marquardt(LM)算法或其改进型算法，迭代求解所述目标定位装置的位置和方向。

根据本发明的第二方面，提供一种计算机可读介质，其上存储有可由处理器执行的指令，所述指令在被处理器执行时，使得处理器执行如本发明第一方面的磁场目标定位计算方法。

本发明利用电磁场定位原理，提供一种磁场目标定位计算方法。本发明的方法尤其可以用于医疗手术中，特别是介入类手术中，能够在不过多占用目标物尺寸的情况下，尽量保证定位的精度。

在磁场定位的应用中，可以将用于产生磁场的励磁线圈与目标定位装置的传感器线圈等效为磁偶极子。这一等效设定在励磁线圈的横截面形状为圆形时可以得到对目标定位的快速而近似的求解。然而，当励磁线圈的横截面形状不是圆形，而是其他形状时，就不能再适用磁偶极子的等效设定了。

本发明针对励磁线圈的各种横截面形状提出了一种更普遍适用的定位计算方法。根据本发明的计算方法，任何横截面形状的励磁线圈都可以被分割为更小的子块，用每一个子块作为电流元，然后再将所有子块在空间中的磁感应强度进行叠加，得到空间中任意一点的磁感应强度，与采集得到的目标传感器线圈上的磁感应强度相对照，就可以列方程求解得到目标定位装置的位置和方向。

此外，对于超定方程组，可以保留计算结果更精确的方程，去除不精确的方程，由此保留适当数量的方程式，进行更精确的定位计算。

类似地，还可以采用Levenberg-Marquardt(LM)算法或其改进型算法，对方程组进行迭代求解。

附图说明

通过以下详细的描述并结合附图将更充分地理解本发明，其中相似的元件以相似的方式编号，其中：

图1是磁场目标定位系统的示意图。

图2是磁场发生装置的示意图。

图3是磁场目标定位原理图。

图4是根据本发明实施例的磁场目标定位计算方法的流程图。

图5是根据本发明实施例的磁感应电动势方程的构建方法的流程图。

图6是根据本发明实施例的采用笛卡尔坐标系励磁线圈的坐标示意图。

图7是根据本发明实施例的励磁线圈截面为圆角矩形的示意图。

图8是根据本发明实施例的励磁线圈截面为圆角三角形的示意图。

图9是根据本发明实施例的励磁线圈截面为另一种圆角矩形的示意图。

图10是根据本发明实施例的目标位置和方向迭代求解流程图。

具体实施方式

下面通过实施例，并结合附图，对本发明的技术方案作进一步详细的说明，但本发明不限于下面的实施例。

下面首先对磁场目标定位的一般原理进行说明。

根据本发明的一个实施例，磁场目标定位计算方法的实现需要首先在空间中生成磁场，然后对目标定位装置的磁感应信号进行采集。

磁场的生成可以分为磁场生成控制和磁场生成两个方面。在本发明的优选实施例中，可以包括使用交变电流、准直流电流或永磁体旋转来驱动生成磁场。磁场可以采用固定在固定装置上的多个磁场生成器来产生的。每个磁场生成器设置在不同位置或不同摆放方向以产生相应的磁场。每个磁场生成器包括励磁线圈。在优选实施例中，磁场生成器可以是至少6个磁场生成器。磁场生成器的形状可根据应用场合做调整，常见的形式为圆柱形、方形或多边形。磁场生成器的相对摆放位置可根据目标物的定位区域范围做调整。磁场生成器的摆放角度可根据目标物采集信号的幅度做调整。

本领域技术人员应当理解，尽管在本发明中使用了“磁场生成器”、“励磁线圈”等术语来描述用于在空间中产生磁场的各个层级的装置，但是也可以使用其他类似的术语，例如磁生成单元、磁发生器、磁生成线圈、定位垫等，来表达相同或相似的含义。

目标定位装置处于励磁线圈所产生的磁场中，产生磁感应信号。根据本发明的实施例，目标定位装置是定位传感器线圈。将目标定位装置安装在需要被定位的目标物上或目标物内。这样，在确定了定位传感器线圈的位置和方向的同时，也就确定了目标物的位置和方向。在一个优选实施例中，目标定位装置位于医疗介入到生物体内的医学装置上。例如，目标物可以是导管，或更具体地，可以是导管上的一个或多个电极；目标定位装置也安装在导管上，靠近电极，用于为导管或电极进行定位。

目标定位装置的位置和方向包括所述目标定位装置的三维坐标、俯仰角和旋转角。更一般地，目标定位装置的方向还可包括滚转角，但在本发明的应用中，并不关注滚转角这个维度。

前文提到，由于磁场的存在，目标定位装置即定位传感器线圈上产生了磁感应信号。采集该磁感应信号以进行分析，从而能够对目标定位装置进行定位。

定位计算的工作就是本发明所要解决的问题。一般来说，定位计算的操作就是利用已经生成的磁场的分布情况，根据采集到的每个励磁线圈作用于目标定位装置的相应磁感应信号，基于毕奥-萨伐尔定律(Biot-Savart Law)建立方程组，求解所述目标定位装置的位置和方向。

本领域技术人员应该理解，定位计算的实现可以通过软件，即完全通过算法编程在通用计算机上完成所述的计算操作。定位计算的实现也可以通过硬件或固件，即通过在诸如现场可编程门阵列(FPGA)、专用集成电路(ASIC)或数字信号处理器(DSP)等专门硬件处理器中进行编程，来完成所述的定位计算操作。

下面通过附图来更具体解释说明根据本发明实施例的磁场目标定位计算方法。

图1是一种磁场目标定位系统的示意图。

如图1中所示，磁场目标定位系统100包括磁场发生装置101。磁场发生装置101包含多个磁场生成器组102A、102B、102C、102D，每个磁场生成器组包括一个或多个磁场生成器。例如，每个磁场生成器组包括3个磁场生成器，用于生成磁场。系统100还包括用于采集所生成的磁场中的调制信号的信号采集模块107、主要功能是调制信号以驱动磁场生成器生成磁场的磁场生成控制模块108、以及用于求解目标物位置和方向的定位计算模块109。如前所述，这里所述的位置可以用三维坐标来表示，而这里所述的方向可以用俯仰角和旋转角来表示。目标定位装置(也称传感器或探测器)103处于目标物上，具有检测磁场的功能，即可以在磁场中通过定位传感器线圈而产生磁感应信号。线缆104的一端连接到目标定位装置103，另一端连接到信号采集模块107。线缆105的一端连接到磁场发生装置101，另一端连接到磁场生成控制模块108。此外，系统100还可以包括显示器106，用于显示信号采集模块107所采集到的磁感应信号或者定位计算模块109所计算得到的定位信息。例如，如图1中所示，在显示器106上显示的是目标定位装置103(也就是目标物)的三维坐标、俯仰角度以及旋转角度的数值。

一般来说，为了建立方程组求解位置与角度，磁场发生装置内应布置至少6个磁场生成器。这里以12个磁场生成器为例描述。磁场生成器的外部形状可设计为圆柱状、方形状或其他各种形状。

此外，本领域技术人员应当理解，尽管在本文中将磁场生成器(乃至磁场发生装置)与磁场生成控制模块作为磁场目标定位系统的两个组成部分来进行描述，然而在很多情况下磁场生成器与磁场生成控制模块可以集成在一起。因此，本发明并不限制磁场生气器与磁场生成控制模块在物理上分离或集成，只是在功能上进行区分。换句话说，在将磁场生成器与磁场生成控制模块集成在一起的实施例中，二者的关系可以看作硬件及其驱动的关系，或者将二者的集成看作是一种固件。

图2是磁场发生装置的示意图。如图2所示，磁场发生装置201包括磁场生成器组202A、202B、202C、202D，每组包括3个彼此正交的磁场生成器。以磁场生成器组202A内部结构图204为例进行说明。磁场生成器组202A包括生成X方向磁场的磁场生成器205，生成Y方向磁场的磁场生成器206，生成Z方向磁场的磁场生成器207。

图3是磁场目标定位原理图。在图3所示的坐标系300中，磁场发生装置包括磁场生成器组301A、301B、301C、301D，每组包括3个磁场生成器。举例来说，磁场生成器302为磁场生成器组301C中的一个磁场生成器，其位置和摆放角度已知P(x_i,y_i,z_i,α_i,β_i)。目标定位装置(定位传感器线圈)303也处于坐标系300中。医疗领域中常见安装有目标定位装置的目标物包括导管、导丝、导引器(鞘管)、探针等，应用领域包括心脏介入治疗导航、肺支气管定位导航、肾动脉消融导航等。目标定位装置303的空间位置和摆放角度P(x,y,z,α,β)为待求解变量。

磁场定位计算方法

在将励磁线圈直接等效为磁偶极子的计算操作中，根据毕奥-萨伐尔定律(Biot-Savart Law)列方程时，由于磁场生成器与目标物相距远大于磁场生成器自身尺寸，故将二者看做磁偶极子，推导出了方程：

Vol_i＝γ*(B_(x,i)*cos(α)*cos(β)+B_(y,i)*cos(α)*sin(β)+B_(z,i)*sin(α))

其中，(x,y,z)为目标物的三维空间位置，(α,β)为传感器线圈的俯仰角度和旋转角度，γ为增益系数，Vol_i为磁感应信号分量，B_(x,i)是第i个磁场生成器在传感器线圈处产生的磁感应强度的x分量，B_(y,i)是第i个磁场生成器在传感器线圈处产生的磁感应强度的y分量，B_(z,i)是第i个磁场生成器在传感器线圈处产生的磁感应强度的z分量。

在该方法中，有一个前提条件是分别将磁场生成器(励磁线圈)和目标定位装置(定位传感器线圈)等效为磁偶极子才能直接应用毕奥-萨伐尔定律。也就是说，在以上的实施例中，磁场生成器的励磁线圈的横截面的形状使得磁场生成器的励磁线圈和目标定位装置可以被等效为磁偶极子，由此才能近似计算目标定位装置的位置和方向。

因此，在实际产品生产的过程中，通常需要将磁场生成器中的励磁线圈的横截面设置为圆形，使得励磁线圈最大程度的接近磁偶极子的结构特征。即，在这样的实施例中，磁场生成器的励磁线圈的横截面为圆形的形状。

这一条件限制了励磁线圈的结构。然而，实际应用过程中，需要将磁场生成器中的线圈的横截面设置为其他形状，例如图7中的圆角矩形。上述励磁线圈的横截面需要设置为圆形的这一结构特征限制了磁导航中励磁线圈结构的设计，不利于励磁线圈的安装和批量生产。

在接下来将要描述的一般情况下的实施例中，打破了将励磁线圈等效为磁偶极子，励磁线圈的横截面为圆形这一限制。即，以下的实施例不仅适用于磁场生成器的励磁线圈的横截面的形状使得所述磁场生成器的励磁线圈可以被等效为磁偶极子的情况，而且也适用于磁场生成器的励磁线圈的横截面的形状使得所述磁场生成器的励磁线圈不能被等效为磁偶极子的情况。在后一种情况中，更具体地说，磁场生成器的励磁线圈的横截面为除圆形以外的形状。

无论是哪种情况，本发明提出了一种磁场目标定位计算方法，使得不同形状横截面的励磁线圈都能进行磁定位计算，并实现目标物的精准定位。

图4是根据本发明实施例的磁场目标定位计算方法的流程图。如图4所示，磁场目标定位计算方法400包括以下步骤：

S410：将用于产生磁场的励磁线圈分割为励磁线圈子块；

S420：将每个励磁线圈子块作为电流元，计算出每个励磁线圈子块在空间中任意一点P的磁感应强度；

S430：将各个励磁线圈子块在空间中任意一点P的磁感应强度进行叠加，得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向(俯仰角度和旋转角度)的关系；

S440：基于所采集的目标定位装置处的磁感应强度信号和步骤S430中得到的励磁线圈在空间中产生的磁感应强度与空间位置和方向(俯仰角度和旋转角度)的关系，求解得到目标定位装置(传感器线圈)的空间位置坐标和方向(俯仰角度和旋转角度)。

步骤S420中，将每个励磁线圈子块作为电流元，计算出每个励磁线圈子块在空间中任意一点P的磁感应强度具体是指，根据作为电流元的每个励磁线圈子块的位置和摆放方向，以及所述励磁线圈子块的电流强度，基于毕奥-萨伐尔定律，得到励磁线圈子块在空间中任意一点P的磁感应强度。

本领域技术人员应当理解，上文以及下文中提到的磁感应强度指的是一种矢量，即磁感应强度矢量或向量信号不仅包括幅度大小，还包括方向。

与现有技术相比，此方法的优点在于，在采用毕奥-萨伐尔定律列在P传感器线圈的磁感应电动势方程时，将励磁线圈截面轮廓的形状作为必要因素考虑进去，而不是将励磁线圈和传感器线圈直接等效为磁偶极子。具体的方法中，将励磁线圈截面轮廓分为各微分段，分别计算各微分段在磁场中的磁感应强度分量，再通过积分进行累加，最后得到励磁线圈整体在空间中任意一点P的磁感应强度的计算公式，从而列出在P传感器线圈的磁感应电动势方程。

图5是根据本发明实施例的磁感应电动势方程的构建方法的流程图。

磁感应电动势方程的构建方法是本发明的重要步骤，磁感应电动势方程的构建方法500的流程图如图5所示，具体包括以下步骤：

S510：将励磁线圈沿轴向分成M段，得到M个子线圈片，将子线圈片等效为子线圈片的轮廓后，再将轮廓分段，该步骤是对图4的步骤S410的进一步扩展；

S520：采用毕奥-萨伐尔定律计算出轮廓各段在磁场中任意一点P的磁感应强度分量，该步骤是对图4的步骤S420的进一步扩展；

S530：将轮廓各段在P的磁感应强度分量进行叠加，得到轮廓在P的磁感应强度；将M个轮廓在P的磁感应强度在轴向上进行叠加，得到励磁线圈在P的磁感应强度，P的磁感应强度的表达式中包含了P的三维空间位置坐标与角度，从而得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，该步骤是对图4的步骤S430的进一步扩展；

S540：根据电磁感应定律列出传感器线圈在P法向向量方向上的磁感应强度，法向向量是指在传感器线圈截面的法向单位向量，并且法向向量由俯仰角度和旋转角度表征；由步骤S530中得到的励磁线圈在P的磁感应强度和本步骤中得到的传感器线圈在P法向向量方向上的磁感应强度相等的原理列出磁感应电动势方程，从而求解得到目标定位装置(传感器线圈)的位置和方向，该步骤是对图4的步骤S440的进一步扩展。

在本发明的优选实施例中，上面所述的轮廓为封闭轮廓。例如，轮廓可以由圆弧段和直线段构成。

在以下将要详细说明的优选实施例中，轮廓的形状是圆角矩形。但是，本领域技术人员应该理解，轮廓的形状也可以是圆形。另外，也可以是圆角三角形等其他的形状。

图6是根据本发明实施例的采用笛卡尔坐标系励磁线圈的坐标示意图。如图9所示，采用笛卡尔坐标系，励磁线圈中心位于坐标原点，轴向向量指向Z方向，截面向量指向X方向和Y方向，沿Z方向，将励磁线圈切片，把长度为H的励磁线圈等效为M个长度为H/M的薄线圈(子线圈片)，其中，励磁线圈的中心位置Z＝0，第i个薄线圈的中心位置为

图7是根据本发明实施例的励磁线圈截面为圆角矩形的示意图。图8是根据本发明实施例的励磁线圈截面为圆角三角形的示意图。图9是根据本发明实施例的励磁线圈截面为另一种圆角矩形的示意图。

在任一个薄线圈(中心位置为[0,0,Z])上，以截面为圆角矩形为例，说明磁感应电动势方程的构建方法。薄线圈等效为圆角矩形，圆角矩形如图7所示，直线段为K1、K2、K3和K4，弧线段为S1、S2、S3和S4，圆角矩形依次由S1、K1、S2、K2、S3、K3、S4和K4拼接为闭合形状，并且图形对称，弧线段为S1、S2、S3和S4可组合为一个圆，也即是说弧线段S1、S2、S3和S4分别为同一个圆的四分之一，将所在的圆均分为四等份。用于计算磁感应强度的薄线圈还可以是其他的形状，例如图8中的圆角三角形，图9中的另一种圆角矩形。图7、图8和图9的共同点在于，轮廓可以被切分成线段和弧线，并可以用积分求得该轮廓各段在空间中某一点的磁感应强度，从而叠加得到该轮廓在空间中某一点的磁感应强度。

下面以圆角矩形为例说明励磁线圈在空间中某一点的磁感应强度的计算方法以及磁场中传感器线圈的三维空间位置和角度的求解方法。

显然，圆角矩形由4个1/4圆弧和4个直线段组成，圆角矩形直线段边长分别为L、W，四角圆弧半径为R，(X，Y，Z)是圆角矩形上的坐标点，这八个线段上分别取任意点：

圆弧S1(圆心[L/2,W/2,Z]，Ψ＝[0,π/2])任一点坐标为：

M1＝[L/2+R*cos(Ψ),W/2+R*sin(Ψ),Z]；

直线段K1([L/2,(W/2+R),Z]到[-L/2,(W/2+R),Z])任一点坐标为：

M2＝[X,(W/2+R),Z]；

圆弧S2(圆心[-L/2,W/2,Z]，Ψ＝[π/2,π])任一点坐标为：

M3＝[-L/2+R*cos(Ψ),W/2+R*sin(Ψ),Z]；

直线段K2([-(L/2+R),W/2,Z]到[-(L/2+R),-W/2,Z])任一点坐标为：

M4＝[-(L/2+R),Y,Z]；

圆弧S3(圆心[-L/2,-W/2,Z]，Ψ＝[π,3π/2])任一点坐标为：

M5＝[-L/2+R*cos(Ψ),-W/2+R*sin(Ψ),Z]；

直线段K3([-L/2,-(W/2+R),Z]到[L/2,-(W/2+R),Z])任一点坐标为：

M6＝[X,-(W/2+R),Z]；

圆弧S4(圆心[L/2,-W/2,Z]，Ψ＝[3π/2,2π])任一点坐标为：

M7＝[L/2+R*cos(Ψ),-W/2+R*sin(Ψ),Z]；

直线段K4([(L/2+R),-W/2,Z]到[(L/2+R),W/2,Z])任一点坐标为：

M8＝[(L/2+R),Y,Z]；

在这些线段上截取任一电流元I(dl)，根据毕奥-萨伐尔(Biot-Savart)定律，该电流元在磁场中产生的磁感应强度为：

多线段产生的磁感应强度B，是由各线段的dB_n积分后叠加得到。

其中，dli是M1～M8的微分：

dl1＝diff(M1,Ψ)；

dl2＝diff(M2,X)；

dl3＝diff(M3,Ψ)；

dl4＝diff(M4,Y)；

dl5＝diff(M5,Ψ)；

dl6＝diff(M6,X)；

dl7＝diff(M7,Ψ)；

dl8＝diff(M8,Y)。

这里需要注意的是，diff是matlab中的微分函数。例如，diff(M1,Ψ)，是对M1以Ψ微分。即：

ai是M1～M8指向磁场空间中某一点P(x,y,z)的向量：

a1＝cp-M1；

a2＝cp-M2；

a3＝cp-M3；

a4＝cp-M4；

a5＝cp-M5；

a6＝cp-M6；

a7＝cp-M7；

a8＝cp-M8。

因为：

而所以：

上述|a_i|^-3难以得到积分解析式，需进行近似计算处理。由于(1+x)^m的泰勒展开式为：

只取第一项，(1+x)^m≈1+m·x。

由此，

对进行积分可得：

b1＝int(dl1×a1*a1^(-3),Ψ,0,π/2)；

b2＝int(dl2×a2*a2^(-3),X,L/2,-L/2)；

b3＝int(dl3×a3*a3^(-3),Ψ,π/2,π)；

b4＝int(dl4×a4*a4^(-3),Y,W/2,-W/2)；

b5＝int(dl5×a5*a5^(-3),Ψ,π,3π/2)；

b6＝int(dl6×a6*a6^(-3),X,-L/2,L/2)；

b7＝int(dl7×a7*a7^(-3),Ψ,3π/2,2π)；

b8＝int(dl8×a8*a8^(-3),Y,-W/2,W/2)。

这里需要注意的是，int是matlab中的积分函数。例如，int(dl1×a1*a1^(-3),Ψ,0,π/2)，是dl1×a1*a1^(-3)以Ψ在区间[0,π/2]积分，b1～b8是圆角矩形分成八段后每段对应的磁感应强度。写为一般的数学式为：

当i＝1,2,3,4,5,6,7,8时，分别记为：

圆角矩形的磁感应强度是各段磁感应强度的矢量积分，用公式可表达为：B＝b1+b2+b3+b4+b5+b6+b6+b8。更通用的表达式为：其中，B_j是第j个轮廓在P的磁感应强度，N是轮廓被分成的段数，b_i是第j个轮廓中第i段在磁场中任意一点P的磁感应强度分量在对应长度或角度范围内的积分。

特别地，当励磁线圈为螺线管线圈时，由于L＝0，W＝0，Z＝0，励磁线圈的截面轮廓为圆形，因此，圆形轮廓在P的磁感应强度的表达式在XYZ坐标系下表示为：

其中，Bx、By和Bz为所述轮廓的磁感应强度在X、Y、Z轴向上的分量，N是励磁线圈匝数，R为四角圆弧半径，μ是磁导率，(x，y，z)是P点的三维坐标。

求得励磁线圈每一个截面轮廓在P点的磁感应强度之后，将M个轮廓在P的磁感应强度在轴向上进行叠加，就得到了整个励磁线圈在P的磁感应强度，所述励磁线圈在P的磁感应强度的表达式为其中，B是所述励磁线圈在P点的磁感应强度，B_j是第j个轮廓在P点的磁感应强度，M是将所述励磁线圈沿轴向切分的段数。

励磁线圈上施加励磁电压U，由可得，励磁电流变化率式中，L’是励磁线圈电感量。设B’＝(Bx’,By’,Bz’)，则μ是磁导率，N是励磁线圈的线圈匝数，U是励磁线圈上施加励磁电压，R为四角圆弧半径，L’是励磁线圈电感量，B’是励磁线圈在P点的磁感应强度B经过坐标转换后的磁感应强度，所述坐标转换是指将以励磁线圈中心点为原点表示的P的磁感应强度B转换到与传感器线圈同一坐标系下的P的磁感应强度B’。传感器线圈所在空间的坐标系并不是以励磁线圈中心点为原点建立的，所以要进行转换，而使得传感器线圈的空间坐标和P的磁感应强度B位于同一坐标系下，转换后的P的磁感应强度表示为B’。

根据电磁感应定律，传感器的感应电动势式中，n是传感器线圈匝数，Φ是穿过传感器线圈的磁通量。而Φ＝B·S，式中，B是励磁线圈产生的磁场在传感器线圈(P)处的磁感应强度，S是传感器线圈截面积，S＝(π·r²)·vp'，式中，r是传感器线圈圆周半径，vp’(xv’,yv’,zv’)是传感器线圈截面的法向单位向量，法向单位向量vp’可以通过俯仰角度和旋转角度表征。

设传感器的磁感应电动势ε＝k·(B'·vp')。将实测的磁感应电动势和计算得到的系数k代入ε＝k·(B'·vp')中，由此可以联立方程组，由于B’中包含了P点的三维空间位置坐标，vp’通过俯仰角度和旋转角度表征，所以通过联立方程组，可以算出传感器线圈的坐标(三维坐标)和姿态(俯仰角度和旋转角度)。优选地，用LM算法解算传感器坐标和姿态。

在以上的方法中，所采集的目标定位装置处的磁感应强度信号包括每个磁场生成器(励磁线圈)作用于所述目标定位装置的相应磁感应强度信号分量。因此，在以上的方法中，在图4的步骤S440，基于所采集的目标定位装置处的磁感应强度信号以及励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，得到所述目标定位装置的位置和方向包括：基于每个磁场生成器作用于所述目标定位装置的相应磁感应强度信号分量以及每个磁场生成器的励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，列出磁感应电动势方程组，从而求解得到所述目标定位装置的位置和方向。

求解超定方程组问题实则是非线性模型求解问题，可以按照一定筛选准则挑选出其中部分(大于等于6个)或全部方程联立求解，常用求解方法是LM(Levenberg-Marquardt)算法或其改进型，本发明的优选实施例采用其改进型，可在3～8次迭代内获得收敛。

上述方程组是根据毕奥-萨伐尔定律，对计算公式进行泰勒展开后，取一次谐波分量而得到的近似计算公式，所以，求解超定方程得出的位置和方向是近似的值。为了提高计算结果的准确性，目标定位装置(目标物)需要放在距离磁场生成器一定的距离以内，求出的数据才比较准确。目标定位装置距离磁场生成器太近，定位传感器线圈通过磁场生成器的励磁信号感应到的响应信号就很强，该响应信号代入方程，计算出来的位置和方向误差较大，不准确；目标定位装置距离磁场生成器太远，定位传感器线圈通过磁场生成器的励磁信号感应到的响应信号就很弱，该响应信号代入方程，计算出来的位置和方向误差也较大，不准确。所以在实际计算过程中，需要剔除距离太近或距离太远的信号分量对应的方程，也可以将距离太近或距离太远的信号分量对应的方程都剔除。使得列出的方程组中各方程的响应信号分量在合理的范围，从而提高计算的准确性。具体步骤包括：

A301、将信号分量均分为多个组，例如N组，并计算各组的信号模量之和；

A302、找到模量之和值最大和/或最小的信号分量组，并将所述模量之和值最大和/或最小的信号分量组对应的方程删除，将剩下的方程组成优化超定方程组，参与最终求解。如果未知量为6个(除三维坐标、俯仰角度和旋转角度之外，还要一个未知数是增益系数)，剔除部分方程后，需要确保剩下的方程个数大于或等于6，通常方程个数为6～12个，作为优选方案，方程个数为6个、9个或12个。比较各组的信号模量之和，得到模量之和值最大和/或最小的信号分量组，并将该信号分量组对应的方程去除；将剩下的方程组成优化超定方程组，求解所述目标定位装置的位置和方向。

作为优选方案，给出了从12个方程中筛选出优化方程组合的方法，具体包括以下步骤：

A3001、分别统计目标定位装置103针对102A、102B、102C、102D磁场生成器组采集的信号模量之和，计算公式为：

其中，Vol_i是第i个磁场生成器生成磁场作用于目标物上产生的信号量，i是各信号分量的编号或索引，Vol^A、Vol^B、Vol^C和Vol^D分别是相邻的三个信号量的模量之和。

A3002、比较Vol^A、Vol^B、Vol^C和Vol^D的大小，从中筛选出值最大的模量之和，并找到值最大的模量之和对应的3个信号量，从12个方程组中剔除该3个信号量对应的方程，剩余的9个方程联立组成优化超定方程组，参与最终求解。

以上的流程可以被总结为：将信号分量均分为多个组，并计算各组的信号模量之和；比较各组的信号模量之和，得到模量之和值最大和/或最小的信号分量组，并将该信号分量组对应的方程去除；将剩下的方程组成优化超定方程组，求解所述目标定位装置的位置和方向。

图10是根据本发明实施例的目标位置和方向迭代求解流程图。如图10所示，目标位置和方向求解流程1000开始于步骤1001，在此步骤，输入各个磁场生成器作用于目标定位装置上所产生的信号分量。接下来，在步骤1002，判别当前输入条件下，是否连续多次未解出结果。若次数超过预设值(如3次)，即步骤1002的“是”分支，则在步骤1004停止本轮求解，输出求解失败。若次数为超过预设值，即步骤1002的“否”分支，则进入步骤1003，判断当前输入对应的目标物是否是第一次求解。如果是第一次求解，即步骤1003的“是”分支，则进入步骤1005，随机生成初始位值作为初始迭代值；如果不是第一次求解，即步骤1003的“否”分支，则进入步骤1006，利用前次求解结果作为迭代求解的初始值。确定初始值之后，在步骤1007，迭代求解目标物坐标和方向，常用方法是LM(Levenberg-Marquardt)算法或其改进型。在步骤1008判断迭代是否收敛。若收敛，即步骤1008的“是”分支，则在步骤1009标记成功求解，然后在步骤1010输出求解结果；若收敛失败，即步骤1008的“否”分支，则返回步骤1002判断失败次数是否超过预设。

也就是说，可以根据Levenberg-Marquardt(LM)算法或其改进型算法，迭代求解所述目标定位装置的位置和方向。

直接磁偶极子等效

另外，下面附上将磁场生成器或励磁线圈直接等效为磁偶极子的算法。该算法可以用于励磁线圈截面为圆形形状的情况。

因磁场生成器与目标物相距远大于磁场生成器自身尺寸，可将二者直接看做磁偶极子。根据毕奥-萨伐尔定律(Biot-Savart Law)，定位原理详述如下：

根据磁场生成器302的位置和摆放角度，可得其归一化的磁场生成器方向向量

Dir_(x,i)＝cos(α_i)*cos(β_i)

Dir_(y,i)＝cos(α_i)*sin(β_i)

Dir_(z,i)＝sin(α_i)

其中(x_i,y_i,z_i)为三维空间位置，(α_i,β_i)为磁场生成器的俯仰角度(极角)和旋转角度(方位角)，其中i表示磁场生成器的编号或索引，例如当存在N个磁场生成器时，i＝1,2,…,N，N≥6。

目标物到磁场生成器距离：

第i个磁场生成器生成磁场作用于目标物上产生的信号量Vol_i，对应于所采集到的磁感应强度信号：

Vol_i＝γ*(B_(x,i)*cos(α)*cos(β)+B_(y,i)*cos(α)*sin(β)+B_(z,i)*sin(α))

其中，(x,y,z)为目标物的三维空间位置，(α,β)为定位传感器线圈的俯仰角度(极角)和旋转角度(方位角)，γ为增益系数，P(x,y,z,α,β,γ)为6个待解未知量。以12个磁场生成器为例，可以获得包含6个未知量的12个方程，联立组成超定方程组。

上述磁场目标定位方法适用于所有需要进行目标定位的领域与应用场景，例如，不仅可以用于医疗应用场景，还可以用于佩戴VR眼镜、AR头盔后确定头部的位置和方向等应用场景。

计算机程序或计算机程序产品以及计算机可读介质

此外，本领域普通技术人员应该认识到，本公开的方法可以实现为计算机程序。如上结合附图所述，通过一个或多个程序执行上述实施例的方法，包括指令来使得计算机或处理器执行结合附图所述的算法。这些程序可以使用各种类型的非瞬时计算机可读介质存储并提供给计算机或处理器。非瞬时计算机可读介质包括各种类型的有形存贮介质。非瞬时计算机可读介质的示例包括磁性记录介质(诸如软盘、磁带和硬盘驱动器)、磁光记录介质(诸如磁光盘)、CD-ROM(紧凑盘只读存储器)、CD-R、CD-R/W以及半导体存储器(诸如ROM、PROM(可编程ROM)、EPROM(可擦写PROM)、闪存ROM和RAM(随机存取存储器))。进一步，这些程序可以通过使用各种类型的瞬时计算机可读介质而提供给计算机。瞬时计算机可读介质的示例包括电信号、光信号和电磁波。瞬时计算机可读介质可以用于通过诸如电线和光纤的有线通信路径或无线通信路径提供程序给计算机。

例如，根据本公开的一个实施例，可以提供一种计算机可读介质，其上存储有可由处理器执行的指令，所述指令在被处理器执行时，使得处理器执行如前所述的磁场目标定位计算方法。

因此，根据本发明公开的内容，还可以提议一种计算机程序或计算机程序产品，当所述计算机程序被执行时，可实现如前所述的磁场目标定位计算方法。

另外，本发明还涉及一种计算装置或计算系统，包括处理器和存储器，所述存储器中存储有计算机程序，当所述计算机程序由所述处理器执行时，可实现如前所述的磁场目标定位计算方法。

本发明的有益效果

综上所述，除了之前已经描述过的效果，本发明的有益效果还可以被总结如下：

1、本发明公开了一种磁场目标位置计算方法，在构建一种磁场目标位置跟踪定位系统后，不再将磁场生成器和定位传感器等效为磁偶极子，而是将磁场生成器的励磁线圈分成励磁线圈子块，将每个励磁线圈子块作为电流元，计算出每个励磁线圈子块在空间中任意一点P的磁感应强度，然后将励磁线圈子块在空间中任意一点P的磁感应强度进行叠加，建立起整个励磁线圈在空间中产生的磁感应强度与空间位置、俯仰角度和旋转角度的关系。这样的好处在于，在磁定位系统中，励磁线圈的横截面的形状有了更多的扩展，不是非得是圆形，还可以是圆角矩形，还可以是圆角三角形以及线段和弧线的组合，给励磁线圈的制作和安装提出了更多的可能性。

2、即使励磁线圈之间距离很近，无法等效为磁偶极子，采用本发明的方法仍然能精确计算出励磁线圈作用于传感器线圈上产生的信号分量，从而实现目标物(传感器线圈)的精确定位。

3、在根据采用多个信号分量联立组成超定方程组后，为了提高计算效率，还提出了从多个方程中筛选出优化方程组合的方法，以减少方程组中方程的个数，提高计算效率。

本发明的实施方式并不限于上述实施例所述，在不偏离本发明的精神和范围的情况下，本领域普通技术人员可以在形式和细节上对本发明做出各种改变和改进，而这些均被认为落入了本发明的保护范围。

Claims (14)

1.一种用于磁场目标定位计算的方法，包括：

将用于产生磁场的励磁线圈分割为励磁线圈子块；所述励磁线圈子块为励磁线圈截面切割后形成的直线段或弧线；

将每个励磁线圈子块作为电流元，计算出每个励磁线圈子块在空间中任意一点P的磁感应强度；

将各个励磁线圈子块在空间中任意一点P的磁感应强度进行叠加，得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系；以及

基于所采集的目标定位装置处的磁感应强度信号以及励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，得到所述目标定位装置的位置和方向，包括：基于每个磁场生成器作用于所述目标定位装置的相应磁感应强度信号分量以及每个磁场生成器的励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，列出磁感应电动势方程组，计算所述信号分量的模量值；

得到模量值最大和/或最小的信号分量，并将该信号分量对应的方程去除；

将剩下的方程组成优化超定方程组，求解所述目标定位装置的位置和方向；

所述将各个励磁线圈子块在空间中任意一点P的磁感应强度进行叠加，得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系进一步包括：将所述直线段或弧线在P的磁感应强度分量进行叠加，得到所述直线段或弧线构成的励磁线圈截面在P的磁感应强度；将若干励磁线圈截面在P的磁感应强度在轴向上进行叠加，得到励磁线圈在P的磁感应强度，从而得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系；

所述磁场是采用固定在固定装置上的多个磁场生成器来产生的，

每个磁场生成器设置在不同位置或不同摆放方向以产生相应的磁场，

每个磁场生成器包括所述励磁线圈。

2.根据权利要求1所述的方法，其中，所述的将每个励磁线圈子块作为电流元，计算出每个励磁线圈子块在空间中任意一点P的磁感应强度，进一步包括：

根据作为电流元的每个励磁线圈子块的位置和摆放方向，以及所述励磁线圈子块的电流强度，基于毕奥-萨伐尔定律，得到励磁线圈子块在空间中任意一点P的磁感应强度。

3.根据权利要求1所述的方法，其中：

所述的将用于产生磁场的励磁线圈分割为励磁线圈子块，进一步包括：将励磁线圈沿轴向分成M段，得到M个子线圈片，将子线圈片等效为子线圈片的轮廓后，再将所述轮廓分段，

所述的将每个励磁线圈子块作为电流元，计算出每个励磁线圈子块在空间中任意一点P的磁感应强度，进一步包括：采用毕奥-萨伐尔定律计算出轮廓各段在磁场中任意一点P的磁感应强度分量，

所述的将各个励磁线圈子块在空间中任意一点P的磁感应强度进行叠加，得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，进一步包括：将轮廓各段在P的磁感应强度分量进行叠加，得到所述轮廓在P的磁感应强度；将M个轮廓在P的磁感应强度在轴向上进行叠加，得到励磁线圈在P的磁感应强度，从而得到整个励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，

所述的基于所采集的目标定位装置处的磁感应强度信号以及励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，得到所述目标定位装置的位置和方向，进一步包括：根据电磁感应定律列出目标定位装置的传感器线圈在P法向向量方向上的磁感应强度，其中，法向向量是指在传感器线圈截面的法向单位向量；根据励磁线圈在P的磁感应强度和传感器线圈在P法向向量上的磁感应强度相等的原理列出磁感应电动势方程，从而求解得到所述目标定位装置的位置和方向。

4.根据权利要求3所述的方法，其中，所述轮廓为封闭轮廓。

5.根据权利要求4所述的方法，其中，所述轮廓由圆弧段和直线段构成。

6.根据权利要求5所述的方法，其中，所述轮廓的形状为圆角矩形。

7.根据权利要求4所述的方法，其中，所述轮廓的形状为圆形。

8.根据权利要求1所述的方法，其中，所述目标定位装置位于医疗介入到生物体内的医学装置上。

9.根据权利要求1所述的方法，其中，所述目标定位装置的位置和方向包括所述目标定位装置的三维坐标、俯仰角和旋转角。

10.根据权利要求1所述的方法，其中，所述多个磁场生成器是至少6个磁场生成器。

11.根据权利要求1所述的方法，其中，所采集的目标定位装置处的磁感应强度信号包括每个磁场生成器作用于所述目标定位装置的相应磁感应强度信号分量。

12.根据权利要求1所述的方法，其中，所述的基于每个磁场生成器作用于所述目标定位装置的相应磁感应强度信号分量以及每个磁场生成器的励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，列出磁感应电动势方程组，从而求解得到所述目标定位装置的位置和方向，进一步包括：

将信号分量均分为多个组，并计算各组的信号模量之和；

比较各组的信号模量之和，得到模量之和值最大和/或最小的信号分量组，并将该信号分量组对应的方程去除；

将剩下的方程组成优化超定方程组，求解所述目标定位装置的位置和方向。

13.根据权利要求1所述的方法，其中，所述的基于每个磁场生成器作用于所述目标定位装置的相应磁感应强度信号分量以及每个磁场生成器的励磁线圈在空间中产生的磁感应强度与空间位置和方向的关系，列出磁感应电动势方程组，从而求解得到所述目标定位装置的位置和方向，进一步包括：

根据Levenberg-Marquardt（LM）算法或其改进型算法，迭代求解所述目标定位装置的位置和方向。

14.一种计算机可读介质，其上存储有可由处理器执行的指令，所述指令在被处理器执行时，使得处理器执行如权利要求1 所述的磁场目标定位计算方法。