How to generate AWGN noise in Matlab/Octave (without using in-built awgn function)

1 Star2 Stars3 Stars4 Stars5 Stars (7 votes, average: 4.57 out of 5)
Loading...

This article is also available in PDF format – Free Download here1357657499_ACP_PDF_3_file_document

Looking for the proper way to generate AWGN noise in Matlab/Octave ? Here you go…

AWGN – the in-built function

Matlab/Octave communication toolbox has an inbuilt function named – \(awgn()\) with which one can add an Additive Gaussian White Noise to obtain the desired \(Signal-to-Noise Ratio\) (SNR). The main usage of this function is to add AWGN to a clean signal (infinite SNR) in order to get a resultant signal with a given SNR (usually specified in dB). This usage is often found in signal processing/digital communication applications. For example, in Monte Carlo simulations involving modeling of modulation/demodulation systems, the modulated symbols at the transmitter are added with a random noise of specific strength, in order to simulate a specific \(E_b/N_0\) or \(E_s/N_0\) point.

The function \(y=awgn(x,SNR,’measured’)\), first measures the power of the signal vector \(x\) and then adds white Gaussian Noise to \(x\) for the given \(SNR\) level in dB. The resulting signal y is guaranteed to have the specified SNR.

Custom function to add AWGN noise

If you do not have the communication toolbox, or if you would like to mimic the in-built AWGN function in any programming language, the following procedure can be used.
1) Assume, you have a vector \(x\) to which an AWGN noise needs to be added for a given \(SNR\) (specified in dB).
2) Measure the power in the vector \(x\)
$$E_s=\frac{1}{L} \sum_{i=0}^{L-1} \left|x[i] \right|^2 ; \;\;\; where \;\;\;\; L=length(x) $$
3) Convert given \(SNR\) in dB to linear scale (\(SNR_{lin}\)) and find the noise vector (from Gaussian distribution of specific noise variance)  using the equations below

$$ noise= \begin{cases} \sqrt{\frac{E_s}{SNR_{lin}}}*randn(1,L) & if \;\; x \;\; is \;\; real \\ \\ \sqrt{\frac{E_s}{2* SNR_{lin}}}*\left [randn(1,L)+j*randn(1,L) \right] & if \;\; x \;\; is \;\; complex
\end{cases} $$

4) Finally add the generated noise vector to the signal \(x\)

$$y=x+noise$$

The custom function – add_awgn_noise():

The custom function written in Matlab, that mimics the \(awgn\) function is given below. It can be easily ported to Octave.

Comparison & Testing:

Let’s cross check the results obtained from the above function with that of the standard in-built \(awgn\) function in Matlab. Testing and comparison is done using two test waveforms – 1) sawtooth waveform (represented by a vector containing only real numbers) , 2) A complex sinusoidal waveform (vector having both real and imaginary part). The testing below is done for SNR=5 dB. Users are advised to test the function for various ranges of SNRs.

Test code and results for sawtooth waveform (the input signal x is a vector of real numbers) is given next. Test of linearity indicates that the results from the custom function matches with that of the in-built function.

Sawtooth waveform with AWGN noise
Sawtooth waveform with AWGN noise generated using custom method and Matlab’s in-built method

Output from add_awgn_function Vs inbuilt awgn function
Output from add_awgn_function Vs inbuilt awgn function

Test code and results for complex sinusoidal waveform (the input signal x is a vector of complex numbers) is given next. Test of linearity indicates that the results from the custom function matches with that of the in-built function.

Complex sinusoidal waveform with AWGN noise generated using custom method and Matlab's in-built method
Complex sinusoidal waveform with AWGN noise generated using custom method and Matlab’s in-built method

Output from add_awgn_function Vs inbuilt awgn function
Output from add_awgn_function Vs inbuilt awgn function

Creative Commons License
How to generate AWGN noise in Matlab/Octave by Mathuranathan Viswanathan is licensed under a Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License.

Recommended Books:

  • Chema

    Hi,

    I have the next question: the awgn object has a parameter called “SamplesPerSymbol”, where this parameter is used in the custom equation ?

    Thanks.

    • The custom function given here, mimics the ‘awgn’ function given here : http://www.mathworks.com/help/comm/ref/awgn.html
      It computes the required noise from the given SNR that is defined as ‘SNR per sample’. So there is no connotation of ‘samplesPersymbol’ in this case.

      Whereas, the comm toolbox object that you are referring to is more sophisticated, that includes different arguments like EbNo, EsNo, SNR, BitsPerSymbol, SignalPower, SamplesPerSymbol, and Variance etc..,

      Of course, the function given here can be extended to include all these parameters above.

  • Alam

    If I have multiple complex signals in parallel , e.g a 6x 5000 matrix. How can I add awgn noise to them using this function?

    • The function given here treats the signal as a 1-dimensional vector The above function needs to be modified for your specific case (for example noise addition in MIMO case). I will post an article on this soon. Thanks for the question.

      • Alam

        Thank you . looking forward for the article

  • Waqas Anjum

    Nice code ….I want to ask about my specific problem…..
    I am building an end to end basic Communication transceiver.My signal at Transmitter before adding AWGN noise is over sampled by factor of 10.I am adding a noise as per your code without considering oversampling factor during noise addition.I have verified the channel block and at receiver processed the channel signal(catering the upsampling factor at receiver) but upon comparing BER curve with theoretical results (berawgn function) ,there is 4dB difference/offset in BER curves(My results are 4dB better than theoretical BER curve as per attachment)…..The only thing I am worried is about upsampling factor incorporation in AWGN block while my transmitter and receiver is bug free…..
    Please tell how to cater this upsampling factor during noise addition.

    • That curve is expected, if you do not incorporate the oversampling factor in the code.
      To incorporate the oversampling factor (OF), change line 17 as in add_awgn_noise function as
      Esym=OF*sum(abs(x).^2)/(L); %Calculate actual symbol energy

      • Waqas Anjum

        Thanks dear for your response.I did as per your suggestion but this time the ber performance becomes 6dB poor as per attachment below.Please suggest Whats wrong …Do I need to use
        Esym=sqrt(OF)*sum(abs(x).^2)/(L); %Calculate actual symbol energy
        instead.Thanks for ur time again.

        • what is the format of the modulated signal from the transmitter ? Are you using a complex number representation to represent inphase and quadrature signal- Complex baseband equivalent model? In this case, you do not need a oversampling factor in the simulation.

          If the modulation signal is not in complex format, you can add noise as follows (I have verified it for all PSK modulations, QAM should be okay too). The only difference here is the noise sigma set to N0/2 (same as in the case of complex signal)

          SNR = 10ˆ(SNR_dB/10); %SNR to linear scale
          Esym=OF*sum(abs(x).^2)/(L); %Calculate actual symbol energy
          N0=Esym/SNR; %Find the noise spectral density
          noiseSigma = sqrt(N0/2) #note the factor N0/2
          n = noiseSigma*randn(1,L);%computed noise
          y = x+n

          Instead, if your simulation is a passband simulation, then the oversampling factor is relevant.

          • Waqas Anjum

            Yes the transmitted signal is a complex signal and noise I added is both to real and imaginary part catering the No/2 factor you mentioned.
            Also,What I foresee about the possible mistake is that I am passing my vector named snrVecdB to berawgn function instead of Eb/No as berawgn function demands Eb/No instead of snr Vector.
            snrVecdB = -4:2:10; % in db
            bertheory= berawgn(snrVecdB(snrCount),’qam’,M);
            Remember it will matter as I am computing Es/No in simulation during adding noise in AWGN channel.

            • Possible. The method given here will give identical results as with commtool box’s awgn function. You could verify your model using the standard function too.

  • RACHIT

    sir,
    how would i generate uniformly distributed two INDEPENDENT random sequences???

  • Tariq

    Dear Sir @Mathuranathan:disqus , How about the Bandwidth of the signal ,I guess we have to adjust the values of the noise power which located only within this banwidth right? … For AWGN case, when adding the noise , its distribution is Flat along the whole frequencies , so how much the noise power within the signal spectrum? our interest should be for the noise power values located within this bandwidth ,isn’t it? . I am very thankful for this great website which benefiting me a lot.

    • dara_parsavand

      If you first generate a real noise signal x, with zero mean and variance sigma^2 (so the power is sigma^2 too) and you use it in a simulation for a particular sample rate Fs, then the spectral density is sigma^2/(Fs/2). For band-limited noise, It’s clearer to most people I think if you first generate the noise power based on the density you want (and thus you need more power the bigger Fs is). Then you filter it with some filter (typically FIR to match some analog filter or the particular digital filter you are using) and choose the filter so that it has unity DC gain (sum of coefficients is 1). Then after you filter the white noise, it will have power = BW*sigma^2/(Fs/2), where BW is the one sided bandwidth of the filter (0 to cutoff freq). Note that one BW = Fs/2 (which is no filtering), you are still at sigma^2 but for any smaller BW, you have reduced the power.

      • Tariq

        Dear @dara_parsavand:disqus , I really appreciate your response . Thanks a lot

  • mayank

    Want to add simple Phase noise in my signal. when i am analyse the signals in frequency domain i do not get desired output. For the ideal signal, it shoud be only single peak whereas after adding noise it shoud be having some fluctuations. What is the mistake i am making?

    sig = cos(2*pi*30*t); %%%% original signal
    sig_fft = fft(sig);
    figure,plot(abs(sig_fft));
    sigma = 2
    sig2 = sig + sigma*randn(1,length(t)); %%%%% adding phase noise
    sig2_fft = fft(sig2);
    figure,plot(abs(sig2_fft));

    • Your code adds noise to the signal amplitude.
      For phase noise, use this equation
      Cos(2pi 30 t + theta)
      Theta is the phase and it should be a random variable

      • mayank

        Thanks Mathuranathan for your reply.

        ahh ok. That means i can use Theta = randn(1,length(t)); right?

        • Usually, the phase noise is considered UNIFORM drawn from the interval [-pi,pi].
          Therefore phase noise needs to be from -pi to +pi and it can be drawn from a uniform distribution. Use rand function for uniform rvs instead of randn function which is for Gaussian rvs.

          The rand function must generate uniform phase from -pi to +pi. Read thorough the matlab doc on how to generate a uniform distribution in the interval [a,b] and substitute a by -pi and b by +pi.

          • mayank

            Thanks for your reply.

            I will go through matlab doc :).

            • mayank

              Hey Mathuranathan,
              As you told me i used rand to generate phase noise. However, when i see the graph in frequency domain the spectrum is not broading. Can you tell me what mistake i am making?

              sigma = 0.5;
              t = 0:0.01:1;
              sig = cos(2*pi*4*t); % original signal
              A = sigma *randn(1,length(t)); % standard deviation of phase noise, increase
              sig2 = (A).*cos(2*pi*4*t + sigma*rand(1,length(t))); % signal with phase noise & amplitud4e noise
              figure,plot(t,sig,’r’,t,sig2, ‘b’);
              figure, plot(abs(fft(sig,10*length(sig))));
              figure, plot(abs(fft(sig2,10*length(sig2))));
              return

              • (A).*cos(2*pi*4*t + sigma*rand(1,length(t)))

                1) Implies you are randomizing the amplitude (A) too.
                2) why there is a sigma for uniform distribution (rand) ? Sigma quantifies a normal distribution (randn)
                function
                The rand function must generate a random number between -pi and pi. The code does not do that.

                • mayank

                  ohh ok. thanks for your comments :)

                • no name

                  method : sum of sinusodial model of channel???
                  write eq plz and cod in matlab ????

          • dara_parsavand

            Phase noise better not be uniform from -pi to pi. That much phase noise will ruin any comm system. If you want to see the PDF of what an example phase noise would look like and you have the Matlab Comm toolbox you can type:

            Nsamp = 1e5;
            Nhist = 1e3;
            hPhNoise = comm.PhaseNoise(‘Level’,[-60 -80], ‘FrequencyOffset’,[20 200], ‘SampleRate’,1024); % an example call from the doc pages
            y = step(hPhNoise, 1j*ones(Nsamp,1)); % instantiates phase noise on a trivial constant phasor which is complex (Matlab will complain if not complex)
            figure, hist((angle(y)-pi/2), Nhist)

            Sure enough the PDF is Gaussian. Why would that be? Because of the Law of Large Numbers – phase noise is the sum of a lot of little things going on with analog hardware that creates an imperfect sinusoid. So in fact mayank is correct about using randn and scaling it by sigma for phase noise. The problem is: though you can adjust a sigma to match an RMS phase noise spec, and then write code to use sigma*randn(1, N) to put in the argument of a cos or exp function, this will give white phase noise (independent with each time sample as is AWGN) when phase noise is never like this – it has a (more difficult to simulate) dBc/Hz spectrum that can only be viewed clearly on a log-log plot. So you have to filter the noise sequence before applying it to the argument. Again, if you have the Comm toolbox, there is a lot of info on doc comm.PhaseNoise. If you don’t have the toolbox (or even if you do), this link looked good: http://www.dsplog.com/2012/09/30/modeling-phase-noise-frequency-domain/

            • That’s not necessarily true. Phase noise assumption depends on the context. For example, the oscillator phase can be made random from one symbol time to another and still a modulation technique like fsk can be detected non-coherently. If the context is not given, I would generally assume noise to be uniform

              • dara_parsavand

                Mathuranathan, Good point about context. I agree that a slowly varying random phase can be tolerated in an FSK system (or any comm system for that matter if it varies slowly enough and can be tracked out). This random phase can walk around all over the place and modulo 2 pi, it will over the very long term be uniform from 0 to 2 pi. My basic understanding of phase noise is this: in the analog domain it is a random process that has all frequency content in it, but has infinite power when integrated from 0 freq to any finite BW. The way I typically model it is to ignore the section from 0 freq up to the PLL BW if a PLL is used (or if it is not used as in a non-coherent FSK system that you mention, some low freq limit that has no effect on BER). And there is a high speed cut off based on filtering that is used. So one only needs to consider the phase noise density from f1 to f2 and typically when integrated between these two frequencies, you might end up with 1 deg RMS or 20 deg RMS (tolerated in BPSK but not 256 QAM) or whatever, but you will have a zero mean random sequence (because you are ignoring the 0 freq to f1 part) which means the random variable will have a PDF that is a certain shape and I maintain at that point, the shape is typically Gaussian and not uniform. I do agree that context matters though and I suppose I need to know a lot more about mayank’s exact simulation task before commenting any further. Did you try my example code that I posted?

                • Dara_parsavand, Thanks for adding your insights and knowledge. Yes, I am aware of the phase noise function in the comm toolbox and I always use it in my simulations wherever required. Yes, as you mentioned, the phase level is all the more important for phase shift keying systems.

  • kool learner

    would the function be the same if things are in Eb/N0 instead? please explain.

    • If the Eb/N0 values are given, the noise should be generated accordingly as follows.

      Actually, the above function calculates the symbol energy (Es) and adds noise to the incoming signal accordingly. The input argument SNR_dB is equivalent to Es/N0 in dB.

      If Eb/N0 is given, it needs to be converted to Es/N0 appropriately. Examples follows

      Es/N0 dB = 10*log10(2) + Eb/N0 dB; %in linear scale for QPSK , 1 symbol = 2bits
      Es/N0 dB = 10*log10(4) + Eb/N0 dB; %in linear scale for 16-QAM , 1 symbol = 4bits

      Generic case:
      Es/N0 dB = 10*log10(k) + Eb/N0 dB; %in linear scale for any modulation where 1 symbol = k bits

      Now call the function

      y = add_awgn_noise(x,EsN0dB); %where x is the variable that contains the modulated samples for which the noise needs to be added

      • kool learner

        thank you for the clarification. much appreciated!

  • mohammad

    hi
    i want to generate dvb-t signal (digital video broadcasting torestrial) in matlab
    plz if any one have this codes give me (my mail:md.morady@gmail.com)
    tnx alot of this site member