Skip to main content

Add AWGN Directly to PSD in MATLAB

 

In general, we compute the power spectral density (PSD) of a noisy periodic signal. However, in this article, you will learn how to add noise directly to the PSD of a signal. This process is approximately equivalent to adding noise to a clean signal and then computing its PSD. Here, I will discuss both the theoretical background and the MATLAB implementation.

Steps

1. First, compute the Fast Fourier Transform (FFT) of the clean signal. Then, calculate the Power Spectral Density (PSD) from the FFT.

2. In our case, ensure that the PSD is in the linear scale. Next, compute the noise power from the given Signal-to-Noise Ratio (SNR) using:

    noise_power = signal power / linear SNR
    

3. Then, generate Additive White Gaussian Noise (AWGN) using the formula:

    AWGN noise = sqrt(noise_power) * randn
    

      where randn generates a Gaussian-distributed signal with a mean of 0 and a variance of 1.

 

MATLAB Code 

clc; clear; close all;

%% Define Parameters
fs = 1000; % Sampling frequency (Hz)
T = 0.2; % Time period of sine wave (s)
A = 1; % Amplitude
N = 1024; % Number of samples
t = linspace(-0.5, 0.5, N); % Time vector
f_sin = 5; % Frequency of sine wave (Hz)

%% Generate Periodic Sine Wave
sine_wave = A * sin(2 * pi * f_sin * t);

%% Compute PSD using FFT
Xf = fftshift(fft(sine_wave)); % Compute FFT and shift
PSD = abs(Xf).^2 / N; % Compute Power Spectral Density

%% Generate AWGN in Frequency Domain (Method 1)
snr_dB = 20; % SNR in dB
snr_linear = 10^(snr_dB/10); % Convert SNR to linear scale
signal_power = mean(PSD); % Approximate power of the original spectrum
noise_power = signal_power / snr_linear; % Compute noise power
noise_spectrum = sqrt(noise_power) .* (randn(size(PSD)) + 1j*randn(size(PSD))); % AWGN

%% Add AWGN Directly to PSD
noisy_PSD = PSD + abs(noise_spectrum).^2; % Add noise power to PSD

%% Generate AWGN in Time Domain (Method 2)
noise_time = sqrt(noise_power) * randn(size(sine_wave)); % AWGN in time domain
noisy_sine = sine_wave + noise_time; % Add noise to signal

%% Compute PSD of Noisy Sine Wave
Xf_noisy = fftshift(fft(noisy_sine)); % Compute FFT of noisy signal
PSD_noisy = abs(Xf_noisy).^2 / N; % Compute Power Spectral Density

%% Plot Results
freq = linspace(-fs/2, fs/2, N); % Frequency axis

figure;

% Plot Time-Domain Sine Wave
subplot(3,1,1);
plot(t, sine_wave, 'b', 'LineWidth', 1.5); hold on;
plot(t, noisy_sine, 'r', 'LineWidth', 1.2);
xlabel('Time (s)');
ylabel('Amplitude');
title('Sine Wave Before and After AWGN');
legend('Original Sine Wave', 'Noisy Sine Wave');
grid on;

% Plot PSD Comparison (Direct AWGN to PSD)
subplot(3,1,2);
plot(freq, 10*log10(PSD + eps), 'b', 'LineWidth', 1.5); hold on;
plot(freq, 10*log10(noisy_PSD + eps), 'r', 'LineWidth', 1.5);
xlabel('Frequency (Hz)');
ylabel('Power Spectral Density (dB)');
title('AWGN Added Directly to PSD');
legend('Original PSD', 'PSD with Direct AWGN');
grid on;

% Plot PSD Comparison (AWGN in Time Domain)
subplot(3,1,3);
plot(freq, 10*log10(PSD + eps), 'b', 'LineWidth', 1.5); hold on;
plot(freq, 10*log10(PSD_noisy + eps), 'g', 'LineWidth', 1.5);
xlabel('Frequency (Hz)');
ylabel('Power Spectral Density (dB)');
title('PSD: Original vs. PSD from Noisy Sine Wave');
legend('Original PSD', 'PSD from Noisy Signal');
grid on;

Output

 





Copy the MATLAB Code from here 

 

Further Reading 

  1. Periodogram in MATLAB

People are good at skipping over material they already know!

View Related Topics to







Admin & Author: Salim

profile

  Website: www.salimwireless.com
  Interests: Signal Processing, Telecommunication, 5G Technology, Present & Future Wireless Technologies, Digital Signal Processing, Computer Networks, Millimeter Wave Band Channel, Web Development
  Seeking an opportunity in the Teaching or Electronics & Telecommunication domains.
  Possess M.Tech in Electronic Communication Systems.


Contact Us

Name

Email *

Message *

Popular Posts

MATLAB code for MSK

 Copy the MATLAB Code from here % The code is developed by SalimWireless.com clc; clear; close all; % Define a bit sequence bitSeq = [0, 1, 0, 0, 1, 1, 1, 0, 0, 1]; % Perform MSK modulation [modSignal, timeVec] = modulateMSK(bitSeq, 10, 10, 10000); % Plot the modulated signal subplot(2,1,1); samples = 1:numel(bitSeq); stem(samples, bitSeq); title('Original message signal'); xlabel('Time (s)'); ylabel('Amplitude'); % Plot the modulated signal subplot(2,1,2); samples = 1:10000; plot(samples / 10000, modSignal(1:10000)); title('MSK modulated signal'); xlabel('Time (s)'); ylabel('Amplitude'); % Perform MSK demodulation demodBits = demodMSK(modSignal, 10, 10, 10000); % Function to perform MSK modulation function [signal, timeVec] = modulateMSK(bits, carrierFreq, baudRate, sampleFreq) % Converts a binary bit sequence into an MSK-modulated signal % Inputs: % bits - Binary input sequence % carrierFreq - Carri...

Difference between AWGN and Rayleigh Fading

Wireless Signal Processing Gaussian and Rayleigh Distribution Difference between AWGN and Rayleigh Fading 1. Introduction Rayleigh fading coefficients and AWGN, or additive white gaussian noise [↗] , are two distinct factors that affect a wireless communication channel. In mathematics, we can express it in that way.  Fig: Rayleigh Fading due to multi-paths Let's explore wireless communication under two common noise scenarios: AWGN (Additive White Gaussian Noise) and Rayleigh fading. y = h*x + n ... (i) Symbol '*' represents convolution. The transmitted signal  x  is multiplied by the channel coefficient or channel impulse response (h)  in the equation above, and the symbol  "n"  stands for the white Gaussian noise that is added to the signal through any type of channel (here, it is a wireless channel or wireless medium). Due to multi-paths the channel impulse response (h) changes. And multi-paths cause Rayleigh fa...

Constellation Diagrams of ASK, PSK, and FSK

BASK (Binary ASK) Modulation: Transmits one of two signals: 0 or -√Eb, where Eb​ is the energy per bit. These signals represent binary 0 and 1.    BFSK (Binary FSK) Modulation: Transmits one of two signals: +√Eb​ ( On the y-axis, the phase shift of 90 degrees with respect to the x-axis, which is also termed phase offset ) or √Eb (on x-axis), where Eb​ is the energy per bit. These signals represent binary 0 and 1.  BPSK (Binary PSK) Modulation: Transmits one of two signals: +√Eb​ or -√Eb (they differ by 180 degree phase shift), where Eb​ is the energy per bit. These signals represent binary 0 and 1.  Key Points For Binary Amplitude Shift Keying (BASK), binary bit '0' can be represented as lower level voltage or no signal and bit '1' as higher level voltage.  For Binary Frequency Shift Keying (BFSK), you can map binary bit '0' to 'j' and bit '1' to '1'. So, signals are in phase.  A phase shift of 0 degrees could represent a binary '1...

BER vs SNR for M-ary QAM, M-ary PSK, QPSK, BPSK, ...

Modulation Constellation Diagrams BER vs. SNR BER vs SNR for M-QAM, M-PSK, QPSk, BPSK, ... What is Bit Error Rate (BER)? The abbreviation BER stands for bit error rate, which indicates how many corrupted bits are received (after the demodulation process) compared to the total number of bits sent in a communication process. It is defined as,  In mathematics, BER = (number of bits received in error / total number of transmitted bits)  On the other hand, SNR refers to the signal-to-noise power ratio. For ease of calculation, we commonly convert it to dB or decibels.   What is Signal the signal-to-noise ratio (SNR)? SNR = signal power/noise power (SNR is a ratio of signal power to noise power) SNR (in dB) = 10*log(signal power / noise power) [base 10] For instance, the SNR for a given communication system is 3dB. So, SNR (in ratio) = 10^{SNR (in dB) / 10} = 2 Therefore, in this instance, the s...

MATLAB code for BER vs SNR for M-QAM, M-PSK, QPSk, BPSK, ...

Modulation Constellation Diagrams BER vs. SNR MATLAB code for BER vs SNR for M-QAM, M-PSK, QPSk, BPSK, ...   MATLAB Script for  BER vs. SNR for M-QAM, M-PSK, QPSk, BPSK %Written by Salim Wireless %Visit www.salimwireless.com for study materials on wireless communication %or, if you want to learn how to code in MATLAB clc; clear; close all; % Parameters num_symbols = 1e5; % Number of symbols snr_db = -20:2:20; % Range of SNR values in dB % PSK and QAM orders to be tested psk_orders = [2, 4, 8, 16, 32]; qam_orders = [4, 16, 64, 256]; % Initialize BER arrays ber_psk_results = zeros(length(psk_orders), length(snr_db)); ber_qam_results = zeros(length(qam_orders), length(snr_db)); % BER calculation for each PSK order and SNR value for i = 1:length(psk_orders) psk_order = psk_orders(i); for j = 1:length(snr_db) % Generate random symbols data_symbols = randi([0, psk_order-1], 1, num_symb...

Comparisons among ASK, PSK, and FSK | And the definitions of each

Modulation ASK, FSK & PSK Constellation MATLAB Simulink MATLAB Code Comparisons among ASK, PSK, and FSK    Comparisons among ASK, PSK, and FSK Comparison among ASK,  FSK, and PSK Performance Comparison: 1. Noise Sensitivity:    - ASK is the most sensitive to noise due to its reliance on amplitude variations.    - PSK is less sensitive to noise compared to ASK.    - FSK is relatively more robust against noise, making it suitable for noisy environments. 2. Bandwidth Efficiency:    - PSK is the most bandwidth-efficient, requiring less bandwidth than FSK for the same data rate.    - FSK requires wider bandwidth compared to PSK.    - ASK's bandwidth efficiency lies between FSK and PSK. Bandwidth Calculator for ASK, FSK, and PSK The baud rate represents the number of symbols transmitted per second Select Modulation Type: ASK...

Fundamentals of Channel Estimation

Channel Estimation Techniques Channel Estimation is an auto-regressive process that may be performed with a number of iterations. There are commonly three types of channel estimation approaches. 1. Pilot estimation  2. Blind estimation  3. Semi-blind estimation. For Channel Estimation,  CIR [↗] is used. The amplitudes of the impulses decrease over time and are not correlated. For example, y(n) = h(n) * x(n) + w(n) where y(n) is the received signal, x(n) is the sent signal, and w(n) is the additive white gaussian noise At the next stage, h(n+1) = a*h(n) + w(n) The channel coefficient will be modified as stated above at the subsequent stage. The scaling factor "a" determines the impulse's amplitude, whereas "h(n+1)" represents the channel coefficient at the following stage. Pilot Estimation Method To understand how a communication medium is currently behaving, a channel estimate is necessary. In order to monitor a channel's behavior in practice communication ...

MATLAB Code for Pulse Amplitude Modulation (PAM) and Demodulation

  Pulse Amplitude Modulation (PAM) & Demodulation MATLAB Script clc; clear all; close all; fm= 10; % frequency of the message signal fc= 100; % frequency of the carrier signal fs=1000*fm; % (=100KHz) sampling frequency (where 1000 is the upsampling factor) t=0:1/fs:1; % sampling rate of (1/fs = 100 kHz) m=1*cos(2*pi*fm*t); % Message signal with period 2*pi*fm (sinusoidal wave signal) c=0.5*square(2*pi*fc*t)+0.5; % square wave with period 2*pi*fc s=m.*c; % modulated signal (multiplication of element by element) subplot(4,1,1); plot(t,m); title('Message signal'); xlabel ('Time'); ylabel('Amplitude'); subplot(4,1,2); plot(t,c); title('Carrier signal'); xlabel('Time'); ylabel('Amplitude'); subplot(4,1,3); plot(t,s); title('Modulated signal'); xlabel('Time'); ylabel('Amplitude'); %demdulated d=s.*c; % At receiver, received signal is multiplied by carrier signal filter=fir1(200,fm/fs,'low'); % low-pass FIR fi...