이전에 아래 포스팅에서 QAM constellation 에 대해 소개를 드렸었고~~

 

http://iamaman.tistory.com/205

 

아래 포스팅에서 16 QAM BER 시뮬레이션에 대해 소개한적이 있는데~

 

http://iamaman.tistory.com/120

 

간만에 64 QAM BER 시뮬레이션에 대해 소개 드리려 합니다.

 

생각도 안하고 있었던 내용인데 홍팡이라는 분이 부탁을 해서 하드를 뒤져 보니 있긴 있더군요.

 

이전에는 그냥 MATLAB Communications System Toolbox 의 이론적 BER 함수들을 사용했었는데 ~

 

오늘은 아래 주소에 있던 이론적 BER(Bit Error Rate) 수식을 사용했습니다.

 

http://www.raymaps.com/index.php/ber-64-qam-awgn/

 

코드는 아래 Zip 파일 다운로드 받아서 BER_simulation_64QAM.m 파일을 돌려 보시기 바랍니다.

 

64qam.zip

 

최신의 MATLAB 을  사용한다면 Randint()함수 관련 경고가 뜰텐데 무시하고 시뮬레이션 하시면 됩니다. 경고가 보기 싫다면 Randi() 함수로 변경해서 사용하면 됩니다.


다음과 같은 Constellation 과 64 QAM BER 결과를 확인 할 수 있습니다. 위 링크에 첨부되어 있는 결과와 마찬가지로 이론적인 성능과 거의 근사치의 시뮬레이션 결과가 나오는 것을 확인 할 수 있습니다.

 

 


MATLAB BER simulation 관련 검색을 하다가 아래 주소의 글을 발견했습니다.

 

http://www.mathworks.com/matlabcentral/fileexchange/22316-communication-systems-reference-curves

 

위 파일에서 설명하는 내용은 아래와 같습니다.

 

- PSK and QAM over AWGN Channel (BER and SER)

- BPSK over Rayleigh fading channel (BER)

- Convolutional Coded BPSK over AWGN (BER)

 

일단 BER 이나 SER 시뮬레이션을 어떻게 해야 할 지 모르는 분들이 보기에 참 좋은 코드라는 생각이 들더군요.

 

위 주소에 들어가서 아래 그림과 같이 Download Submission 버튼을 누르면 파일을 다운로드 할 수 있습니다.

 

위 코드가 실행되기 위해서는 MATLAB 뿐만 아니라 Simulink, Communications Toolbox, Communications Blockset 이 필요 합니다.

 

다운로드 받은 Submit_v2.zip 파일을 압축해제 하고~ RunMe.m 파일을 실행시키면 아래 그림과 같이 웹 화면이 나옵니다.

 

웹 화면에서 뒤로 가기를 누르면 정상적으로 실행 안 되는 경우가 있으므로 이전 페이지의 내용을 다시 돌려 보고 싶다는 분들은 RunMe.m 파일을 다시 실행하고 들어가서 실행시키면 됩니다.

 

다양한 내용들이 있는데~ 저는 그 중에서 Rayleigh Flat Fading Channel Bit Error Rate curves 만 한번 실행해 봤습니다.

 

아래 그림과 같이 Simulink 파일이 실행되면서 시뮬레이션이 진행 됩니다.

 

시뮬레이션이 진행 될 때마다~ 아래 그림과 같이 BER 그래프에 결과가 하나씩 추가되며~

 

맨 마지막에는 다음과 같은 결과를 확인 할 수 있습니다.

 

BER 또는 SER 시뮬레이션을 해 보려고 하시는 엔지니어 분들이나 디지털 통신을 공부하시는 학생 분들에게 대단히 좋은 예가 될 것 같네요~


아래 포스팅에서 BPSK, QPSK BER(Bit error rate) simulation in AWGN channel 에 대해 설명 드렸었는데~

   

2011/03/27 - [programming language/MATLAB] - MATLAB QPSK BER simulation in AWGN channel


2011/03/19 - [programming language/MATLAB] - MATLAB BPSK BER simulation in AWGN channel

 

8 PSK(Phase-shift keying) 시뮬레이션에 대해 질문하신 분이 있어서 답변 드립니다.

 

먼저 제 경험을 말씀 드리면, 8 PSK 나 16 PSK 의 경우 책에서 공부한적은 있지만 실제 통신 시스템에서 사용하는 것을 본적은 없는 것 같습니다. 사용하지 않는 이유에 대해서는 아래 글을 참조해 보시기 바랍니다. 아래 글은 16 QAM 과 16 PSK 의 SER (Symbol Error Rate)성능 비교에 대한 내용인데 16 PSK 가 16 QAM 대비하여 확실히 SER 성능이 안 좋다는 것을 확인 할 수 있습니다.

 

http://www.dsplog.com/2008/03/29/comparing-16psk-vs-16qam-for-symbol-error-rate/

 

일단 QPSK, BPSK 의 경우 각 constellation 지점들이 90 도 또는 180 도 위상차이를 보이기 때문에 심볼에 대해 위상을 구할 필요가 없이 BPSK 의 경우 간단하게 x=0 지점을, QPSK 의 경우 x=0, y=0 지점을 기준으로 심볼 demapping 을 할 수 있었습니다.

 

하지만 8 PSK 의 constellation 은 아래 그림과 같이 8개 지점(45도 간격)으로 구성되어 있습니다. 이런 경우에는 당연히 45/2 = 22.5 도 간격으로 symbol demapping 을 해야 하므로 수신 심볼의 위상을 구해야 할 것입니다.

 

이미지 출처 : http://commons.wikimedia.org/wiki/File:8PSK_Gray_Coded.svg

 

PSK 시뮬레이션을 위해서는 PSK 의 Gray code 에 대해 공부를 하셔야 하고~ 아래 주소에 굉장히 좋은 설명이 있더군요. MATLAB/OCTAVE 코드도 있으니 실습해 보시면서 천천히 공부해 보시기 바랍니다.

 

http://www.dsplog.com/2008/05/12/gray-code-to-binary-conversion-for-psk-pam/

 

http://www.dsplog.com/2008/05/11/binary-to-gray-code-conversion-psk-pam/

 

다음으로 질문하신 8 PSK BER 시뮬레이션에 대해 소개하겠습니다.

 

아래 주소에 16 PSK BER 시뮬레이션에 대해 소개하는데~

 

http://www.dsplog.com/2008/05/18/bit-error-rate-for-16psk-modulation-using-gray-mapping/

 

코드는 아래 주소에서 다운로드 받을 수 있습니다.

 

http://www.dsplog.com/db-install/wp-content/uploads/2008/05/script_16psk_gray_mapping_bit_error_rate.m

 

위 코드를 실행해보면~ 16 PSK 뿐만 아니라 8 PSK, QPSK도 시뮬레이션 해 볼 수 있습니다.

 

위 주소에서 설명한 16 PSK 부터 시뮬레이션을 해보죠~ 아래 코드의 모든 권리는 http://www.dsplog.com 사이트의 Krishna Pillai 에게 있습니다. 재 배포를 하실 경우 원저자를 명시 하시기 바랍니다.

 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% All rights reserved by Krishna Pillai, http://www.dsplog.com

% The file may not be re-distributed without explicit authorization

% from Krishna Pillai.

% Checked for proper operation with Octave Version 3.0.0

% Author    : Krishna Pillai

% Email        : krishna@dsplog.com

% Version    : 1.0

% Date        : 17 May 2008

% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 

% Bit Error Rate for 16-PSK modulation using Gray modulation mapping

 

clear

N = 10^5; % number of symbols

M = 16; % constellation size

k = log2(M); % bits per symbol

 

 

thetaMpsk = [0:M-1]*2*pi/M; % reference phase values

 

Eb_N0_dB = [0:25]; % multiple Es/N0 values

Es_N0_dB = Eb_N0_dB + 10*log10(k);

 

% Mapping for binary <--> Gray code conversion

ref = [0:M-1];

map = bitxor(ref,floor(ref/2));

[tt ind] = sort(map);

 

ipPhaseHat = zeros(1,N);

for ii = 1:length(Eb_N0_dB)

 

% symbol generation

% ------------------

ipBit = rand(1,N*k,1)>0.5; % random 1's and 0's

bin2DecMatrix = ones(N,1)*(2.^[(k-1):-1:0]) ; % conversion from binary to decimal

ipBitReshape = reshape(ipBit,k,N).'; % grouping to N symbols having k bits each

ipGray = [sum(ipBitReshape.*bin2DecMatrix,2)].'; % decimal to binary

 

% Gray coded constellation mapping

ipDec = ind(ipGray+1)-1; % bit group to constellation point

ipPhase = ipDec*2*pi/M; % conversion to phase

ip = exp(j*ipPhase); % modulation

s = ip;

 

% noise

% -----

n = 1/sqrt(2)*[randn(1,N) + j*randn(1,N)]; % white guassian noise, 0dB variance

 

y = s + 10^(-Es_N0_dB(ii)/20)*n; % additive white gaussian noise

 

% demodulation

% ------------

% finding the phase from [-pi to +pi]

opPhase = angle(y);

% unwrapping the phase i.e. phase less than 0 are

% added 2pi

opPhase(find(opPhase<0)) = opPhase(find(opPhase<0)) + 2*pi;

 

% rounding the received phase to the closest constellation

ipPhaseHat = 2*pi/M*round(opPhase/(2*pi/M))    ;

% as there is phase ambiguity for phase = 0 and 2*pi,

% changing all phases reported as 2*pi to 0.

% this is to enable comparison with the transmitted phase

ipPhaseHat(find(ipPhaseHat==2*pi)) = 0;

ipDecHat = round(ipPhaseHat*M/(2*pi));

 

% Decimal to Gray code conversion

ipGrayHat = map(ipDecHat+1); % converting to decimal

ipBinHat = dec2bin(ipGrayHat,k) ; % decimal to binary

 

% converting binary string to number

ipBinHat = ipBinHat.';

ipBinHat = ipBinHat(1:end).';

ipBinHat = str2num(ipBinHat).' ;

 

% counting errors

nBitErr(ii) = size(find([ipBit- ipBinHat]),2); % couting the number of errors

 

end

simBer = nBitErr/(N*k);

theoryBer = (1/k)*erfc(sqrt(k*10.^(Eb_N0_dB/10))*sin(pi/M)); % PSK 의 이론적 BER

 

close all; figure

semilogy(Eb_N0_dB,theoryBer,'bs-','LineWidth',2);

hold on

semilogy(Eb_N0_dB,simBer,'mx-','LineWidth',2);

axis([0 20 10^-5 1])

grid on

legend('theory', 'simulation');

xlabel('Eb/No, dB')

ylabel('Bit Error Rate')

title('Bit error probability curve for 16-PSK modulation')

 

위 코드를 돌려보면~ 아래 그림과 같이 16 PSK 결과가 나옵니다. 아래 BER 그래프에서 16 dB 이상에서 이론 값과 약간 달라 보이는 것은 시뮬레이션 횟수가 부족해서 BER 이 0 으로 나오는 경우가 생겨서 그런 것입니다. 고 dB 에서도 이론값과 같은 결과를 보이고 싶다면 심볼수인 N = 10^5 값을 좀더 늘려 보시기 바랍니다. 

 

위 코드를 약간 수정해서~질문하신 8 PSK 시뮬레이션을 해보죠~ 위 코드에서 빨간색 표시한 부분들만 아래와 같이 변경합니다.

 

M = 8; % constellation size

Eb_N0_dB = [0:15]; % multiple Es/N0 values

axis([0 15 10^-5 1])

title('Bit error probability curve for 8-PSK modulation')

 

시뮬레이션을 해보시면~ 다음과 같이 8 PSK 결과를 확인 할 수 있습니다.

 

다음과 같이 변경하면 QPSK 시뮬레이션이 되는거죠~

 

M = 4; % constellation size

Eb_N0_dB = [0:10]; % multiple Es/N0 values

axis([0 10 10^-5 1])

title('Bit error probability curve for QPSK modulation')

 


  1. 2014.04.17 15:33

    비밀댓글입니다

    • 남성 2014.04.17 18:23 신고

      네~ 방문해 주셔서 감사합니다. ^^ 열공하세요~

  2. 곰비 2014.05.12 00:34

    QPSK의 BER이론치를 구하는 공식을 알수잇을까요.??ㅠㅠㅠ

    • 남성 2014.05.12 08:02 신고

      PSK 의 이론적 BER 이라고 되 있는 부분의 코드를 보시면 수식을 아실 수 있을 겁니다.

  3. yj 2017.09.25 20:31

    BPSK는 해당 이론식에 적용되지 않는건가요?

    • 남성 2017.09.26 04:21 신고

      아니요 위 포스팅에서는 그냥 시뮬레이션 안 한 것뿐입니다. 똑같이 적용 됩니다

  4. 궁금해 2018.11.10 00:54

    이론값과 실제값이 차이나는게 심볼 수 때문이라고 하는데 왜 그런 차이가 나는건가요??
    오차가 생기는 원인이 궁금합니다ㅠㅠ

    • 남성 2018.11.10 21:43 신고

      확률적인 시뮬레이션이니까 시뮬레이션 샘플수가 많을수록 이론적인 값과 근사하게 되겠죠.

오늘은 MATLAB 을 이용한 16 QAM BER Simulation코드에 대해 소개해 보려 합니다.

 

예전에 제가 작성했던 코드들도 있지만~

 

아래 주소를 보니 16 QAM 의 Gray coding 부터 이론적인 BER 성능 까지 자세히 설명되어 있더군요.

   

http://www.dsplog.com/2008/06/05/16qam-bit-error-gray-mapping/

   

코드는 위 블로그의 약간 아래 쪽을 보시면~ 링크가 되어 있습니다.

 

못찾으실 분들을 위해 링크를 걸죠, 아래 주소를 오른쪽 클릭한 후에 파일로 다운로드 받거나 그냥 클릭하고 들어가서 전체 선택후에 m 파일에 붙여넣기 해도 됩니다.

 

http://www.dsplog.com/db-install/wp-content/uploads/2008/06/script_16qam_gray_mapping_bit_error_rate.m

 

혹시 못 찾으실 분들을 위해 코드도 넣습니다. 다시 한번 밝히지만 아래 코드는 Krishna Pillai, http://www.dsplog.com 님이 작성하신 겁니다.

 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% All rights reserved by Krishna Pillai, http://www.dsplog.com

% The file may not be re-distributed without explicit authorization

% from Krishna Pillai.

% Checked for proper operation with Octave Version 3.0.0

% Author    : Krishna Pillai

% Email        : krishna@dsplog.com

% Version    : 1.0

% Date        : 05 June 2008

% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 

% Bit Error Rate for 16-QAM modulation using Gray modulation mapping

 

clear

N = 10^5; % number of symbols

M = 16; % constellation size

k = log2(M); % bits per symbol

 

% defining the real and imaginary PAM constellation

% for 16-QAM

alphaRe = [-(2*sqrt(M)/2-1):2:-1 1:2:2*sqrt(M)/2-1];

alphaIm = [-(2*sqrt(M)/2-1):2:-1 1:2:2*sqrt(M)/2-1];

k_16QAM = 1/sqrt(10);

 

Eb_N0_dB = [0:15]; % multiple Es/N0 values

Es_N0_dB = Eb_N0_dB + 10*log10(k);

 

% Mapping for binary <--> Gray code conversion

ref = [0:k-1];

map = bitxor(ref,floor(ref/2));

[tt ind] = sort(map);

 

for ii = 1:length(Eb_N0_dB)

 

% symbol generation

% ------------------

ipBit = rand(1,N*k,1)>0.5; % random 1's and 0's

ipBitReshape = reshape(ipBit,k,N).';

bin2DecMatrix = ones(N,1)*(2.^[(k/2-1):-1:0]) ; % conversion from binary to decimal

% real

ipBitRe = ipBitReshape(:,[1:k/2]);

ipDecRe = sum(ipBitRe.*bin2DecMatrix,2);

ipGrayDecRe = bitxor(ipDecRe,floor(ipDecRe/2));

% imaginary

ipBitIm = ipBitReshape(:,[k/2+1:k]);

ipDecIm = sum(ipBitIm.*bin2DecMatrix,2);

ipGrayDecIm = bitxor(ipDecIm,floor(ipDecIm/2));

% mapping the Gray coded symbols into constellation

modRe = alphaRe(ipGrayDecRe+1);

modIm = alphaIm(ipGrayDecIm+1);

% complex constellation

mod = modRe + j*modIm;

s = k_16QAM*mod; % normalization of transmit power to one

 

% noise

% -----

n = 1/sqrt(2)*[randn(1,N) + j*randn(1,N)]; % white guassian noise, 0dB variance

 

y = s + 10^(-Es_N0_dB(ii)/20)*n; % additive white gaussian noise

 

% demodulation

% ------------

y_re = real(y)/k_16QAM; % real part

y_im = imag(y)/k_16QAM; % imaginary part

 

% rounding to the nearest alphabet

ipHatRe = 2*floor(y_re/2)+1;

ipHatRe(find(ipHatRe>max(alphaRe))) = max(alphaRe);

ipHatRe(find(ipHatRe<min(alphaRe))) = min(alphaRe);

ipHatIm = 2*floor(y_im/2)+1;

ipHatIm(find(ipHatIm>max(alphaIm))) = max(alphaIm);

ipHatIm(find(ipHatIm<min(alphaIm))) = min(alphaIm);

 

% Constellation to Decimal conversion

ipDecHatRe = ind(floor((ipHatRe+4)/2+1))-1; % LUT based

ipDecHatIm = ind(floor((ipHatIm+4)/2+1))-1; % LUT based

 

% converting to binary string

ipBinHatRe = dec2bin(ipDecHatRe,k/2);

ipBinHatIm = dec2bin(ipDecHatIm,k/2);

 

% converting binary string to number

ipBinHatRe = ipBinHatRe.';

ipBinHatRe = ipBinHatRe(1:end).';

ipBinHatRe = reshape(str2num(ipBinHatRe).',k/2,N).' ;

 

ipBinHatIm = ipBinHatIm.';

ipBinHatIm = ipBinHatIm(1:end).';

ipBinHatIm = reshape(str2num(ipBinHatIm).',k/2,N).' ;

 

% counting errors for real and imaginary

nBitErr(ii) = size(find([ipBitRe- ipBinHatRe]),1) + size(find([ipBitIm - ipBinHatIm]),1) ;

 

end

simBer = nBitErr/(N*k);

theoryBer = (1/k)*3/2*erfc(sqrt(k*0.1*(10.^(Eb_N0_dB/10))));

 

close all; figure

semilogy(Eb_N0_dB,theoryBer,'bs-','LineWidth',2);

hold on

semilogy(Eb_N0_dB,simBer,'mx-','LineWidth',2);

axis([0 15 10^-5 1])

grid on

legend('theory', 'simulation');

xlabel('Eb/No, dB')

ylabel('Bit Error Rate')

title('Bit error probability curve for 16-QAM modulation')

 

실행 시키고 조금만 기다려 보면~~ Krishna Pillai의 블로그에 있는것과 같이 다음과 같은 시뮬레이션 결과를 확인 할 수 있습니다.

 

 

위 코드를 배포 하실 때는 원저자를 명시하셔야 합니다.


  1. 홍팡 2015.07.31 17:20

    안녕하세요
    매틀랩 공부를 하다가 우연히 들어오게 되었습니다!
    혹시 64QAM도 시뮬레이션 돌리려면 어디어디 수정해야하는지 알수있을까요?ㅜ

    • 남성 2015.08.01 13:45 신고

      간단하게 말씀 드리면 mapping 과 demapping 부분 그리고 파워 계산 부분이 바뀌어야하는데 위 코드에서 바꾸려면 대부분 다 바뀌어야 할 것 같네요.

  2. 홍팡 2015.08.03 15:49

    죄송합니다만 조금만 더 알려주실 수 있을까요??
    16QAM 에 대해 도움을 많이 얻어서 감사합니다

    • 남성 2015.08.03 19:39 신고

      for 문 내의 시뮬레이션 과정은 반복적으로 비트를 생성해서 16qam 으로 변조하고 잡음을 넣어주고 복조하는 과정인데 결국에는 이 내부가 대부분 바뀌어야 할거라는 말씀 입니다.

      아래 주소의 글을 보시면 constellation 을 확인 할 수 있는데 64 qam constellation 대로 비트를 심볼로 맵핑하고 수신단에서는 심볼을 다시 비트로 demapping 하는 과정이 필요 할 겁니다.

      http://iamaman.tistory.com/205

      이 과정을 생각하면 위 코드의 많은 부분이 바뀌게 되어 어디 한군데만 바뀌면 된다라고 하기가 힘들것 같습니다.

    • 홍팡 2015.08.04 16:04

      그렇군요.. 감사합니다.

      마지막으로 죄송합니다만 혹시 64QAM 코드도 가지고 계시다면
      공유 부탁드려도 될까요?? 제 비루한 코딩실력으로 변경하기가 쉽지가 않아서요.
      감사합니다.

    • 남성 2015.08.04 19:56 신고

      아래 주소에 포스팅 했습니다.

      http://iamaman.tistory.com/1631

      자주 방문해 주세요. ^^

  3. 홍팡 2015.08.09 16:50

    아 정말 감사합니다 많은 도움이 되었습니다!

MATLAB 을 이용한 통신 시뮬레이션 코드를 찾다 보니 아래 주소의 블로그를 발견했다.

 

http://www.dsplog.com/

 

대단히 자세한 이론 설명과 훌륭한 MATLAB 예제 코드들을 보면 디지털 통신을 공부하는데 정말 부족함이 없어 보인다.

 

학부 시절에 FSK(Frequency Shift Keying)를 공부하면서~ BER(Bit Error Rate) 시뮬레이션 까지는 못해봤던 것 같은데~

 

아래 주소에 Frequency Shift Keying BER 시뮬레이션 코드가 있었다.

 

http://www.dsplog.com/2007/08/30/bit-error-rate-for-frequency-shift-keying-with-coherent-demodulation/

 

이론적 설명도 자세하고 세부 코드까지 포함되어 있어서 굉장히 도움이 많이 되는 포스팅이라는 생각이 든다.

 

아래 코드는 저자의 시뮬레이션 코드에 스펙트럼 확인을 하기 위해 몇 가지 코드를 추가한 예이다. 내가 추가한 부분은 빨간색으로 표시했다.

 

아래 코드에서 T 변수의 주석에 symbol duration 이라고 되어 있어서 약간 헷갈렸지만~ 변수 T 가 Sampling Rate(Fs)이다.

 

암튼 저자의 FSK simulation 코드를 보면 Coherent FSK 복조까지 명확하게 알 수 있어서 FSK 를 공부하는 분들이 참조하기에 대단히 좋은 코드라는 생각이 든다.

 

아래 코드에서 내가 스펙트럼을 확인하기 위해 사용한 pwelch() 함수를 사용하기 위해서는 MATLAB Signal Processing Toolbox 가 있어야 한다.

 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% All rights reserved by Krishna Pillai, http://www.dsplog.com

% The file may not be re-distributed without explicit authorization

% from Krishna Pillai.

% Checked for proper operation with Octave Version 3.0.0

% Author    : Krishna Pillai

% Email        : krishna@dsplog.com

% Version    : 1.0

% Date        : 05 June 2008

% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 

% Simple Matlab example of binary Frequency Shift Keying using

% coherent demodulation

 

clear

N = 10^5 % number of bits or symbols

T = 8; % symbol duration

Fs=T; % Sampling Rate

t = [0:1/T:0.99]; % sampling instants

tR = kron(ones(1,N),t); % repeating the sampling instants

 

Eb_N0_dB = [0:11]; % multiple Eb/N0 values

 

for ii = 1:length(Eb_N0_dB)

% generating the bits

ip = rand(1,N)>0.5; % generating 0,1 with equal probability

freqM = ip+1; % converting the bits into frequency, bit0 -> frequency of 1, bit1 -> frequency of 2

freqR = kron(freqM,ones(1,T)); % repeating

x = (sqrt(2)/sqrt(T))*cos(2*pi*freqR.*tR); %generating the FSK modulated signal

 

% noise

n = 1/sqrt(2)*[randn(1,N*T) + j*randn(1,N*T)]; % white gaussian noise, 0dB variance

 

% coherent receiver

y = x + 10^(-Eb_N0_dB(ii)/20)*n; % additive white gaussian noise

op1 = conv(y, sqrt(2/T)*cos(2*pi*1*t)); % correlating with frequency 1

op2 = conv(y, sqrt(2/T)*cos(2*pi*2*t)); % correlating with frequency 2

 

% demodulation

ipHat = [real(op1(T+1:T:end)) < real(op2(T+1:T:end))]; %

nErr(ii) = size(find([ip - ipHat]),2); % counting the number of errors

 

end

simBer = nErr/N;

theoryBer = 0.5*erfc(sqrt((10.^(Eb_N0_dB/10))/2)); %theoretical BER

 

close all

figure

semilogy(Eb_N0_dB,theoryBer,'b-');

hold on

semilogy(Eb_N0_dB,simBer,'mx-');

axis([0 11 10^-4 0.5])

grid on

legend('theory:fsk-coh', 'sim:fsk-coh');

xlabel('Eb/No, dB')

ylabel('Bit Error Rate')

title('Bit error probability curve')

 

 

%%

 

figure

pwelch(y,1024,512,1024,Fs)

 

위 MATLAB 코드를 실행시켜보면 다음과 같이 Coherent Frequency Shift Keying BER 과 FSK Power Spectral Density 그래프를 확인 할 수 있다.


< Coherent Frequency Shift Keying BER >

< Coherent FSK Power Spectral Density >

 

정확히 주파수 1, 2 Hz 에 톤이 뜨는 것을 확인 할 수 있다.


  1. Genius 2014.02.03 13:53

    좋은 자료 감사합니다.
    담아갈께요. ^^

    • 남성 2014.02.04 10:38 신고

      방문 감사합니다. 공지에 밝혔다 시피 링크만 허용합니다.

  2. ryu 2014.06.19 17:29

    conv에서 t 설정은 어떻게 해주는건지 잘 모르겠어요ㅠㅠ

    • 남성 2014.06.20 19:58 신고

      질문이 무슨 뜻인지 잘 몰겠네여... 위에 보시다시피 t = [0:1/T:0.99];

      t 는 1/T 간격으로 해 준건데요....

  3. chobo 2016.06.10 21:35

    안녕하세요 좋은자료 잘봣습니다.
    혹지 저기서 fs 는 뭘 말하는 건가요?
    FSK에서의 sampling rate가 의미하는걸 잘 모르겟어서..

    • 남성 2016.06.12 21:26 신고

      FSK 에서의 의미라기 보다는 디지털에서의 의미이라고 볼 수 있을것 같네요. 초당 Fs 샘플이 있다는 거죠.

이번 포스팅에서는 AWGN (Additive white gaussian noise) 채널 환경에서의 QPSK Bits Error Rate(BER)에 성능 검증 simulation 에 대해 설명한다.

 

QPSK 및 잡음의 분산(σ2) 과 No(noise power spectral density) 의 관계에 대한 내용은 다음 포스팅을 참조 하기 바란다.




2011/03/19 - [MATLAB] - MATLAB BPSK BER simulation in AWGN channel


2011/03/08 - [통신] - [디지털 통신] BPSK, QPSK

 

QPSK 에서는 아래 식의 M =4 이다.

 

 

따라서 다음 식과 같이 정리 된다.

 

 

QPSK 시뮬레이션 시의 주의할 점신호의 파워를 1 로 정규화 해 줘야 하며 잡음 역시 real , imaginary 에 대해 각각 독립적이어야 한다.

 

이제 QPSK 시뮬레이션을 해 보자.

 

  

Command 창의 결과는 다음과 같다.

Eb/No= 0 BER: 0.078801

Eb/No= 2 BER: 0.0374294

Eb/No= 4 BER: 0.0124703

Eb/No= 6 BER: 0.00238721

Eb/No= 8 BER: 0.0001905

Eb/No= 10 BER: 3.07143e-006

 

 

위 그래프를 보면 이론적인 결과와 근사적으로 일치 함을 확인 할 수 있다.


  1. 남성 2012.11.25 18:47 신고

    이렇게 바꿔보세요.

    Symbols=((bits_v(:,1)*(-2)+1) + i*(bits_v(:,2)*(-2)+1))./sqrt(2);

  2. 과객 2013.06.18 23:52

    BPSK QPSK 정말 감사하게 잘 실행되었습니다.
    다만 QPSK에서는 비트 2개를 하나는 리얼로 하나는 이미지로 생각해서 처리하셨는데 8PSK에서는 비트 3개라 어떻게 해야 할지 모르겠네요 8PSK에서는 어떻게 해야 할까요?

    • 남성 2013.06.19 06:45 신고

      8 psk 에서도 3 비트를 한 점의 symbol 로 mapping 하시면 됩니다.

      아래 주소의 그림처럼 말이죠.

      http://commons.wikimedia.org/wiki/File:8PSK_Gray_Coded.svg


      방문해 주셔서 감사합니다. ^^

    • 과객 2013.06.20 20:11

      감사합니다.

  3. 답변부탁드릴게요 2013.11.17 23:06

    도움 감사합니다 하지만 실행이 되지 않아서 질문드립니다
    답변해주시면 정말 감사하겠습니다.

    clc
    clear
    close all

    N_bits=70000000; %생성할 비트 수 / 2
    Eb_No_dB=0:2:10; %Eb_No dB scale
    BER_buffer=zeros(size(Eb_No_dB)); %BER 저장할 버퍼


    for n=1:length(Eb_No_dB)
    Eb_No_ral_scale=10^(Eb_No_dB(n)/10); %Eb/No real scale 값 반환
    sigma_v=sqrt(1./(4*Eb_No_ral_scale)); %잡음의 표준 편차 값 계산

    bits_v=randi([0 1],N_bits,2); %비트 생성

    Symbols=(bits_v(:,1)*(-2)+1 + i *(bits_v(:,2)*(-2)+1))./sqrt(2); %symbol mapping & normalizat

    noise_v=(randn(N_bits,1)*sigma_v)+ i * (randn(N_bits,1)*sigma_v); % complex AWGN 생성

    tx_signal = Symbols * noise_v; % 신호*잡음

    demapped_bits=[real(tx_signal) imag(tx_signal)] < 0; %수신단 demapping

    BER_buffer(n)=sum(bits_v(:) ~= demapped_bits(:)) ./ (N_bits*2); %BER calculation

    fprintf('Eb/No=%g BER: %g\n' ,Eb_No_dB(n) ,BER_buffer(n))
    end

    %%이론적인 QPSK BER 계산 Using communication toolbox

    Eb_No_theory = 0:0.05:10;
    Ber_value_theory = berawgn(Eb_No_theory,'psk',4,'nondiff');
    %% 그래프그리기

    figure(1),semilogy(Eb_No_dB,BER_buffer,'b+',Eb_No_theory,Ber_value_theory,'r--'),x;xlabel('Eb/No [dB]'),ylabel('BER')
    grid.axis([0 10 0.000001 1]), legend('Simulated result','Theoretical BPSK result')


    제시해주신 코드와 똑같이 작성하려고 노력하였습니다,
    하지만 아래와 같은 오류가 나와서 질문드립니다.

    ??? Error using ==> randn
    Out of memory. Type HELP MEMORY for your options.

    Error in ==> QPSK_BER at 18
    noise_v=(randn(N_bits,1)*sigma_v)+ i * (randn(N_bits,1)*sigma_v); % complex
    AWGN 생성

    어떤 문제인지 알려주시면 정말 감사하겠습니다!

    • 남성 2013.11.17 23:59 신고

      randn() 에서 메모리 에러 나는것은 N_bits 를 너무 크게 줘서 그렇습니다.

      아래 세 부분이 문제더군요.

      다음과 같이 바꾸세요.

      N_bits=7000000; %생성할 비트 수 / 2, 너무 크게 주면 컴퓨터뻗어요.

      tx_signal = Symbols + noise_v; % 신호 + 잡음, 곱하기가 아닙니다.


      % 아래 부분에서는 단순 오타가 좀 있더군요.
      figure(1),semilogy(Eb_No_dB,BER_buffer,'b+',Eb_No_theory,Ber_value_theory,'r--');xlabel('Eb/No [dB]'),ylabel('BER')
      grid on , axis([0 10 0.000001 1]), legend('Simulated result','Theoretical BPSK result')

    • 답변부탁드릴게요 2013.11.18 01:06

      감사합니다! 정상작동됩니다!

    • 남성 2013.11.18 01:39 신고

      네~ 열심히 공부하세욧 !! 방문해 주셔서 감사합니다. ^^

  4. 답변좀 부탁드려요 2013.11.20 17:06

    안녕하세요~ 과제 때문에 들렀다가 오류가 떠서 질문드립니다..ㅠ
    Eb_N0_theory=0:0.05:10 ; 이 부분에서요
    다음과 같은 오류가 뜹니다..
    The expression to the left of the equals sign is not a valid target for an assignment.
    Eb_N0 대신 혹시나 해서 Eb_No로 바꿔서 입력해보아도 같은 오류 메시지가 뜨는데 어떻게 해야 하나요?

    • 남성 2013.11.20 18:30 신고

      Eb_N0_theory=0:0.05:10 ; 부분에 코드상의 문제는 없습니다.
      다른 부분에 문제가 있는것 같은데요.

  5. 119 2014.04.01 17:59

    음.. ㅠㅠ 저기 tx_signal 차원이 [N_bits,1]맞죠?? symbol하구 noise 왜 더하는데 둘이 자꾸 차원이 안맞다구 나오는지 ㅠㅠ

  6. 119 2014.04.01 19:58

    저는 sigma(n)으로 고쳐서 곱해주니까다시됏어요.. 그리구 demappping이랑 BER계산 어떻게 하신건가요ㅠㅠ??

    • 남성 2014.04.01 21:23 신고

      demapping은 0 보다 작으면 1 그렇지 않으면 0 으로 demapping 한 거구요. 이론적 ber 게산은 berawgn 이라는 통신 툴박스의 함수를 사용했습니다. 위에 보시다 시피 코드는 저게 다에요.

  7. 2014.04.16 23:14

    비밀댓글입니다

    • 남성 2014.04.17 01:06 신고

      아래 포스팅에 답변 합니다.

      http://iamaman.tistory.com/1251

  8. 공돌공돌 2014.06.06 22:06

    clc
    clear all
    close all

    N_bits=7000000; %생성할 비트 수/2
    Eb_No_dB=0:2:10; %Eb/No dB scale
    BER_buffer=zeros(size(Eb_No_dB)); %BER 저장할 버퍼

    for n=1:length(Eb_No_dB)
    Eb_No_ral_scale=10^(Eb_No_dB(n)/10); %Eb/No real scale 값 변환
    sigma_v=sqrt(1./(4*Eb_No_ral_scale)); %잡음의 표준 편차 값 계산

    bits_v=randi([0 1],N_bits,2); %비트 생성

    Symbols=(bits_v(:,1)*(-2)+1 +i*(bits_v(:,2)*(-2)+1))./sqrt(2); %symbol mapping & normalization

    noise_v=(randn(N_bits,1)*sigma_v) +i*(randn(N_bits,1)*sigma_v); %complex AWGN 생성

    tx_signal=Symbols+noise_v; %신호 + 잡음

    demapped_bits=[real(tx_signal) imag(tx_signal)]<0; %수신단 demapping

    BER_buffer(n)=sum(bits_v(:)~=demapped_bits(:))./(N_bits*2); %BER calculation

    fprint('Eb/No=%g BER:%g\n',Eb_No-dB(n), BER_buffer(n))
    end

    % 이론적인 QPSK BER 계산
    Eb_No_theory=0:0.05:10;
    Ber_value_theory=berawgn(Eb_No_theory,'psk',4,'nondiff');

    % 그래프 그리기
    figure(1),semilogy(Eb_No_dB,BER_buffer,'b+',Eb_No_theory,Ber_value_theory,'r--');xlabel('Eb/No[dB]'), ylabel('BER')
    grid on, axis([0 10 0.000001 1]), legend('Simulated result','Theoritical QPSK result')

    이렇게 했는데 왜 BERTool에서 run 누르면 invalid simulation M-file or model 이라는 에러가 나오면서 그림이 안그려지죠? ㅜㅜ 제발 알려주세요 ㅜㅜㅜㅜ

  9. 남성 2014.06.07 06:07 신고

    위 코드에서는 25 번째 줄에 오타가 몇개 있더군요. 아래와 같이 하면 정상적으로 동작 합니다.

    fprintf('Eb/No=%g BER:%g\n',Eb_No_dB(n), BER_buffer(n))

    위 질문에서 bertool 이라고 얘기하시는데... 위 코드를 m 파일로 작성하시고 bertool 에서 monte carlo simulation 하셨다는 건가요?

    • 남성 2014.06.07 10:34 신고

      먼가 잘 못 이해하고 계신거 같은데... 위 코드는 그냥 m 파일로 돌리면 되는건데요? bertool 안 써도 되는겁니다.

    • 공돌공돌 2014.06.07 10:35

      25번째줄 말씀하신대로 수정하고 저장한 후에 bertool 켜서 monte carlo simulation 했는데 또 invalid simulation M-file or model 이라는 에러가 뜨네요 ㅜㅜ 대체 왜 이러는지 모르겠습니다 ㅜ 위에 님이 포스팅하신 것처럼 이론적인 QPSK 그래프와 코드 시뮬레이션 한 값과 비교하는 그래프를 봐야하는데 그게 안뜹니다 ㅜㅜ

    • 2014.06.07 10:38

      비밀댓글입니다

    • 남성 2014.06.07 11:51 신고

      bertool을 사용하신다면 변수명 등을 bertool 에 맞게 써 줘야 하는것 같더군요.

  10. 공돌이 2014.06.08 15:36

    안녕하세요 자료를 찾아보다가 오게 되어서 질문드립니다 ㅠ

    제가 지금 설계프로젝트를 하고 있는데요 32QAM을 써서 통신시스템 설계를 진행하고 있는데 awgn만 첨가했을 때

    bertool이랑은 전혀 다른 결과가 나와서 문의드립니다.

    들어가는 비트 Eb - 채널코딩 한 후 Ec - 변조 한 후 Es로 해서

    N0=10.^(EcN0dB(snrindex)/10);
    gain=sqrt(N0/2);
    noise=randn(200,1)*gain

    으로 했는데 bertool보다 훨씬 더 좋게 나오네요,,, 이러면 안되는데,,, awgn 채널을 잘못구성한건가요??

    • 남성 2014.06.08 17:08 신고

      채널 코딩도 들어갔으면 파워 계산 과정에서 채널 코드 레이트에 대한 수치가 반영 되어야 합니다. 채널 코딩시 encoder 의 입력 비트에 비해 출력비트가 더 많이 나오므로 출력의 각 비트당 파워는 줄어들게 되겠죠. 이를 정확히 반영하셔야 정확한 시뮬레이션이 될 겁니다. 이를 제대로 반영안했다면 이론 값보다 당연히 좋게 나오겠죠.

  11. 멍청이 2015.07.22 20:32

    질문이 있는데 qpsk 에서 위에 코드에서 Symbols=(bits_v(:,1)*(-2)+1 + i *(bits_v(:,2)*(-2)+1))./sqrt(2); 루트 2 이걸로 나눠주는데 이건 왜 나눠 주는거죠???

    • 남성 2015.07.23 04:18 신고

      위의 포스팅에서 적어 놨지만 신호의 파워, 여기서는 심볼의 파워를 1로 맞춰주는겁니다.
      Symbols = ±1/sqrt(2) ± j · 1/sqrt(2) 이므로 절대값이 1 이라는것을 확인 할 수 있습니다. 심볼 파워를 1로 정규화한 상태에서 잡음의 파워 즉 분산을 계산해서 넣어주는 겁니다.

오늘은 AWGN (Additive white gaussian noise) 채널 환경에서의 BPSK Bits Error Rate(BER)에 대해 설명한다. 일단 BPSK 가 뭔지 잘 모르는 사람은 다음 페이지를 참조 하기 바란다.


2011/03/08 - [통신] - [디지털 통신] BPSK, QPSK

 

위 페이지를 참조하면 bit 를 BPSK 심볼로 mapping 하는 방식에 대해 알았을 테고, 이제 잡음을 생성하는 방법을 알아 보자. 잡음을 생성하기 위해서는 몇 지 통신 관련 수식에 대해 알아야 한다. 보통 BER 그래프를 그릴 때는 SNR 축에 따른 그래프가 아닌 Eb/No 에 따른 그래프를 그리곤 한다. 따라서 BER 시뮬레이션을 정확히 하려면 Eb/No 에 따른 잡음을 생성할 수 있어야 하겠다.

 

잡음의 분산(σ2) 과 No(noise power spectral density) 는 다음과 같은 관계가 있다.

 

 

주의 : 본 식에서는 oversample 이나 channel coding 에 대해 고려 하지 않은 것이다.

 

 

 

위 식에서 각 기호의 의미는 아래와 같다.

 

  • Es : Signal energy

 

  • Eb : Bit energy

 

  • N0 : Noise power spectral density

 

  • M : PSK modulation order

 

 

시뮬레이션의 Es 는 1 로 정규화 한다. 그럼 잡음의 분산 값은 다음과 같은 식이 된다.

 

 

BPSK 의 경우 성좌도의 경우의 수가 2 이므로 M = 2 이다. 그럼 이제 최종적으로 BPSK 시뮬레이션시의 잡음 분산은 다음과 같이 된다.

 

 

이제 BPSK 의 BER 시뮬레이션을 해 보자.

 

MATLAB 코드 내용은 다음과 같다.


2014-07-03 추가 내용 

캡쳐 화면을 올렸더니 오타로 고생하시는 분들이 너무 많더군요. Github 에 파일 올립니다. 

파일 다운로드 링크 입니다. 

https://github.com/ssgkd/MATLAB-BPSK-AWGNSimulation/archive/master.zip

 

위 코드를 실행 했을때의 command 창 결과는

 

Eb/No= 0 BER: 0.0784773

Eb/No= 2 BER: 0.0374897

Eb/No= 4 BER: 0.0125513

Eb/No= 6 BER: 0.00238729

Eb/No= 8 BER: 0.000199714

Eb/No= 10 BER: 4.85714e-006

 

이며 그래프는 다음과 같이 나온다.

 

 


 

이론값과 시뮬레이션 값이 거의 일치하는 것을 확인 할 수 있고, 10 dB 일때 약간 값이 다른 것은 시뮬레이션 bits 수를 늘리면 해결 된다. 그리고 위 코드에서 이론적인 BER 값을 계산할 때는 communication toolbox berawgn() 함수를 이용하였다.



  1. 이전 댓글 더보기
  2. 흐미 2014.04.13 14:56

    좋은자료 감사합니다. 그런데 왜 저는 simulated된 BER값이 일정하게 나올까요? 어떤게 잘못됬을까요?

    • 남성 2014.04.13 15:11 신고

      방문해 주셔서 감사합니다. ^^
      일정하게 나오다는건 0.5 정도로 나온다는 거죠? 0.5 로 나온다는건 시뮬레이션이 잘못 됐다는겁니다. 오타일 확률이 클것 같네요.

    • 흐미 2014.04.13 20:57

      계속 찾아보려고 해도 오타가 보이지 않네요ㅜ
      한번 봐주실 수 있으세요?
      clc
      clear
      close all

      N_bits=7000000;
      Eb_No_dB=0:1:10;
      BER_buffer=zeros(size(Eb_No_dB));


      for n=1:length(Eb_No_dB)
      Eb_No_real_scale =10^(Eb_No_dB(n)/10);
      sigma_v=sqrt(1./(2*Eb_No_real_scale));

      bits_v=randi([0,1],N_bits,1);
      Symbols=bits_v*(-2)*1;

      noise_v=randn(N_bits,1)*sigma_v;

      tx_signal=Symbols + noise_v;

      demapped_bits= tx_signal < 0;

      BER_buffer(n)=sum(bits_v~= demapped_bits)./N_bits;

      fprintf('Eb/No = %g BER: %g\n', Eb_No_dB(n), BER_buffer(n))
      end

      %% Theory BPSK with BER
      Eb_NO_theory = 0:0.05:10;
      Ber_value_theory=berawgn(Eb_NO_theory,'psk',2,'nondiff');

      %% 그래프 그리그

      figure(1),semilogy(Eb_No_dB,BER_buffer,'b+',Eb_NO_theory,Ber_value_theory, 'r--');xlabel('Eb/No[dB]'),ylabel('BER')
      grid on, axis([0 10 0.000001 1]), legend('Simulated result','Theoretical BPSK result')

      덕분에 좋은 자료로 디지털 통신 이해가 되고 있습니다 감사합니다

    • 남성 2014.04.13 22:32 신고

      다 하셨는데 딱 한군데 틀리셨네요.
      Symbols=bits_v*(-2)*1; 를 Symbols=bits_v*(-2)+1;
      로 바꾸세요. 곱하기 1을 더하기 1로 바꾸시면 됩니다.
      바꾸고 돌려보니 잘 되네요.

    • 흐미 2014.04.13 22:34

      정말 감사합니다 ^^ 진짜 많이 배우고 갑니다.
      좋은 자료 정말 감사합니다!

    • 남성 2014.04.13 22:37 신고

      학교 숙젠가봐여? 도움 되셨다니 다행이네요.

    • 2014.04.13 23:07

      비밀댓글입니다

  3. 급해요 2014.04.15 13:29

    close all

    N_bits=200; %비트생성
    Eb_No_dB=0:1:10; %0부터 10까지 1씩지점을 체크하겠다.dB스케일로
    BER_buffer=zeros(size(Eb_No_dB)); %저장할 버퍼

    for n=1:length(Eb_No_dB)
    Eb_No_ral_scale=10^(Eb_No_dB(n)/10); %Eb/No값스케일변환
    sigma_v=sqrt(1./(2+Eb_No_ral_scale)); %잡음분산하려고 표준편차값 계산
    bits_v=randi([0 1],N_bits,1); %비트생성
    Symbols=bits_v*(-2)+1; %심볼 매핑
    noise_v=randn(N_bits,1)*sigma_v; %AWGN생성

    tx_signal=Symbols + noise_v; %원신호에 잡음더해주기
    res_bits=tx_signal < 0 %수신된것
    BER_buffer(n)=sum(bits_v ~= res_bits) ./ N_bits; %버퍼값 계산
    fprintf('Eb/No= %g BER: %g\n',Eb_no_dB(n), BER_buffer(n))
    end

    Eb_No_theory=0:0.05:10;
    Ber_value_theory=berawgn(Eb_No_theory,'psk',2,'mpmdoff');

    figure(1),semilogy(Eb_No_dB,BER_buffer,'b+',Eb_No_theory,Ber_value_theory,'r--'),xlabel('Eb/No [dB]'),ylabel('BER')
    grid,axis([0 10 0.000001 1]), legend('계산결과','이론적인 BPSK결과') %이론적값은 berawgn()함수사용

  4. 급해요 2014.04.15 13:30

    위에처럼 코딩했는데요
    Run하면 커맨드창에

    Undefined function 'Eb_no_dB' for input arguments of type 'double'.

    Error in Untitled2 (line 17)
    fprintf('Eb/No= %g BER: %g\n',Eb_no_dB(n), BER_buffer(n))

    나옵니다

  5. 급해요 2014.04.15 13:30

    그리고 그래프나오게 어떻게하죠??

  6. 부탁드려요 2014.05.08 11:09

    암만봐도어느부분이틀린건지알수가없네요 ㅠㅠ 한번만 확인해주시면 안될까요
    clc
    clear
    close all


    N_bits=7000000; % 생성할 비트 수
    Eb_No_dB=0:2:10; % Eb/N0 dB scale
    BER_buffer=zeros(size(Eb_No_dB)); % BER 저장할 버퍼


    for n=1:length(Eb_No_dB)
    Eb_No_ral_scale=10^(Eb_No_dB(n)/10); %Eb/N0 real scale 값 변환
    sigma_v=sqrt(1./(2+Eb_No_ral_scale)); %잡음 표준편차 값 계산

    bits_v=randi([0 1],N_bits,1); % 비트 생성
    Symbols = bits_v*(-2)+1; % symbol mapping

    noise_v=randn(N_bits,1)*sigma_v; %AWGN 생성

    tx_signal=Symbols + noise_v; % 신호 + 잡음

    demapped_bits=tx_signal < 0; % 수신단 demapping

    BER_buffer(n)=sum(bits_v ~= demapped_bits) ./ N_bits; %BER calculation

    fprintf('Eb/No= %g BER: %g\n', Eb_No_dB(n), BER_buffer(n))
    end

    %% 이론적인 BPSK BER 계산 Using Communication Toolbox

    Eb_No_theory=0:0.05:10;
    Ber_value_theory=berawgn(Eb_No_theory,'psk',2,'nondiff');

    %% 그래프 그리기

    figure(1), semilogy(Eb_No_dB,BER_buffer,'b+',Eb_N0_theory,Ber_value_theory,'r--'),xlabel('Eb/No [dB]'),ylabel('BER')
    grid on,axis([0 10 0.000001 1]), legend('Simulated result','Theoretical BPSK result')

    에러가 이렇게 뜨는데 봐도봐도 안보이네요 ㅠㅠ
    ??? Undefined function or method 'berawgn' for input
    arguments of type 'double'.

    Error in ==> BPSK at 32
    Ber_value_theory=berawgn(Eb_No_theory,'psk',2,'nondiff');

    공부하는데 정말 큰 도움이 되고 있습니다 고맙습니다~

    • 남성 2014.05.08 20:01 신고

      두군데 고치시면 되겠네요. 두군데 다 오타네요.


      아래 부분에서는 더하기를 곱하기로 고침

      sigma_v=sqrt(1./(2 * Eb_No_ral_scale)); %잡음 표준편차 값 계산



      아래 부분에서는 Eb_N0_theory 를 Eb_No_theory 로 고침, No의 O 를 소문자로

      figure(1), semilogy(Eb_No_dB,BER_buffer,'b+',Eb_No_theory,Ber_value_theory,'r--'),xlabel('Eb/No [dB]'),ylabel('BER')



      위 코드들은 수정한 코드이니 그대로 사용하시면 됩니다.

  7. ㅠㅠ 주무시겟죠 2014.05.12 01:16

    demapped_bits=tx_signal<0; % 수신된 비트신호
    이렇게 생각했는데 이게 왜 수신 신호인지 모르겟어요 ㅠㅠ

    • 남성 2014.05.12 08:01 신고

      mapping 을 할때 비트 0 을 1로 비트 1을 -1 로 mapping 을 하기 때문에 수신 신호가 0 보다 크면 비트 0으로 0보다 작으면 비트 1로 mapping 하는 겁니다.

  8. 2014.11.11 01:46

    비밀댓글입니다

    • 남성 2014.11.11 13:11 신고

      방문해 주셔서 감사합니다. ^^ 수정했습니다.

  9. 디시설 2014.11.16 23:32

    방금 위에서 물어 봣는데 그럼 혹시 저코드를 사진을 이용해 영상처리를 해서 잡음이 많은 거랑 적은걸 표현을 어떻게 하나요 ㅜㅜ

    • 남성 2014.11.17 01:25 신고

      사용하는 영상이 몇비트로 구성되어 있는지는 모르겟지만 영상의 값들을 비트로 표현한 후에 잡음을 섞어서 imshow() 와 같은 함수로 이미지를 보여 주면 잡음이 섞인 정도에 따라 비교가 될것 같은데요.

  10. 급한질문 2014.11.29 12:16

    혹시 이 코드에서 엔코더와 디코더를 넣을려면 어떻게해야하는지 알수잇을까요?

    • 남성 2014.11.29 18:22 신고

      어떤 인코더와 디코더를 말씀 하시는지... 위 코드의 시작과 끝이 비트 니까 각각의 외부에 넣으면 되겠네요.

  11. 서영화 2014.12.03 18:42

    송신에 Raised cosine 수신에 matched filter(square root RC로)추가하려면 어디에 너야할까요??

    • 남성 2014.12.03 21:29 신고

      송신 필터는 송신단 끝단에 넣으면 되고 수신 필터는 채널 통과후 수신 첫단에 넣으면 됩니다.

  12. 2014.12.10 18:17

    비밀댓글입니다

    • 남성 2014.12.11 18:12 신고

      아래 주소 참조하시면 쉽게 할 수 있을것 같네요.

      http://www.mathworks.com/matlabcentral/fileexchange/37455-bpsk-modulatian-demodulation-by-using-simulink

      Communication system toolbox 가 있어야 할 겁니다.

    • 2014.12.12 00:58

      비밀댓글입니다

    • 남성 2014.12.12 08:49 신고

      위에 그래프 그리는 코드도 들어가 있는데요.

    • 2014.12.13 01:40

      비밀댓글입니다

    • 남성 2014.12.13 02:30 신고

      아래 글에 작성한 바와 같이 simulink 에는 to workspace 라는 블록이 있습니다.

      http://iamaman.tistory.com/485

      simulink 에서 BER 계산 결과와 Eb/N0 값을 to workspace 블록을 이용하여 MATLAB workspace 에 저장한 후에 그래프를 그리면 될것 같네요.

    • 도와주세요 ㅠㅠ 2014.12.16 12:57

      감사합니다. 별거 아닌데 너무 어렵게 생각했네요. 많이배워갑니다

    • 남성 2014.12.16 13:21 신고

      도움되었다니 다행이네요. 방문해 주셔서 감사합니다. ^^

  13. 디시설 2014.12.20 15:37

    clc
    clear
    close all


    N_bits=7000000; % 생성할 비트 수
    Eb_No_dB=0:2:10; % Eb/N0 dB scale
    BER_buffer=zeros(size(Eb_No_dB)); % BER 저장할 버퍼
    double berawgn;


    for n=1:length(Eb_No_dB)
    Eb_No_ral_scale= 10 ^(Eb_No_dB(n)/10); % Eb/N0 real scale 값 변환
    sigma_v=sqrt(1./(2*Eb_No_ral_scale)); % 잡음의 표준 편차 값 계산

    bits_v=randi([0 1],N_bits,1); % 비트 생성
    Symbols=bits_v*(-2)+1; % symbol mapping

    noise_v=randn(N_bits,1)*sigma_v; % AWGN 생성

    tx_signal=Symbols + noise_v ; % 신호 + 잡음

    demapped_bits= tx_signal < 0 ; % 수신단 demapping

    BER_buffer(n)=sum(bits_v ~= demapped_bits) ./ N_bits; % BER calculation

    fprintf('Eb/No= %g BER: %g\n',Eb_No_dB(n) ,BER_buffer(n))
    end

    %% 이론적인 BPSK BER 계산 Using communication toolbox

    Eb_N0_theory=0:0.05:10 ;
    Ber_value_theory=berawgn(Eb_N0_theory,'psk',2,'nondiff');
    %% 그래프 그리기

    figure(1),semilogy(Eb_No_dB,BER_buffer,'b*',Eb_N0_theory,Ber_value_theory,'r--'), xlabel('Eb/No [dB]'), ylabel('BER')
    grid,axis([0 10 0.000001 1]), legend('Simulated result','Theoretical BPSK result')

    'berawgn'은(는) 유형 'double'의 입력 인수에 대해 정의되지 않은 함수입니다.

    오류 발생: Untitled (line 33)
    Ber_value_theory=berawgn(Eb_N0_theory,'psk',2,'nondiff');

    이렇게오류가 나는데 무슨 문제인가요 ㅜㅜ

    • 남성 2014.11.16 20:02 신고

      Communications System Toolbox 가 있는건지 확인해 보세요.

  14. dduss 2015.01.14 13:42

    디지털공부하는데 정말 많은 도움이 되고 있습니다! 이글 보고 BER 코드를 짜고있는데 송수신단에 각각 다른 필터를 사용해 성능비교를 하려고 하는데 raised cosine filter를 추가할때 oversample을 해서인지 demapped bit 와 bits_v가 사이즈가 맞지않아 연산이 수행되지 않는것 같습니다.. 어떻게 해결해야할까요? ㅠㅠ

    • 남성 2015.01.14 23:03 신고

      송신단에서 oversample 을 했으면 수신단 필터링 후에 downsample 을 해야 겠죠. 그리고 필터 delay 때문에 송신단에서는 0 을 delay 만큼 추가를 해 줘야 송신한 모든 데이터를 수신 할 수 있을 겁니다.

    • dduss 2015.01.15 13:29

      그럼 oversampling만큼 upsampling 하고 필터링 후에 downsampling하면 되는 건가요?? 필터를 추가하면 noise를 더할 때도 비트수가 안맞는 것 같아요.. 공부를 더 많이 해야겠네요 ㅠㅠ

    • 남성 2015.01.15 20:58 신고

      upsampling 하고 필터링 하는게 oversampling 이라고 하고요. oversampling 한 다음에 잡음을 섞어야 하고 잡음을 섞을때는 sample 에 섞어야 하니까 sample 에 따른 잡음 분산을 다시 계산해야겠죠. 그리고 수신단에서 필터링 한다음에 downsampling 을 하는 겁니다.

    • dduss 2015.01.15 22:18

      개념을 제가 잘못이해하고있었네요! 정말 감사합니다. 본문도 댓글도 정말많은 도움이 되었어요~

    • 남성 2015.01.16 08:19 신고

      도움되었다니 다행이네요. ^^

  15. 온지 2016.04.19 20:33

    통신이론에 대해 공부하는 학생인데요 BER그래프를 그릴때 X축을 SNR이 아닌 Eb/No으로 해야되는건지 정확한 이유좀 알수있을까요?ㅠㅠㅠ

    • 남성 2016.04.19 21:30 신고

      SNR 축으로 그릴 수는 있습니다. 어떤걸 확인하는게 목적이냐에 따라 다르겠죠. 보통 BER 은 비트의 에너지에 따른 BER 성능을 구하려고 하는거죠. BER 을 구한다 함은 변조에 따른 심볼당 비트수 차에 의한 신호의 파워 차이를 없애고 순수 비트 에너지에 따른 BER 성능을 구하는게 타당 할 겁니다. BPSK 와 QPSK 를 Eb/No 를 x 축으로 해서 그리면 성능이 같지만 SNR 을 축으로 해서 BER 을 구한다면 성능이 다르겠죠? 이때 Eb/No 를 기준으로 한 BER 을 보고 BPSK 와 QPSK 는 BER 성능이 같다고 해석하는게 맞을까요? 아니면 SNR 축에 따라 그려진 BER 을 보고 BPSK 가 QPSK 보다 BER 성능이 좋다라고 해석하는게 맞을까요?

  16. 임명호 2016.11.28 20:22

    안녕하세요. 통신을 전혀 모르는 직장인인데, 아는 지인 분께 통신 관련 논문으로 BPSK랑 QPSK BER을 시뮬레이션 해달라고 하셔서 이리저리 헤메다가 질문 올려봅니다. (통신에 대한 지식이 전무합니다 ㅠㅠ...엉뚱한 소리라도 이해좀 부탁드립니다)

    1. Eb/N0 가 왜 Matlab 소스에서 0에서 10까지로 2간격으로 되어있는지요? 단순히 축을 설정하기 위한 것인가요?
    단순히 축 설정을 위한 것이라면 X축이 같으므로 Y 값이 다 똑같은 결과가 나오는 것이 아닙니까?
    제가 참조하는 논문에는 Eb/N0 = (S/Rb) / (N/W) 라고 되어있습니다. (S: 신호전력 Rb:전송속도 N:잡음전력 W:대역폭)

    밑에는 제가 올린 지식인 질문입니다. 해당 내용에 논문에 대한 참조 데이터가 있습니다. (원하신다면 논문을 보내드리겠습니다.)
    http://kin.naver.com/qna/detail.nhn?d1id=1&dirId=108&docId=265206138

    • 남성 2016.11.28 21:21 신고

      ber 시뮬레이션은 eb/n0 에 따른 비트 에러율을 구하는 겁니다 따라서 eb/n0 값을 시뮬레이션 때마다 조절해줘야 합니다. 즉 S값을 달리하면서 시물레이션 하는 거죠

    • 2016.11.28 23:45

      비밀댓글입니다

    • 남성 2016.11.29 00:13 신고

      제가 질문자님의 의도를 제대로 이해한건지 모르겠는데... 비콘과 발리스의 비트레이트가 다르다고해서 비콘은 bpsk로 변조하고 발리스는 qpsk로 변조해서 에러 레이트 비교한 다음에 둘중에 어느하나가 더 좋다라고 결론내려는건가요?

  17. 임명호 2016.11.30 16:02

    안녕하세요 우선 답변주셔서 감사드립니다. 제가 글을 명확히 안 써서 혼돈이 온 점 죄송스럽게 생각합니다.

    비콘
    발리스에
    대해 각각
    Bpsk ber 시뮬레이션
    Qpsk ber 시뮬레이션
    을 하여 성능을 비교해보려합니다.

    1. Bpsk ber
    비콘 vs 발리스

    2. Qpsk ber
    비콘 vs 발리스

    이런 시뮬레이션을 하려합니다.

    • 남성 2016.11.30 16:13 신고

      해당 내용에 대해 위 포스팅 대로 시뮬레이션을 하면 성능차이가 없을거 같네요

  18. 임명호 2016.12.28 10:19

    남성님 안녕하세요. 1달만에 다시 질문 드립니다.
    현재 BPSK와 QPSK 시뮬레이션을 하려하는데
    현재 SNR 값을 구하였습니다. 제가 조사한 바에 의하면 SNR은 아날로그 통신에서 주로 사용하고
    디지털통신에서는 SNR[DB]을 정규화한 Eb/N0[DB] 로 변환하는걸로 알고있습니다.
    최종적으로 정리한 공식을 찾아보니 아래와 같습니다.

    SNR=10log(S/N) [DB]
    Eb/N0 = (S/N) * (W/R)

    여기서,
    S: 신호전력 [V/m]
    N: 노이즈전력[V/m]
    W : 대역폭 [HZ]
    R : 비트전송속도 [bits/s]

    현재 저는 S,N,R 값을 알고있는데 W 값을 모르겠습니다.
    제가 가지고 있는 정보는

    샘플링 시간 : 4.4325 * 10^(-8)
    심볼 시간 :1.773 * 10^(-6)
    전송 주파수 4.23 MHZ 입니다.

    그리고 BPSK와 QPSK로 시뮬레이션을 하려하는데 제가 W 값을 위 정보를 바탕으로 구할 수 있는지요?

    • 남성 2016.12.28 10:32 신고

      아래 주소 보시면 심볼레이트와 대역폭 사이의 관계를 알 수 있습니다.

      http://www.ktword.co.kr/abbr_view.php?m_temp1=4185

  19. zzzxxxccc 2017.11.25 14:21

    BPSK에서 unipolar signal과 bipolar signal 을 동시에 그래프상에서 보여주고 싶은데 어떻게하나요??

    • 남성 2017.11.25 14:31 신고

      두 시그널에 따른 ber 성능을 비교하고 싶다는 건가요? 아니면 그래프 두개를 동시에 그리는 방법을 물어보시는 건가요?

  20. here 2018.11.25 12:46

    안녕하세요! AWGN에서 BER값을 내기 때문에 BER값을 구할때 계속 다른값이 나오잖아요!
    이걸 백번돌려서 나온 백개들의 평균을 구해서 그래프를 그릴려고 합니다! 이 때 이걸 백번돌리는 함수는 어떻게 만들면 될까요?

    • 남성 2018.11.25 15:27 신고

      백번 돌려서 평균을 구하려면 당연히 버퍼를 만들어서 저장을 하는 형태로 만들어야 되겠죠.

  21. 안녕하세요 2018.12.11 16:35

    제가 QPSK BER curve를 plot하려고 하는데요
    QPSK지만 1bit에 대해서만 ber curve를 그리려고 합니다. 그래서 matlab코드를 이렇게 짜봣는데요
    clc;
    clear all;
    %QPSK BER
    %1->√(2Eb/Tb)*cos(2πfct)이고 0->-√(2Eb/Tb)*cos(2πfct)로 mapping 된다.
    Eb=1;
    sample_num=100000;%각 Eb/No당 샘플갯수
    Eb_over_No_dB=-2.5:0.1:7.5;
    Eb_over_No=10.^(Eb_over_No_dB./10);
    No=1./Eb_over_No.*Eb;
    Bit_error=zeros(1,length(Eb_over_No_dB));%각 Eb/No값마다 bit error값을 저장할 행렬
    for k=1:length(Eb_over_No_dB)% Noise variance를 바꿈.즉, Eb/No값을 -2.5dB에서 7.5dB까지 바꿈. 각각의 Noise variance에 대하여 BER을 찾는다.
    A=randi([0,1],1,sample_num);%랜덤한 0,1 sequence를 sample number만큼 생성.즉, binary데이터 열을 만든것.
    No_over_2=No(1,k)./2;%Noise Variance No/2
    %QPSK oddbit에 대한 receiver 출력구하기.평균이 ±√Eb variance가 No/2인
    %Random Varable(=Output)
    QPSK_odd=(A*2-1)*sqrt(Eb);
    Noise=sqrt(No_over_2)*(randn(1,sample_num));
    Output=QPSK_odd+Noise;


    %BER계산
    for i=1:sample_num
    if(A(1,i)==1) %송신된 Binary data가 1일때
    if( Output(1,i)<0)%receiver 출력이 0보다 작으면 -----
    Bit_error(1,k)=Bit_error(1,k)+1;%-->error
    end
    else if(A(1,i)==0) %송신된 Binary data가 1일때
    if( Output(1,i)>0)%receiver 출력이 0보다 크면 -----
    Bit_error(1,k)=Bit_error(1,k)+1;%error
    end
    end
    end
    end
    end

    BER=Bit_error./sample_num.*100;
    BER_Theory=1/2*erfc(sqrt(Eb_over_No));
    semilogy(Eb_over_No_dB,BER);

    방법에 대해 설명해보자면 binary data와 receive된 data를 비교해서 1을 보냈을때 receive signal이 0보다 작으면 error이고
    0을 보냈을때 receive signal이 0보다 크면 error가 나오는 이런 식으로 짰는데 BER이 너무 높게 나옵니다.
    이론값과 100배이상 나는데 어디가 잘못된 걸까요?

    • 남성 2018.12.11 17:42 신고

      위 코드 링크 넣어놨는데 다운로드 받아서 돌려 보시면서 비교해 보는게 좋을것 같네요

+ Recent posts