雙Johnson SB 分布湍流擾動OAM 光通信接收信號生成研究

陳純毅，熊琴琴，于海洋，姚海峰，宋佳雪，婁巖

（1.長春理工大學計算機科學技術學院，吉林 長春 130022；2.長春理工大學空地激光通信技術國防重點學科實驗室，吉林 長春 130022）

1 引言

光波的軌道角動量（OAM,orbital angular momentum）橫向空間模為經典和量子光信息科學中的信息編碼提供了一種額外的自由度，近年來受到國內外研究者的廣泛關注[1-3]。在經典光通信中，人們使用OAM 來實現空間模分復用，以便提高光通信系統容量[1-2]。在光量子密鑰分發領域，OAM被用于實現高維量子態，以便把超過1 bit 的信息編碼在一個光子中，從而增加每個光子的密鑰生成率[3]。從是否需要光纜的角度可將光通信分成光纖通信和自由空間光通信2 種形式。在這2 種形式中，前者穩定性好，但需要預先在通信雙方之間鋪設光纜，實施成本高；后者不需要鋪設光纜，靈活機動，實施成本低，但容易受到外界環境的干擾[4]。綜合比較來說，對于難以鋪設光纜或者對靈活性要求很高的應用場景，自由空間光通信是一個較好的方案。

大氣湍流是自由空間光通信中最重要的外界環境干擾因素之一[4-5]。為了理解自由空間光通信中大氣湍流的影響，國內外學者對受大氣湍流擾動的光信號的統計特性做了大量研究。目前，學者已經提出了多種描述受大氣湍流擾動的光信號起伏概率分布模型，例如對數正態分布、K分布和Gamma-Gamma 分布等[6-7]。然而，這些傳統分布模型在建立時沒有考慮OAM 橫向空間模接收的特殊性。大氣湍流對OAM 光通信的影響表現為信號光功率在不同OAM 橫向空間模之間的隨機再分布。Chen 等[8]發現，上述光信號起伏概率分布模型不太適合用來描述因大氣湍流擾動而導致的OAM接收光信號起伏，因此提出了雙Johnson SB分布模型。模擬實驗表明該分布模型能夠很好地描述受大氣湍流擾動的 OAM 接收光信號起伏統計特性。

在用蒙特卡洛方法來分析 OAM 光通信系統性能時，需要模擬生成受大氣湍流擾動的OAM 接收光信號的時域隨機樣本。受大氣湍流擾動的光信號起伏具有時間相關性，屬于時域自相關信號。模擬生成受大氣湍流擾動的OAM 接收光信號的關鍵是，要使生成的隨機光信號時域樣本在統計上服從雙Johnson SB分布，同時具有指定的時間自相關。本文將滿足這2 個條件的隨機光信號稱為雙Johnson SB分布湍流擾動OAM 光通信自相關接收信號。文獻[9-10]基于隨機微分方程（SDE,stochastic differential equation）分別生成了具有指數自相關函數的Gamma-Gamma 分布和Johnson SB分布隨機數序列。但目前還沒有生成雙Johnson SB分布湍流擾動OAM 光通信自相關接收信號時域樣本的報道。文獻[11]利用自回歸（AR,autoregression）模型生成了具有指定時間自相關的復高斯隨機樣本。文獻[12]利用逆變換法將自相關高斯隨機樣本變換為服從特定統計分布的隨機樣本。然而，采用逆變換法將自相關高斯樣本轉換為雙 Johnson SB分布樣本的過程較復雜。文獻[13]提出的秩匹配技術通過排序操作將自相關瑞利序列的時間自相關傳遞給服從Nakagami 分布的白噪聲隨機數序列，避免了逆變換法中對目標分布的求逆運算，對生成服從雙Johnson SB分布的時間自相關隨機光信號的時域樣本具有參考價值。

本文基于AR 模型和秩匹配法，提出一種模擬生成雙Johnson SB分布湍流擾動OAM 光通信自相關接收信號的方法。所提方法分為2 個步驟：首先用AR 模型生成具有指定時間自相關的高斯隨機樣本；然后生成一個服從雙Johnson SB分布的白噪聲隨機數序列，并對該白噪聲隨機數序列進行適當的重新排列，以匹配指定的時間自相關。其中，重排的順序由自相關高斯過程的統計秩決定。根據指定的指數形式的時間自相關函數，本文利用所提方法生成了雙Johnson SB分布湍流擾動OAM 光通信接收信號，并將生成信號的自相關系數與理論自相關系數的均方誤差作為生成信號的自相關性匹配精度，由此分析了不同的雙Johnson SB分布模型參數對OAM接收光信號時域樣本的自相關性匹配精度的影響。

2 OAM 光通信湍流擾動接收信號模擬模型

生成自相關高斯隨機數序列的常用方法為采用線性濾波，高斯白噪聲隨機數序列經過線性濾波后可獲得指定的自相關，并且服從其原有的高斯分布[11]。當目標分布為非高斯分布時，僅進行線性濾波不再可行，例如，對指數序列進行線性濾波時，輸出的自相關隨機數序列將不再服從指數分布。對于雙Johnson SB分布來說，這個結論同樣適用。為了對雙Johnson SB分布湍流擾動OAM 光通信自相關接收信號進行蒙特卡洛模擬，本文提出如圖1 所示的原理框架。首先使用AR 模型生成具有指定自相關的高斯隨機信號，將其作為輔助集；然后利用蒙特卡洛接受?拒絕采樣法[15]生成服從雙Johnson SB分布的白噪聲隨機信號，即樣本集；最后通過秩匹配將樣本集按照輔助集的統計秩進行重新排列，在不改變分布特性的情況下使樣本集具有指定的自相關。

2.1 基于AR 模型的自相關高斯過程

AR 模型本質上是一種線性濾波器，能生成具有指定自相關的高斯隨機過程。由文獻[9]可知，指數形式的歸一化時間自相關函數如式(1)所示。

其中，τ表示時延，τc表示時間自相關長度。基于AR模型生成的自相關高斯過程稱為AR高斯過程，一個p階AR高斯過程通過時域遞歸[11]表示為

其中，ak(k=1,2,…,p)表示AR 模型系數，p表示AR 模型的階數；w(n)表示高斯白噪聲信號，其均值為0、方差為。AR 模型的參數包括系數ak和高斯白噪聲的方差，它們和指定的自相關序列R(n)間的非線性關系用Yule-Walker方程[14]可表示為

其中，R表示由R(n)構成的自相關矩陣，R(n)表示式(1)對應的自相關值序列，具體為

其中，Δt表示采樣時間間隔；對于本文的p階AR高斯過程，n滿足?p≤n≤p的條件。

本文利用Levinson-Durbin 遞推算法來快速地求解式(3)，得到AR 模型的參數。該算法的初始值[14]為

對于k=2,3,…,p，遞推寫為[14]

其中，

其中，ρk為k階預測誤差方差。遞推完成得到參數集{a1(1),ρ1},{a2(1),a2(2),ρ2},…,{ap(1),ap(2),…,ap(p),ρp}。最后的p階參數集就是Yule-Walker方程的解，即ap(i)=ai,i=1,2,…,p，。對式(4)的自相關序列進行一次Levinson-Durbin 遞推后，式(9)中的預測誤差方差便不再隨p的增加而變化，由此確定AR 模型的階數p為1。

由式(2)可以看出，若要得到x(n)，除了AR 模型的參數外，還需要計算前p個平穩輸入序列元素的初始值[11]。對于本文的一階AR 模型，其初始條件x(0)是均值為0、方差為1 的高斯變量。最后將AR 模型的參數和初始條件代入式(2)，就得到自相關高斯過程。

圖1 本文的原理框架

2.2 自相關雙Johnson SB 分布接收光信號

由自相關高斯過程生成自相關非高斯過程的標準方法是非線性逆變換法[12]。在這種方法中，目標分布的累積分布函數（CDF,cumulative distribution function）非常重要，因為在構建合適的非線性變換和高斯過程的自相關函數時都需要目標CDF的反函數。對雙Johnson SB分布來說，其概率密度函數比較復雜，不便于獲得其CDF 的反函數，因此本文不考慮使用該方法而是利用統計秩來匹配指定自相關[13]。

統計秩是按指定（遞減或遞增）順序排列的樣本列表中一個值的序數[13]，本文選擇升序排序。假設樣本長度為5，對應的統計秩向量如圖2 所示，其中，第i個元素表示排列中第i個原始樣本的統計秩。對于平穩隨機過程，樣本的自相關函數近似于其統計秩的自相關函數[13]。因此，AR 高斯過程和雙Johnson SB分布隨機過程具有式(10)所示關系。

圖2 統計秩向量原理

其中，RGau和分別表示AR 高斯過程及其統計秩的自相關函數；RJoh和分別表示雙Johnson SB分布隨機過程及其統計秩的自相關函數。若雙Johnson SB分布隨機過程的統計秩與AR 高斯過程的統計秩相同，則二者統計秩的自相關函數相同，為

由式(10)和式(11)可得

由此可知，將服從雙Johnson SB分布的樣本按照AR 高斯過程的統計秩進行重排，可以得到與AR高斯過程秩匹配的自相關雙Johnson SB分布隨機過程，實現過程如圖3 所示，具體步驟如下。

步驟1生成服從雙Johnson SB分布的白噪聲隨機信號。雙Johnson SB分布的概率密度函數（PDF,probability density function）為[8]

其中，4 個模型參數滿足條件δ1＞0，δ2＞0，γ1≤γ2；pz(·)表示Johnson SB分布的PDF。為了避免數值擬合雙Johnson SB分布CDF 的反函數，本文采用蒙特卡洛接受?拒絕采樣法[15]生成服從雙Johnson SB分布的白噪聲隨機信號，計算方法 如下。

1)分別生成在[c,d]和[0,1]范圍內均勻分布的隨機變量u和v。

圖3 雙Johnson SB 分布湍流擾動OAM 光通信自相關接收信號的秩匹配流程

3)重復1)和2)，直到生成N個被接受樣本，其中，N表示樣本長度。

步驟2找到自相關高斯隨機信號x(n)的統計秩向量。由圖2 可以看出，原始樣本的統計秩向量與其進行遞增排序后的索引向量相同。因此，對x(n)升序排序，得到x(Jn),n=0,…,N?1，其中索引Jn使x(J0)≤x(J1)≤…≤x(JN?1)，Jn即為x(n)的統計秩向量。

步驟3根據Jn對雙Johnson SB分布白噪聲隨機信號z(n)進行重新排列，使z(n)的最小值放置在與x(n) 最小值相同的位置，z(n)的第二小值位于x(n) 第二小值的位置，以此類推。首先將z(n)升序排序，使z(I0)≤z(I1)≤…≤z(IN?1)，有序值記為z(In),n=0,…,N?1；然后重排z(n)，即

最后得到的重排信號z(n) 就具有指定的時間自相關并且服從雙Johnson SB分布。

3 仿真實驗與討論

通過模擬仿真實驗來驗證圖1 所示的原理框架，所有仿真的樣本長度N=105、采樣時間間隔Δt=0.1 ms。為了更準確地描述雙Johnson SB分布OAM 光通信接收信號的自相關性，所有仿真實驗重復100 次，并根據這100 次仿真得到的數據來計算自相關系數。

圖4(a)為τc=10 ms、δ1=3、γ1=1.1、δ2=5、γ2=2.4的雙Johnson SB分布湍流擾動OAM 光通信自相關接收信號，其中，雙Johnson SB分布的4 個模型參數的取值源自文獻[8]。圖4(b)為對應仿真信號的自相關系數隨時延變化的曲線，與式(1)的理論歸一化自相關函數基本一致。為了進一步分析仿真信號的自相關性是否與理論值匹配，本文以仿真信號的自相關系數與理論自相關系數的均方誤差（MSE,mean square error）作為自相關性匹配精度，即

圖4 雙Johnson SB 分布OAM 光通信接收信號及其自相關函數示例

為了驗證所提方法的有效性，根據文獻[8]提供的查找表，另外設置了6 組雙Johnson SB分布的模型參數進行仿真實驗，如表1 所示。當τc=10ms 時，根據表1 的參數仿真生成信號的概率密度如圖5所示。OAM 功率傳輸系數為接收的特定OAM 模中的功率與特定OAM 模的發射功率之比[8]。本文通過生成的雙Johnson SB分布湍流擾動OAM 光通信自相關信號對應的功率傳輸系數擾動來表征OAM光信號功率接收時域擾動。

表1 6 組雙Johnson SB 分布模型參數值

圖5 根據表1 的參數仿真生成信號的樣本概率密度與理論概率密度的對比(τ c=10ms)

圖6 表示當τc為10 ms 和20 ms 時，根據表1的參數仿真生成信號的自相關系數與理論自相關系數。從圖6 可以看出，仿真信號的自相關性基本匹配理論自相關性。當雙Johnson SB分布的形狀與高斯分布比較相似時，誤差會減小。通過式(16)的自相關性匹配精度對此進行詳細分析，M的取值如圖4 所示。

表1 中模型參數對應的匹配精度如表2 所示。從表2 可以看出，τc=10 ms 時的匹配精度與τc=20 ms 時的精度基本上接近。其中第4 組和第6 組參數的匹配精度最好，結合圖5 的概率密度曲線可以看出，這是由于第4 組和第6 組模型參數對應的概率密度曲線形狀類似于高斯分布。值得注意的是，AR 模型的輸出是均值為0、方差為1 的自相關高斯過程。進一步分析可知，當進行秩匹配的2 種分布相似時，本文所提方法的仿真精度更高。

為分析雙Johnson SB分布的模型參數值對仿真信號的自相關性的影響，對表2 進一步研究：將第1 組參數與第6 組參數的自相關性匹配精度進行比較可以看出，隨著δ1和δ2的增加，接近或者大于2 時，自相關性匹配更好；根據第6 組參數的匹配精度優于第4 組參數的精度和第2 組參數的匹配精度優于第5 組參數的精度可以發現，當越小時，自相關性匹配越好；通過對比第3 組與第5 組的匹配精度發現，當的比率接近于1 時，匹配精度更高。對于表2 中的6 組雙JohnsonSB分布的模型參數，對應的自相關性匹配精度基本保持在10?3以下，再次說明了所提的雙Johnson SB分布湍流擾動OAM 光通信自相關接收信號蒙特卡洛模擬方法的有效性。該方法生成的信號可用于研究大氣湍流對OAM 光通信的影響。

表2 表1 的6 組模型參數的自相關性匹配精度

4 結束語

本文提出一種基于AR 模型和秩匹配的雙Johnson SB分布湍流擾動OAM 光通信自相關接收信號蒙特卡洛模擬方法。由于AR 模型的參數計算簡單，秩匹配法僅涉及排序，因此本文所提方法具有復雜度低、便于工程應用的優點。通過模擬仿真實驗發現，模擬生成的信號在統計上滿足雙Johnson SB分布且具有指定的時間自相關，進而證明了本文所提方法的可靠性。本文所提方法生成的雙Johnson SB分布湍流擾動OAM 光通信自相關接收信號，可以為用蒙特卡洛方法分析大氣OAM 光通信系統性能提供支持。

圖6 不同條件下的仿真生成信號自相關函數與理論自相關函數的對比