Camera Calibration Using Cross-Section Waterline Orientation for Video-based Flow Measurement
-
摘要: 针对现有基于直接线性变化法(DLT)的图像法测流技术依赖于地面控制点,存在效率低、风险高、宽断面天然河流操作难度大等问题,该文提出一种基于断面水边线定向的摄像机姿态角标定方法(CSWO)。该方法将标定过程分解为实验室内参标定和现场外参标定两步,其中后者又被划分为摄像机定位和定向两个环节。定向环节中首先在摄像机安装时将光轴与断面方向对齐,使方位角置零。然后利用无畸变图像中人工标注的平直断面水边线的斜率计算出横滚角。接下来通过计算水边线与图像纵轴交点作为其亚像素像方坐标。最后依据透视投影成像模型,联合水位与插值断面高程求得的物方坐标解算出俯仰角。该方法应用于时空图像测速法(STIV),实现了200 m宽河流的免像控表面流速测量。结果表明:起点距的最大绝对误差为0.59 m,最大相对误差为0.45%,表面流速的最大相对误差小于6.3%。Abstract: The image-based water flow measurement technology based on the Direct Linear Transformation (DLT) method relies on ground control points, and has problems such as low efficiency, high risk, and difficult operation in natural rivers with wide sections. A camera attitude angle calibration method based on Cross-Section Waterline Orientation (CSWO) is proposed. In this method, the calibration process is divided into two steps: laboratory calibration and field calibration, and the latter is divided into camera positioning and orientation. In the orientation step, the optical axis is required to be aligned with the cross-section when the camera is installed, so that the azimuth angle is set to zero. The roll angle is calculated by the slope of the straight waterline marked manually in the undistorted image. Then the sub-pixel image coordinate of the intersection point between the waterline and the vertical axis of the image is calculated. Finally, the pitch angle is calculated according to the perspective projection imaging model combined with the object coordinates obtained by the water level and the interpolated elevation of cross-section. This method has been applied to Space-Time Image Velocimetry (STIV) to measure the image-free surface velocity of a river with a width of 200 m. The results show that the maximum absolute error of starting distance is 0.59 m, the maximum relative error is 0.45%, and the maximum relative error of surface velocity is less than 6.3%.
-
1. 引言
能量收集(Energy Harvesting, EH)技术能够将不同类型的绿色能量转化为电能,在通信系统中采用EH技术有利于环境的保护和节点的部署。对于依靠收集的能量运转的无线通信系统,能量供应具有随机波动性,无线信道随机时变,如何对发送功率和传输速率进行动态的控制,实现对能量和信道的高效利用是一个重要课题[1]。EH通信系统中的功率控制算法包括离线功率控制和在线功率控制两类。离线功率控制算法应用于事先知道整个通信过程中收集的能量、信道衰落和数据到达等信息的情况下,离线注水算法是常见的离线算法[2]。实际系统中,能量到达和信道衰落无法预知,不依赖于未来系统状态信息的在线功率控制算法才具有实用性。在线算法根据能量到达、信道衰落等的统计信息,以及当前和过去的系统状态做出决策。将功率控制过程建模为马尔可夫决策过程(Markov Decision Process, MDP),并采用动态规划进行求解是常用的在线功率控制算法[3]。MDP加动态规划的方法依赖于状态转移概率等统计信息,在系统状态空间较大时求解复杂度很高。Lyapunov优化技术在控制论中应用广泛,它不依赖于系统的统计信息和未来的状态信息,而是根据系统当前和过去的状态做出决策,优化长期平均性能。保持队列长期时间意义上稳定是Lyapunov优化的基本特征[4]。Lyapunov优化框架的基本方法是把优化目标作为惩罚项,将队列长度漂移与惩罚项的和作为目标函数,通过最小化该目标函数来实现队列的稳定和性能的优化。对于包含约束条件的优化问题,通过构造虚队列并使其稳定来满足约束,可显著降低优化问题求解的复杂度。与MDP方法相比,Lyapunov优化不依赖于状态转移概率等统计信息,复杂度更低。近年来,有不少文献研究利用Lyapunov优化框架来求解EH通信系统中的功率控制问题。文献[5]针对发送端由EH装置供电的点对点通信系统平均速率最大化问题,利用Lyapunov优化框架将约束转换为虚队列稳定性要求,将长期平均优化问题转化为仅依赖于当前电池状态和信道状态的单时隙在线优化问题。文献[6]考虑两跳解码转发中继系统,针对源节点平均能耗最小化和中继传输速率最大化为目标的源节点和中继节点发送功率的联合优化问题,利用Lyapunov优化框架将优化问题转换为数据和能量队列稳定下的两个优化子问题,并求解。
研究EH系统中功率控制的文献,包括以上介绍的文献[2-6],多以优化系统理论上的最大传输速率为目标。在实际的系统中,需要采用适当的信道编码和调制方式才能高效率地利用能量和信道的传输能力。目前有少量文献研究了EH通信系统中的发送功率和调制方式的控制问题[7-10]。文献[7]针对点对点EH无线通信系统,提出了一种基于接收端反馈的1 bit信道状态信息的后向迭代算法来获得最优调制方式与发送功率组合。文献[8]针对点对点EH无线通信系统,将系统误码性能约束下的发送功率和调制阶数的控制过程建模为一个状态转移概率未知的MDP,并使用深度Q网络算法在连续状态空间中求解。文献[9]研究EH无线通信系统中传输速率最大化问题,在误比特率约束和有限种调制方式下,提出一种基于电池电量和信道状态的发送功率与调制方式的联合优化策略。文献[10]考虑点对点EH无线通信系统,以最大化系统长期平均实际可实现传输速率为目标构造发送功率、调制方式和帧长的联合优化问题,然后利用Lyapunov优化框架对优化问题进行转换并求解。以上文献中虽然在优化发送功率的同时优化了调制方式,但都没有考虑信道编码。在实际的通信系统中,要实现接近于信道容量的传输速率,信道编码是必不可少的。
本文针对能量收集无线通信系统中实际可达的信息传输速率最大化问题,研究发送功率、调制方式、信道编码码率的联合优化策略。系统的发送端配备有能量收集装置和容量有限的可充电电池,发送信号的能量来自从周围环境中收集的能量。首先构造可用能量和电池使用约束下的长期优化问题,然后利用Lyapunov优化框架对优化问题进行转换,再给出求解方法。与现有相关文献不同,本文的优化目标不是理论上传输速率的极限,即信道容量,而是实际系统可实现的传输速率,对发送功率、调制方式、信道编码参数进行联合优化。仿真的结果表明,本文所提算法实际可达的信息传输速率要高于以理论极限传输速率为优化目标的算法。
2. 系统模型
本文的系统模型如图1所示。发送端配备能量收集设备和容量有限的可充电电池,能量收集设备从周围环境中收集能量并转换为电能后存储在电池中,为发送端提供能量。信源产生的信息经过信道编码、调制后生成符号周期为Ts的符号送入信道传输。记时隙的长度为Tts,
Tts≫Ts 。发送端根据瞬时的信道状态以及能量收集情况,动态地调整发送功率、调制方式以及编码参数,最大化长期时间平均信息传输速率。信源生成的信息序列以k bit为一组,经二进制信道编码后得到n bit的码字,再经过M元调制后得到发送符号。设信道带宽为B Hz,则符号速率Rs=1/Ts=B Baud。信息的发送速率为
Rbt=B⋅log2M⋅kn (bit/s) (1) 记x(t)为发送端调制后得到的单位功率信号。假设信道为平坦衰落信道,信道系数为h(t),在1个时隙内保持不变。发送信号经过信道传输,接收端接收到的信号为
y(t)=√PT(t)h(t)x(t)+n(t) (2) 其中,PT(t)为信号发送功率,n(t)为均值为0、功率谱密度为N0/2的加性高斯白噪声。假设可选用的调制方式的集合为Ω={BPSK, QPSK, 16QAM, 64QAM, 256QAM}。接收端解调后的误比特率(Bit Error Rate, BER)为
Peb={Q(√2EbavN0),M=2,41−[1−2(1−1√M)Q(√3Ebavlog2M(M−1)N0)]2log2M,M=16,64,256 (3) 其中,
Q(x)=1√2π∫∞xe−t22dt 为高斯Q函数,Ebav 为接收信号的平均比特能量Ebav=|h(t)|2⋅Ts⋅PT(t)log2M (4) 解调后的二进制码元序列以码字为单位送入信道编码的译码器进行译码,正确译码的概率为[11]
Pcf=r∑i=0CinPieb(1−Peb)n−i (5) 其中,r为编码的纠错能力。译码错误的码字将被丢弃。本文以单位带宽上正确传输信息的速率作为衡量系统性能的指标,其表达式为
Rb=RbtPcf/B=log2M⋅k⋅r∑i=0CinPieb(1−Peb)n−i/n (6) 信息传输速率与调制方式、编码码率和误比特率有关,而误比特率则与发送功率和调制方式有关。调制阶数越大,虽然每个符号可携带更多的信息比特,但错误概率也越大。
发送端的功耗P(t)除了发送信号的功耗PT(t)外,还包括电路功耗Pcc(t),由编码电路功耗Pc(t)、调制电路功耗Pm(t)及其他电路功耗PA等几部分组成,即
P(t)=Pc(t)+Pm(t)+PA+PT(t)=Pcc(t)+PT(t) (7) 编码器的功耗与信息速率有关,可以表示为
Pc(t)=σ1⋅Rbt ,其中,σ1为一个常系数,与具体采用的编码方式、编码电路、电路参数等有关;调制器的功耗与每秒钟调制的二进制码元数量成正比,可表示为Pm(t)=σ2⋅Blog2M ,其中,σ2为一个常系数,具体数值与电路或信号处理器件有关。发送端在时隙t的总功耗P(t) 受到电池存储电量Eb(t)和最大放电功率Pd,max的限制0≤P(t)≤min(Eb(t)Tts,Pd,max) (8) 设电池最大容量为
EB ,可充电电池的最大充电功率为Pc,max,时隙t收集的电能为Ea(t) ,则时隙t存储进电池的能量EH(t) 可以表示为0≤EH(t)≤min(Ea(t),TtsPc,max,EB−(Eb(t)−TtsP(t))) (9) 其中,最小值函数最后一项为电池剩余存储容量。经过1个时隙的充电和放电后,电池中的电量更新为
Eb(t+1)=Eb(t)−Tts⋅P(t)+EH(t) (10) 3. 优化问题的构造与求解
3.1 优化问题的构造
发送功率越大,当前时隙的信息传输速率越高,但发送端的可用能量是受限的,需要在可用能量的约束下,根据信道状态合理地在时隙间分配和使用能量,才能使平均传输速率更高。如式(1)所示,M, k越大,发送端发送的信息速率越高,但是错误概率也越高,而M, k还要影响编码和调制电路的功耗。本文的优化问题就是在收集能量的约束下,对每个时隙的发送功率、调制阶数M和信息位长度k进行联合优化,最大化系统的长期平均信息传输速率
P1:max{PT(t),M,k} limT→∞1TT−1∑t=0E[Rb(t)],s.t.式(8)和式(10),M∈Ω,k∈K (11) 其中,
E[⋅] 表示期望运算,K为一个码字中可选信息位长度的集合,与选用的信道编码有关。将式(10)改写为
Eb(t+1)−Eb(t)=EH(t)−Tts⋅P(t) (12) 对式(12)中时隙t从0~T–1的各式两端进行叠加,并求期望可得
E[Eb(T)]−E[Eb(0)]=T−1∑t=0E[EH(t)−Tts⋅P(t)] (13) 式(13)左右两端同时除以T,并求T→∞的极限,得到
ˉEH−Tts⋅ˉP=0 (14) 其中,
ˉEH = limT→∞1TT−1∑t=0E[EH(t)] ,ˉP = limT→∞1TT−1∑t=0E[P(t)] ,分别表示每个时隙平均收集的能量和平均功耗。式(14)表示从长期来看,能量收集与使用间应保持平衡。将P1中的单时隙电池电量约束式(10)放松为长期时间电量约束式(14),优化问题放松为P2:max{PT(t),M,k} limT→∞1TT−1∑t=0E[Rb(t)],s.t.式(8)和式(14),M∈Ω,k∈K (15) 3.2 优化问题的转换
优化问题P2的目标函数和约束条件涉及长期平均,需要进行整个传输过程的联合优化,求解的复杂度很高。本文采用Lyapunov优化框架将长期平均约束转化为保持虚队列稳定,并将长期平均优化通过单时隙优化来逼近。首先构造发送端的电池能量虚队列
X(t)=Eb(t)−A (16) 其中,A为偏移量,用于控制电池的电量水平。根据式(10)易得能量虚队列的更新公式为
X(t+1)=X(t)−Tts⋅P(t)+EH(t) (17) 定义2次Lyapunov函数
L(X(t))≜ (18) Lyapunov漂移定义为
\Delta (X(t)) \triangleq {\rm{E}}\left[ {L(X(t + 1)) - L(X(t))\left| {X(t)} \right.} \right] (19) 漂移越小,电池电量越稳定。若能保持能量虚队列稳定,即电池的电量在一个有限的范围内波动,则从长期来看收集的能量与消耗的能量就是相等的,约束式(14)就能满足。将待最大化的信息传输速率Rb的负值作为惩罚项,与Lyapunov漂移一起构造“漂移加惩罚”
\Delta (X(t)) - V{\rm{E}}[{R_{\text{b}}}(t)\left| {X(t)} \right.] (20) 其中,V是漂移和惩罚项之间的权重。若能使“漂移加惩罚”最小化,则就在保持虚队列(即电池电量)稳定的同时,最大化信息的传输速率。由式(17)和式(18)易得漂移的一个上界为
\Delta X(t) \le {{C}} + X(t){\rm{E}}[{E_{\text{H}}}(t) - {T_{{\text{ts}}}} \cdot P(t)\left| {X(t)} \right.] (21) 其中,C为一个有限非负参数,满足
{{C}} \ge \frac{1}{2}{\rm{E}}[{({E_{\text{b}}}(t) - {T_{{\text{ts}}}} \cdot P(t))^2}\left| {X(t)} \right.] (22) 将式(20)中的漂移用上式代替,可得“漂移加惩罚”的上界为
\begin{split} &\Delta X(t) - V{\rm{E}}[{R_{\text{b}}}(t)\left| {X(t)} \right.] \le {{C}} + X(t){\rm{E}}[{E_{\text{H}}}(t) \\ & \quad - {T_{{\text{ts}}}} \cdot P(t)\left| {X(t)} \right.] - V{\rm{E}}[{R_{\text{b}}}(t)\left| {X(t)} \right.] \end{split} (23) 优化问题P2可以用去除了约束式(14)的“漂移加惩罚”上界的最小化来逼近。进一步去掉上界中的均值运算,并将速率的表达式(6)代入,优化问题P2可以转化为单时隙的优化问题
\begin{split} & \text{P}3:\underset{{P}_{\text{T}}(t),M,k}{\mathrm{max}}P(t)X(t)+V{\mathrm{log}}_{2}M\cdot k\\ & \qquad \cdot {\displaystyle \sum _{i=0}^{r}{C}_{n}^{i}{P}_{\text{eb}}{}^{i}{(1-{P}_{\text{eb}})}^{n-i}}/n,\\ & \qquad \text{s}\text{.t}\text{.}\;式(8),\;\;\, M\in {\boldsymbol{\varOmega}},\;\;\, k\in {\boldsymbol{K}} \end{split} (24) 3.3 优化问题的求解
记优化问题P3中的优化目标函数为
F({P_{\text{T}}}(t),M,k) ,展开其中总功率P(t)的表达式,得到\begin{split} F({P_{\text{T}}}(t),M,k) = & \left[ {{P_{\text{T}}}(t) + {\sigma _1} \cdot {R_{{\text{bt}}}} + {\sigma _2} \cdot B{{\log }_2}M} \right] \\ & \cdot X(t) + V{\log _2}M \cdot k \\ &\cdot \sum\limits_{i = 0}^r {C_n^i{P_{{\text{eb}}}}^i{{(1 - {P_{{\text{eb}}}})}^{n - i}}} /n \\[-21pt] \end{split} (25) 目标函数中包含3个优化变量,其中调制阶数M和信息位长度k为离散取值,且取值数量有限,而发送功率PT(t)可连续取值,因此这是一个混合整数优化问题。直接联合优化这3个变量十分困难,本文采用搜索的方式获得M和k的优化解,PT(t)则通过优化方法求解,即遍历所有M和k可能的组合,优化每个组合下的PT(t),选择使目标函数最大的M, k以及对应的PT(t)为优化问题的最优解。由于M和k的取值集合不大,遍历的复杂度可以接受。下面讨论给定M, k下最优PT(t)的求解。
当前时隙的传输速率随PT(t)增大而增大,因此目标函数式(25)中第2项的值也增大。易知当
X(t) \ge 0 时,目标函数(25)是信号发送功率PT(t)的单调增函数,为使目标函数最大化,PT(t)应取在电池放电功率和存储能量支持下的最大值{P_{{\text{T,max}}}}(t) = \min \left( {{P_{{\text{d,max}}}},\frac{{{E_{\text{b}}}(t)}}{{{T_{{\text{ts}}}}}}} \right) - {P_{{\text{cc}}}}(t) (26) 从实际物理意义上看,
X(t) \ge 0 表示电池中有足够的电量,可以最大功率发送信号。当X(t)<0时,式(25)中的第1项随PT(t)的增大而减小,而第2项增大,因此不能直接判断优化目标函数的单调性,需要进一步分析。求目标函数对发送功率PT(t)的偏导,得到\begin{split} \frac{{\partial F({P_{\text{T}}}(t),M,k)}}{{\partial {P_{\text{T}}}(t)}} =& X(t) + \frac{{V{{\log }_2}M \cdot k}}{n} \\ & \cdot \sum\limits_{i = 0}^r {C_n^i{P_{{\text{eb}}}}^{i - 1}(1 - {P_{{\text{eb}}}}){}^{n - i - 1} } \\ & \cdot\frac{{\partial {P_{{\text{eb}}}}}}{{\partial {P_{\text{T}}}(t)}} \cdot (i - n \cdot {P_{{\text{eb}}}})\\[-21pt] \end{split} (27) 若
\dfrac{{\partial F({P_{\text{T}}}(t),M,k)}}{{\partial {P_{\text{T}}}(t)}}{\text{ < 0}} ,目标函数式(25)单调递减,发送功率应取为0;否则需要找到\left[ {0,{P_{{\text{T,max}}}}(t)} \right] 范围内目标函数的极大值点,即满足\dfrac{{\partial F({P_{\text{T}}}(t),M,k)}}{{\partial {P_{\text{T}}}(t)}}{\text{ = 0}} 的解。观察偏导的表达式(27)和BER的表达式(3),易知该方程是一个包含积分和指数函数的、关于PT(t)的高次非线性方程,无法获得该方程的解析解,因此本文采用数值方法进行求解。由于极值点可能不止一个,先将发送功率范围[{P_{{\text{T,min}}}},{P_{{\text{T,max}}}}] ({P_{{\text{T,min}}}} 为发送功率的最小值,后面介绍如何确定)等分为数个大小为δ1的小区间,判断每一个小区间中是否包含极大值点,若该区间内包含极大值点,再采用黄金分割法找到该极大值点。记Pi为第i个小区间的左端点功率值。若在该区间内存在极大值点,则目标函数在Pi的右邻域是递增的,而在Pi+δ1的左邻域是递减的。因此,如果F({P_i}) < F({P_i} + {\delta _2}) 且F({P_i} + {\delta _1} - {\delta _2}) > F({P_i} + {\delta _1}) ,则该区间存在一个极大值点,这里{\delta _2} < {\delta _1} 为一个很小的正常数。在找到所有的极值点后,比较各极值点和边界点{P_{{\text{T,min}}}} ,{P_{{\text{T,max}}}} 处的目标函数值,选择使目标函数最大的点为最优功率点。搜索的流程见图2。任何一种信道编码都有这样一个特性,当进入译码器的序列错误概率超过某个阈值时,译码错误概率会很高,而只有在输入序列错误概率低于该阈值时,编码的纠错能力才能发挥出来,输出序列的错误概率会随输入错误概率的下降而快速下降。根据译码器的输入错误概率门限Peb,max的限制,相应有一个最小接收信号功率要求,与当前信道状态一起决定了最小发送功率要求PT,min。
优化问题P3的求解算法总结在算法1中,其中步骤(1)—步骤(16)遍历调制方式集合,步骤(3)—步骤(15)遍历编码信息位长度集合,步骤(4)—步骤(14)是在给定的调制方式、信息位长度下优化PT(t)。
算法1 P3的求解算法 设定参数:V, A, Peb,max; 输入:{\boldsymbol{\varOmega}} , K, δ1, δ2; 输出:M, k, PT(t), Rb(t); 在时隙t: (1) for M\in {\boldsymbol{\varOmega}} do (2) 将Peb,max代入式(3)求得PT,min; (3) for k\in {\boldsymbol{K}} do (4) Pcc(t)=Pc(t)+Pm(t)+PA; (5) 由式(26)计算得到PT,max; (6) if PT,min<PT,max (7) if X(t)>0 (8) PT(t)=PT,max; (9) else (10) 利用算法1搜索最优发送功率PT(t); (11) end if (12) else (13) PT(t)=0; (14) end if (15) end for (16) end for (17) 选择最大目标函数对应的PT(t), M, k ; (18) return PT(t), M, k, Rb; 4. 权重和偏移量的自适应调整
在优化问题P3中有两个参数,分别为“漂移加惩罚”中的权重V和电池电量虚队列中的偏移量A。权重反映对目标优化的重视程度,权重值越大,则倾向于消耗更多的电量使得当前时隙获得更高的传输速率,但同时会增大电量波动性,并降低电池的平均电量;权重值越小,虽然单时隙的传输速率优化强度降低,但电池的平均电量升高,在时隙间有更大的能量调度空间。偏移量用于控制电池的电量水平,偏移量越大,相同的电池电量下虚队列的值越小,优化中会更偏向于选择更小的发送功率,电池的平均电量越高。这两个参数过大或过小都会对系统的长期性能产生不利的影响,其最优取值与信道和能量的随机特性等有关,无法通过数学分析的方法得到。Lyapunov方法的特点之一就是不依赖系统的统计信息和未来的状态信息自适应地对系统进行控制,因此这两个参数也应能自适应地进行调整。
权重V和虚队列偏移量A都对电池电量波动和平均电量水平有直接的影响,所以可根据过去一段时间的电池电量状态对其进行调整。聚类可将数据分为不同的类,使得同类中的元素尽可能相似,不同类中的元素差别尽可能大,K-means算法简单且时间复杂度与数据集大小呈线性关系。本文利用基于滑动窗口的K-means聚类方法将一段时间内的电池电量点分为不同的簇,利用簇内和簇间电量关系反映电池电量的波动和平均特征。记滑动窗口的大小为l个时隙。自适应调整算法思路如下:对当前滑动窗口的电池电量进行K-means聚类,分析聚类后的均值电量
{\bar E_{\text{b}}}(t) 和加权均值电量{\bar E_{{\text{b,w}}}}(t) ,以{\bar E_{{\text{b,w}}}}(t) 是否处在允许的波动范围内为依据判断电量是否稳定,若不稳定,则对权重进行调整,而以{\bar E_{{\text{b,w}}}}(t) 是否超出上下阈值为依据进行偏移量的调整。将时间和电量作为2维欧氏空间中的坐标,时间和该时间的电池电量确定一个欧氏空间中的点,计算各点之间的欧氏距离,根据欧氏距离将滑动窗口内的点分为m个簇,聚类簇中心的电量表示为
{E_{{\text{c}},i}},i = 1,2,\cdots,m ,i越大,聚类簇中心越靠近当前时隙。滑动窗口内的均值电量与加权均值电量分别为{\bar E_{\text{b}}}(t) = 1/m \cdot \displaystyle\sum\nolimits_{i = 1}^m {{E_{{\text{c}},i}}} 和{\bar E_{{\text{b,w}}}}(t) = \displaystyle\sum\nolimits_{i = 1}^m {{w_i}{E_{{\text{c}},i}}} ,其中wi为第i簇的权重,满足\displaystyle\sum\nolimits_{i = 1}^m {{w_i} = 1} ,时间上距离当前时隙越近的簇中心的权重越大,即{w_1} < {w_2} < \cdots < {w_m} 。用加权均值电量与均值电量的比值{q_1}(t) = {\bar E_{{\text{b,w}}}}(t)/{\bar E_{\text{b}}}(t) 衡量时隙t滑动窗口内电池电量的稳定性。q1(t)>1表示当前滑动窗口内电量有增加的趋势,q1(t)<1则表示有减少的趋势。每隔d个时隙更新一次权重和偏移量的值,若q1(t)处于[{\gamma _{{\text{min}}}},{\gamma _{{\text{max}}}}] 范围内,则认为电量是稳定的,不对权重进行调整;若{q_1}(t) > {\gamma _{\max }} ,表明电池电量增加迅速,应适当增大权重;若{\text{ }}{q_1}(t) < {\gamma _{\min }} ,表明电池电量减少过快,应适当减小权重。权重调整的步长与滑动窗口内的加权均值电量和均值电量的差值大小成正比,该差值反映了当前滑动窗口下电量增长或降低的程度。权重调整的公式为V(t) = \left\{ \begin{aligned} & V(t - d) + {{{b}}_1} \cdot \left| {{{\bar E}_{{\text{b,w}}}}(t) - {{\bar E}_{\text{b}}}(t)} \right|,{\text{ }} {{q_1}(t) > {\gamma _{\max }}} \\ & V(t - d),\qquad\qquad\qquad\qquad\qquad\;\; {{\gamma _{\min }} \le {q_1}(t) \le {\gamma _{\max }}} \\ & V(t - d) - {{{b}}_1} \cdot \left| {{{\bar E}_{{\text{b,w}}}}(t) - {{\bar E}_{\text{b}}}(t)} \right|,\; {{q_1}(t) < {\gamma _{\min }}} \end{aligned} \right. (28) 其中,b1为一个常数,
V(t - d) 为上一个参数调整时刻权重的值。偏移量A直接影响电池的均值电量。定义滑动窗口中的均值电量与电池容量的比值为相对均值电量
{q_{\text{2}}}(t) = {\bar E_{\text{b}}}(t)/{E_{\text{B}}} ,当它处于上下阈值ηmax和ηmin之间时,表明电池电量维持在合适的范围内,偏移量不需要调整,否则需要对偏移量进行调整。调整的公式为A(t) = \left\{ \begin{aligned} & A(t - d) + {{{b}}_2} \cdot ({\eta _{\min }} - {q_2}(t)),{\text{ }}{{q_2}(t) < {\eta _{\min }}} \\ & A(t - d),\qquad\qquad\qquad\qquad\quad {{\eta _{\min }} \le {q_2}(t) \le {\eta _{\max }}}\\ & A(t - d) - {{{b}}_2} \cdot ({q_2}(t) - {\eta _{\max }}), {{q_2}(t) > {\eta _{\max }}} \end{aligned} \right. (29) 其中,b2为固定的常数。
5. 仿真结果与分析
本节通过仿真来验证所提出算法的性能。仿真中,信道系数h(t)服从均值为0、方差为
2 \times {10^{ - 9}} 的复高斯分布,即信道衰减为87 dB,信道系数在一个时隙内保持不变,在时隙间独立随机变化;发送端能量到达过程服从复合均匀分布的泊松过程,到达率为λ=0.5 单位/时隙,每个单位能量服从[0, 0.2] J之间的均匀分布;电池容量EB=50 J,电池初始电量为15 J,最大放电和充电功率为Pc,max=0.6 W和Pd,max=1.2 W;时隙长度Tts=1 s,仿真时长T=1×104 s;信道编码译码器限制的最大BER为Peb,max=2.5×10–2;噪声功率谱密度N0=1×10–17 W/Hz;带宽B=1×106 Hz,符号速率Rs=1×106 Baud。算法1中的{\delta _{\text{1}}} = ({P_{{\text{T,max}}}} - {P_{{\text{T,min}}}})/20 ,{\delta _{\text{2}}} = ({P_{{\text{T,max}}}} - {P_{{\text{T,min}}}})/1000 ;参数调整中的滑动窗口长度l=100 s,聚类数量m=4,聚类后中心电量点的权重{w_i},i = 1,2,\cdots,m 集合为{0.1, 0.2, 0.3, 0.4},权重和偏移量的调整间隔d=10 s,阈值为γmin=0.6, γmax=1.4, ηmax=0.8,ηmin=0.1;式(28)、式(29)中的b1和b2分别为0.025和0.01。如无特殊说明,电池电量虚队列的偏移量A的初始值为42,惩罚权重V的初始值为2.5。选用(n, k, r) Bose-Chaudhuri-Hocquenghem(BCH)码作为信道编码,码长为n=255,信息位长度k的取值集合为K={247, 239, 231, 223, 215, 207, 199, 191, 187, 179, 171, 163, 155, 147, 139, 131}。编码电路采用反馈移位寄存器电路[11],包括k个移位寄存器和平均约k/2个异或器(加法器)。假设采用互补金属氧化物半导体(Complementary Metal Oxide Semiconductor, CMOS)电路,忽略电路的静态功耗,1个移位寄存器或1个异或器在每个时钟周期的能耗为
{E_{{\text{co}}}} = \alpha \cdot {C_{\text{L}}} \cdot V_{{\text{dd}}}^{\text{2}} [12],其中,α为活动因子,为1个时钟周期内电路状态翻转的平均次数,这里考虑最大电路功耗,故设为1;CL=15 pF为负载电容,Vdd=1.0 V为工作电压。编码电路在1个时钟周期完成1次移位运算,时钟频率等于编码器输出码元的速率{R_{\text{c}}} = B \cdot {\log _2}M ,因此编码器的功耗为{P_{\text{c}}}(t) = \dfrac{{(k + 0.5k){E_{{\text{co}}}}{R_{\text{c}}}}}{{{T_{{\text{ts}}}}}} = \dfrac{{3k{E_{{\text{co}}}}B{{\log }_2}M}}{{2n \cdot {T_{{\text{ts}}}}}} 。调制采用数字信号处理器来完成,以TI的C55x系列处理器为例,其功耗为0.15 mW/MHz。进行M阶调制方式时,每秒钟需要进行Blog2M次运算。假设每两个时钟完成1次运算,则调制的功耗为{P_{\text{m}}}(t) = 3 \times {10^{ - 10}} \cdot B{\log _2}M 。其他电路功耗设为PA=0.011 W。先对权重V和偏移量A自适应调整算法的有效性进行验证。首先通过计算机用2维搜索的方法找到设定的仿真条件下,能获得最高长期平均信息传输速率的权重和偏移量的最优值分别为Vopt=2.5, Aopt=42。将本文第3节优化算法中的权重和偏移量固定为该最优值,仿真得到长期平均速率作为权重和偏移量自适应调整性能评估的基准。图3给出了将权重V和偏移量A固定为最优值,以及给定不同的初始值自适应调整时,从仿真开始到当前时隙的平均信息传输速率随时间变化的轨迹。可以看到,权重和偏移量自适应调整的算法中,初始取值不同时,仿真初段的传输速率值会有所不同,但随着仿真的进行,不同初始值下的平均信息传输速率会逐渐趋于一致,稍低于权重和偏移量固定为最优值的平均传输速率。这是因为权重和偏移量的自适应调整过程是一个逐步收敛的过程,仿真开始时由于参数偏离最优值较多,会有性能的损失。随着仿真的进行,参数逐渐收敛到最优值,性能也逐步提高并逼近最优。仿真结果表明本文给出的参数自适应调整算法是有效的。
为验证本文算法的性能,将其与离线注水算法、半功率算法以及固定比例算法进行性能比较。 (1) 离线注水算法:传输前已经获得整个传输过程中信道的变化情况和能量收集的情况,不考虑电量使用的因果性和电池电量的溢出,在收集总能量的约束下,采用交替迭代注水算法得到平均信道容量最大化的最优发送功率。(2) 文献[5]的在线功率控制算法:该文献的系统模型与本文类似,提出了一种以平均信道容量最大化为目标的在线功率控制算法,同样采用Lyapunov优化框架求解。仿真中权重和偏移量分别设置为2.5和40,此时算法性能最佳。(3) 半功率算法(Half Power Algorithm, HPA):每个时隙消耗的总能量为电池中现有电量的1/2,即
{P_{{\text{HPA}}}}(t) = {E_{\text{b}}}(t)/(2 \cdot {T_{{\text{ts}}}}) 。(4) 固定比例算法(Fixed Fraction Algorithm, FFA)[13]:每时隙消耗的总能量为固定比例部分的电池存储能量,即{P_{{\text{FFA}}}}(t) = \dfrac{{2{{\bar E}_{\text{H}}}(t)}}{{{E_{\text{B}}}}} \cdot \dfrac{{{E_{\text{b}}}(t)}}{{{T_{{\text{ts}}}}}} 。4种对比算法的调制阶数和信息位长度在发送功率确定后采用搜索的方法得到。HPA和FFA中,如果根据电池存储电量确定的总功率超过了电池的放电最大功率{P_{{\text{d,max}}}} ,则设为{P_{{\text{d,max}}}} 。图4给出了不同算法平均信息传输速率随着时间的变化轨迹,每时隙的平均信息传输速率为从仿真开始到当前时隙的信息传输速率的平均值。HPA和FFA都只根据电池电量计算得到发送功率,而没有考虑信道状态变化带来的影响,因此性能要明显低于根据信道状态和能量约束进行了优化的本文算法。离线注水算法以及文献[5]算法的优化目标为最大化理论上的长期平均信息传输速率,因此这两种算法在平均信道容量上要优于本文算法,相较于本文算法分别有1.6%,1.2%的优势。本文算法则是以实际可达的传输速率最大化为目标进行这3个变量的联合优化,因此尽管信道容量相对较低,但实际可达的信息传输速率却要高于这两种算法,分别有9.4%,10.3%的优势。图5给出了4种在线算法电池电量随时间的变化轨迹,离线注水算法不需要考虑能量使用的因果性和电池溢出情况,电池电量变化轨迹没有意义,因此没有给出。仿真结果显示,本文所提出的算法及文献[5]算法的电池电量能在一定水平上下波动,每时隙电池中都有足够的电量供信息传输,也为收集的能量留有足够的存储空间。而HPA在很短时间内消耗完事先存储的电量,随后电量稳定在一个很低的水平上;FFA在消耗电量时相对保守,仅使用现有电池中存储电量的固定比例部分,没有根据信道状态进行调整,因此性能较低。图6给出了长期平均信息传输速率随能量到达率λ变化的仿真结果。随着能量到达率的增加,可用能量增加,信号的平均发送功率增大,因此所有算法的平均信息传输速率都增大。相较于对比算法,本文算法由于是以实际可达的传输速率最大化对发送功率、调制方式和编码参数进行联合优化,因此能获得最高的信息传输速率。
为验证本文算法的适用性,还分别在随机风力发电模型、伯努利能量到达模型下进行了平均传输速率的仿真。风力发电模型下,风力发电模块输出电功率为
{P_{\text{a}}} = 0.5\rho S{C_{\text{p}}}{v_{\text{T}}}^3{\eta _1}{\eta _2} [14],其中:ρ=1.22 kg/m3为空气密度;vT为随机变化的风速,采用4分量随机模型仿真得到[15];Cp=0.593为风能利用系数;S = \pi {R^{\text{2}}} 为扇叶的扫掠面积,R为扇叶半径;η1=0.94为齿轮箱的能量转换效率;η2=0.52为电机的能量转换效率。伯努利能量到达模型下,每时隙到达的能量为独立同分布的伯努利随机变量,即{E_{\text{a}}}(t) = \left\{ \begin{aligned} & 0.5,p \\ & 0, \;\; {1 - p} \end{aligned} \right. ,其中概率p为能量到达的概率。图7给出风力发电模型下,不同扇叶半径R时长期平均传输速率的仿真结果,图8则给出伯努利能量到达模型下,不同能量到达概率p时长期平均传输速率的仿真结果。随着扇叶半径的增加,风力发电机的扫掠面积增大,输出的电功率增大,或者随着能量到达概率的增加,平均到达的能量增加,所有算法都能获得更高的信息传输速率。与图6类似,仍然是本文算法能获得的传输速率最高,其次是离线注水算法,HPA算法最低。6. 结束语
本文针对发送端由能量收集设备供电的无线通信系统,在能量收集和信道状态先验信息未知的条件下,研究以最大化实际可达传输速率的发送功率、调制方式和信道编码码率(信息位长度)的联合优化问题。基于Lyapunov优化框架,将能量使用的长期约束转换为能量虚队列的稳定性要求,将能量使用约束下的长期时间平均实际可达传输速率最大化问题转化为单时隙的、仅依赖于当前信道状态和电池状态的“漂移加惩罚”项上界的最小化问题。由于优化变量中的调制阶数和信息位长度为离散取值,且组合数量不大,故将原优化问题转化为调制阶数和信息位长度组合下求解发送功率的单变量优化问题,再遍历各组合获得原问题的解。求解发送功率时,首先根据电池电量虚队列的状态对目标函数的单调性进行分析,对于单调增或单调减的情况直接获得发送功率的最优值;对于存在极值点的情况,由于不能通过解析的方法获得最优功率的解,故给出了一种高效的数值求解方法。另外还给出了基于滑动窗口的K-means聚类方法的“漂移加惩罚”中权重和电池电量虚队列偏移量两个参数的自适应调整算法。在不同能量到达随机模型下与对比算法进行了性能的仿真对比,结果表明,本文所提算法在各种能量到达模型下都能获得更高的长期实际平均传输速率。另外,通过与参数固定为最优情况下算法性能的对比,证明算法对参数自适应调整正确、有效。本文算法不依赖于能量到达和信道状态的非因果信息,也不依赖于它们的统计信息,仅依赖于当前的电池状态和信道状态信息,计算复杂度低,且优化的是在实际可用的编码和调制方式的可达信息传输速率,因此是一种实用的算法。本算法中译码后错误概率的表达式适用于信道编码采用硬判决译码的情况,若信道编码采用软判决译码,该表达式不再适用,因此本文提出的算法不能直接使用。信道编码采用软判决译码时,发送功率、调制方式和信道编码的码长、码率等参数的联合优化将是下一步研究中需要解决的问题。
-
表 1 内参矩阵和畸变参数
m (pixel) n (pixel) {f_x} (pixel) {f_y} (pixel) {C_x} (pixel) {C_y} (pixel) {k_1} {k_2} {p_1} {p_2} 标定结果 3840 2160 2876.507 2884.631 1947.382 1043.743 –0.407 0.001 0.001 –0.001 表 2 各水位级下俯仰角、横滚角标定结果对比
测次 日期 时间 场景 水位级 水位(m) R\left(x_{1}, y_{1}\right) R_{2}\left(x_{2}, y_{2}\right) 俯仰角(°) 横滚角(°) 1 210407 15:00 晴 低 985.34 (990, 443) (2071, 435) 20.16 0.42 2 210204 10:30 晴 986.19 (982, 445) (2125, 441) 19.70 0.20 3 210416 20:05 夜 986.34 (984, 439) (2127, 433) 19.79 0.30 4 210416 20:10 夜 986.34 (984, 439) (2125, 431) 19.81 0.40 5 201118 10:30 阴 986.37 (960, 437) (2151, 429) 19.84 0.39 6 201118 10:35 阴 986.37 (966, 437) (2143, 431) 19.81 0.29 7 210415 14:18 晴 986.37 (984, 437) (2141, 427) 19.87 0.50 8 210318 9:00 阴 986.40 (978, 437) (2027, 427) 19.88 0.55 9 210414 13:00 晴 986.42 (984, 437) (2141, 431) 19.79 0.30 10 210408 14:00 晴 986.42 (1006, 439) (2079, 429) 19.82 0.53 11 210408 16:00 晴 986.62 (974, 437) (2187, 431) 19.79 0.28 12 210316 9:00 阴 986.90 (978, 427) (2143, 419) 19.81 0.39 13 210223 11:00 晴 987.03 (982, 423) (2095, 415) 19.84 0.41 14 220704 9:35 阴 中 993.80 (1100, 257) (2057, 255) 19.91 0.12 15 220704 9:20 阴 993.84 (1096, 257) (2073, 255) 19.89 0.12 16 220704 9:26 阴 993.85 (1092, 257) (2071, 255) 19.89 0.12 17 200821 10:00 阴 高 997.96 (1008, 206) (1824, 202) 19.88 0.28 表 3 DSO法和CSWO法实验结果
测次 水位(m) 起点距真值(m) 绝对误差(m) 相对误差(%) DSO CSWO DSO CSWO 1 985.34 55.00 0.29 0.02 0.22 0.02 2 986.19 90.00 0.64 0.12 0.51 0.09 3 986.34 65.00 0.24 0.18 0.18 0.14 4 986.34 90.00 0.81 0.02 0.62 0.01 5 986.37 55.00 0.13 0.17 0.10 0.13 6 986.37 65.00 0.64 0.23 0.50 0.18 7 986.37 90.00 0.79 0.01 0.61 0.01 8 986.40 120.00 1.64 0.27 1.27 0.21 9 986.42 65.00 0.34 0.09 0.26 0.07 10 986.42 105.00 0.94 0.13 0.72 0.10 11 986.62 135.00 1.18 0.59 0.90 0.45 12 986.90 65.00 0.26 0.28 0.19 0.21 13 987.03 90.00 0.72 0.19 0.55 0.15 14 993.80 90.00 0.97 0.29 0.57 0.17 15 993.84 55.00 0.42 0.05 0.25 0.03 16 993.85 65.00 0.66 0.00 0.39 0.00 17 997.96 55.00 0.31 0.28 0.17 0.15 -
[1] LE COZ J, HAUET A, PIERREFEU G, et al. Performance of Image-based Velocimetry (LSPIV) applied to flash-flood discharge measurements in mediterranean rivers[J]. Journal of Hydrology, 2010, 394(1/2): 42–52. doi: 10.1016/j.jhydrol.2010.05.049. [2] FUJITA I, WATANABE H, and TSUBAKI R. Development of a non-intrusive and efficient flow monitoring technique: The Space-Time Image Velocimetry (STIV)[J]. International Journal of River Basin Management, 2007, 5(2): 105–114. doi: 10.1080/15715124.2007.9635310. [3] TAURO F, TOSI F, MATTOCCIA S, et al. Optical Tracking Velocimetry (OTV): Leveraging optical flow and trajectory-based filtering for surface streamflow observations[J]. Remote Sensing, 2018, 10(12): 2010. doi: 10.3390/rs10122010. [4] 张振, 徐枫, 王鑫, 等. 河流水面成像测速研究进展[J]. 仪器仪表学报, 2015, 36(7): 1441–1450. doi: 10.19650/j.cnki.cjsi.2015.07.001.ZHANG Zhen, XU Feng, WANG Xin, et al. Research progress on river surface imaging velocimetry[J]. Chinese Journal of Scientific Instrument, 2015, 36(7): 1441–1450. doi: 10.19650/j.cnki.cjsi.2015.07.001. [5] 蒋海军. 大视场条件下的摄像机标定方法研究[D]. [硕士论文], 国防科技大学, 2017. doi: 10.27052/d.cnki.gzjgu.2017.000625.JIANG Haijun. Research on camera calibration method in large field of view[D]. [Master dissertation], National University of Defense Technology, 2017. doi: 10.27052/d.cnki.gzjgu.2017.000625. [6] 张振, 徐枫, 沈洁, 等. 基于变高单应的单目视觉平面测量方法[J]. 仪器仪表学报, 2014, 35(8): 1860–1868. doi: 10.19650/j.cnki.cjsi.2014.08.025.ZHANG Zhen, XU Feng, SHEN Jie, et al. Plane measurement method with monocular vision based on variable-height homography[J]. Chinese Journal of Scientific Instrument, 2014, 35(8): 1860–1868. doi: 10.19650/j.cnki.cjsi.2014.08.025. [7] 杨聃, 邵广俊, 胡伟飞, 等. 基于图像的河流表面测速研究综述[J]. 浙江大学学报:工学版, 2021, 55(9): 1752–1763. doi: 10.3785/j.issn.1008-973X.2021.09.017.YANG Dan, SHAO Guangjun, HU Weifei, et al. Review of image-based river surface velocimetry research[J]. Journal of Zhejiang University:Engineering Science, 2021, 55(9): 1752–1763. doi: 10.3785/j.issn.1008-973X.2021.09.017. [8] AL-MAMARI M M, KANTOUSH S A, KOBAYASHI S, et al. Real-time measurement of flash-flood in a Wadi area by LSPIV and STIV[J]. Hydrology, 2019, 6(1): 27. doi: 10.3390/hydrology6010027. [9] 曹列凯, DETERT M, 李丹勋. 基于无人机的长河段表面流场测量系统与应用[J]. 清华大学学报:自然科学版, 2022, 62(12): 1922–1929. doi: 10.16511/j.cnki.qhdxxb.2022.21.012.CAO Liekai, DETERT M, and LI Danxun. Airborne image velocimetry system and its application to measure the surface flow fields of long river reaches[J]. Journal of Tsinghua University:Science and Technology, 2022, 62(12): 1922–1929. doi: 10.16511/j.cnki.qhdxxb.2022.21.012. [10] TSAI R. A versatile camera calibration technique for high-accuracy 3D machine vision metrology using off-the-shelf TV cameras and lenses[J]. IEEE Journal on Robotics and Automation, 1987, 3(4): 323–344. doi: 10.1109/JRA.1987.1087109. [11] HOLLAND K T, HOLMAN R A, LIPPMANN T C, et al. Practical use of video imagery in nearshore oceanographic field studies[J]. IEEE Journal of Oceanic Engineering, 1997, 22(1): 81–92. doi: 10.1109/48.557542. [12] BECHLE A J, WU C H, LIU Wencheng, et al. Development and application of an automated river-estuary discharge imaging system[J]. Journal of Hydraulic Engineering, 2012, 138(4): 327–339. doi: 10.1061/(ASCE)HY.1943-7900.0000521. [13] HUANG Weiche, YOUNG C C, and LIU Wencheng. Application of an automated discharge imaging system and LSPIV during typhoon events in Taiwan[J]. Water, 2018, 10(3): 280. doi: 10.20944/preprints201802.0089.v1. [14] COSTA F A L and MITISHITA E A. An approach to improve direct sensor orientation using the integration of photogrammetric and lidar datasets[J]. International Journal of Remote Sensing, 2019, 40(14): 5651–5672. doi: 10.1080/01431161.2019.1580794. [15] ZHANG Zhen, ZHAO Lijun, LIU Boyuan, et al. Free-surface velocity measurement using direct sensor orientation-based STIV[J]. Micromachines, 2022, 13(8): 1167. doi: 10.3390/mi13081167. [16] LIU Wencheng and HUANG Weiche. Development of a three-axis accelerometer and Large-Scale Particle Image Velocimetry (LSPIV) to enhance surface velocity measurements in rivers[J]. Computers & Geosciences, 2021, 155: 104866. doi: 10.1016/j.cageo.2021.104866. [17] 张振, 吕莉, 石爱业, 等. 基于物像尺度变换的河流水面流场定标方法[J]. 仪器仪表学报, 2017, 38(9): 2273–2281. doi: 10.19650/j.cnki.cjsi.2017.09.023.ZHANG Zhen, LV Li, SHI Aiye, et al. River surface flow field calibration method based on object-image scaling[J]. Chinese Journal of Scientific Instrument, 2017, 38(9): 2273–2281. doi: 10.19650/j.cnki.cjsi.2017.09.023. [18] TSUBAKI R. On the texture angle detection used in Space-Time Image Velocimetry (STIV)[J]. Water Resources Research, 2017, 53(12): 10908–10914. doi: 10.1002/2017WR021913. [19] FUJITA I, KOSAKA Y, and YOROZUYA A. Tracking of river surface features by space time imageing[C]. The 15th International Symposium on Flow Visualization, Minsk, Belarus, 2012: 25–28. [20] ZHANG Zhen, ZHOU Yang, LIU Haiyun, et al. Visual measurement of water level under complex illumination conditions[J]. Sensors, 2019, 19(19): 4141. doi: 10.3390/s19194141. [21] 王慧斌, 董伟, 张振, 等. 基于时空图像频谱的时均流场重建方法[J]. 仪器仪表学报, 2015, 36(3): 623–631. doi: 10.19650/j.cnki.cjsi.2015.03.018.WANG Huibin, DONG Wei, ZHANG Zhen, et al. Time-averaged flow field reconstruction method based on spectrum of spatio-temporal image[J]. Chinese Journal of Scientific Instrument, 2015, 36(3): 623–631. doi: 10.19650/j.cnki.cjsi.2015.03.018. [22] ZHANG Zhengyou. A flexible new technique for camera calibration[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2000, 22(11): 1330–1334. doi: 10.1109/34.888718. [23] FUJITA I, NOTOYA Y, TANI K, et al. Efficient and accurate estimation of water surface velocity in STIV[J]. Environmental Fluid Mechanics, 2019, 19(5): 1363–1378. doi: 10.1007/s10652-018-9651-3. [24] 赵浩源, 陈华, 刘维高, 等. 基于河流表面时空图像识别的测流方法[J]. 水资源研究, 2020, 9(1): 1–11. doi: 10.12677/JWRR.2020.91001.ZHAO Haoyuan, CHEN Hua, LIU Weigao, et al. Application of flow measurement method based on space-time image recognition of river surface[J]. Journal of Water Resources Research, 2020, 9(1): 1–11. doi: 10.12677/JWRR.2020.91001. 期刊类型引用(0)
其他类型引用(1)
-