Effects of stochastic models on BDS DCB estimation and ionospheric modeling
-
摘要: 差分码偏差(DCB)作为电离层建模和导航定位中一项重要的误差源，对其进行估计求解至关重要. 为提高北斗卫星导航系统(BDS) DCB估计和电离层建模精度，提出了一种综合高度角、卫地距和测站纬度多因素的随机模型，并对比分析了不同随机模型对BDS DCB估计和电离层垂直总电子含量(VTEC)建模精度的影响. 结果表明：不同随机模型对卫星端DCB解算产生约0.2 ns差异. 相较于高度角随机模型，采用高度角、卫地距组合模型测站DCB估计精度平均提高0.13 ns，电离层建模精度提高了约0.2 TECU. 新提出的随机模型，在低纬度测站DCB解算精度上差于高度角模型和高度角、卫地距组合模型，但在高纬度测站DCB解算结果上更优，且对电离层VTEC建模精度提升效果明显，与前两种随机模型相比分别提升了0.88 TECU和0.68 TECU.
-
关键词:
- 差分码偏差(DCB) /
- 随机模型 /
- 北斗卫星导航系统(BDS) /
- 电离层模型 /
- 总电子含量(TEC)
Abstract: As an important error source in ionospheric modeling and navigation, so it is very important to estimate and solve it. In order to improve the accuracy of differential code bias (DCB) estimation and ionospheric modeling, this paper proposes a multi-factor stochastic model based on altitude angle, the distance between the station and the satellite and station latitude, and analyzes the influence of different stochastic models on the accuracy of DCB estimation and ionospheric vertical total electron content (VTEC) modeling. The results indicated that different stochastic models may produce about 0.2 ns difference in satellite DCB. Compared with altitude angle weighting method, the estimation accuracy of station DCB was improved by 0.13 ns on average, and the ionospheric modeling accuracy was improved by about 0.2 TECU when the combine model of the elevation angel and the distance from the satellite to the ground station was used. Under the new stochastic model, the DCB solution of low latitude stations was worse than the first two methods, but the DCB solution of high latitude stations was better, and the modeling accuracy of ionospheric VTEC was improved significantly, which was 0.88 TECU and 0.68 TECU respectively, compared with the first two stochastic models. -
0. 引 言
同一时刻不同频率或同一频率上不同测距信号在经过发射链路和接收链路中所产生时间延迟的差值称为差分码偏差(DCB)，是影响电离层总电子含量(TEC)监测和建模的主要误差源，也是精密单点定位(PPP)中必须考虑的误差项[1-2]. 目前常用的DCB估计方法主要有两种：1)与电离层参数同步估计[3]；2)先用经验模型改正电离层TEC，再解算DCB[4]. LI等[5]提出的IGG DCB方法克服了前者对测站过分依赖的缺点，WANG等[6]将该方法进一步拓展至多系统多频DCB上，并通过中国科学院评估研究中心（CAS）向用户发布. 考虑到利用载波相位平滑伪距提取观测量，易受平滑弧段长度和多路径效应影响，文献[7-8]提出了采用非差非组合PPP提取电离层观测值，进而解算DCB参数，实验验证了非组合PPP方法更利于电离层延迟提取、建模和预报. 随着低轨卫星的升空作业，采用低轨卫星估计DCB的方法也不断被提出[9].
北斗卫星导航系统(BDS)星座包含有倾斜地球同步轨道(IGSO)卫星、地球静止轨道(GEO)卫星、中圆地球轨道(MEO)卫星三种不同类型的卫星，且不同类型卫星的轨道高度、轨道精度和DCB解算残差存在差异. 另外，由于观测站所处纬度不同，其天顶(U)方向电离层变化和活跃程度也不尽相同. 目前关于DCB解算随机模型的研究相对较少，现有研究也大多基于随机模型或高度角模型展开[10-12]. 为了提高BDS DCB估计和电离层建模精度，提出了一种综合高度角、卫地距和测站纬度的随机模型，并对比分析了不同随机模型对BDS DCB估计以及电离层模型精度的影响.
1. DCB与电离层参数估计方法
通过双频无几何组合，可得到伪距和载波相位无几何观测量，考虑到伪距观测量噪声较大，本文采用载波相位平滑伪距的方式得到较高精度的电离层观测量，计算公式为
$$ \left\{ \begin{gathered} \mathop {\tilde P}\nolimits_4 = \frac{{\mathop {40.28(f}\nolimits_2^2 {{ - }}\mathop f\nolimits_1^2 )}}{{\mathop f\nolimits_1^2 \mathop f\nolimits_2^2 }}\mathop S\nolimits_{{\text{TEC}}} + c(\mathop D\nolimits^{\text{s}} + \mathop D\nolimits_{\text{r}} ) \\ \mathop S\nolimits_{{\text{TEC}}} = \mathop V\nolimits_{{\text{TEC}}} /\mathop M\nolimits_f (z) \\ \mathop M\nolimits_f (z){\text{ = }}\cos \left[ {\arcsin \left( {\frac{R}{{R + H}}\sin (\alpha z)} \right)} \right] \\ \end{gathered} \right.. $$ (1) 式中：
$ \mathop {\tilde P}\nolimits_4 $ 为经载波相位平滑后的伪距无几何观测量；$\mathop S\nolimits_{{\text{TEC}}}$ 、$\mathop V\nolimits_{{\text{TEC}}}$ 分别为倾斜方向TEC和U方向TEC；$ \mathop D\nolimits^{\text{s}} $ 、$\mathop D\nolimits_{\text{r}} $ 为卫星端$ {\text{DCB}} $ 和接收机端$ {\text{DCB}} $ ；s、r分别为卫星、接收机；$c$ 为光速；$ \mathop M\nolimits_f (z) $ 为投影函数，$z$ 为卫星高度角的余角；$R$ 为地球平均半径；$H$ 为电离层薄层高度，取值为506.7 km；$\alpha $ 为比例因子，取值为0.978 2.本文采用低阶球谐函数对U方向电离层TEC进行建模，表达式为
$$ \begin{split} \mathop V\nolimits_{{\text{TEC}}} \left( {\varphi ,s} \right) =& \sum\limits_{n = 0}^{\mathop n\nolimits_{\max } } {\sum\limits_{m = 0}^n {\mathop {\overline P }\nolimits_{nm} } } \left( {\sin\;\varphi } \right) \\& \left[ {\mathop a\nolimits_{nm} \cos \left( {ms} \right) + \mathop b\nolimits_{nm} \sin \left( {ms} \right)} \right] \end{split} . $$ (2) 式中：
$ \mathop n\nolimits_{\max } $ 为球谐函数最大展开度数；$\mathop {\overline P }\nolimits_{nm} \left( {\sin\;\varphi } \right)$ 为$ n $ 度$ m $ 阶的归化勒让德函数；$\varphi $ 为穿刺点的地理纬度；$s = \lambda - \mathop \lambda \nolimits_0 $ 为穿刺点的日固经度，其中$\lambda $ 为穿刺点的地理经度，$\mathop \lambda \nolimits_0 $ 为太阳的地理经度；$\mathop a\nolimits_{nm} $ 和$\mathop b\nolimits_{nm} $ 为球谐函数模型中的待求系数. 联立式(1)和式(2)，改写矩阵形式可得$$ \left\{ \begin{gathered} {\boldsymbol{L}} = {\boldsymbol{B}}\widehat {\boldsymbol{X}} \\ \widehat {\boldsymbol{X}} = \mathop {\left[ {\mathop {\boldsymbol{K}}\nolimits^{\text{T}} ,\mathop {\boldsymbol{R}}\nolimits^{\text{T}} ,\mathop {\boldsymbol{S}}\nolimits^{\text{T}} } \right]}\nolimits^{\text{T}} \end{gathered} \right. . $$ (3) 式中：
$ {\boldsymbol{L}} $ 为平滑的无几何观测量；$ {\boldsymbol{B}} $ 为参数向量的系数矩阵；$\widehat {\boldsymbol{X}}$ 为参数向量，其中包含有球谐函数的待求系数向量${\boldsymbol{K}}$ 和接收机端DCB参数向量$ {\boldsymbol{R}} $ 和卫星端DCB参数向量$ {\boldsymbol{S}} $ . 由于卫星DCB和接收机DCB参数列相关，需要引入一个零均值基准约束条件，为了便于与CAS DCB产品进行比较，本文采用的零均值基准约束与CAS解算DCB采用的基准约束保持一致，即$ {\boldsymbol{H}}\widehat {\boldsymbol{X}} = 0 $ ，$ {\boldsymbol{H}} $ 为约束向量，其中与卫星DCB对应元素为1，其他元素为0. 基于附有限制条件的间接平差，可求解得到未知参数$$ \widehat {\boldsymbol{X}} = \mathop {\left( {\mathop {\boldsymbol{B}}\nolimits^{\text{T}} {\boldsymbol{PB}} + \mathop {\boldsymbol{H}}\nolimits^{\text{T}} {\boldsymbol{H}}} \right)}\nolimits^{ - 1} \mathop {\boldsymbol{B}}\nolimits^{\text{T}} {\boldsymbol{PL}} . $$ (4) 式中，
${\boldsymbol{P}}$ 为观测量权阵. 将求解得到的电离层参数${\boldsymbol{K}}$ 代入式(2)中可得到相应的电离层垂直总电子含量(VTEC)值.2. 随机模型
由式(4)可知，电离层观测值的权阵对解算DCB、电离层球谐参数有着重要作用. DCB估计中几种常见的随机模型有：
1) 等权随机模型
等权随机模型假设各颗卫星在不同历元的电离层观测值误差符合相同的统计分布，即所占权重是大小相等且互不相关. 等权模型实际上并没有进行加权，即权矩阵等于单位阵
${\boldsymbol{E}}$ ，表达式可写为：${\boldsymbol{P}} = {\boldsymbol{E}}$ .2) 高度角随机模型
一般地，当卫星高度角较大时，受到大气延迟等误差的影响越小，且测站接收机受到多路径效应越弱，将会大大降低观测量中包含的观测噪声. 基于卫星高度角的权重计算公式[11]可写为
$$ \left\{ \begin{gathered} \mathop P\nolimits_e = \frac{{\mathop \sigma \nolimits_0^2 }}{{\mathop \sigma \nolimits^2 }} \\ \mathop \sigma \nolimits_0^2 = \mathop {(f + d)}\nolimits^2 \\ \mathop \sigma \nolimits^2 = \mathop {\left[ {0.05 + \frac{{0.02}}{{\mathop {\sin (z)}\nolimits^2 }}} \right]}\nolimits^2 \end{gathered} \right. . $$ (5) 式中：
$ \mathop P\nolimits_e $ 为基于高度角定权得到的权重；z为某一观测值对应的卫星高度角；$ \mathop \sigma \nolimits^2 $ 为该观测值的方差；$ \mathop \sigma \nolimits_0^2 $ 为单位权方差，其中，$ f $ 、$ d $ 为两个经验值常量，分别为5 cm和2 cm.3) 高度角、卫地距组合模型
BDS卫星星座包含IGSO、GEO、MEO三种不同类型的卫星，三种卫星所处的轨道高度、轨道精度、DCB解算残差都不尽相同. 为了更好区分不同类别卫星间的差异，基于卫星高度角和卫地距综合确权的权重计算公式为
$$ \mathop P\nolimits_{e{{ - }}\rho } = \frac{{\mathop P\nolimits_e }}{{\mathop \rho \nolimits^2 }} . $$ (6) 式中：
$ \mathop P\nolimits_{e - \rho } $ 为基于高度角卫地距定权得到的权重；$ \mathop P\nolimits_e $ 为按高度角随机模型计算得到的权重；$\;\mathop \rho \nolimits^2$ 为卫星和测站间的卫地距.在综合多个测站对参数进行估计时，由于各个测站所处纬度不同，U方向上电离层活跃程度在一天中不尽相同. 其中，低纬度地区电离层电子密度最高且变化最为剧烈. 比如赤道附近在当地时间正午时，上空的电离层变得最为活跃，其TEC含量达到最大值，且此时电离层建模较为复杂，存在的模型误差较大，应赋予更小的权重. 考虑到纬度因素，文中提出了一种综合卫星高度角、卫地距、测站纬度三因素的随机模型，具体权重公式为
$$ \mathop P\nolimits_{e - \rho - b} = \frac{{\mathop {\sin (B)P}\nolimits_e }}{{\mathop \rho \nolimits^2 }} . $$ (7) 式中：
$ \mathop P\nolimits_{e - \rho - b} $ 为基于高度角卫地距纬度定权方式获得的权重；$ B $ 为当前测站所处纬度；其他参数同上式.3. 数据处理及案例分析
3.1 数据处理
实验数据取自2021年年积日第175—184下MGEX (Multi-GNSS Experiment)测站网BDS C2I、C6I类型观测数据，在中国及周边区域选取了13个测站，具体测站位置分布如图1所示. 数据采样间隔为30 s，卫星截止高度角设置为10°. 本文采用4阶球谐函数模型模拟电离层VTEC的时空变化，并将DCB在一天内视为常数，和时间分辨率为2 h的球谐函数模型参数进行同步估计. 为分析不同随机模型对DCB估计、电离层建模的影响，实验分3种方案进行：方案1，基于高度角随机模型解算；方案2，基于高度角、卫地距组合模型解算；方案3，基于高度角、卫地距和纬度三因素模型解算.
3.2 DCB精度分析
在解算C2I-C6I类型DCB过程中，除C15、C17、C18和C31卫星不再提供服务外，其余卫星参加了解算，由于卫星数目众多，表1中仅列举了GEO卫星、MEO卫星、IGSO卫星中各5颗卫星在年积日第175—184下的卫星DCB平均值. 图2展示了在不同定权方式下解算出卫星DCB的均方根误差(RMSE)，本文以CAS发布的DCB产品作为真值进行比较.
表 1 三种方案下北斗C2I-C6I类型卫星DCB平均估值ns 卫星 PRN 方案1 方案2 方案3 GEO C01 −1.37 −1.34 −1.23 C02 4.25 4.22 4.36 C03 5.62 5.62 5.59 C04 8.82 8.85 9.06 C05 2.73 2.74 3.05 MEO C11 7.25 7.23 7.26 C12 6.71 6.67 6.73 C20 21.48 21.66 21.71 C24 7.98 7.94 7.81 C30 −8.58 −8.55 −8.55 IGSO C06 11.73 11.74 11.64 C07 18.10 18.10 18.02 C10 9.44 9.37 9.28 C38 4.63 4.63 4.70 C40 1.35 1.33 1.57 由表1可知，采用方案1和方案2解算出的DCB平均估值相差较为细微，除C20卫星两者差值在0.18 ns外，其他卫星二者差值皆小于0.1 ns；方案3和两者的差值较大，但也都在0.2 ns左右. 由图2可知，方案1和方案2下卫星DCB的RMSE值在除C20卫星存在较大差异外，其他卫星处两者RMSE波动较小，而方案3下统计的卫星DCB RMSE值和前两者之间的差距较大. 这说明额外添加卫地距定权，对卫星端DCB估计有着一定影响，但影响效果不太明显；采用方案3随机模型时，卫星DCB的平均估值和RMSE与前两个方案可能会存在着约0.2 ns的偏差.
表2中统计了年积日第175—184下各测站DCB估值与CAS产品真值的平均差值，图3展示了各测站DCB的RMSE情况，由于LCK3、SIN1、USUD测站在实验时段存在着CAS发布真值丢失情况，下文统计结果并未包含该三个测站. 由表2和图2可以看出，在方案2下测站DCB的平均差值和RMSE普遍优于方案1，GUAM测站除外. 对比方案1，采用方案2平均每个测站DCB估计精度提高了0.13 ns. 这说明相较于单独的高度角定权，基于高度角-卫地距的定权方式解算出的接收端DCB更加接近理论真值. 对比方案3和前两种方案可以发现，低纬度测站DCB结果不论是在平均差值还是RMSE上都会更差，但高纬度测站DCB结果会得到较好地提升，即基于高度角-卫地距-纬度定权方式对于较高纬度测站求解DCB优于方案1和方案2，对低纬度测站求解DCB结果要差于前两者.
表 2 各测站DCB与CAS DCB的平均差值ns 测站 方案1 方案2 方案3 SGOC 0.62 0.51 1.12 IISC 0.59 0.59 1.44 GUAM 0.83 1.03 1.22 CUSV 0.75 0.71 0.71 PTGG 1.76 1.60 1.88 WUH2 0.92 0.61 0.39 GAMG 0.60 0.56 0.66 POL2 2.51 2.39 2.15 URUM 3.27 3.04 2.80 ULAB 3.38 2.90 2.69 3.3 电离层VTEC精度分析
通过上文解算出的球谐函数模型参数，可实现对电离层VTEC的解算，生成区域电离层产品文件. 本文选取的格网数据空间覆盖纬度范围为5°N~45°N，分辨率为2.5°，经度范围为80°E~140°E，分辨率为5°. 图4展示了本文选取格网点具体分布情况. 为分析不同随机模型对电离层模型精度的影响，以国际GNSS服务(IGS)研究中心发布的IGSG电离层格网产品作为真值，此处采用VTEC模型值与真值的RMSE来反映三种方案下模型解算的精度状况，计算公式为
$$ \mathop R\nolimits_j = \sqrt {\frac{{\mathop {\displaystyle\sum\limits_{i = 1}^n {\left( {\mathop v\nolimits_{i,j} - \mathop {\tilde v}\nolimits_{i,j} } \right)} }\nolimits^2 }}{{n - 1}}} . $$ (8) 式中：
$ \mathop v\nolimits_{i,j} $ 是由$ j $ 时刻电离层参数计算得到的某一格网点VTEC值；$ \mathop {\tilde v}\nolimits_{i,j} $ 是$ j $ 时刻该格网点对应的IGSG发布的真值；$ n $ 为所选用格网点总数.图5中统计了在年积日第175—179中电离层VTEC在对应时段的RMSE. 从图5中可以明显看出，方案2和方案3解算出的RMSE皆小于方案1，其中方案1和方案2的RMSE较为接近，两者在5天内的平均RMSE分别为4.54 TECU、4.35 TECU，方案3和方案1、方案2的RMSE偏差较大，其平均RMSE为3.66 TECU，且在一天中起始时刻和结束时刻三者RMSE差距较小，中间时段内存在明显差距. 这说明了相较于高度角随机模型，基于高度角-卫地距，高度角-卫地距-纬度这两种随机模型可以提高电离层VTEC建模精度，其中综合高度角、卫地距、测站纬度的随机模型最适合电离层建模.
为了能够直观反映三种方案下VTEC建模的准确性，将得到的区域VTEC产品与IGS发布的电离层VTEC产品作差取绝对值，得到的电离层残差分布图如图6所示. 图6展示了三种方案下5°N~45°N，80°E~140°E区域一天内VTEC残差分布状况，自左向右依次为方案1、方案2、方案3下残差分布. 限于篇幅原因，在此只展示了4个时段的VTEC分布图. 从图6中可以看出三种模型下，电离层残差较大的地方主要位于建模的边缘地带，可能是由于边缘地带穿刺点分布较少，观测量稀疏导致的；从协调世界时(UTC) 00：00和22：00时刻对应的分布图中，可以看出三者残差分布差距较小，不易发现三者差距，这印证了图5中年积日175起始时刻和结束时刻三者RMSE差距较小. 从UTC 08：00和16：00时刻的分布图来看，方案2与方案1得到的电离层残差分布相对接近，两者在偏差较小区域差距不大，但方案1残差较大区域所占空间比方案2更大，即相较方案2，方案1存在更大偏差. 方案3得到的残差分布图中靠近残差为0的区域比方案1、方案2更多，这说明方案3得到的VTEC更为接近IGS提供的电离层产品.
4. 结 论
在电离层反演与导航定位中，DCB是必须要考虑的误差项. 现阶段关于DCB估计的方法大多基于函数模型进行研究，随机模型研究较少. 本文提出了一种基于高度角、卫地距、测站纬度三种因素的定权方式，并对比分析了三种不同定权方案对DCB估计以及电离层建模的影响. 得到的结论如下：
1)不同随机模型下卫星DCB解算结果存在差异，差值一般约在0.2 ns.
2)和高度角模型相比，高度角、卫地距组合随机模型得到的结果更接近IGS产品，解算精度平均提升了0.13 ns. 在新提出的随机模型下，低纬度测站DCB解算上差于前两种模型，但在高纬度测站DCB解算结果上更优.
3)相较于高度角模型，高度角、卫地距组合模型对电离层建模有一定的提升，格网点RMSE降低了约0.2 TECU；新提出的模型对电离层建模精度提升尤为显著，格网点RMSE降低了0.88 TECU.
-
表 1 三种方案下北斗C2I-C6I类型卫星DCB平均估值
ns 卫星 PRN 方案1 方案2 方案3 GEO C01 −1.37 −1.34 −1.23 C02 4.25 4.22 4.36 C03 5.62 5.62 5.59 C04 8.82 8.85 9.06 C05 2.73 2.74 3.05 MEO C11 7.25 7.23 7.26 C12 6.71 6.67 6.73 C20 21.48 21.66 21.71 C24 7.98 7.94 7.81 C30 −8.58 −8.55 −8.55 IGSO C06 11.73 11.74 11.64 C07 18.10 18.10 18.02 C10 9.44 9.37 9.28 C38 4.63 4.63 4.70 C40 1.35 1.33 1.57 
下载: 导出CSV
表 2 各测站DCB与CAS DCB的平均差值
ns 测站 方案1 方案2 方案3 SGOC 0.62 0.51 1.12 IISC 0.59 0.59 1.44 GUAM 0.83 1.03 1.22 CUSV 0.75 0.71 0.71 PTGG 1.76 1.60 1.88 WUH2 0.92 0.61 0.39 GAMG 0.60 0.56 0.66 POL2 2.51 2.39 2.15 URUM 3.27 3.04 2.80 ULAB 3.38 2.90 2.69
下载: 导出CSV
-
[1] 邓远帆, 郭斐, 张小红, 等. 北斗三号卫星多频多通道差分码偏差估计与分析[J]. 测绘学报, 2021, 50(4): 448-456. [2] 刘姣, 陈俊平, 王彬. 估计卫星频间偏差的GLONASS伪距定位改进模型及其验证[J]. 大地测量与地球动力学, 2021, 41(1): 27-33. DOI: 10.14075/j.jgg.2021.01.006 [3] 王健, 党亚民, 王虎, 等. 融合BDS/GPS/GLONASS解算卫星差分码偏差及精度分析[J]. 测绘通报, 2018(10): 6-9,21. DOI: 10.13474/j.cnki.11-2246.2018.0304 [4] 姚宜斌, 刘磊, 孔建, 等. GIM和不同约束条件相结合的BDS差分码偏差估计[J]. 测绘学报, 2017, 46(2): 135-143. DOI: 10.11947/j.AGCS.2017.20160375 [5] LI Z S, YUAN Y B, LI H, et al. Two-step method for the determination of the differential of the differential code biases of COMPASS satellites[J]. Journal of geodesy, 2012, 86(11): 1059-1076. DOI: 10.1007/s00190-012-0565-4
[6] WANG N B, YUAN Y B, LI Z S, et al. Determination of differential code biases with multi-GNSS observations[J]. Journal of geodesy, 2016, 90(3): 209-228. DOI: 10.1007/s00190-015-0867-4
[7] 张宝成, 欧吉坤, 袁运斌, 等. 利用非组合精密单点定位技术确定斜向电离层总电子含量和站星差分码偏差[J]. 测绘学报, 2011, 40(4): 447-453. [8] WANG J, HUANG G W, ZHOU P Y, et al. Advantages of uncombined precise point positioning with fixed ambiguity resolution for slant total electron content (STEC) and differential code bias (DCB) estimation[J]. Remote sensing 2020, 12(2): 304. DOI: 10.3390/rs12020304
[9] LI X X, MA T Z, XIE W L, et al. FY-3D and FY-3C onboard observations for differential code biases estimation[J]. GPS solutions, 2019, 23(2): 57. DOI: 10.1007/s10291-019-0850-2
[10] 徐丹丹, 赵兴旺, 徐跃, 等. GPS双频观测值差分码偏差估计及精度分析[J]. 全球定位系统, 2018, 43(5): 17-22. DOI: 10.13442/j.gnss.1008-9268.2018.05.003 [11] ELGHAZOULY A, DOMA M, SEDEEK A. Estimating satellite and receiver differential code bias using a relative Global Positioning System network[J]. Annales geophysicae, 2019, 37(6): 1039-1047. DOI: 10.5194/ANGEO-37-1039-2019
[12] 刘乾坤, 隋立芬, 肖国锐, 等. 北斗系统差分码偏差解算中一种新的定权方法[J]. 测绘科学技术学报, 2015, 32(5): 473-478. DOI: 10.3969/j.issn.1673-6338.2015.05.008 -
期刊类型引用(0)
其他类型引用(1)
计量
- 文章访问数: 273
- HTML全文浏览量: 88
- PDF下载量: 23
- 被引次数: 1
