室内定位的方法、系统、电子设备和存储介质

1.本技术属于室内定位领域，尤其涉及室内定位的方法、系统、电子设备和存储介质。


背景技术：

2.当前，基于视觉和激光雷达传感器的高精度跟踪是室内定位研究的主流方向。研究者们通过利用摄像头、激光雷达等传感器融合imu实现了高精度、高保真度的定位以及室内环境重建。然而视觉传感器受光照、环境纹理等因素影响严重，激光雷达价格昂贵并且受天气、空气湿度等影响严重，导致现有的同步定位与地图构建(simultaneous localization and mapping，slam)算法无法广泛应用于复杂室内环境中。
3.针对激光雷达和视觉传感器的局限性，研究人员在当前同步定位与地图构建slam框架中加入一些低成本且不受上述环境条件约束的传感器。射频(radio-frequency，rf)传感器对动态照明条件和结构化的室内环境的鲁棒性，即使在黑暗，弱纹理和烟雾环境下也能正常工作。其中wifi信号相比于其他无线信号在部署区域和覆盖范围上有巨大优势。目前的一些工作提出了wifi-slam系统。使用wifi信号的信号强度信息(received signal strength，rssi)进行slam。但rssi测量在大多数室内场景中并不可靠的，因为室内环境的多径衰落会严重干扰无线信号强度。随着无线技术的发展，使得商用网卡能够获取wifi正交频分复用信号的信道状态信息(channel sate information，csi)。与信号强度不同，csi提供了信号的相位信息和多个子载波的信道响应，从而大大降低了定位模糊度并提高了室内的多径分辨率。一些工作利用csi估计移动端的位置、速度以及ap(access point，ap)位置的系统。研究了通过wifi信号赋能自主导航的可行性。但这些系统只能工作在无多径衰落的场景，如空旷的室外环境。对于室内等遭受严重多径衰落的环境下并不适用。
4.综上所述，现有基于视觉和激光雷达的slam系统受到光照、环境纹理，算力、载重以及成本等方面的限制，无法广泛应用。无线信号具有丰富的环境感知能力，有着巨大的潜力为手持移动端提供更鲁棒的状态估计和环境感知。然而，目前基于rssi进行的slam算法由于rssi在室内场景下的随机性，导致效果并不理想。而利用csi的slam系统中，缺乏有效特征参数的提取，在室内环境下鲁棒性较差。


技术实现要素：

5.本发明实施例的主要目的在于提供室内定位的方法、系统、电子设备和存储介质，使得多径aoa可以克服传统wifi定位算法中需要获得ap位置先验信息的问题，而且通过多径aoa也解决了在nlos或多径干扰严重场景下无法进行定位的问题。
6.第一方面，提供了室内定位的方法，所述方法包括：
7.通过待定位设备获取室内的信号状态信息csi和传感信息；
8.根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实时粗略位置信息；
9.通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。
10.在一个可能的实现方式中，所述根据所述csi获取多径信道到达角aoa，包括：
11.将所述待定位设备的行动轨迹虚拟为大孔径天线，根据合成孔径雷达sar的原理通过所述大孔径天线获取所述多径aoa。
12.在另一个可能的实现方式中，所述将所述待定位设备的行动轨迹虚拟为大孔径天线，根据合成孔径雷达sar的原理通过所述大孔径天线获取所述多径aoa，包括：
13.根据预设的信号接收公式获取所述待定位设备接收到的信号，记录预设时间段内《h
ij
(t),pi(t)》的数据对，将所述数据对设置为多天线阵列，其中，h
ij
(t)为t时刻待接收设备接收的信号，pi(t)为t时刻待接收设备的位置，d
ij
(t)为t时刻待接收设备与信号发射端的距离；
14.根据预设的运动轨迹公式p
ij
(t)＝pi(t)-pj(t)获取所述待定位设备的相对运动轨迹；
15.根据预设的转向向量公式和所述相对运动轨迹获取转向向量其中，(ρk,ξk,ψk),(ρ
l
,ξ
l
,ψ
l
)为k、l时刻待定位设备在世界坐标系下的坐标；
16.根据预设的能量谱公式获取能量谱，进而获取所述f
ij
(θ,φ)的峰值对应的(θ,φ)，并将所述峰值对应的(θ,φ)设置为多径aoa。
17.在另一个可能的实现方式中，所述通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化，包括：
18.根据预设的融合优化算法将所述实时粗略位置信息、历史位置信息、真实ap和虚拟ap的估计位置信息融合优化，其中，ρ(
·
,c)为huber代价函数，c为经验阈值。
19.第二方面，提供了一种室内定位的系统，所述系统包括：
20.信息获取模块，用于通过待定位设备获取室内的信号状态信息csi和传感信息；
21.初略位置信息获取模块，用于根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实时粗略位置信息；
22.估计及优化模块，用于通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。
23.在一个可能的实现方式中，所述根据所述csi获取多径信道到达角aoa，包括：
24.将所述待定位设备的行动轨迹虚拟为大孔径天线，根据合成孔径雷达sar的原理通过所述大孔径天线获取所述多径aoa。
25.在另一个可能的实现方式中，所述初略位置信息获取模块，包括：
26.多天线阵列设置单元，用于根据预设的信号接收公式获取所述待定位设备接收到的信号，记录预设时间段内《h
ij
(t),pi(t)》的数据对，将所述数据对设置为多天线阵列，其中，h
ij
(t)为t时刻待接收设备接收的信号，pi(t)t时刻待接收设备的位置，d
ij
(t)为t时刻待接收设备与信号发射端的距离；
27.相对运动轨迹获取单元，用于根据预设的运动轨迹公式获取所述待定位设备的相对运动轨迹；
28.转向向量获取单元，用于根据预设的转向向量公式和所述相对运动轨迹获取转向向量其中，(ρk,ξk,ψk),(ρ
l
,ξ
l
,ψ
l
)为k、l时刻待定位设备在世界坐标系下的坐标；
29.多径aoa设置单元，用于根据预设的能量谱公式获取能量谱，进而获取所述f
ij
(θ,φ)的峰值对应的(θ,φ)，并将所述峰值对应的(θ,φ)设置为多径aoa。
30.在另一个可能的实现方式中，通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化，包括：
31.根据预设的融合优化算法将所述实时粗略位置信息、历史位置信息、真实ap和虚拟ap的估计位置信息融合优化，其中，ρ(
·
,c)为huber代价函数，c为经验阈值。
32.第三方面，提供了一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，处理器执行程序时实现如第一方面提供的室内定位的方法。
33.第四方面，提供了一种非暂态计算机可读存储介质，其上存储有计算机程序，计算机程序被处理器执行时实现如第一方面提供的室内定位的方法。
附图说明
34.为了更清楚地说明本技术实施例中的技术方案，下面将对本技术实施例描述中所需要使用的附图作简单地介绍。
35.图1为本发明一个实施例提供的室内定位的方法的流程图；
36.图2为本发明再一个实施例提供的室内定位的方法的流程图；
37.图3为本发明一个实施例提供的室内定位的系统的结构图；
38.图4为本发明再一个实施例提供的室内定位的系统的结构图；
39.图5为本发明一种电子设备的实体结构示意图。
40.具体实现方式
41.下面详细描述本技术的实施例，所述实施例的示例在附图中示出，其中自始至终相同或类似的标号表示相同或类似的模块或具有相同或类似功能的模块。下面通过参考附图描述的实施例是示例性的，仅用于解释本技术，而不能解释为对本发明的限制。
42.本技术领域技术人员可以理解，除非特意声明，这里使用的单数形式“一”、“一个”、“所述”和“该”也可包括复数形式。应该进一步理解的是，本技术的说明书中使用的措辞“包括”是指存在所述特征、整数、步骤、操作、模块和/或组件，但是并不排除存在或添加一个或多个其他特征、整数、步骤、操作、模块、组件和/或它们的组。应该理解，当我们称模块被“连接”或“耦接”到另一模块时，它可以直接连接或耦接到其他模块，或者也可以存在中间模块。此外，这里使用的“连接”或“耦接”可以包括无线连接或无线耦接。这里使用的措辞“和/或”包括一个或更多个相关联的列出项的全部或任一模块和全部组合。
43.为使本技术的目的、技术方案和优点更加清楚，下面将结合附图对本技术实现方式作进一步地详细描述。
44.下面以具体地实施例对本技术的技术方案以及本技术的技术方案如和解决上述技术问题进行详细说明。下面这几个具体的实施例可以相互结合，对于相同或相似的概念或过程可能在某些实施例中不再赘述。下面将结合附图，对本技术的实施例进行描述。
45.如图1所示为本发明一个实施例提供的室内定位的方法的流程图，所述方法包括：
46.步骤101，通过待定位设备获取室内的信号状态信息csi和传感信息；
47.步骤102，根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实时粗略位置信息；
48.步骤103，通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。
49.在本发明实施例中，待定位设备为具有wifi信号接收模块以及惯性测量单元的设备，通常为手机、平板电脑等移动终端。待定位设备通过wifi信号接收模块获取室内信道的csi(channel sate information，信道状态信息)，通过惯性测量单元获取包含待定位设备的加速度、待定位设备的陀螺仪信息等传感信息。获取的csi通过计算可以获取多径aoa(angle of arrival，信道到达角)，而传感信息可以用于获取待定位设备在室内的实时大致位置信息，即实时粗略位置信息，对于实时粗略位置信息由于其只是通过惯性测量单元获取的大致信息，准确度不够，而且随着时间的增加，惯性信息获取的位置信息会越来越不准确，因此，在实时粗略位置信息的基础上，引入多径aoa和融合优化算法，通过融合优化算
法对多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息待定位设备的当前位置信息进行估计，并对待优化参数进行实时优化，获取位置更为准确的优化后的位置信息。
50.其中，所述根据所述csi获取多径信道到达角aoa，包括：
51.将所述待定位设备的行动轨迹虚拟为大孔径天线，根据合成孔径雷达sar的原理通过所述大孔径天线获取所述多径aoa。
52.在本发明实施例中，由于待定位设备通常只拥有1-2根天线，该天线数量并不足以通过多天线阵列的方式估计出多径aoa，因此，需要获取待定位设备的行动轨迹，将行动轨迹虚拟为大孔径天线，根据sar的原理通过该大孔径雷达获取多径aoa信息。
53.惯性导航系统具有抗干扰性强，不需要借助外部设备等优点，目前大多智能手机都配备着惯性测量单元，imu由三个正交陀螺仪和三个正交加速度计组成，分别用于测量系统的加速度和比力。通过积分可以从加速度和角速度中得到位置、速度和姿态(也称为旋转)三个量，分别用p、v、q表示。因此ins是获取载体状态信息的有效方案。此外，imu测量信息频率相对于大多数传感器更高。
54.在惯性导航中，对于姿态表述常用到的方法有四元数法、欧拉角法等。欧拉角使用俯仰角、横滚角、航向角来描述姿态发生的变化，其优点在于更符合人们的思维习惯，更加直观的描述了角度的变化。然而其有一个重大缺陷，在俯仰角绝对值为90
°
时，翻滚角与俯仰角无法进行区分，出现“锁死”现象。四元数可用于描述载体的旋转，它将三维空间中的旋转用一个单位长度的四维向量来表示，一个四维向量包含一个实部和三个虚部，其优势在于能够无歧义的描述姿态。本发明采用四元素法来对姿态进行描述：
55.q＝q0+q1i+q2j+q3k
56.其中i,j,k为三个虚部。i2＝-1,j2＝-1,k2＝-1，且满足此外，需要对载体坐标系和世界坐标系进行区分，用b表示系统自身坐标系，w表示世界坐标系。是用四元数表示的从b系相对于w系的旋转。此外可以通过四元数转化得到旋转矩阵
57.通过对imu的测量数据进行积分，可以推导出从第m时刻到n时刻的位置速度旋转
[0058][0059][0060][0061]
其中，bg为陀螺仪的偏置，ba为加速度计的偏置，gw为重力加速度，ng陀螺仪的噪声，na为加速度计的噪声，ω(
·
)为四元数的右乘矩阵。
[0062]
其中，所述通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化，包括：
[0063]
根据预设的融合优化算法
将所述实时粗略位置信息、历史位置信息、真实ap和虚拟ap的估计位置信息融合优化，其中，ρ(
·
,c)为huber代价函数，c为经验阈值。
[0064]
在本发明实施例中，ap并不总是与用户放置在相同的高度，而是被放置在天花板附近，ap的高度也是未知的，这增加了优化的复杂性。然而，ap高度的变化对ap到用户的aoa影响很小，使得最终对定位性能的影响可以忽略，因此为简化系统复杂性本发明对室内场所用户的二维位置和方向进行优化，这也包含了目前室内场景大多数的应用需求。
[0065]
用户手机中的加速度计和陀螺仪传感器提供跨时间的相对位姿测量，这里将手机与用户近似为同一参考系。令状态空间s为用户在p个时间步长的一组姿态，则状态空间可表示为：
[0066][0067]
其中每个pi属于特殊欧几里得群。令pi与p
i+1
时刻的相对位姿测量为zi，同时其对角协方差矩阵为∑i。该协方差矩阵可以在imu预积分阶段通过对误差项求积分得到。用来表示残差的置信程度，协方差矩阵越大，表明当前的测量数据的可信度越低，其在优化过程中对结果的影响程度越小。
[0068]
考虑两个连续的位姿pi与p
i+1
。则相对位姿的测量模型可定义为：
[0069][0070]
其中r(
·
)∈so(3)是旋转矩阵。通过该预测模型，可以将两个连续的用户位姿测量和预测之间的误差函数定义为：
[0071][0072]
在考虑imu的误差模型后，需要继续对ap与vap的测量模型与误差模型进行定义。首先，由于本设计考虑无法获得ap与vap位置及方向先验信息的情况，因此需要在状态空间s中附加ap与vap的位置。为方便讨论，下面以仅使用单ap的情况进行分析，但该方法易扩展为多ap场景。因此，将状态空间s改写为：
[0073][0074]
其中a＝(a
x
,ay)∈se(2)表示ap的位置。表示第j个vap的位置。通过这种重新定义，本发明将ap的位姿纳入了优化的范围。
[0075]
下面对真实ap到用户间aoa测量进行建模：
[0076][0077]
其中旋转使ap到用户的aoa转化到用户坐标系下。由于用户手机天线数量限制无法进行aoa参数估计，通过上述sar方法可对ap与用户间aoa进行估计，选择每接收到m次csi数据进行一次sar参数估计，并取水平角作为每一次观测量。则ap与用户直射路径的
aoa测量和预测之间的误差函数定义为：
[0078][0079]
其中是ap到用户在第i个位置上的直射路径aoa测量值。需要说明的是，由于采用sar技术，因此获得aoa参数估计的速率相比于有m根天线的接收机下降了m倍数。
[0080]
同理对vap到用户间aoa测量进行建模：
[0081][0082]
vap与用户间aoa测量和预测之间的误差函数定义为：
[0083][0084]
需要说明的是，由于vap的信号本质上为真实ap的反射信号，因此用户无法区分解算的多径aoa分别来自哪个vap或真实ap。然而由于用户行走速度的限制，相邻两个位置间距离较小，因此对于同一个vap或真实ap在下一个时刻的aoa会受到当前的aoa数值以及两个时刻间的相对运动约束。由于两个时刻的间隔时间较短，可以通过imu在短时间内的积分获得。
[0085]
令k到k+1时刻手机的相对旋转表示为两个时刻的相对位移表示为以给定直射路径θk为例，通过第k个位置与第k+1个位置间的imu测量对θ
k+1
进行预测。
[0086]
由于ap相对于手机的aoa主要由ap与手机的相对位置和手机的姿态决定。因此该预测过程主要包含两部分的预测：一部分是由于手机位移导致角度变化α，另一部分是由于姿态变化导致的角度变化β。因此对于k+1时刻的预测值可以表示为
[0087]
β可以由旋转矩阵获得。但α在求解过程中存在需要取近似的地方，这是由于ap与手机距离l未知，但对于k到k+1时刻弧长求解过程中需要利用该参数。弧长的计算公式为：
[0088][0089]
因此，假设k和k+1时刻相对位移较小，弧长可以用位移近似代替。因此近似的弧长计算为：
[0090][0091][0092]
其中a
t
表示imu的线加速度，表示k时刻相对初始时刻的旋转。δt表示k到k+1之间的时间差。初始距离l可以选取经验值。则，初始距离对最终优化效果影响不大，因为本系统的wifi-imu融合优化模块可以不断的优化用户的位置和方向以及ap，vap的位置，从而更新了l。因此通过上述方法可以通过θk预测同理，该方法可以对其他vap下一时刻的aoa进行预测。因此，对于k+1时刻利用sar方法获得的n+1个多径aoa可以分别与预测的aoa进行对照，将与θ
k+1
最接近的匹配为同一vap或ap。
[0093]
最后，为了将用户转化在一个固定的全局参考系上，需要提供了一个初始先验位置综上所述，完整的优化函数可以定义为：
[0094][0095]
其中ρ(
·
,c)采用huber代价函数，c为经验阈值，能够进一步加强对异常值的优化。
[0096][0097]
此外，σ
ap
和σ
vap
分别表示ap与vap到用户的aoa方差，其数值可以通过记录一段时间数据，结合统计学公式获得。值得说明的一点是，本发明创新性地提出将vap误差加入总体的因子图中进行优化。之前大多融合算法使用的扩展卡尔曼滤波算法(extend kalman filter，ekf)且仅仅利用当前wifi的直射aoa进行融合，而将vap的aoa加入总体的因子图中，可以有效克服nlos情况，且能够利用当前与历史定位结果之间的联系，进一步降低定位误差。
[0098]
构建完上述的目标函数后，采用levenberg-marquadt方法求解该优化问题。通过以上方式，在无先验ap位置信息和nlos环境下，可以实现高精度和高鲁棒性的室内定位。
[0099]
本发明实施例，通过待定位设备获取室内的信号状态信息csi和传感信息；根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实时粗略位置信息；通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。使得多径aoa可以克服传统wifi定位算法中需要获得ap位置先验信息的问题，而且通过多径aoa也解决了在nlos或多径干扰严重场景下无法进行定位的问题。
[0100]
如图2所示为本发明再一个实施例提供的室内定位的方法的流程图，所述将所述待定位设备的行动轨迹虚拟为大孔径天线，根据合成孔径雷达sar的原理通过所述大孔径天线获取所述多径aoa，包括：
[0101]
步骤201，根据预设的信号接收公式获取所述待定位设备接收到的信号，记录预设时间段内《h
ij
(t),pi(t)》的数据对，将所述数据对设置为多天线阵列，其中，h
ij
(t)为t时刻待接收设备接收的信号，pi(t)t时刻待接收设备的位置，d
ij
(t)为t时刻待接收设备与信号发射端的距离；
[0102]
步骤202，根据预设的运动轨迹公式p
ij
(t)＝pi(t)-pj(t)获取所述待定位设备的相对运动轨迹；
[0103]
步骤203，根据预设的转向向量公式和所述相对运动轨迹获取转向向量
其中，(ρk,ξk,ψk),(ρ
l
,ξ
l
,ψ
l
)为k、l时刻待定位设备在世界坐标系下的坐标；
[0104]
步骤204，根据预设的能量谱公式获取能量谱，进而获取所述f
ij
(θ,φ)的峰值对应的(θ,φ)，并将所述峰值对应的(θ,φ)设置为多径aoa。
[0105]
在本发明实施例中，在的三维极坐标系下，为方便讨论令ap为发射端j且在坐标原点，令手机为接收端i。考虑到该过程在较短时间内完成，因此接收端pi(k)到pi(l)相对运动轨迹可由imu进行估计。接收端运动轨迹为一条由m组坐标构成的三维轨迹[pi(k)，...，pi(l)]，l＝k+m。接收信号可表示为：
[0106][0107]
这里d
ij
(t)代表发送端与接收端之间的距离。用pi(t)，pj(t)代表t时刻接收端与发送端的位置，那么，可以在接收端记录一段时间里的《h
ij
(t)，pi(t)》数据对，并将每个这样的数据对看作三维空间中的一根虚拟天线。这样，通过采集一段时间窗口内的《h
ij
(t)，pi(t)》数据对，便可以得到一个虚拟的三维空间上的多天线阵列。
[0108]
值得注意的是，信道状态信息采集h
ij
(t)所需的带宽开销是很小的。以wifi信号为例，只需要5kb/s的通信速率就能完成h
ij
(t)的采集，不会影响收发端之间的正常通信。
[0109]
有了虚拟多天线阵列，就可以仿照之前一维线性天线阵列，进行到达角估计。首先可以得到sar的信号接收模型：
[0110]yij
(t)＝[y(tk),...,y(t
l
)]
t
＝a(θ,φ)(t)h
ij
(t)+n(t)
[0111]
这里y
ij
(t)代表t时刻时间窗口[tk，...，t
l
]上产生的接收信号。由于天线阵列拓展到了三维极坐标系下，转向向量是θ，φ的函数。
[0112]
这样，同样可以利用bartlett波束赋形法来估计到达角，信号能量谱可以写作如下形式：
[0113][0114]fij
(θ，φ)是一个二维矩阵，它的每一个元素可以看作其代表的信号方向是真实的信号到达角的可能性。得到能量谱f
ij
(θ，φ)后，只需要知道其峰值对应的(θ，φ)便可以得到信号到达角。
[0115]
在分析一维线性天线阵列的信号接收模型时，引入了转向向量来帮助分析各天线上的相位关系。在给定接收端移动轨迹条件下，只需导出其对应的转向向量便可以计算到达角了。
[0116]
由于需要的只是手机与ap的相对方位，与手机和ap收发系统所处的绝对位置无关。因此对于发送端位置并无要求，接收端相对运动轨迹可以表示为：
[0117]
p
ij
(t)＝pi(t)-pj(t)
[0118]
利用相对轨迹，导出极坐标系下的三维天线阵列的转向向量：
[0119][0120]
这里(ρk，ξk，ψk)，(ρ
l
，ξ
l
，ψ
l
)代表k，l时刻接收端在世界坐标系下的坐标。
[0121]
得到了转向向量以及接收端在每个虚拟天线位置上收集到的信号h
ij
(t)，便可以通过bartlett波束赋形到达角解算法得到能量谱f
ij
(θ，φ)：
[0122][0123]fij
(θ，φ)的各个峰值对应的(θ，φ)即为三维空间中的信号多径到达角，在los情况下最高峰为直射路径到达角，其余峰值为各主反射面反射信号到达角度，即各vap到达角度，峰值大小可用于表征其置信度。
[0124]
然而要实现到达角估计，最重要的一步就是获取信道信息h
ij
(t)，即csi。对于一对收发端来说，它们使用的射频电路的晶振并不相同，这会导致收发设备之间存在载波频率的差异，这种差异称为载波频偏(carrier frequency offset，cfo)。cfo会导致接收机采集到的csi包含噪声，从而影响到达角预测精度。
[0125]
一种消除cfo的常用手段是信道互惠。信道互惠利用了前向信道和反向信道之间的关系，来消除cfo。因为不仅可以在接收端获得csi，在发射端也能接收到来自接收端信号的csi(反向信道)。接收端每收到一个来自发射端的数据包，便立即发回一个ack信号，这样发射端也能采集到一帧csi。
[0126]
收发端采集的csi可分别表达如下：
[0127][0128][0129]
和分别为前向信道和反向信道的csi，δf为cfo，τ为信号传播时间。这里假设在接收端接收到前向信道csi到发射机接收到反向信道csi的过程中接收端位移忽略不计，故上两式中信号传播时间τ相同。
[0130]
得到了双向信道信息后，可以由下式计算出不含cfo的原始信道：
[0131][0132]
因此可以通过和得到不包含cfo的原始信道h
ij
(t)。这里的κh是一个不随时间变化的常数。能量谱函数：
[0133][0134]
通过上式计算能量谱并取各个峰值对应的(θ，φ)，便可得到不受cfo影响的多径
到达角估计值。
[0135]
如图3所示为本发明一个实施例提供的室内定位的系统的结构图，所述系统包括：
[0136]
信息获取模块301，用于通过待定位设备获取室内的信号状态信息csi和传感信息；
[0137]
实时粗略位置信息获取模块302，用于根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实时粗略位置信息；
[0138]
估计及优化模块303，用于通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。
[0139]
在本发明实施例中，待定位设备为具有wifi信号接收模块以及惯性测量单元的设备，通常为手机、平板电脑等移动终端。待定位设备通过wifi信号接收模块获取室内信道的csi(channel sate information，信道状态信息)，通过惯性测量单元获取包含待定位设备的加速度、待定位设备的陀螺仪信息等传感信息。获取的csi通过计算可以获取多径aoa(angle of arrival，信道到达角)，而传感信息可以用于获取待定位设备在室内的实时大致位置信息，即实时粗略位置信息，对于实时粗略位置信息由于其只是通过惯性测量单元获取的大致信息，准确度不够，而且随着时间的增加，惯性信息获取的位置信息会越来越不准确，因此，在实时粗略位置信息的基础上，引入多径aoa和融合优化算法，通过融合优化算法对多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息待定位设备的当前位置信息进行估计，并对待优化参数进行实时优化，获取位置更为准确的优化后的位置信息。
[0140]
其中，所述根据所述csi获取多径信道到达角aoa，包括：
[0141]
将所述待定位设备的行动轨迹虚拟为大孔径天线，根据合成孔径雷达sar的原理通过所述大孔径天线获取所述多径aoa。
[0142]
在本发明实施例中，由于待定位设备通常只拥有1-2根天线，该天线数量并不足以通过多天线阵列的方式估计出多径aoa，因此，需要获取待定位设备的行动轨迹，将行动轨迹虚拟为大孔径天线，根据sar的原理通过该大孔径雷达获取多径aoa信息。
[0143]
惯性导航系统具有抗干扰性强，不需要借助外部设备等优点，目前大多智能手机都配备着惯性测量单元，imu由三个正交陀螺仪和三个正交加速度计组成，分别用于测量系统的加速度和比力。通过积分可以从加速度和角速度中得到位置、速度和姿态(也称为旋转)三个量，分别用p、v、q表示。因此ins是获取载体状态信息的有效方案。此外，imu测量信息频率相对于大多数传感器更高。
[0144]
在惯性导航中，对于姿态表述常用到的方法有四元数法、欧拉角法等。欧拉角使用俯仰角、横滚角、航向角来描述姿态发生的变化，其优点在于更符合人们的思维习惯，更加直观的描述了角度的变化。然而其有一个重大缺陷，在俯仰角绝对值为90
°
时，翻滚角与俯仰角无法进行区分，出现“锁死”现象。四元数可用于描述载体的旋转，它将三维空间中的旋转用一个单位长度的四维向量来表示，一个四维向量包含一个实部和三个虚部，其优势在于能够无歧义的描述姿态。本发明采用四元素法来对姿态进行描述：
[0145]
q＝q0+q1i+q2j+q3k
[0146]
其中i,j,k为三个虚部。i2＝-1,j2＝-1,k2＝-1，且满足此外，
需要对载体坐标系和世界坐标系进行区分，用b表示系统自身坐标系，w表示世界坐标系。是用四元数表示的从b系相对于w系的旋转。此外可以通过四元数转化得到旋转矩阵
[0147]
通过对imu的测量数据进行积分，可以推导出从第m时刻到n时刻的位置速度旋转
[0148][0149][0150][0151]
其中，bg为陀螺仪的偏置，ba为加速度计的偏置，gw为重力加速度，ng陀螺仪的噪声，na为加速度计的噪声，ω(
·
)为四元数的右乘矩阵。
[0152]
其中，所述通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化，包括：
[0153]
根据预设的融合优化算法将所述实时粗略位置信息、历史位置信息、真实ap和虚拟ap的估计位置信息融合优化，其中，ρ(
·
,c)为huber代价函数，c为经验阈值。
[0154]
在本发明实施例中，ap并不总是与用户放置在相同的高度，而是被放置在天花板附近，ap的高度也是未知的，这增加了优化的复杂性。然而，ap高度的变化对ap到用户的aoa影响很小，使得最终对定位性能的影响可以忽略，因此为简化系统复杂性本发明对室内场所用户的二维位置和方向进行优化，这也包含了目前室内场景大多数的应用需求。
[0155]
用户手机中的加速度计和陀螺仪传感器提供跨时间的相对位姿测量，这里将手机与用户近似为同一参考系。令状态空间s为用户在p个时间步长的一组姿态，则状态空间可表示为：
[0156][0157]
其中每个pi属于特殊欧几里得群。令pi与p
i+1
时刻的相对位姿测量为zi，同时其对角协方差矩阵为∑i。该协方差矩阵可以在imu预积分阶段通过对误差项求积分得到。用来表示残差的置信程度，协方差矩阵越大，表明当前的测量数据的可信度越低，其在优化过程中对结果的影响程度越小。
[0158]
考虑两个连续的位姿pi与p
i+1
。则相对位姿的测量模型可定义为：
[0159][0160]
其中r(
·
)∈so(3)是旋转矩阵。通过该预测模型，可以将两个连续的用户位姿测
量和预测之间的误差函数定义为：
[0161][0162]
在考虑imu的误差模型后，需要继续对ap与vap的测量模型与误差模型进行定义。首先，由于本设计考虑无法获得ap与vap位置及方向先验信息的情况，因此需要在状态空间s中附加ap与vap的位置。为方便讨论，下面以仅使用单ap的情况进行分析，但该方法易扩展为多ap场景。因此，将状态空间s改写为：
[0163][0164]
其中a＝(a
x
,ay)∈se(2)表示ap的位置。表示第j个vap的位置。通过这种重新定义，本发明将ap的位姿纳入了优化的范围。
[0165]
下面对真实ap到用户间aoa测量进行建模：
[0166][0167]
其中旋转使ap到用户的aoa转化到用户坐标系下。由于用户手机天线数量限制无法进行aoa参数估计，通过上述sar方法可对ap与用户间aoa进行估计，选择每接收到m次csi数据进行一次sar参数估计，并取水平角作为每一次观测量。则ap与用户直射路径的aoa测量和预测之间的误差函数定义为：
[0168][0169]
其中是ap到用户在第i个位置上的直射路径aoa测量值。需要说明的是，由于采用sar技术，因此获得aoa参数估计的速率相比于有m根天线的接收机下降了m倍数。
[0170]
同理对vap到用户间aoa测量进行建模：
[0171][0172]
vap与用户间aoa测量和预测之间的误差函数定义为：
[0173][0174]
需要说明的是，由于vap的信号本质上为真实ap的反射信号，因此用户无法区分解算的多径aoa分别来自哪个vap或真实ap。然而由于用户行走速度的限制，相邻两个位置间距离较小，因此对于同一个vap或真实ap在下一个时刻的aoa会受到当前的aoa数值以及两个时刻间的相对运动约束。由于两个时刻的间隔时间较短，可以通过imu在短时间内的积分获得。
[0175]
令k到k+1时刻手机的相对旋转表示为两个时刻的相对位移表示为以给定直射路径θk为例，通过第k个位置与第k+1个位置间的imu测量对θ
k+1
进行预测。
[0176]
由于ap相对于手机的aoa主要由ap与手机的相对位置和手机的姿态决定。因此该预测过程主要包含两部分的预测：一部分是由于手机位移导致角度变化α，另一部分是由于姿态变化导致的角度变化β。因此对于k+1时刻的预测值可以表示为
[0177]
β可以由旋转矩阵获得。但α在求解过程中存在需要取近似的地方，这是由于ap
与手机距离l未知，但对于k到k+1时刻弧长求解过程中需要利用该参数。弧长的计算公式为：
[0178][0179]
因此，假设k和k+1时刻相对位移较小，弧长可以用位移近似代替。因此近似的弧长计算为：
[0180][0181][0182]
其中a
t
表示imu的线加速度，表示k时刻相对初始时刻的旋转。δt表示k到k+1之间的时间差。初始距离l可以选取经验值。则，初始距离对最终优化效果影响不大，因为本系统的wifi-imu融合优化模块可以不断的优化用户的位置和方向以及ap，vap的位置，从而更新了l。因此通过上述方法可以通过θk预测同理，该方法可以对其他vap下一时刻的aoa进行预测。因此，对于k+1时刻利用sar方法获得的n+1个多径aoa可以分别与预测的aoa进行对照，将与θ
k+1
最接近的匹配为同一vap或ap。
[0183]
最后，为了将用户转化在一个固定的全局参考系上，需要提供了一个初始先验位置如默认进入大门时刻的位置为初始位置。综上所述，完整的优化函数可以定义为：
[0184][0185]
其中ρ(
·
,c)采用huber代价函数，c为经验阈值，能够进一步加强对异常值的优化。
[0186][0187]
此外，σ
ap
和σ
vap
分别表示ap与vap到用户的aoa方差，其数值可以通过记录一段时间数据，结合统计学公式获得。值得说明的一点是，本发明创新性地提出将vap误差加入总体的因子图中进行优化。之前大多融合算法使用的扩展卡尔曼滤波算法(extend kalman filter，ekf)且仅仅利用当前wifi的直射aoa进行融合，而将vap的aoa加入总体的因子图中，可以有效克服nlos情况，且能够利用当前与历史定位结果之间的联系，进一步降低定位误差。
[0188]
构建完上述的目标函数后，采用levenberg-marquadt方法求解该优化问题。通过以上方式，在无先验ap位置信息和nlos环境下，可以实现高精度和高鲁棒性的室内定位。
[0189]
本发明实施例，通过待定位设备获取室内的信号状态信息csi和传感信息；根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实
时粗略位置信息；通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。使得多径aoa可以克服传统wifi定位算法中需要获得ap位置先验信息的问题，而且通过多径aoa也解决了在nlos或多径干扰严重场景下无法进行定位的问题。
[0190]
如图4所示为本发明再一个实施例提供的室内定位的系统的结构图，所述实时粗略位置信息获取模块302，包括：
[0191]
多天线阵列设置单元401，用于根据预设的信号接收公式获取所述待定位设备接收到的信号，记录预设时间段内《h
ij
(t),pi(t)》的数据对，将所述数据对设置为多天线阵列，其中，h
ij
(t)为t时刻待接收设备接收的信号，pi(t)t时刻待接收设备的位置，d
ij
(t)为t时刻待接收设备与信号发射端的距离；
[0192]
相对运动轨迹获取单元402，用于根据预设的运动轨迹公式p
ij
(t)＝pi(t)-pj(t)获取所述待定位设备的相对运动轨迹；
[0193]
转向向量获取单元403，用于根据预设的转向向量公式和所述相对运动轨迹获取转向向量其中，(ρk,ξk,ψk),(ρ
l
,ξ
l
,ψ
l
)为k、l时刻待定位设备在世界坐标系下的坐标；
[0194]
多径aoa设置单元404，用于根据预设的能量谱公式获取能量谱，进而获取所述f
ij
(θ,φ)的峰值对应的(θ,φ)，并将所述峰值对应的(θ,φ)设置为多径aoa。
[0195]
在本发明实施例中，在的三维极坐标系下，为方便讨论令ap为发射端j且在坐标原点，令手机为接收端i。考虑到该过程在较短时间内完成，因此接收端pi(k)到pi(l)相对运动轨迹可由imu进行估计。接收端运动轨迹为一条由m组坐标构成的三维轨迹[pi(k)，...，pi(l)]，l＝k+m。接收信号可表示为：
[0196][0197]
这里d
ij
(t)代表发送端与接收端之间的距离。用pi(t)，pj(t)代表t时刻接收端与发送端的位置，那么，可以在接收端记录一段时间里的《h
ij
(t)，pi(t)》数据对，并将每个这样的数据对看作三维空间中的一根虚拟天线。这样，通过采集一段时间窗口内的《h
ij
(t)，pi(t)》数据对，便可以得到一个虚拟的三维空间上的多天线阵列。
[0198]
值得注意的是，信道状态信息采集h
ij
(t)所需的带宽开销是很小的。以wifi信号为例，只需要5kb/s的通信速率就能完成h
ij
(t)的采集，不会影响收发端之间的正常通信。
[0199]
有了虚拟多天线阵列，就可以仿照之前一维线性天线阵列，进行到达角估计。首先
可以得到sar的信号接收模型：
[0200]yij
(t)＝[y(tk),...,y(t
l
)]
t
＝a(θ,φ)(t)h
ij
(t)+n(t)
[0201]
这里y
ij
(t)代表t时刻时间窗口[tk，...，t
l
]上产生的接收信号。由于天线阵列拓展到了三维极坐标系下，转向向量是θ，φ的函数。
[0202]
这样，同样可以利用bartlett波束赋形法来估计到达角，信号能量谱可以写作如下形式：
[0203][0204]fij
(θ，φ)是一个二维矩阵，它的每一个元素可以看作其代表的信号方向是真实的信号到达角的可能性。得到能量谱f
ij
(θ，φ)后，只需要知道其峰值对应的(θ，φ)便可以得到信号到达角。
[0205]
在分析一维线性天线阵列的信号接收模型时，引入了转向向量来帮助分析各天线上的相位关系。在给定接收端移动轨迹条件下，只需导出其对应的转向向量便可以计算到达角了。
[0206]
由于需要的只是手机与ap的相对方位，与手机和ap收发系统所处的绝对位置无关。因此对于发送端位置并无要求，接收端相对运动轨迹可以表示为：
[0207]
p
ij
(t)＝pi(t)-pj(t)
[0208]
利用相对轨迹，导出极坐标系下的三维天线阵列的转向向量：
[0209][0210]
这里(ρk，ξk，ψk)，(ρ
l
，ξ
l
，ψ
l
)代表k，l时刻接收端在世界坐标系下的坐标。
[0211]
得到了转向向量以及接收端在每个虚拟天线位置上收集到的信号h
ij
(t)，便可以通过bartlett波束赋形到达角解算法得到能量谱f
ij
(θ，φ)：
[0212][0213]fij
(θ，φ)的各个峰值对应的(θ，φ)即为三维空间中的信号多径到达角，在los情况下最高峰为直射路径到达角，其余峰值为各主反射面反射信号到达角度，即各vap到达角度，峰值大小可用于表征其置信度。
[0214]
然而要实现到达角估计，最重要的一步就是获取信道信息h
ij
(t)，即csi。对于一对收发端来说，它们使用的射频电路的晶振并不相同，这会导致收发设备之间存在载波频率的差异，这种差异称为载波频偏(carrier frequency offset，cfo)。cfo会导致接收机采集到的csi包含噪声，从而影响到达角预测精度。
[0215]
一种消除cfo的常用手段是信道互惠。信道互惠利用了前向信道和反向信道之间的关系，来消除cfo。因为不仅可以在接收端获得csi，在发射端也能接收到来自接收端信号的csi(反向信道)。接收端每收到一个来自发射端的数据包，便立即发回一个ack信号，这样
发射端也能采集到一帧csi。
[0216]
收发端采集的csi可分别表达如下：
[0217][0218][0219]
和分别为前向信道和反向信道的csi，δf为cfo，τ为信号传播时间。这里假设在接收端接收到前向信道csi到发射机接收到反向信道csi的过程中接收端位移忽略不计，故上两式中信号传播时间τ相同。
[0220]
得到了双向信道信息后，可以由下式计算出不含cfo的原始信道：
[0221][0222]
因此可以通过和得到不包含cfo的原始信道h
ij
(t)。这里的κh是一个不随时间变化的常数。能量谱函数：
[0223][0224]
通过上式计算能量谱并取各个峰值对应的(θ，φ)，便可得到不受cfo影响的多径到达角估计值。
[0225]
图5示例了一种电子设备的实体结构示意图，如图5所示，该电子设备可以包括：处理器(processor)501、通信接口(communications interface)502、存储器(memory)503和通信总线504，其中，处理器，通信接口，存储器通过通信总线完成相互间的通信。处理器可以调用存储器中的逻辑指令，以执行室内定位的方法，该方法包括：通过待定位设备获取室内的信号状态信息csi和传感信息；根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实时粗略位置信息；通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。
[0226]
此外，上述的存储器中的逻辑指令可以通过软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：u盘、移动硬盘、只读存储器(rom，read-only memory)、随机存取存储器(ram，random access memory)、磁碟或者光盘等各种可以存储程序代码的介质。
[0227]
另一方面，本发明实施例还提供一种计算机程序产品，所述计算机程序产品包括存储在非暂态计算机可读存储介质上的计算机程序，所述计算机程序包括程序指令，当所述程序指令被计算机执行时，计算机能够执行上述各方法实施例所提供的室内定位的方法，该方法包括：通过待定位设备获取室内的信号状态信息csi和传感信息；根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实时粗略
位置信息；通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。
[0228]
又一方面，本发明实施例还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现以执行上述各实施例提供的室内定位的方法，该方法包括：通过待定位设备获取室内的信号状态信息csi和传感信息；根据所述csi获取多径信道到达角aoa，并根据所述传感信息获取所述待定位设备在所述室内的实时粗略位置信息；通过预设的融合优化算法根据所述实时粗略位置信息、多径aoa、历史位置信息、真实ap和虚拟ap的估计位置信息对所述待定位设备的当前位置进行估计，并对待优化参数进行实时优化。
[0229]
应该理解的是，虽然附图的流程图中的各个步骤按照箭头的指示依次显示，但是这些步骤并不是必然按照箭头指示的顺序依次执行。除非本文中有明确的说明，这些步骤的执行并没有严格的顺序限制，其可以以其他的顺序执行。而且，附图的流程图中的至少一部分步骤可以包括多个子步骤或者多个阶段，这些子步骤或者阶段并不必然是在同一时刻执行完成，而是可以在不同的时刻执行，其执行顺序也不必然是依次进行，而是可以与其他步骤或者其他步骤的子步骤或者阶段的至少一部分轮流或者交替地执行。
[0230]
以上所述仅是本发明的部分实现方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。