一种基于LINE-MOD模板匹配的增强现实三维跟踪注册方法与流程

文档序号:17445992发布日期:2019-04-17 05:37阅读:1773来源:国知局
一种基于LINE-MOD模板匹配的增强现实三维跟踪注册方法与流程

本发明属于增强现实技术领域,具体涉及一种基于line-mod模板匹配的增强现实三维跟踪注册方法。



背景技术:

刚性物体检测和位姿估计是近年来增强现实技术中研究的热点也是研究难点。根据物体表观特征,可以将其划分为无纹理刚性物体检测与位姿估计和有纹理刚性物体检测与位姿估计。当刚性物体表面纹理丰富时,能够提取到数量稳定的特征点和特征描述符,基于这些稳定的特征点和特征描述符,物体的位姿就能被准确的检测和描述。其中sift(scaleinvariantfeaturetransform)、pca-sift、surf(speededuprobustfeatures)等特征点和特征描述符对刚性纹理物体的检测和位姿估计方法的发展起到了积极的推动作用。sift将多尺度空间上提取的高斯差分图极值点作为特征点,并将每个特征点邻域内计算获得的梯度方向直方图作为特征描述子。由于sift具有良好的尺度、旋转以及一定视角和光照不变性,并且具有很强的可区分性,因此自被提出以来,在物体检测与位姿估计领域获得到了广泛应用。surf是sift的改进算法,它利用haar小波替换sift方法中的梯度直方图计算,利用积分图技术对运算过程进行加速,使surf的运算速度较sift提高了3-7倍,同时在性能上维持与sift相当,因此它在很多应用中得到了推广,尤其是在对运行速度要求高的场景。brief(binaryrobustindependentelementfeature)利用图像邻域内随机点对的灰度大小差异测试,将图像块表示为二进制字符串,此描述符的构造和匹配不仅速度快,而且具有更高的识别率,因此在运算能力有限的手机等移动设备上具有很好的应用前景。除了brief特征描述子,近几年还提出了许多二进制特征描述子,例如orb、brisk、freak。gordon采用sift特征描述符以及从运动中恢复结构的理论框架,并利用已知对象模型完成物体位姿估计。collet等利用sift特征描述符在离线训练阶段构建的三维局部描述符,在线运行时,通过sift特征提取与匹配,对物体姿态进行估计。lowe等采用sift特征进行物体检测,该方法通过对待检测物体提取特征点,并利用快速邻域法(fastnearest-neighboralgorithm)与已知物体数据库进行特征匹配。然后用hough变换鉴别物体的类属,最后通过姿态参数一致性进行物体确认。但是由于机械零件表面缺少足够的纹理特征,导致很难提取到稳定且判别力强的特征点,因此上述方法描述的基于特征点的刚性物体检测和跟踪方法难以直接适用于机械产品装配现场。

模板匹配近年来在刚性物体检测和位姿估计中发挥了重要作用,该类方法操作过程简单、具有强大的处理不同对象的能力。它既不需要一个庞大的训练集,也不需要一个耗时的训练阶段,并且对于纹理刚性物体和无纹理刚性物体都能处理。然而,模板匹配鲁棒性的增加往往伴随着计算负荷的增大,导致模板匹配难以满足实时性要求。到目前为止,许多学者试图降低这种复杂性。早期的模板匹配方法和它的扩展使用chamfer距离度量模板和输入图像轮廓之间的差异。例如gavrila介绍了一个在形状和参数空间上的由粗到精的二值边缘图像chamfer距离度量。chamfer匹配能最大限度地减少了两组边缘点之间的广义距离,实时性较好,但这种方法对遮挡等造成的外点极为敏感。另一种二值图像的距离度量方式是hausdorff距离,它测量从图像中的每个边缘点到模板中最近邻的所有距离的最大值。但是它易受遮挡和复杂背景的影响。huttenlocher等试图引入广义hausdorff距离来克服这一缺陷,该方法在一定程度上克服了遮挡和杂乱背景的影响,但是该方法需要预先对背景杂波进行估计,此外,当模板较多时,计算代价很高。上述这些方法使用的二值图像大多是通过类似canny的边缘提取算法得到的,因此它们对光照变化、噪声和模糊极为敏感。为了避免上述算法的缺陷,hinterstoisser等提出使用图像梯度而非图像轮廓作为匹配特征的line-mod方法,方法采用梯度方向特征的二进制方式表示,巧妙地运用了现代计算机的超高速缓存进行并行处理。该方法能够在背景复杂环境下实时对多类无纹理刚性物体进行检测。然而,该方法只能在固定尺度对刚性物体进行检测与跟踪,同时该方法对物体位姿估计精度较低。

综上所述,基于特征点的刚性物体检测和位姿估计方法需要物体表面具有丰富的纹理特征,对于表面光滑的机械零部件该方法难以直接适用。模板匹配是解决无纹理刚性物体检测和位姿估计的有效方法,该类方法实时性好,不需要物体表面有足够的纹理特征。其中,line-mod方法是近年来最先进的模板匹配方法之一,本文的刚性物体检测与位姿估计是在line-mod的启发下开展的,但是该方法只能在固定尺度对刚性物体进行检测与位姿估计,同时该方法对物体位姿估计精度较低。如何提高算法多尺度检测与位姿估计精度,仍有待研究。



技术实现要素:

本发明旨在解决以上现有技术的问题。提出了一种提高算法多尺度检测与位姿估计精度的一种基于line-mod模板匹配的增强现实三维跟踪注册方法。本发明的技术方案如下:

一种基于line-mod模板匹配的增强现实三维跟踪注册方法,其包括离线阶段和在线阶段,具体包括以下步骤:

设计改进的line-mod模板匹配,并以改进的line-mod模板匹配为基础,离线阶段在cad环境中对待装配零部件进行训练,采集多视角下的参考图像,并在采集过程中记录相机的采集位置和姿态,然后对采集的参考图像进行描述符提取并存储到xml文件中;在在线阶段,采集真实物体视频帧,并对每一帧图像进行物体梯度方向提取,然后将在线提取的梯度方向描述符与离线阶段提取的梯度方向描述符进行匹配,获取与当前输入视频帧最相似的关键帧和关键帧对应的相机位姿,最后利用关键帧的空间拓扑关系,采用插值的方法对待装配对象的姿态角和位置进行估计,完成对待装配对象的跟踪和位姿估计。

进一步的,所述离线阶段在cad环境中对待装配零部件进行训练,采集多视角下的参考图像,并在采集过程中记录相机的采集位置和姿态,具体包括:首先在cad环境中,通过虚拟相机采集装配对象不同视角下的投影视图,每一幅视图中包含了装配对象轮廓的不同形状及每一幅视图对应相机采集位姿,采用在正二十面体表面进行采样,将正二十面体的每一个面平均分成四份并进行迭代划分,迭代2次,使每一个面上形成16个等边三角形,在cad环境中,将物体cad模型放置在正二十面体中心,虚拟相机光轴始终通过中心。采样视点位于每一个三角形顶点处。

进一步的,当机械零件存在大量的对称结构时,对于这类特殊零件,降低采样数量;对于轴套类的某些零件,只需在经度0°,纬度0~90°区间进行采样;对于轮盘类的某些零件,在经度上的采样区间为0~90°,纬度上的采样区间为0~90°;对于叉架类、盖板类和箱壳类的某些零件,在经度上的采样区间为0~180°,纬度上的采样区间为0~90°;对于其他非对称零件,在经度上的采样区间为0~360°,纬度上的采样区间为0~90°;对于每一个视点,在-45°~45°旋转区域内形成四个旋转角度,用以表示同一位置不同视角的情况。最后提取所有视图中零件轮廓的形状特征描述符,保存所有零件轮廓的形状特征描述符和对应的位姿参数到xml文件中。

进一步的,所述在线阶段,采集真实物体视频帧,并对每一帧图像进行物体梯度方向提取,然后将在线提取的梯度方向描述符与离线阶段提取的梯度方向描述符进行匹配,具体包括:

对于每一幅采样得到的视图,计算目标物体轮廓每一个位置{r,g,b}三个颜色通道上的梯度向量,并将三个通道中的梯度向量最大值作为该点处的梯度方向,一幅输入图像i,在位置x处的梯度方向为:

其中:

式中:r,g,b是视图i中位置x处的三个颜色通道。c是r,g,b三个颜色通道上的某一颜色通道。ori()表示梯度方向;将梯度向量方向共线的作为一个方向所有的梯度向量方向被统一到0~180°区间范围。然后将0~180°的梯度向量方向空间平均分成8份,将每一个位置x的3×3像素邻域内出现次数最多且梯度向量的模大于给定阈值的方向作为该点处的梯度方向。

进一步的,在离线阶段提取物体的轮廓梯度后,对每一个梯度方向进行发散处理。

进一步的,在线阶段待装配对象在线识别具体包括:

为了完成待装配物体在任意尺度下的识别,将深度信息引入相似性评价方法,使其具有尺度不变性,改进后的相似性评价公式为:

式中:ori(o,so(c′o,r′))—在参考图像o位置点处的梯度方向的弧度表示;ori(i,t)—在输入图像i位置t点处的梯度方向的弧度表示;

表示以为中心,τ为邻域的区域;p—是位置r的列表;t=(o,p)—表示物体的模板;d(c′x)—表示在点c′x处的深度值;d(c'o)—表示离线训练时参考图像o相机到正多面体中心的距离。

进一步的,所述.姿态角估计包括:在模板训练时,相机光轴始终通过零件坐标系原点,但在线匹配时,相机光轴可能出现在空间任意位置,(1)首先对相机相对关键帧的姿态角进行求解;(2)为了简化相机姿态角计算过程,假设一个虚拟成像平面,c-x1y1z1是当前相机的坐标系,z1为相机的光轴,与虚拟像平面1相交于点o1,c-x2y2z2为零件训练时相机坐标系,z2为零件训练时相机的光轴,它穿过零件坐标系原点,与像平面1和像平面2分别交于点o1'和o2,z2与x1z1平面夹角为θ2,在x1z1平面的投影为直线cp,cp与z1轴的夹角为θ1;(3)变换相机坐标系(4)除了相机光轴偏转带来的姿态角变化外,还有零件绕自身坐标系旋转带来的相机相对姿态角的变化;(5)对零件绕自身坐标系旋转带来的姿态角变化进行计算。

进一步的,所述位置估计包括:

为了更加准确的获得相机位置,采用插值方法进行精确计算:(1)首先采用姿态角中方法将“轴侧”视图变成“轴向”视图并提取形状描述符d',然后将其与数据库中的所有训练视图的形状描述符n(d,k)进行匹配,其中,dk(k=1,2,3...n)表示训练图像形状描述符集合,k表示采样视点之间的拓扑关系;(2)从而获得与d'最相似的对应网格上的图相描述符di及其相邻节点形状描述符dj;相机位置目标优化函数为:(3)这里通过计算可以得到权因子ωi及ωij的值,优化计算出来的权因子ωi及ωij是表示离线零件形状描述符库n中与当前输入轴向视图描述符d'最接近的图像描述符di的重心坐标数据,再根据权因子ωi及ωij的值来完成相机位置的插值计算,ωi和ωij是加权因子的形状描述符,因此公式能被写成:这里可以得到加权因子ωi和ωij;(4)最后,根据权因子ωi及ωij的值来完成相机在正二十面体表面位置的插值计算;

(5)摄像机位置通过将经过正二十面体表面位置0的点进行平移获得,其平移位置可以通过深度相机获得,因此摄像机的最终位置和姿态可以用公式表示:p=fm(do,r,t,mr)式中:do——深度传感器获取的o点的深度信息,r,t——深度传感器和颜色传感器的相对位置,可以通过标定获得,mr——是相机相对零件的姿态,可以通过姿态角估计获得。

本发明的优点及有益效果如下:

本发明设计了一种基于无纹理目标物体跟踪与位姿估计的刚性待装配对象与虚拟装配引导信息的虚实遮挡处理方法。首先对算法的基本思想进行了介绍,概述了本章方法的整体思路。引入深度信息,使line-mod模板匹配方法能完成任意尺度下的物体识别与匹配。以改进的line-mod模板匹配为基础,通过当前视角图像和离线不同视角参考视图的匹配,获取与当前视角最接近的关键帧及关键帧对应位姿。通过参考视图之间的拓扑关系,利用插值的方法对零部件位姿进行精确估计。基于跟踪与位姿估计结果在待装配零部件表面叠加相应透明几何模型,通过模型之间空间位置关系的判断,完成遮挡关系处理。实验结果表明,本方法优于hinterstoisser在accv上提出的无纹理目标物体跟踪与位姿估计精度,且算法实时性好(≥15帧/s),可较好地满足ar装配系统待装配对象和引导信息之间的虚实遮挡应用的需求。

附图说明

图1是本发明提供优选实施例工作流程示意图;

图2是待装配零件位姿估计示意图;

图3是相机小孔成像模型;

图4是相机自身坐标系变换示意图;

图5是当前视角下物体在虚拟平面1上的投影。

具体实施方式

下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、详细地描述。所描述的实施例仅仅是本发明的一部分实施例。

本发明解决上述技术问题的技术方案是:

1.算法基本思想和工作流程

系统总体方案流程图如图1所示。首先离线阶段在cad环境中对待装配零部件进行训练,采集多视角下的参考图像,并在采集过程中记录相机的采集位置和姿态。然后对采集的参考图像进行描述符提取并存储到xml文件中。在在线阶段,采集真实物体视频帧,并对每一帧图像进行物体梯度方向提取,然后将在线提取的梯度方向描述符与离线阶段提取的梯度方向描述符进行匹配,获取与当前输入视频帧最相似的关键帧和关键帧对应的相机位姿。最后利用关键帧的空间拓扑关系,采用插值的方法对待装配对象的姿态角和位置进行估计,完成对待装配对象的跟踪和位姿估计。最后基于跟踪结果,在待装配对象表面叠加相应的三维几何模型,通过模型之间空间位置判断,完成装配对象和虚拟装配引导信息之间的虚实遮挡关系处理。

2.多视角零件图像采集

目标物体的训练主要有在线和离线两种方式,但由于在线学习需要人工采集样本训练集,过程复杂,且在线学习采用的贪心法很难在速度和鲁棒性之间得到很好的权衡。因此,本文采用离线的方式,利用物体的cad模型生成不同视角下的参考图像。

首先在cad环境中,通过虚拟相机采集装配对象不同视角下的投影视图,每一幅视图中包含了装配对象轮廓的不同形状及每一幅视图对应相机采集位姿。为了实现采样视角的均匀分布,避免出现半球表面采样极点处采样密度过大导致的视图冗余问题,本文采用在正二十面体表面进行采样。将正二十面体的每一个面平均分成四份并进行迭代划分。为了权衡运算速度和精度,最终选择迭代2次,使每一个面上形成16个等边三角形。在cad环境中,将物体cad模型放置在正二十面体中心,虚拟相机光轴始终通过中心。采样视点位于每一个三角形顶点处。

由于机械零件存在大量的对称结构,对于这类特殊零件,可降低采样数量,以提高识别和跟踪速度。对于轴套类的某些零件,只需在经度0°,纬度0~90°区间进行采样。对于轮盘类的某些零件,在经度上的采样区间为0~90°,纬度上的采样区间为0~90°。对于叉架类、盖板类和箱壳类的某些零件(e)),在经度上的采样区间为0~180°,纬度上的采样区间为0~90°。对于其他非对称零件,在经度上的采样区间为0~360°,纬度上的采样区间为0~90°。对于每一个视点,在-45°~45°旋转区域内形成四个旋转角度,用以表示同一位置不同视角的情况。最后提取所有视图中零件轮廓的形状特征描述符,保存所有零件轮廓的形状特征描述符和对应的位姿参数到xml文件中。

3.梯度向量描述符提取

在计算梯度向量之前,首先简单介绍使用梯度向量和梯度向量作为描述符的优势。相比于其他形式的局部描述符,图像梯度具有良好的可识别性,并且对光照强度变化和噪音具有良好的鲁棒性。特别是针表面光滑且无纹理的机械零部件时,图像梯度通常是唯一可靠的特征,因此本文将物体的梯度向量作为特征描述符。

对于每一幅采样得到的视图,计算目标物体轮廓每一个位置{r,g,b}三个颜色通道上的梯度向量,并将三个通道中的梯度向量最大值作为该点处的梯度方向。例如,一幅输入图像i,在位置x处的梯度方向为:

其中:

式中:

r,g,b是视图i中位置x处的三个颜色通道。c是r,g,b三个颜色通道上的某一颜色通道。ori()表示梯度方向。

由于背景、光照等因素的影响,物体轮廓像素点的梯度方向可能发生反转。为了量化这些梯度向量,同时提高对背景、光照等因素的鲁棒性,将梯度向量方向共线的作为一个方向,因此所有的梯度向量方向被统一到0~180°区间范围。然后将0~180°的梯度向量方向空间平均分成8份。同时为了提高对噪声的鲁棒性,将每一个位置x的3×3像素邻域内出现次数最多且梯度向量的模大于给定阈值的方向作为该点处的梯度方向。

4.梯度方向发散

在线匹配过程中可能会由于输入图像的微小形变和位置偏移,导致提取的物体轮廓梯度方向与离线提取的物体轮廓梯度方向不匹配。为了提高匹配过程的鲁棒性,在离线阶段提取物体的轮廓梯度后,对每一个梯度方向进行发散处理,。例如将零部件表面某一位置x的梯度方向发散,就是以当前位置为中心的3×3像素坐标邻域内,都采用位置x的梯度方向进行表示。梯度方向发散后形成模板。为在线过程提取到的同一物体同一视角的轮廓梯度方向。可以看出,由于透视形变或视角微小偏差,提取到的轮廓梯度方向与离线阶段提取到的物体轮廓梯度方向不一致。但由于采用了梯度方向发散的方法,当输入图像在模板图像上进行滑动时,仍能与模板图形获得良好的匹配效果,从而提高了方法对微小形变和位置偏移的鲁棒性。

5.待装配对象在线识别

为了将在线阶段提取的物体表面梯度向量与离线阶段产生的梯度向量模板进行匹配,hinterstoisser等提出采用式(3)对其相似性进行度量(后文称line-mod):

式中:

ε——表示待识别零部件与离线训练阶段获得的零部件训练模板之间的相似性。

ori(o,r)—在参考图像o位置r点处的梯度方向的弧度表示;

ori(i,t)—在输入图像i位置t点处的梯度方向的弧度表示;

—表示c+r为中心,τ为邻域的区域。

p—是位置r的列表;

t=(o,p)—表示物体的模板。

式3表示,在模板匹配过程中,包含零部件训练模板的滑动窗口在输入图像上进行滑动,并统计零部件表面轮廓处每一个位置与训练模板对应位置的相似性。当所有位置的梯度方向相似性总和大于阈值δ,则待识别零部件在图像中的初始位置以滑动窗口中心为中心,以零部件训练模板最小包围盒为边界框进行确定。在搜索结束后,对所有具有重叠阈值为0.5的潜在边界框执行非极大值抑制,剩余边界框即为输入图像中零部件的位置,从而完成物体的识别。

虽然上述方法能很好地在具有微小位置移动和形变的图像中识别出待装配对象,但是该方法只能在固定尺度对待装配对象进行识别。例如,如果离线阶段生成模板的图像采集时远离相机,而在线运行时采集的图像距离相机较近,那么即使在线采集的图像中物体与模板中物体具有相同的梯度方向,上述方法仍然会显示不匹配。

为了完成待装配物体在任意尺度下的识别,将深度信息引入相似性评价方法,使其具有尺度不变性。改进后的相似性评价公式为:

式中:

ori(o,so(c′o,r′))—在参考图像o位置点处的梯度方向的弧度表示;

ori(i,t)—在输入图像i位置t点处的梯度方向的弧度表示;

表示以为中心,τ为邻域的区域。

p—是位置r的列表;

t=(o,p)—表示物体的模板。

d(c′x)—表示在点c′x处的深度值。

d(c'o)—表示离线训练时参考图像o相机到正多面体中心的距离。

为了避免直接进行max运算符的运算,本文采用与line-mod相同的方法,对max操作进行简化。同时预先计算响应图,建立查阅表,并将内存线性化以便并行处理,从而实现对相似性评价过程的加速。

从式4中不难看出,由于在离线训练时物体距离相机每一个采样位置都已知,所以根据在线深度信息就能获取到模板的缩放尺度。

采用原始line-mod的相似性评价方法无法完成零件识别。采用本文改进的相似性评价方法,首先利用深度信息获取模板的缩放尺度,然后通过滑动窗口在输入图像上进行滑动。在滑动窗口内,将位置列表p处的所有相似度相加,从而确定出输入图像和参考图像之间的相似度,完成物体的识别过程。

6.待装配对象位姿精确计算

待装配零件位姿估计示意图如图2所示。通过当前视角采集图像与离线训练的物体模板库的匹配,便可识别出与当前视角采集的图像最相似的参考视图(后文称之为关键帧)和其对应的训练时标记的相机3d位姿。然后通过关键帧周边投影视图的空间拓扑关系,确定出相机在待装配对象坐标系下的位置和姿态。最后通过坐标系变换,求解待装配对象在相机坐标系下的位姿。本文将对相机相对零件的姿态和位置分别求解。

7.姿态角估计

在模板训练时,相机光轴始终通过零件坐标系原点。但在线匹配时,相机光轴可能出现在空间任意位置。本部分将首先对相机相对关键帧的姿态角进行求解。图3为相机小孔成像模型,光线从p(x,y,z)经过透镜c在像平面呈倒立像p(x,y),为了简化相机姿态角计算过程,我们假设一个虚拟成像平面。

图4为在线匹配过程中待装配零件与相机坐标系的关系示意图,图中c-x1y1z1当前相机的坐标系,z1为相机的光轴,与虚拟像平面1相交于点o1,c-x2y2z2为零件训练时相机坐标系,z2为零件训练时相机的光轴,它穿过零件坐标系原点,与像平面1和像平面2分别交于点o1'和o1',z2与x1z1平面夹角为θ2,在x1z1平面的投影为直线cp,cp与z1轴的夹角为θ1。变换相机坐标系使其光轴从z2旋转到z1,过程如下:相机绕x2旋转θ2使得光轴cz2与cp重合,再绕y2旋转θ1使得光轴z2与z1重。相机的变换矩阵r为:

为了简化θ1和θ2的计算,光轴cz1和光轴cz1与虚拟像平面1的交点关系如图5所示,其中xc表示在o1x方向上的距离,yc表示在o1y方向上的距离,则θ1和θ2分别为:

除了相机光轴偏转带来的姿态角变化外,还有零件绕自身坐标系旋转带来的相机相对姿态角的变化。本文采用对零件绕自身坐标系旋转带来的姿态角变化进行计算。

8.位置估计

上文不仅能识别出待装配物体,还能识别出与当前视角输入图像接近的训练图像对应的相机位置和姿态。但是受采集投影图像时视点密度的影响,相机位置并不一定出现在交点位置,更多的情况是出现在空间任意位置。为了更加准确的获得相机位置,本文采用插值方法进行精确计算。

首先采用方法将“轴侧”视图变成“轴向”视图并提取形状描述符d',然后将其与数据库中的所有训练视图的形状描述符n(d,k)进行匹配,其中,dk(k=1,2,3...n)表示训练图像形状描述符集合,k表示采样视点之间的拓扑关系。从而获得与d'最相似的对应网格上的图相描述符di及其相邻节点形状描述符dj。位置0是当前相机视角所在位置,点4位置是与当前视角图像形状描述符最相似离线训练视角位置,其他位置点为相邻节点,则相机位置目标优化函数为:

这里通过计算可以得到权因子ωi及ωij的值,优化计算出来的权因子ωi及ωij是表示离线零件形状描述符库n中与当前输入轴向视图描述符d'最接近的图像描述符di的重心坐标数据,我们再根据权因子ωi及ωij的值来完成相机位置的插值计算,ωi和ωij是加权因子的形状描述符,因此公式(8)能被写成:

这里可以得到加权因子ωi和ωij。最后,根据权因子ωi及ωij的值来完成相机在正二十面体表面位置的插值计算。

在前面的步骤中,我们始终假设摄像机在正二十面体表面进行移动,但在实际中摄像机可能会出现在任意位置。

此时摄像机位置可以通过将经过正二十面体表面位置0的点进行平移获得,其平移位置可以通过深度相机获得。因此摄像机的最终位置和姿态可以用公式(10)表示:

p=fm(do,r,t,mr)(12)

式中:

do——深度传感器获取的o点的深度信息,

r,t——深度传感器和颜色传感器的相对位置,可以通过标定获得。

mr——是相机相对零件的姿态。

以上这些实施例应理解为仅用于说明本发明而不用于限制本发明的保护范围。在阅读了本发明的记载的内容之后,技术人员可以对本发明作各种改动或修改,这些等效变化和修饰同样落入本发明权利要求所限定的范围。

当前第1页1 2 
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1