0
$\begingroup$

I'm trying to implement the Rayleigh flat fading channel in MATLAB using SRRC (Square-root raised cosine) pulse shape. The modulation scheme is BPSK. Here is my code:

EbNo_vec = -4:2:15; % Eb/N0 vector
number_of_packets = 3e3; % Number of packets for each Eb/N0
packet_len = 501; % Length of each packet
span = 30; % Span of SRRC pulse
beta = 1; % Roll-off factor of SRRC pulse
resample_rate = 10; % Number of samples per symbols
g_square_root_rc = rcosdesign(beta , span , resample_rate); % SRRC pulse
ber_vec_no_eq = ones(length(EbNo_vec) , number_of_packets); % Bit error rate matrix
for i=1:length(EbNo_vec)
    for j=1:number_of_packets
    d_k = randi(2 , [1 , packet_len]); % index of data symbols
    u_k = d_k;
    u_k(d_k == 2) = -1; % BPSK symbols
    v_n = upsample(u_k , resample_rate); % Upsample BPSK symbols
    s_n = conv(v_n,g_square_root_rc , "same"); % Transmitted signal
    signal_power = sum(abs(s_n).^2)/length(s_n); % Calculate power of signal
    s_n = s_n/sqrt(signal_power); % Normalize the power of signal
    rayleigh_samples = (1/sqrt(2))*(randn(1,packet_len)+1j*randn(1,packet_len));
    rayleigh_channel = repelem(rayleigh_samples,resample_rate); % Constant channel for each symbol
    faded_signal = s_n.*rayleigh_channel;
    current_EbNo = EbNo_vec(i); % Select Eb/N0
    current_SNR = current_EbNo - 10*log10(resample_rate); % Calculate SNR
    noise_power_dB = -current_SNR; % Calculate power of noise
    noise_power_linear = 10^(noise_power_dB/10); % Convert to linear scale
    w_n = sqrt(noise_power_linear/2)*randn(1,length(faded_signal)); % Create noise vector
    %w_n = 0;
    r_n = w_n + faded_signal; % Add noise
    v_hat_n = conv(r_n , g_square_root_rc , "same"); % Matched filtering
    uu_hat_k = v_hat_n./rayleigh_channel; % Equalization
    u_hat_k = downsample(uu_hat_k , resample_rate); % Downsample
    u_hat_k = sign(real(u_hat_k)); % Detection
    ber_vec_no_eq(i , j) = mean(u_k ~= u_hat_k); % Calculate BER
    end
    disp(i)
end
average_ber = mean(ber_vec_no_eq , 2);
EbNo_linear_vec = 10.^(EbNo_vec/10);
berTheory = 0.5.*(1-sqrt(EbNo_linear_vec./(EbNo_linear_vec+1)));
figure;
semilogy(EbNo_vec, average_ber,'o', "MarkerSize",5);  
hold on;
semilogy(EbNo_vec, berTheory, "LineWidth",2);
grid on;
title('Bit error ratio, packet length = ' + string(packet_len))
xlabel('E_b/N_0(dB)');
ylabel('BER')
legend('Simulated BER','Theory BER')

When I set the additive noise to zero and the number of samples per symbol to one, the resulting BER is zero but increasing the number of samples per symbol resulted in increase of BER. When I interchanged the matched filtering step and the equalization step, BER stayed zero even for greater than one samples per symbol. I don't know why this happens. Also when the samples per symbol was set to one, the resulting BER curve differed from the theoretical curve:enter image description here I think my simulation model for the Rayleigh flat fading channel is wrong but I don't know where is the problem. If we don't want to simulate the pulse shaping step, it's straightforward to write the program which generates a BER curve that agrees completely with the theoretical BER curve. For example, take a look at Krishna Pillar’s code from https://dsplog.com/2008/08/10/ber-bpsk-rayleigh-channel/ :

for ii = 1:length(Eb_N0_dB)
   n = 1/sqrt(2)*[randn(1,N) + j*randn(1,N)]; % white gaussian noise, 0dB variance 
   h = 1/sqrt(2)*[randn(1,N) + j*randn(1,N)]; % Rayleigh channel
   % Channel and noise addition
   y = h.*s + 10^(-Eb_N0_dB(ii)/20)*n; 
   % equalization
   yHat = y./h;
   % receiver - hard decision decoding
   ipHat = real(yHat)>0;
   % counting the errors
   nErr(ii) = size(find([ip- ipHat]),2);
end 

Edit: In order to ensure that the power of signal and the noise is selected correctly (in the absence of fading), I used the following code which is based on the @MBaz's clever idea:

EbNo_vec = -4:2:10; % Eb/N0 vector
number_of_packets = 1e4; % Number of packets for each Eb/N0
packet_len = 700; % Length of each packet
span = 10; % Span of SRRC pulse
beta = 1; % Roll-off factor of SRRC pulse
resample_rate = 10; % Number of samples per symbols
g_square_root_rc = rcosdesign(beta , span , resample_rate); % SRRC pulse
energy_signal = zeros(length(EbNo_vec) , number_of_packets);
energy_noise = zeros(length(EbNo_vec) , number_of_packets);
for i=1:length(EbNo_vec)
    for j=1:number_of_packets
    d_k = randi(2 , [1 , packet_len]); % index of data symbols
    u_k = d_k;
    u_k(d_k == 2) = -1; % BPSK symbols
    v_n = upsample(u_k , resample_rate); % Upsample BPSK symbols
    s_n = conv(v_n,g_square_root_rc , "same"); % Transmitted signal
    signal_power = sum(abs(s_n).^2)/length(s_n); % Calculate power of signal
    s_n = (s_n)/(sqrt(signal_power)); % Normalize the power of signal
    current_EbNo = EbNo_vec(i); % Select Eb/N0
    current_SNR = current_EbNo  - 10*log10(resample_rate); % Calculate SNR
    noise_power_dB = -current_SNR; % Calculate power of noise
    noise_power_linear = 10^(noise_power_dB/10); % Convert to linear scale
    w_n = sqrt(noise_power_linear/2)*randn(1,length(s_n)); % Create noise vector
    v_hat_n1 = conv(s_n , g_square_root_rc , "same"); % Matched filtering (signal)
    v_hat_n2 = conv(w_n , g_square_root_rc , "same"); % Matched filtering (noise)
    u_hat_k1 = downsample(v_hat_n1 , resample_rate); % Downsample (signal)
    u_hat_k2 = downsample(v_hat_n2 , resample_rate); % Downsample (noise)
    energy_signal(i , j) = mean(abs(u_hat_k1).^2); % Eb
    energy_noise(i,j) = mean(abs(u_hat_k2).^2); % N0/2
    end
    disp(i)
end
Eb_signal = mean(energy_signal , 2);
N0_half_noise = mean(energy_noise , 2);
EbNo_linear_vec = 10.^(EbNo_vec/10);
disp(Eb_signal./(2*N0_half_noise));
disp(EbNo_linear_vec);

The obtained result was close to the true Eb/N0:

Simulated Eb/N0:[0.3983,0.6308,1.0004,1.5854,2.5149,3.9829,6.3138,9.9988]

True Eb/N0: [0.3981,0.6310,1.0000,1.5849,2.5119,3.9811,6.3096,10.0000]

$\endgroup$
8
  • $\begingroup$ @DanBoschen Thanks for pointing out, I forgot to add the source. What do you mean by "not selecting a new gain from the Rayleigh distribution"? $\endgroup$
    – S.H.W
    Jul 6, 2023 at 19:30
  • $\begingroup$ @DanBoschen Thanks. I added the code which implements your suggestion. Sorry for possible mistakes. $\endgroup$
    – S.H.W
    Jul 6, 2023 at 23:47
  • $\begingroup$ If you remove both fading and noise from your simulation do you see zero BER? If you then add just noise, do you get the expected BER? $\endgroup$
    – MBaz
    Jul 7, 2023 at 14:06
  • 1
    $\begingroup$ @MBaz Yes, after removing fading and noise the BER is zero. Also with additive noise (and no fading), the BER curve matches the theoretical BER curve for the AWGN channel. I think the main issue is simulating the Rayleigh flat fading channel. $\endgroup$
    – S.H.W
    Jul 7, 2023 at 14:24
  • 2
    $\begingroup$ Yeah, there may be different issues. One thing you may try: remove all noise and find Eb (average energy of u_hat_k before quantization). Then, remove the signal and leave only noise, and repeat to find N0. Verify that the ratio is what you expect. $\endgroup$
    – MBaz
    Jul 7, 2023 at 19:52

1 Answer 1

1
$\begingroup$

The code the OP references from dsplog.com is really only useful for predicting BER but not channel simulation for receiver testing. This is because the simulation is changing the channel magnitude and phase for every symbol according to samples drawn from a complex Gaussian Process (which itself is a Rayleigh fading channel), but immediately equalizes each symbol to correct for that distortion, with what is equivalent to a single tap equalizer in that the gain and phase coefficient is inverted in the presence of noise. An actual receiver would take several symbols to determine what that perfect equalization needs to be, so that simulation is essentially showing what the average BER performance would be after successful equalization- and importantly over many experimental runs (or over a long time duration with changing fading conditions). How would any receiver be able to correct a phase change that is independent and uniformly distributed from sample to sample without any other information about that change?

This is not a flaw with the code provided by dsplog.com but a clarification on its purpose. The phase and amplitude changes from sample to sample, and is then corrected (together with the modifying the AWGN that is added) instantly. This allows us to simulate the BER temporarily over $N$ samples, just the same as if we ran $N$ receivers simultaneously to determine the average bit error rate across all of them on one single symbol (ergodic process, meaning we can simulate one changing experiment in time or simulate multiple experiments in parallel, it will be the same result).

Please refer to this post where I detail further simulation approaches for Rayleigh channels for cases of fast and slow fading, and frequency selective or flat fading. The possible situations are fast with either frequency selective or flat, or slow with either frequency selective or flat. Fading simulations for mobile applications will also include Doppler spread which is the frequency domain corollary to delay spread. Delay spread if wide enough results in frequency selective which is deep nulls in the frequency domain within the channel bandwidth. Similarly Doppler spread if wide enough in the frequency domain can cause nulls in the time domain within a symbol duration.

Simulating BER is straightforward using the approach given by dsplog.com, but simulating an actual receiver operation with a Rayleigh channel would be more challenging as we would require the channel to be time varying (if we simulate a single channel; in order to get the average result over the range of phase and amplitudes the statistics will supply), and we would need to model the complete receiver including timing recovery and channel equalization. If we are designing an actual receiver, then this could be a good test of the algorithms for that receiver implementation, and practical if done at the lower SNRs (such as 15 dB SNR and lower where close to 0.01 BER and higher is expected).

$\endgroup$
2
  • $\begingroup$ I removed that part from the question. Practical matters aside, is it possible to simulate the Rayleigh flat fading with pulse shaping such that the resulting BER curve matches the theoretical curve? Also, as I mentioned in the question, I couldn't resolve the issue with the sample per symbol greater than one. $\endgroup$
    – S.H.W
    Jul 8, 2023 at 1:01
  • 1
    $\begingroup$ Yes the same way as long as your simulation works in that you get the expected BER without the Rayleigh channel. For that ensure proper matched filter and symbol timing in the receiver. If you mean modeling a Rayleigh channel to test an actual receiver algorithm (vs just confirming what is given by the tail probabilities in the decision process), that is a lot more challenging as I elaborated in my answer. $\endgroup$ Jul 8, 2023 at 3:23

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.