QAM modulation, pulse shaping and raised cosine filter BER

Hi! I am trying to simulate a communication channel at 100Mb/s data rate, using ideal rectangular pulses for sampling, with 16 QAM modulation and raised cosine filters for pulse shaping, adding some white gaussian noise and then receive it in a filter, demodulate and find out the bit error ratio between the original signal and the one at the end, to then plot it. I have managed to do this process without the cosine filters for pulse shaping and the results turned out well, but when I add the filters there is definitely something wrong in the resulting plot and I can't tell what. I also have doubts regarding other factors, like

1. if the data rate is correct and how exactly do I define the data rate with these parameters?
2. is it correct to use the rectpulse AND raised cosine filters at the same time?
3. how do I account for the delay of the raised cosine filters? That might be why the result is not as the theory? Are the filters matched?

M = 16;  % Alphabet size, the number of different symbols available
k = log2(M) % Number of bits per symbol
NSymbols = 5000;  % Number of Symbols in the test message
Ts = 1e-5;  % Symbol duration (s)
Nsamp = 10;  %  Oversampling rate - the number of time samples in each symbol

x = randi([0 M-1], NSymbols, 1)  % Test message as a set of random integer symbols from 0 to M-1


y = qammod(x, M)  % Use built-in QAM modulation function

w = rectpulse(y, Nsamp)
calculated_ber = zeros(1, 17) % array to store the bit error ratio results
z = zeros(1,17) % number of bit errors. not using this
i = 1
for EbN0 = 4:20
txfilter = comm.RaisedCosineTransmitFilter('OutputSamplesPerSymbol',10,'RolloffFactor',0.3, 'FilterSpanInSymbols', 10);
rxfilter = comm.RaisedCosineReceiveFilter('InputSamplesPerSymbol',10, ...
'DecimationFactor',10,'RolloffFactor',0.3, 'FilterSpanInSymbols', 10);
snr = EbN0+10*log10(k)-10*log10(Nsamp); % calculating snr from EbN0
Nw = txfilter(w); % applying raised cosine transmit filter
Nw = awgn(Nw, snr,'measured'); % adding white gaussian noise
Nw = rxfilter(Nw);    % applying receive filter                           
received_w = qamdemod (Nw,M); % signal w after Demodulating
received_w = downsample(received_w,10); % downsampling to reach size of original signal and be able to compare
[z, calculated_ber(i)] = biterr(x,received_w); % getting the bit error rate
i = i + 1 % incrementing index for the calculated_ber array
end

%plotting 
EbN0 = 4:20
berTheory = berawgn(EbN0,'qam',16)

figure;
semilogy(EbN0, calculated_ber,'o');  
hold on;
semilogy(EbN0, berTheory);
title('Bit error ratio')
xlabel('E_b/N_0(dB)');
ylabel('BER')
legend('Simulated BER','Theory BER')
axis([4 22 10e-50 1]);

[Rena Berman]

(Answers Dev) Restored edit

Answer 1

[Santosh Fatale]

Hi Dragos-Valentine,

I investigated the code shared by you and made the following observations and incorporated answers to your queries into the same observations.

1. if the data rate is correct and how exactly do I define the data rate with these parameters?

I did not get what you meant by defining the data rate. I understand that you are comparing the theoretical BER with the simulation results. To calculate BER, why do you need to define the data rate? Note that the BER is the ratio of the number of bits received erroneously to the number of bits transmitted.

2. is it correct to use the rectpulse AND raised cosine filters at the same time?

It is worth noting that both raised cosine filters and rectangular filters are pulse shaping filters. You should use either of these at the transmitter and receiver to reverse the operation at the transmitter. You do not need to use both filters at the same time.

https://in.mathworks.com/help/comm/ug/modulation-with-pulse-shaping-and-filtering-examples.html

3. how do I account for the delay of the raised cosine filters? That might be why the result is not as the theory? Are the filters matched?

The raised cosine filter creates a total symbol delay of txfilter.FilterSpanInSymbols/2 + rxfilter.FilterDelayInSymbols/2. To convert this symbol delay into a bit delay, multiply the symbol delay by the number of bits per symbol. You need to consider this delay while calculating BER.

https://in.mathworks.com/help/comm/ug/raised-cosine-filtering.html

The modified code is as follows:

clearvars;
M = 16;  % Alphabet size, the number of different symbols available
bitPerSymbol = log2(M); % Number of bits per symbol
NSymbols = 1e6;  % Number of Symbols in the test message
Ts = 1e-5;  % Symbol duration (s)
% Nsamp = 1;  %  Oversampling rate - the number of time samples in each symbol

x = randi([0 1], NSymbols*bitPerSymbol, 1);  % Test message as a set of random integer symbols from 0 to M-1

y = qammod(x, M, "InputType","bit","UnitAveragePower",true);  % Use built-in QAM modulation function

txfilter = comm.RaisedCosineTransmitFilter;
rxfilter = comm.RaisedCosineReceiveFilter;
% w = rectpulse(y, Nsamp);
EbN0 = 4:20;
calculated_ber = zeros(1, length(EbN0)); % array to store the bit error ratio results
z = zeros(1,17); % number of bit errors. not using this

for idx = 1 : length(EbN0)
    snr = EbN0(idx) + 10*log10(bitPerSymbol) - 10*log10(txfilter.OutputSamplesPerSymbol); % calculating snr from EbN0
    txOut = txfilter(y); % applying raised cosine transmit filter
    chanOut = awgn(txOut, snr,"measured"); % adding white gaussian noise
    rxOut = rxfilter(chanOut);    % applying receive filter                           
    received_w = qamdemod (rxOut,M, "OutputType","bit","UnitAveragePower", true); % signal w after Demodulating
    
    delayInSymbols = txfilter.FilterSpanInSymbols/2 + rxfilter.FilterSpanInSymbols/2;
    delayInBits = delayInSymbols * bitPerSymbol;
    
    received_wAlligned = received_w(1+delayInBits:end);
    xAlligned = x(1:end-delayInBits);
    %     received_w = downsample(received_w,Nsamp); % downsampling to reach size of original signal and be able to compare
    calculated_ber(idx) = nnz(xAlligned ~= received_wAlligned)/ length(xAlligned); % getting the bit error rate
    %     i = i + 1; % incrementing index for the calculated_ber array
end

%% plotting 
% EbN0 = 4:20;
berTheory = berawgn(EbN0,'qam',16);

figure;
semilogy(EbN0, calculated_ber,'o', "MarkerSize",10);  
hold on;
semilogy(EbN0, berTheory, "LineWidth",2);
grid on;
title('Bit error ratio')
xlabel('E_b/N_0(dB)');
ylabel('BER')
legend('Simulated BER','Theory BER')
axis([4 22 10e-20 1]);

Note that for high values of SNR, it is not possible to calculate the BER through simulation with current code due to constraints on the system on which you are running the simulation. For a high SNR, you need to consider a high value for the number of symbols transmitted to capture the erroneous bits in transmission. In the following figure, you can see the simulated BER until "EbN0 = 14" matching the theoretical BER. For higher values, simulated BER is zero and is not plotted on the logarithmic Y-axis.