この記事のまとめ：
- 判定の信頼性とは
- BPSKでの対数尤度比 (L値) の計算方法
- 多値変調での対数尤度比 (L値) の計算方法
- MatlabでのQPSK/16QAMの軟判定復調サンプルコードの紹介
背景
以前、ターボ符号器を紹介して、復号器も紹介するかもと書きましたが、ターボ復号を理解するためには対数尤度比による軟判定を理解しておく必要があるので、今回は軟判定復調 (Soft demapping)における対数尤度比の計算方法等を紹介しようと思います。
判定の信頼性
受信した変調信号 をビット情報に戻す際、硬判定した場合、0/1でしか表現できないため、その判定に対する信頼性は失われてしまいます。それに対して、軟判定であればどのくらい信頼のある判定なのか表現できます。その指標として、判定閾値と受信信号 との距離を使うことができます。等価と復号のパフォーマンスは軟判定によって飛躍的に向上できます。
信頼性の計測は対数尤度比 (Log-likelihood ratio: LLR) もしくは L値 (L-value) と呼ばれます。L値はバイナリ値 (つまり、0/1) をアナログ値 (”ソフトビット”、”ソフト値”) に拡張します。
それでは、以下に様々な変調信号のL値の計算方法について紹介します。
なお、この記事において、 は確率関数を、 は連続変数に対する確率密度関数を表します。（参考：確率の記法＠朱鷲の杜Wiki）
BPSKのL値
まず、前提として次のような通信路を想定します。
は送信信号、 は受信信号であり、AWGN通信路を想定しています。
BPSKにおいては、信号点が下図のように配置され、受信信号を硬判定する場合には、非常にシンプルで、受信信号が 0 を中心にそれより大きいのか小さいのかで判定すればよいです。
例えば、受信信号レベルが 0.01 だったときにその信号を 1 とビット判定するときと、受信信号レベルが 0.99 だった時にその信号を 1 とビット判定するのでは、どちらが信頼性が高いビット判定だったのかが硬判定では失われてしまいます。上記繰り返しになりますが、L値を使ってソフトビットとして軟判定することでどのくらい信頼のある判定なのかを表現できます。それでは、BPSKにおけるL値についてみていきます。
L値はBPSKだと次のように定義されます。
L値を求めるために、まず をベイズの定理で変形します。
そうするとL値は次のように変換できます。
なお、ランダムビット系列の場合には、 となるため、 となります。
また、AWGN通信路におけるノイズの分散を とした場合には、 は に従う正規分布となるため、 となります。このことから、冒頭で述べた「判定閾値（ここでは ）と受信信号 との距離」を表していることがわかります。
多値変調時のL値
シンボルあたり ビットにマッピングする場合、シンボルに対して 個の個別のL値を求めなければなりません。
従って、 は次のように表されます。
なお、 は ( ビット目が ) を有する 個 ( ビット目が で固定されているため、 になる) のシンボル に対する の集合です。つまり、それぞれ次の式で表されます。
また、 内の総和が意味することは、入力ビット系列 について、 である ビット目が であるすべての が対象となります。例えば、入力ビット長 としたときに対象となるビット系列は の4つです。 の場合、 となるのは であり、 のときのみということです。なお、入力ビット系列が のときには空和となります。従って、次のようになります。
これらをまとめると、 ビット目のL値は次のように表せます。
例：4-PAMのL値
4-PAM変調信号の最初のビット のL値について求めてみます。コンスタレーションと入力ビット と のマッピングは次の通りです。
と において、最初のビットが になり、 と において になります。従って、最初のビットのL値は次のように計算できます。
なぜこうなるのかもう少し詳しく見ていきます。
まず、 は次のように分解できます。
について見てみると、自然対数関数 の中の総和の対象となる入力ビット系列 は次のように表すことができます。
,
これをもとに、 の自然対数関数 の中の分子を見ると次のようになります。
同様に分母は次の通りです。
従って、これらをまとめると次のようになります。
なお、 です。
以上のことから、最初の計算結果になることがわかります。
QSPK、QAMへの拡張
基本的には4-PAMと大きな違いはなく、違いは虚数を扱うようになるだけです。具体的には、QPSKの場合を考えると、信号点配置が下図のように のように虚数が含まれるようになります。
これに基づいて、 を求めると、下記のようになります。
虚数を扱うことに伴って は次のように絶対値を以て、計算する必要があります。
違いは以上で、あとはPAMのときと同様に計算すればよいです。
Matlabサンプルコード
QSPK/16QAM変調のSoft demappingのサンプルコードを載せておきます。
上記の式に基づいてL値を計算したものと、Matlabの関数を使ってL値を計算したものの比較しています。なお、L値を計算した後で硬判定しているので、これだけではL値を計算した意味はないです。
また、Matlabの関数を使ったL値の計算の場合、分母が 1 、分子が 0 としており、上記で紹介した計算とは分母分子が逆となっていましたので、サンプルコードではMatlabに合わせる形に変更しました。
% Simulation parameters
snrRange = 0:2:16; % dB range SNR
numSimulation = 100; % number of simulations
K = 512; % number of bits being transmitted
modOrder = 16; % modulation order: 4 for QPSK, 16 for 16QAM
symMap = get_symbol_map(modOrder); % symbol mapping in integer
% Initialize variables for Error count
err1 = zeros(length(snrRange),1);
err2 = zeros(length(snrRange),1);
for snrIndex = 1:length(snrRange)
snr = snrRange(snrIndex);
fprintf('\nSNR = %d dB ', snr);
pause(0.001)
for simIndex = 1:numSimulation
fprintf('\b\b\b\b\b\b\b%3d/%3d', simIndex, numSimulation);
pause(0.001)
%% Transmitter
s = randi([0 1], K, 1);
x = qammod(s, modOrder, symMap, ...
'InputType', 'bit');
%% Channel
y = awgn(x, snr, 'measured'); %add noise according to SNR
%% Receiver
% Calculate noise variance
Eb = log2(modOrder);
N0 = Eb/(10^(snr/10));
sigma = N0/2; % Calculating SNR=1/sigma^2
% Matlab function
sEst = qamdemod(y, modOrder, symMap, ...
'NoiseVariance', sigma, ...
'OutputType', 'llr');
% Original implementation
sEst2 = demapper(y, modOrder, symMap, sigma);
%% Hard Detection
sEstHard = sEst; sEstHard(sEst>0)=0; sEstHard(sEst<=0)=1;
sEstHard2 = sEst; sEstHard2(sEst2>0)=0; sEstHard2(sEst2<=0)=1;
%% Error Count
err1(snrIndex) = err1(snrIndex) + sum(abs(sEstHard - s));
err2(snrIndex) = err2(snrIndex) + sum(abs(sEstHard2 - s));
end
end
ber1 = err1/K/numSimulation; %bit error
ber2 = err2/K/numSimulation; %bit error
figure(2)
% semilogy(snr,berapp),grid
semilogy(snrRange,ber1,snrRange,ber2),grid
title('Performance');xlabel('Eb/N0(dB)');ylabel('Bit Error Rate');
legend("Matlab func", "Implement")
%% Sub-fuctions
function sEst = demapper(y, modOrder, symMap, sigma)
binaryMap = de2bi(symMap);
tmp = binaryMap.';
constellation = qammod(tmp(:), modOrder, symMap, ...
'InputType', 'bit');
L = zeros(length(y),log2(modOrder)); % storage for L-value
for m = 1:log2(modOrder)
Sm1 = constellation(binaryMap(:,m)==1,:); % X_{m=1:Qm,1}
Sm0 = constellation(binaryMap(:,m)==0,:); % X_{m=1:Qm,0}
% Check for positions of Xm1 and Xm0
% hold on
% plot(Sm1, 'ro')
% plot(Sm0, 'bo')
% L(:,m) = log( sum(exp(-abs(y-Xm1.').^2/(2*sigma)), 2) ./ sum(exp(-abs(y-Xm0.').^2/(2*sigma)), 2) );
L(:,m) = log( sum(exp(-abs(y-Sm0.').^2/sigma), 2) ./ sum(exp(-abs(y-Sm1.').^2/sigma), 2) ); % MATLAB definition (see: https://jp.mathworks.com/help/comm/ug/digital-modulation.html#brc6yjx)
end
tmp = L.';
sEst = tmp(:);
end
function symMap = get_symbol_map(modulationOrder)
if modulationOrder==4
symMap = [ 2 3 ...
0 1 ];
elseif modulationOrder==16
symMap = [11 10 14 15 ...
9 8 12 13 ...
1 0 4 5 ...
3 2 6 7 ];
end
end
ちなみに実行結果は次のようになります。
前振りだけしていますが、気が向いたら畳み込み復号のアルゴリズムの記事を書こうと思います。
参考
本記事をまとめるにあたり下記ページを参考にさせていただきました。
今回は以上です。
最後まで読んでいただき、ありがとうございます。
関連記事
コメント
コメントを投稿