首页 / 专利库 / 门,大门和窗户 / 框架 / 一种带双轴转位机构捷联惯导的三位置自对准方法

一种带双轴转位机构捷联惯导的三位置自对准方法

阅读:741发布:2023-02-06

专利汇可以提供一种带双轴转位机构捷联惯导的三位置自对准方法专利检索,专利查询,专利分析的服务。并且本 发明 公开了一种带双轴转位机构捷联惯导的三 位置 自对准方法,包括以下步骤:1)建立改进的卡尔曼滤波模型;2)可观测性分析;3)可观测组合的估计收敛速度分析;4)三位置自对准;本发明通过改进模型可以使捷联惯导系统的等效北向陀螺漂移和等效天向 加速 度计 的零偏快速收敛;提出的带双轴转位机构捷联惯导三位置自对准方法,适用于捷联惯导系统在初始任意 姿态 角 条件下的快速自对准;本发明可以有效缩短对准时间,提高带双轴转位机构的捷联惯导系统的自对准 精度 。,下面是一种带双轴转位机构捷联惯导的三位置自对准方法专利的具体信息内容。

1.一种带双轴转位机构捷联惯导的三位置自对准方法,其特征在于,包括以下步骤:
步骤1,建立改进的卡尔曼滤波模型:
以东北天坐标系作为导航坐标系,其中导航坐标系用n表示、机体坐标系用b表示、惯性坐标系用i表示;ωie中的下标表示地球坐标系相对于惯性坐标系的旋转速度,εb表示等效陀螺漂移ε在b系内的投影, 表示等效加速度计的零偏 在b系内的投影;速度误差向量δVn、姿态误差向量φ、导航坐标系下的陀螺常值漂移εn和加速度计零偏 其形式如下所示:
在自对准过程中,载体无移动,捷联惯导系统的速度误差方程为:
其中, Fvφ=(fn×),“×”表示叉乘运算;
捷联惯导系统的姿态误差方程为:
其中, RM、RN分别为经纬圈和卯酉圈半
径,L、h分别为纬度和高度;
将等效陀螺的漂移εb和等效加速度计零偏 视为随机常值误差,选取状态变量为:
同时以速度误差δVn作为外界量测,则系统模型表示如下:
其中 H=[I 0 0 0],I为3阶单位矩阵,0为3阶零矩阵;V为由干扰
引起的量测噪声,可以近似看作白噪声序列,假设其方差为R;
步骤2,可观测性分析:
依据连续时间线性定常系统可观测性分析相关理论,在单一位置对式(3)所表示的系统进行可观测性分析,构造可观测性矩阵Q=[HT (HF)T … (HF11)T]T,令Q1=[HT (HF)T … (HF4)T]T;
其中
令 其中N21=-Fvv,N31=-Fφv,
对式(4)左乘矩阵N进行行变换,得:
由式N52HF+N54HF3+HF4=0可得,HF4由HF和HF3线性表示,因此HFk均由HF和HF3线性表示,其中k=4,5…11;因此rank(Q)=rank(Q1)=rank(NQ1);
由式(5)可得,rank(NQ1)=9,因此得出结论:在任一固定位置,式(3)表示的系统可观测矩阵的秩为9,不满足完全可观测的条件;但是存在可观测组合,现加以分析:
记 则Y = Q 1 X ,令 其 中
其中
则式(3)表示的系统的可观测组合状态如下:
步骤3,可观测组合的估计收敛速度分析:
由式(7)对式(8)表示的可观测组合的收敛速度进行分析;根据式(7)可知,MNY的前9阶对应可观测组合的量测,因此,仅需对MNY矩阵的前9行进行分析,过程如下:
[MN]1:9Y=[MNQ1]1:9X=PX  (9)
式中,下标1:9表示矩阵的前9行;
由式(10)可知,式(9)的前3行对应状态量可由Z直接得到,第4~6行对应的可观测状态组合可由Z和Z一阶导得到,第7~9行对应的可观测状态组合分析如下:
NY=NQ1X  (11)
由式(12)可得:
即 方程两边同时左乘矩阵
得 取其前两行记为:
取式(15)的前两行,记作:
第7~8行对应的可观测状态组合由Z、 和 得到,第9行对应的可观测状态组合由Z、和Z(3)得到;需要量测量阶次信息的阶次越高对应状态所需估计的时间就越长,因此得出如下结论:除速度误差外状态组合外,(PX)4:6的收敛速度最快,(PX)7:8次之,而(PX)9最慢;
步骤4,三位置自对准:
由式(8)得,若εE、 和 能够快速准确辨识出,即可求解初始失准角φE,φN,φU;双轴捷联惯导系统内框架能够沿惯导系统的Z轴转动,外框架沿惯导系统的Y轴转动;在初始对准时,首先控制内框架绕IMU的Z轴转动-90°,然后控制外框架绕IMU的X轴转动-90°,作为初始对准的第一个位置;按照式(9)在此位置进行改进模型的卡尔曼滤波精对准;
其中θ,γ和ψ分别为载体的俯仰横滚和航向角;εx,εy和εz分别为x,y,z方向的等效陀螺漂移, 和 分别为x,y,z方向的等效加速度计零偏;
控制外框架绕IMU的X轴转动90°,作为初始对准的第二个位置;按照式(9)在此位置进行改进模型的卡尔曼滤波精对准;
控制外框架绕IMU的Z轴转动90°,作为初始对准的第三个位置;按照式(9)在此位置进行改进模型的卡尔曼滤波精对准;
在三个位置分别对等效北向陀螺漂移 及 和等效天向加速度计零偏
及 进行辨识,根据式(24)和式(25),得:
其中
若矩阵A满秩,此时εx、εy和εz由式(28)唯一确定:
此时, 由式(29)得到:
其中C=[cosγcosψ-sinγsinθsinψ -cosθsinψ sinγcosψ+cosγsinθsinψ];
若矩阵B满秩,此时 和 由式(30)唯一确定:
此时, 和 由式(31)得到:
其中
在第三位置进行对准,按照估计出来的 和 代入式(9)估计的PX中,按照式(32)进行计算,得到对准结束时刻的失准角φE、φN和φU;
2.根据权利要求1所述的带双轴转位机构捷联惯导的三位置自对准方法,其特征在于,所述步骤2中,若采用单位置对准方法时,φE的对准误差为 φN的对准误差为 φU的对准误差为 根据步骤4估计出的 和 对φE、φN和φU的
结果进行补偿,计算出对准结束时刻的失准角φE、φN和φU,然后根据式(33)对第三位置的姿态矩阵进行更新
即为经过误差补偿后的最终姿态矩阵。

说明书全文

一种带双轴转位机构捷联惯导的三位置自对准方法

【技术领域】

[0001] 本发明属于捷联惯性导航领域,提供一种带双轴转位机构捷联惯导的三位置自对准方法。【背景技术】
[0002] 现代战争对战略战术导弹提出了越来越高的要求,为提高导弹武器系统的反应速度和命中精度,对惯性导航系统提出了更高的要求。提高惯性传感器精度成本高、周期长、难度大,受国内基础工业薄弱和国外对华武器禁运的影响,短期内难以取得突破性进展。而通过对模型的改进和算法的合理设计可以有效的提高惯导系统的精度。
[0003] 与平台式惯性导航系统相比,捷联惯导系统具有反应快、成本低、信息全、可靠性高的优点,但同时捷联惯导系统的精度,为提高捷联惯导系统的精度,避免反复将捷联惯导系统从载体上拆下进行参数标定,进行系统分解维护,节约使用成本。捷联惯导增加双轴转位机构,使捷联惯导系统的免拆弹自检、标定和维护,有利于提高导弹武器装备的可操作性和机动性,对我国武器装备发展具有重要的使用价值和意义。
[0004] 捷联惯导系统自对准,由惯性敏感元件的输出信息,经过导航计算机按照对准程序实时运算,不断地将数学平台变换到能精确描述理想载体坐标系到地理坐标系的方向余弦阵。传统的捷联式惯导系统由于无转位机构,自对准为单位置自对准,其本质是将重加速度的方向做为平对准的基准,其对准精度主要取决于两个加速度计的精度。用于方位对准本质是将地球速率的北向分量方向做为方位对准的基准,其对准精度主要取决于东向陀螺漂移。带有单轴转位机构的捷联惯导系统能够在姿态接近水平的条件下,对等效东向陀螺的漂移进行辨识,进而提高自对准的精度。对双轴转位机构的捷联惯导系统,由于存在转位机构,可以通过旋转转位机构进行多位置对准,通过对自对准误差模型进行改进,设计三位置的快速自对准方法,可有效提高带双轴转位机构的捷联惯导系统的自对准精度。【发明内容】
[0005] 本发明的目的在于克服现有单位置自对准和双位置自对准技术的不足,提供一种双轴转位机构捷联惯导的三位置自对准方法,该方法通过改进自对准误差模型,设置双轴转位机构的旋转次序和位置,以及设计三位置自对准算法,有效提高了自对准精度。
[0006] 为达到上述目的,本发明采用以下技术方案予以实现:
[0007] 一种带双轴转位机构捷联惯导的三位置自对准方法,包括以下步骤:
[0008] 1)对以速度误差为量测量的卡尔曼滤波模型进行可观测性分析,对滤波模型进行改进,变换出新的可观测状态组合,并对其可观测性和收敛速度进行分析;
[0009] 2)在初始位置按照转位机构先绕内框架转-90°,再绕转位机构的外框架转动-90°,作为捷联惯导系统自对准的第一位置,采集惯导系统的陀螺和加速度计数据。根据记录的 到 时刻陀螺和加速度计的输出,进行惯性系粗对准和改进模型的精对准,完成等效北向陀螺漂移 和天向加速度零位 的测漂;
[0010] 3)绕外框架转90°得到第二位置,在第二位置静止采集陀螺和加速度计的输出。根据记录的 到 时刻陀螺和加速度计的输出,进行惯性系粗对准和改进模型的精对准,完成等效北向陀螺漂移 和天向加速度零位 的测漂。
[0011] 4)绕内框架转90°得到第三位置,在第三位置静止采集陀螺和加速度计的输出。根据记录的 到 时刻陀螺和加速度计的输出,进行惯性系粗对准和改进模型的精对准,完成等效北向陀螺漂移 和天向加速度零位 的测漂。
[0012] 5)依据三个位置姿态矩阵和 和 测漂结果进行反解得到εx、εy和εz。利用第三位置的姿态矩阵和εx、εy和εz计算 依据三个位置姿态矩阵和 和 测漂结果进行反解得到▽x、▽y和▽z。利用第三位置的姿态矩阵和▽x、▽y和▽z计算 和[0013] 6)依据改进模型计算对准结束时刻的失准角φE、φN和φU。进行一次修正后得到最终的精对准结果。
[0014] 与现有技术相比,本发明具有以下有益效果:
[0015] 本发明通过改进模型可以使捷联惯导系统的等效北向陀螺漂移和等效天向加速度计的零偏快速收敛;提出的带双轴转位机构捷联惯导三位置自对准方法,适用于捷联惯导系统在初始任意姿态角条件下的快速自对准;本发明可以有效缩短对准时间,提高带双轴转位机构的捷联惯导系统的自对准精度。【附图说明】
[0016] 图1为零位定时坐标系与两框架轴编排示意图;
[0017] 图2为初始位置及转位过程示意图;
[0018] 图3为第一位置等效北向陀螺漂移估计曲线;
[0019] 图4为第二位置等效北向陀螺漂移估计曲线;
[0020] 图5为第三位置等效北向陀螺漂移估计曲线;
[0021] 图6为两种方法俯仰角误差对比曲线;
[0022] 图7为两种方法滚转角误差对比曲线;
[0023] 图8为两种方法航向角误差对比曲线。【具体实施方式】
[0024] 下面结合附图对本发明做进一步详细描述:
[0025] 参见图1-图8,本发明带双轴转位机构捷联惯导的三位置自对准方法,包括以下步骤:
[0026] 步骤1,建立改进的卡尔曼滤波模型
[0027] 以东北天坐标系作为导航坐标系,其中导航坐标系用n表示、机体坐标系用b表示、惯性坐标系用i表示。ωie中的下标表示地球坐标系相对于惯性坐标系的旋转角速度,εb表示等效陀螺漂移ε在b系内的投影,▽b表示等效加速度计的零偏▽在b系内的投影。速度误差向量δVn、姿态误差向量φ、导航坐标系下的陀螺常值漂移εn和加速度计零偏▽n,其形式如下所示:
[0028]
[0029] 在自对准过程中,载体无移动,捷联惯导系统的速度误差方程为:
[0030]
[0031] 其中, Fvφ=(fn×),“×”表示叉乘运算。
[0032] 捷联惯导系统的姿态误差方程为:
[0033]
[0034] 其中, RM、RN分别为经纬圈和卯酉圈半径,L、h分别为纬度和高度。
[0035] 将等效陀螺的漂移εb和等效加速度计零偏▽b视为随机常值误差,选取状态变量为:X=[δVn φ ▽n εn]Τ,同时以速度误差δVn作为外界量测,则系统模型可表示如下:
[0036]
[0037] 其中 H=[I 0 0 0],I为3阶单位矩阵,0为3阶零矩阵。V为由晃动等干扰引起的量测噪声,可以近似看作白噪声序列,假设其方差为R。
[0038] 步骤2,可观测性分析
[0039] 依据连续时间线性定常系统可观测性分析相关理论,在单一位置对式(3)所表示的系统进行可观测性分析,构造可观测性矩阵Q=[HΤ (HF)Τ …  (HF11)Τ]Τ,令Q1=[HΤ (HF)Τ … (HF4)Τ]Τ。
[0040]
[0041] 其中
[0042] 令 其中N21=-Fvv,N31=-Fφv,
[0043] 可以证明,当L≠0°或90°时A2可逆,当L=0°或90°时用类似的方法也能得到式(6)同样的结果。
[0044] 对式(4)左乘矩阵N进行行变换,得:
[0045]
[0046] 由式(5)中第5行N52HF+N54HF3+HF4=0可得,HF4可由HF和HF3线性表示,因此HFk(k=4,5…11)均可由HF和HF3线性表示。因此rank(Q)=rank(Q1)=rank(NQ1)。
[0047] 由式(5)可得,rank(NQ1)=9,因此可以得出结论:在任一固定位置,式(3)表示的系统可观测矩阵的秩为9,不满足完全可观测的条件。但是存在可观测组合,现加以分析。
[0048] 记 则Y=Q1X,令 其中
[0049]
[0050] 其中
[0051]
[0052] 则式(3)表示的系统的可观测组合状态如下:
[0053]
[0054] 至此可观测性分析完毕,从分析结果可见,除速度误差状态和等效天向加速度计零偏直接可观测以外,其余的状态均为各种形式的组合,其中水平姿态误差角和等效东向及等效北向加速度计零偏为一可观测状态组合,航向误差角和等效东向陀螺漂移以及等效东向加速度计零偏为一可观测状态组合,等效北向陀螺漂移和等效北向加速度计零偏为一可观测状态组合,等效天向陀螺漂移和等效北向加速度计零偏为一可观测状态组合。可以看出,采用单位置对准方法时,φE的对准误差为 φN的对准误差为 φU的对准误差为 由此可见,除等效天向加速度计的零偏以外其余的状态估计精度均受制于相应的惯性器件的精度。
[0055] 步骤3,可观测组合的估计收敛速度分析
[0056] 由式(7)可以对式(8)表示的可观测组合的收敛速度进行分析。根据式(7)可知,MNY的前9阶对应可观测组合的量测,因此,仅需对MNY矩阵的前9行进行分析,过程如下:
[0057] [MN]1:9Y=[MNQ1]1:9X  (9)
[0058] 式中,下标1:9表示矩阵的前9行。
[0059]
[0060] 由式(10)可知,式(9)的前3行对应状态量可由Z直接得到,第4~6行对应的可观测状态组合可由Z和Z一阶导得到,第7~9行对应的可观测状态组合分析如下
[0061] NY=NQ1X  (11)
[0062]
[0063] 由式(12)可得
[0064]
[0065] 即 方 程 两 边 同 时 左 乘 矩 阵可得 取其前两行记为
[0066]
[0067]
[0068] 取式(15)的前两行,记作
[0069]
[0070] 第7~8行对应的可观测状态组合可由Z、 和 得到,第9行对应的可观测状态组合可由Z、 和Z(3)得到。需要量测量阶次信息的阶次越高对应状态所需估计的时间就越长。因此可以得出如下结论:除速度误差外状态组合外,(PX)4:6的收敛速度最快,(PX)7:8次之,而(PX)9最慢,至此收敛速度分析完毕。
[0071] 步骤4,三位置自对准方法
[0072] 由式(8)可得,若εE、▽E和▽N能够快速准确辨识出,即可求解初始失准角φE,φN,φU。因此,自对准的问题转化为如何快速准确辨识εE、▽E和▽N。在实际对准过程中,与▽E,▽N相比,εE的影响较大。为减小εE的的影响,双位置中提出绕方位轴正向转90°再次进行对准。在该位置东向等效陀螺漂移的计算值 其中 为第二个位置处等效的北向陀螺漂移。但存在两个问题,第一:若初始对准时载体不在水平面,即存在一定的俯仰角或者横滚角,如导弹发射车处于斜坡或者一侧轮胎处于台阶上,此时第二位置的等效北向陀螺漂移不等于第一位置的等效东向陀螺漂移,此时,第一位置的等效东向陀螺漂移不可测;第二,两位置测漂法中,无法对与▽E和▽N进行辨识,水平对准的精度无法进一步提高。基于以上两个问题,本专利在改进对准误差模型的基础上,提出了适用于带双轴转位机构捷联惯导系统的三位置对准法。
[0073] 如图1所示,双轴捷联惯导系统内框架可以沿惯导系统的Z轴转动,外框架沿惯导系统的Y轴转动。在初始对准时,首先控制内框架绕IMU的Z轴转动-90°,然后控制外框架绕IMU的X轴转动-90°,作为初始对准的第一个位置。按照式(9)在此位置进行改进模型的卡尔曼滤波精对准。
[0074]
[0075]
[0076]
[0077] 其中θ,γ和ψ分别为载体的俯仰,横滚和航向角。εx,εy和εz分别为x,y,z方向的等效陀螺漂移,▽x,▽y和▽z分别为x,y,z方向的等效加速度计零偏。
[0078] 控制外框架绕IMU的X轴转动90°,作为初始对准的第二个位置。按照式(9)在此位置进行改进模型的卡尔曼滤波精对准。
[0079]
[0080]
[0081]
[0082]
[0083] 控制外框架绕IMU的Z轴转动90°,作为初始对准的第三个位置。按照式(9)在此位置进行改进模型的卡尔曼滤波精对准。
[0084]
[0085]
[0086]
[0087]
[0088] 在三个位置处均可对等效北向陀螺漂移εN和等效天向加速度计零偏▽U进行辨识,即 和 ▽x、▽y、和▽z可以得到。根据式(24),可得:
[0089]
[0090]
[0091] 其中
[0092]
[0093]
[0094] 若矩阵A满秩,此时εx、εy和εz可由式(28)唯一确定。
[0095]
[0096] 此时, 可由式(29)得到。
[0097]
[0098] 其中C=[cosγcosψ-sinγsinθsinψ -cosθsinψ sinγcosψ+cosγsinθsinψ]。
[0099] 若矩阵B满秩,此时▽x、▽y和▽z可由式(30)唯一确定。
[0100]
[0101] 此时, 和 可由式(31)得到。
[0102]
[0103] 其中
[0104] 在第三位置进行对准,按照估计出来的 和 对估计结果进行补偿,即可计算出对准结束时刻的失准角φE、φN和φU。
[0105] 步骤5,殊情况的讨论
[0106] εx、εy、εz、▽x、▽y、和▽z的可观测性问题等价于式(26)和式(27)的求解问题,为此,先引入两个定理。
[0107] 定理1:齐次线性方程组的系数矩阵秩rank(A)=n,方程组有唯一解。齐次线性方程组的系数矩阵秩rank(A)<n,方程组有无数多解。
[0108] 定理2:A是n×n矩阵,则
[0109] a)rank(A)=n的充要条件为|A|≠0,称A为满秩矩阵;
[0110] b)rank(A)<n的充要条件为|A|=0。
[0111] 式(26)能够求解陀螺漂移的前提是矩阵A满秩,因此,求解使得|A|=0时的初始姿态角θ0、γ0和ψ0,出当初始姿态角为θ0、γ0和ψ0时,εx、εy和εz不完全可观测。同理可以求出▽x、▽y、和▽z不可观测的姿态角初值。
[0112] 在导弹武器系统的实际的应用过程中,导弹发射车的初始姿态角受到实际应用环境的限制,俯仰角、横滚角在一定范围内变化,根据实际使用情况,取-45°≤θ≤45°,-45°≤γ≤45°,-180°<ψ≤180°。求解|A|=0得:当θ0=0°,ψ0=0°(或ψ0=180°)时,|A|=0,rank(A)<3。此时εx、εy和εz不完全可观测。当θ0=0°,ψ0=0°时,初始位置时Y向指北,第一位置和第二位置时X向指南,第三位置时Y向指北,此时,X向和Y向的陀螺漂移εx和εy均可通过等效北向陀螺测漂得到,Z向陀螺的漂移εz不可观测。当θ0=0°,ψ0=180°时,初始位置时Y向指南,第一位置和第二位置时X向指北,第三位置时Y向指南,此时,X向和Y向的陀螺漂移εx和εy均可通过等效北向陀螺测漂得到,Z向陀螺的漂移εz不可观测。此时,由分析可知,陀螺漂移εx和εy可以观测,εz不可观测,只能设置为εz=0。
[0113] 当θ0≠0°或ψ0≠0°和180°时,|A|≠0,rank(A)=3。此时εx、εy和εz可由式(28)唯一确定。
[0114] 求解|B|=0得:当θ0=0°和γ0=0°同时成立时,|B|=0,rank(B)<3。此时▽x、▽y和▽z不完全可观测。当θ0=0°,γ0=0°时,初始位置时Z向指天,第一位置Y向指地,第二位置和第三位置时Z向指天,此时,Y向和Z向的加速度计零偏▽y、和▽z均可通过等效天向加速度计测漂得到,X向加速度计的零偏▽x不可观测。此时,由分析可知,加速度计零偏▽y和▽z可以观测,▽x不可观测,只能设置为▽x=0。
[0115] 当θ0≠0°或γ0≠0°时,|B|≠0|,rank(B)=3。此时▽x、▽y和▽z可由式(30)唯一确定。
[0116] 【验证】
[0117] 为验证分析的正确性,进行了仿真验证,陀螺仪的等效漂移取为[0.03 0.03 0.03]°/h,加速度计等效零偏取为[50 50 50]μg,采样频率为200Hz,捷联惯导系统的初始姿态为att_0=[15° 10° 30°]。开始对准时捷联惯导系统的台体先绕内框Z轴转动-90°,再绕外框X轴转动-90°作为对准的第一位置,转动过程10s,在第一位置采集数据90s。然后绕外框X轴转动90°度作为对准的第二位置,转动过程10s,在第二位置采集数据90s。再绕内框Z轴转动90°作为对准的第三位置,转动过程10s,在第三位置采集数据90s,对准时间共计
300s。分别采用单位置卡尔曼滤波对准法、双位置卡尔曼滤波对准法和本专利提出的三位置自对准方法进行对准,图3~图5为三个位置上的等效北向陀螺漂移和等效加速度计天向零偏的估计曲线。表1为利用不同方法进行对准的试验结果,图6~图8为不同方法进行姿态估计的曲线。由图中可以看出,本专利提出的三位置自对准方法航向角收敛速度和精度明显优于另两种方法。
[0118] 表1 不同对准方法的结果
[0119]
[0120] 通过由仿真结果可以看出,在相同的时间内,本专利提出的可以有效的提高捷联惯导系统的初始对准精度,是双位置对准方法误差的10%,是单位置对准方法误差的5%。使用本专利提出的可以使捷联惯导系统的航向对准误差有效降低。
[0121] 以上内容仅为说明本发明的技术思想,不能以此限定本发明的保护范围,凡是按照本发明提出的技术思想,在技术方案基础上所做的任何改动,均落入本发明权利要求书的保护范围之内。
高效检索全球专利

专利汇是专利免费检索,专利查询,专利分析-国家发明专利查询检索分析平台,是提供专利分析,专利查询,专利检索等数据服务功能的知识产权数据服务商。

我们的产品包含105个国家的1.26亿组数据,免费查、免费专利分析。

申请试用

分析报告

专利汇分析报告产品可以对行业情报数据进行梳理分析,涉及维度包括行业专利基本状况分析、地域分析、技术分析、发明人分析、申请人分析、专利权人分析、失效分析、核心专利分析、法律分析、研发重点分析、企业专利处境分析、技术处境分析、专利寿命分析、企业定位分析、引证分析等超过60个分析角度,系统通过AI智能系统对图表进行解读,只需1分钟,一键生成行业专利分析报告。

申请试用

QQ群二维码
意见反馈