This example shows how to measure the bit error rate (BER) of the end-to-end chain of the Consultative Committee for Space Data Systems (CCSDS) telemetry (TM) system. The simulation chain follows the coding and modulation schemes that are specified by these two standards:
Data from various instruments is generated onboard the satellite. This data collectively is called as TM data. CCSDS specifies the coding and modulation schemes for the transmission of TM data from the satellite to an Earth station. This example shows the end-to-end simulation of the satellite to the Earth station communication link. The example shows how to generate a complex baseband CCSDS TM waveform from the randomly generated transfer frames (TFs), introduce radio frequency (RF) impairments to the baseband signal, and add additive white gaussian noise (AWGN) to the impaired signal. Then, the example shows the synchronization, demodulation, and decoding of this impaired noisy signal to get the final bits in the form of TFs. The example also shows how to measure the BER with respect to the signal to noise ratio (SNR) for one configuration of the CCSDS TM signal. This figure shows the end-to-end simulation chain.
This example models these RF impairments:
Carrier frequency offset (CFO)
Carrier phase offset (CPO)
Symbol timing offset (STO)
This figure shows the receiver side operations.
The frame synchronization and channel decoding processing step performs these three tasks.
Perform phase ambiguity resolution
Correctly synchronize the frame to the starting of the attached sync marker (ASM)
Perform channel decoding of the synchronized frame to get the recovered TF
This figure shows these three tasks. To start, form a bank of ASM sequences. Each sequence corresponds to the original ASM value in which phase ambiguity is introduced. Correlate each of these sequences with the demodulated symbols. Choose the phase ambiguity value that has the highest correlation peak. Perform the frame synchronization with this correlation process. The process to perform correlation is illustrated in section 9.3.7 in . This example adopts the simplified Massey algorithm for frame synchronization. Resolve the phase ambiguity on the complete set of demodulated symbols after the frame synchronization process is complete. Finally, perform channel decoding on these symbols to obtain the recovered TFs.
This example uses a QPSK modulation scheme for signal generation and reception and rate 1/2 convolutional coding scheme for channel coding. The end-to-end chain this example shows can also be used for the channel coding schemes that are specified in : Reed-Solomon (RS) codes and concatenated codes. For convolutional and concatenated codes, this example supports rates 1/2 and 2/3 along with pulse code modulation (PCM)-format of non return to zero-line (NRZ-L). The supported modulation schemes in this example are BPSK and QPSK.
seeConstellation = true; % Flag to toggle the visualization of constellation channelCoding = "convolutional"; % Channel coding scheme transferFrameLength = 1115; % In bytes corresponding to 223*5 modScheme = "QPSK"; % Modulation scheme alpha = 0.35; % Root raised cosine filter roll-off factor sps = 8; % Samples per symbol
Set the frequencies that are used for signal generation and the RF impairment values.
fSym = 2e6; % Symbol rate or Baud rate - this is useful while modeling RF impairments cfo = 2e5; % Percentage of symbol rate as CFO
Initialize the energy per bit to noise power ratio (Eb/N0), which is used to calculate the SNR using system parameters.
EbN0 = 10; % To see a BER curve, run the simulation for 3.2:0.2:5
Initialize the parameters to terminate the simulation. The parameters are set to small values in this example to get quick results. Increase these parameters value to get a smoother BER curve.
maxNumErrors = 1e2; % Simulation stops after maxNumErrors bit errors maxNumBits = 1e5; % Simulation stops after processing maxNumBits irrespective of number of errors % Set maxNumBits = 1e8 for a smoother BER curve maxFramesLost = 1e2; % Simulation stops after maxFramesLost frames are lost
Initialize all of the objects that are required for the proper functioning of the end-to-end chain.
Create a CCSDS TM waveform generator with these parameters by using the
ccsdsTMWaveformGenerator System object™. Display the properties of the object.
tmWaveGen = ccsdsTMWaveformGenerator("ChannelCoding",channelCoding, ... "NumBytesInTransferFrame",transferFrameLength, ... "Modulation",modScheme, ... "RolloffFactor",alpha, ... "SamplesPerSymbol",sps); disp(tmWaveGen)
ccsdsTMWaveformGenerator with properties: WaveformSource: "synchronization and channel coding" NumBytesInTransferFrame: 1115 HasRandomizer: true HasASM: true PCMFormat: "NRZ-L" Channel coding properties: ChannelCoding: "convolutional" ConvolutionalCodeRate: "1/2" Digital modulation and filter properties: Modulation: "QPSK" PulseShapingFilter: "root raised cosine" RolloffFactor: 0.3500 FilterSpanInSymbols: 10 SamplesPerSymbol: 8 Use get to show all properties
Calculate the SNR from the Eb/N0 and initialize the parameters related to the calculation of BER.
rate = tmWaveGen.info.ActualCodeRate; M = tmWaveGen.info.NumBitsPerSymbol; numBitsInTF = tmWaveGen.NumInputBits; snr = EbN0 + 10*log10(rate) + 10*log10(M) - 10*log10(sps); % As signal power is scaled to one while introducing noise, snr value should be reduced by a factor of SPS numSNR = length(snr); ber = zeros(numSNR,1); % Initialize the BER parameter bercalc = comm.ErrorRate;
Create a receive filter object by using the
comm.RaisedCosineReceiveFilter System object.
b = rcosdesign(alpha,tmWaveGen.FilterSpanInSymbols,sps); % |H(f)| = 1 for |f| < fN(1-alpha) - Annex 1 in Section 2.4.17A in  Gain = sum(b); rxFilterDecimationFactor = sps/2; rxfilter = comm.RaisedCosineReceiveFilter("DecimationFactor",rxFilterDecimationFactor, ... "InputSamplesPerSymbol",sps, ... "RolloffFactor",alpha, ... "Gain",Gain);
Model frequency and phase offsets by using the
comm.PhaseFrequencyOffset System object. Compensate for frequency and phase offset at the receiver in two steps.
phaseOffset = pi/8; fqyoffsetobj = comm.PhaseFrequencyOffset( ... "FrequencyOffset",cfo, ... "PhaseOffset",phaseOffset, ... "SampleRate",sps*fSym); coarseFreqSync = comm.CoarseFrequencyCompensator( ... "Modulation",modScheme, ... "FrequencyResolution",100, ... "SampleRate",sps*fSym); fineFreqSync = comm.CarrierSynchronizer("DampingFactor",1/sqrt(2), ... "NormalizedLoopBandwidth",0.0007, ... "SamplesPerSymbol",1, ... "Modulation",modScheme);
Create a variable fractional delay object by using the dsp.VariableFractionalDelay System object, which introduces the fractional delay in the transmitted waveform. Create a symbol synchronization object by using the
varDelay = dsp.VariableFractionalDelay("InterpolationMethod","Farrow"); fixedDelayVal = 10.2; Kp = 1/(pi*(1-((alpha^2)/4)))*cos(pi*alpha/2); symsyncobj = comm.SymbolSynchronizer( ... "DampingFactor",1/sqrt(2), ... "DetectorGain",Kp, ... "TimingErrorDetector","Gardner (non-data-aided)", ... "Modulation","PAM/PSK/QAM", ... "NormalizedLoopBandwidth",0.0001, ... "SamplesPerSymbol",sps/rxFilterDecimationFactor);
Demodulate and decode the received signal by using the
HelperCCSDSTMDecoder helper files, respectively. Display the properties of the resulting objects.
demodobj = HelperCCSDSTMDemodulator("Modulation",modScheme,"ChannelCoding",channelCoding); decoderobj = HelperCCSDSTMDecoder("ChannelCoding",channelCoding, ... "NumBytesInTransferFrame",transferFrameLength, ... "Modulation",modScheme); disp(demodobj)
HelperCCSDSTMDemodulator with properties: Modulation: "QPSK" PCMFormat: "NRZ-L" ChannelCoding: "convolutional"
HelperCCSDSTMDecoder with properties: ChannelCoding: "convolutional" HasRandomizer: true HasASM: true DisableFrameSynchronization: 0 DisablePhaseAmbiguityResolution: 0 NumBytesInTransferFrame: 1115 ConvolutionalCodeRate: "1/2" ViterbiTraceBackDepth: 60 ViterbiTrellis: [1x1 struct] ViterbiWordLength: 8 Modulation: "QPSK" PCMFormat: "NRZ-L"
Initialize the constellation diagram object by using the
costellationobj = comm.ConstellationDiagram; % Default view is for QPSK if strcmp(modScheme,'BPSK') costellationobj.ReferenceConstellation = [1, -1] end
To simulate the end-to-end chain and measure the BER of the CCSDS TM system, follow these steps.
Generate random bits to form a TF.
Generate the TM waveform by passing the TF through the
ccsdsTMWaveformGenerator System object.
Introduce RF impairments, such as CFO and symbol delay.
Add AWGN to the RF impaired signal. This noisy signal is considered the received signal.
Pass the received signal through coarse frequency correction, which performs the initial coarse carrier frequency synchronization. Coarse frequency estimation is done using the "FFT-based" algorithm.
Use a matched filter (root raised cosine filter) with the same configuration that is applied at the transmitter end. Because the symbol timing synchronization module works at a sampling rate that is higher than the symbol rate, the complex baseband samples are not down sampled to the symbol rate after filtering. It is down sampled such that at least 2 samples per symbol exist.
Perform symbol timing synchronization by using the Gardner timing error detector (TED) to remove the timing offset that is present in the signal.
Perform carrier frequency and phase tracking by using the
Visualize the constellation after symbol timing and carrier frequency synchronization is complete. Observe how the constellation evolves over multiple iterations.
Demodulate the received signal and verify that the signal is at the symbol rate (that is, the samples per symbol is 1).
Perform frame synchronization and channel decoding to resolve the phase ambiguity, synchronize the frame to the start of the ASM, and then decode the synchronized frame to recover the TF.
numBitsForBER = 8; % For detecting which frame is synchronized numMessagesInBlock = 2^numBitsForBER; for isnr = 1:numSNR rng default; % Reset for every SNR value to get repeatable results reset(bercalc); % Reset the comm.ErrorRate object for every SNR point berinfo = bercalc(int8(1), int8(1)); % Initialize berinfo to be used before BER is calculated tfidx = 1; numFramesLost = 0; prevdectfidx = 0; inputBuffer = zeros(numBitsInTF, 256,"int8"); while((berinfo(2) < maxNumErrors) && (berinfo(3) < maxNumBits) && (numFramesLost < maxFramesLost)) seed = randi([0 2^32-1],1,1); % Generate seed and use in AWGN channel for repeatable simulation % Transmitter side processing bits = int8(randi([0 1],numBitsInTF-numBitsForBER,1)); % The first 8 bits correspond to the TF index modulo 256. When % synchronization modules are included, there can be few frames % where synchronization is lost temporarily and then locks again. % In such cases, to calculate the BER, these 8 bits aid in % identifying which TF is decoded. If an error in these 8 bits % exists, then this error is detected by looking at the difference % between consecutive decoded bits. If an error is detected, then % that frame is considered lost. Even though the data link layer is % out of scope of this example, the data link layer has a similar % mechanism. In this example, only for calculating the BER, this % mechanism is adopted. The mechanism that is adopted in this % example is not as specified in the data link layer of the CCSDS % standard. And this mechanism is not specified in the physical % layer of the CCSDS standard. msg = [de2bi(mod(tfidx-1,numMessagesInBlock),numBitsForBER,"left-msb").';bits]; inputBuffer(:,mod(tfidx-1,numMessagesInBlock)+1) = msg; tx = tmWaveGen(msg); % Introduce RF impairments cfoInroduced = fqyoffsetobj(tx); % Introduce CFO delayed = varDelay(cfoInroduced,fixedDelayVal); % Introduce timing offset rx = awgn(delayed, snr(isnr),'measured',seed); % Add AWGN % Receiver side processing coarseSynced = coarseFreqSync(rx); % Apply coarse frequency synchronization filtered = rxfilter(coarseSynced); % Filter received samples through RRC filter TimeSynced = symsyncobj(filtered); % Apply symbol timing synchronization fineSynced = fineFreqSync(TimeSynced); % Track Frequency and phase % Visualize constellation if seeConstellation % Plot constellation of first 1000 symbols in a TF so % that variable size of fineSynced does not impede the % requirement of constant input size for the % comm.ConstellationDiagram System object. costellationobj(fineSynced(1:1000)); end demodData = demodobj(fineSynced); % Demodulate decoded = decoderobj(demodData); % Perform phase ambiguity resolution, frame synchronization and channel decoding % Calculate BER and adjust all buffers accordingly dectfidx = bi2de(double(decoded(1:8).'),"left-msb")+1; % See the value of first 8 bits if tfidx > 30 % Consider to calculate BER only after 30 TFs are processed % As the value of first 8 bits is increased by one in each % iteration, if the difference between the current decoded % value of decimal value of first 8 bits is not equal to the % previous decoded one, then it indicates a frame loss. if dectfidx - prevdectfidx ~= 1 numFramesLost = numFramesLost + 1; disp(['Frame lost at tfidx: ' num2str(tfidx) '. Total frames lost: ' num2str(numFramesLost)]); else berinfo = bercalc(inputBuffer(:,dectfidx),decoded); if nnz(inputBuffer(:,dectfidx)-decoded) disp(['Errors occurred at tfidx: ' num2str(tfidx) '. Num errors: ' num2str(nnz(inputBuffer(:,dectfidx) - decoded))]) end end end prevdectfidx = dectfidx; % Update tfidx tfidx = tfidx + 1; end fprintf("\n"); currentBer = berinfo(1); ber(isnr) = currentBer; disp(['Eb/N0: ' num2str(EbN0(isnr)) '. BER: ' num2str(currentBer) '. Num frames lost: ' num2str(numFramesLost)]); % Reset objects reset(tmWaveGen); reset(fqyoffsetobj); reset(varDelay); reset(coarseFreqSync); reset(rxfilter); reset(symsyncobj); reset(fineFreqSync); reset(demodobj); reset(decoderobj); end
Eb/N0: 10. BER: 0. Num frames lost: 0
This example demonstrates BER simulation for convolutional codes with QPSK modulation in the presence of several RF impairments. To observe the end-to-end simulation chain for different scenarios, change the properties related to the channel coding and modulation schemes. The modulation schemes that are supported by the receiver in this example are BPSK and QPSK. The channel coding schemes that are supported by the receiver in this example are none (that is no channel coding), RS, convolutional and concatenated codes.
Run a full BER simulation by setting Eb/N0 value to 3.2:0.2:5 and observe the BER by setting
maxNumBits to 1e8. Use this code to plot the BER results.
% semilogy(EbN0,ber); % grid on; % xlabel('E_b/N_0 (dB)'); % ylabel('BER'); % title('BER plot');
Always reserve the initial few TFs for the synchronizers (symbol timing and carrier frequency synchronizers) to lock. This example discards the first 30 TFs. This number can vary based on the SNR at which the receiver is operating and the parameters of the synchronization loops, such as loop bandwidth and damping factor. If you operate the receiver at low SNR and observe large errors in the initial values of
tfidx, then the synchronizers are not yet locked. For the given simulation parameters, discard the initial TFs as appropriate. The second output arguments of
The example uses the following helper files:
 TM Synchronization and Channel Coding. Recommendation for Space Data System Standards, CCSDS 131.0-B-3. Blue Book. Issue 3. Washington, D.C.: CCSDS, September 2017.
 Radio Frequency and Modulation Systems--Part 1: Earth Stations and Spacecraft. Recommendation for Space Data System Standards, CCSDS 401.0-B-30. Blue Book. Issue 30. Washington, D.C.: CCSDS, February 2020.
 TM Synchronization and Channel Coding - Summary of Concept and Rationale. Report Concerning Space Data System Standards, CCSDS 130.1-G-3. Green Book. Issue 3. Washington, D.C.: CCSDS, June 2020.