How to remove the AWGN noise from data?

Let's say we have a vector data Uo=[20 30 40 50];

If we add AWGN noise to it, it becomes data U i.e.,

U=awgn(Uo,30);

Now we see that U and Uo are different. We want to get back our original data Uo. I tried as:

Uo=U-awgn(30);

But it gives error as:

Error using awgn (line 63)

Not enough input arguments.

What commands should I use to get back my vector Uo?

Answer 1

[Chunru]

Uo=[20 30 40 50];          % signal
U=awgn(Uo,30)              % add noise to signal

U = 1×4

19.9837 29.9773 40.0475 50.0504

% Uo=U-awgn(30);

% if you don't know the noise (which is U-Uo), you are not able to fully
% recover the signal.
% However, it is possible to "filter" out noise if signal and noise are
% different in certain characteristics.
% For your case, signal is a low-pass signal and lowpass filter can help to
% remove some noise.

Uo = [20:1:100]/30;
U = awgn(Uo, 10);

b = fir1(20, 0.02);  % low pass filter
Ufiltered = filtfilt(b, 1, U);  
t = 0:length(Uo)-1;
plot(t, Uo, 'r', t, U, 'b', t, Ufiltered, 'k')
legend("Signal", "Signal+noise", "Filtered")

[Sadiq Akbar]

Thanks a lot dear Chunru for your kind response. Actually I want to get back my my original vector Uo if not exactly then aprroximately equal.I don't need the plot. In your piece of code, you have used a low pass filter. But I don't know why did you use it? 2ndly you have not get back the original Uo. Can you do some filtering so that we can get exactly or approximately Uo. 3rdly my original vector consists of only four values but you have taken 81 values i.e., Uo = [20:1:100]/30; Why is this?

[Chunru]

The plot is for u to see how good the filtered result as an approximation. You can always comment it out when you are satisfied the results.

For filter to work. You need sufficient data. If you have only 4 points, you may not be able to do much with the added noise.

Answer 2

[Star Strider]

I would simply use a moving average filter, for example movmean, or choose one of the options in smoothdata, to eliminate normally-distributed additive noise.

[Chunru]

Perhaps that is a demonstration of the differences between data processing and signal processing :-)

[Star Strider]

I’m just suggesting a more straightforward approach. I frequently use (and recommend) smoothdata — particularly the 'sgolay' option — to eliminate broadband noise.

[Chunru]

I am a signal processing person and tend to think of filtering though I am also aware the smoothdata,movmean. Indeed, movmean is more STRAIGHTFORWARD approach (which is essentially a filter with uniform weight).

[Star Strider]

I have a reasonable background in signal processing. I generally favour the sgolayfilt function for broadband noise problems, and wavelets are also appropriate.

[Sadiq Akbar]

Thanks a lot for your valuable discussion dear Star Strider and Chunru. Some one told me to use a matched filter for this but didn't tell me how to use it? Can you guide me what is a matched filter and how to use it for this?

[Star Strider]

See: How to implement a matched filter for one approach. (That example use the filter function, although filtfilt could likely provide a better result, so try both and see what works best for you.)

[Sadiq Akbar]

Thanks a lot for your guidance dear Star Strider. I copied that code and made some changes in it and ran it but still it doesn't give me my original vector temp. I re-produce that changed code here:

clear all
clc
% A template is given
temp = [20 30 40 50]';

% Create a matched filter based on the template
b = flipud(temp(:));%My comment:Just reverse the vector temp up side down

% For testing the matched filter, create a random signal which
% contains a match for the template at some time index
x = [randn(200,1); temp(:); randn(300,1)];
n = 1:length(x);

% Process the signal with the matched filter
y = filter(b,1,x);

% Set a detection threshold (exmaple used is 90% of template)
thresh = 0.9

% Compute normalizing factor
u = temp.'*temp;

% Find matches
matches = n(y>thresh*u);

% Plot the results
plot(n,y,'b', n(matches), y(matches), 'ro');

% Print the results to the console
display(matches);

Can you guide me further what should I change so that I get back my original vector?

[Star Strider]

Getting your original vector may not be possible. It may be necessary to accept the result being ‘close enough’.

% clear all
% clc
% A template is given
temp = [20 30 40 50]';

% Create a matched filter based on the template
b = flipud(temp(:));%My comment:Just reverse the vector temp up side down

% For testing the matched filter, create a random signal which
% contains a match for the template at some time index
x = [randn(200,1); temp(:); randn(300,1)];
n = 1:length(x);

% Process the signal with the matched filter
y = filter(b,1,x);

% Set a detection threshold (exmaple used is 90% of template)
thresh = 0.9

thresh = 0.9000

% Compute normalizing factor
u = temp.'*temp;

% Find matches
matches = n(y>thresh*u);

% Plot the results
plot(n,y,'b', n(matches), y(matches), 'ro')
grid

% Print the results to the console
display(matches);

matches = 204

Check = y(y>500)

Check = 7×1

1.0e+03 * 0.9831 2.2866 3.7918 5.4000 3.8339 2.3067 1.0961

Also, the filter may need to be ‘normalised’, for example:

y = filter(b,sum(b),x);

I will let you experiment further with this.

.

[Sadiq Akbar]

Thanks a lot for your guidance dear Star Strider. If I look at the graph, I see only one peak which is at about 200 and there is no such point in temp. Further, if I look at the Check values, still these are quite different from temp. How can we get original vector approximately if not exactly? Moreover how will we decide range in x i.e.,

x = [randn(200,1); temp(:); randn(300,1)];

I mean how will we choose randn(200,1) and randn(300,1) for our temp?

[Star Strider]

My pleasure!

Yopu are filtering a signal padded on both ends by normally-distributed random numbers. It is unlikely that you will be able to recover your exact signal from it. (I changed the filter call to experiment with a normalised version. Change it back if you want to.)

% clear all
% clc
% A template is given
temp = [20 30 40 50]';

% Create a matched filter based on the template
b = flipud(temp(:));%My comment:Just reverse the vector temp up side down

% For testing the matched filter, create a random signal which
% contains a match for the template at some time index
x = [randn(200,1); temp(:); randn(300,1)];
n = 1:length(x);

% Process the signal with the matched filter
y = filter(b,sum(b),x);                             % Changed

% Set a detection threshold (exmaple used is 90% of template)
thresh = 0.9

thresh = 0.9000

% Compute normalizing factor
u = temp.'*temp;

u = 5400

% Find matches
matches = n(y>thresh*u);

% Plot the results
plot(n,y,'b', n(matches), y(matches), 'ro')
grid

% Print the results to the console
display(matches);

matches =

  1×0 empty double row vector

Check = y(y>mean(y(1:100))+10*std(y(1:100)))        % Recover Peak Data

Check = 7×1

7.4727 16.5828 27.1003 38.5714 27.9014 16.9928 7.4240

.

[Sadiq Akbar]

Thanks a lot for your guidance dear Star Strider. But I didn't get my answer yet. I asked the following adding more too:

1-how will we choose randn(200,1) and randn(300,1) for our temp. I mean why 200 and 300 why not some other?

2-How can we get our original signal approximatey

3-My original vector is 4x1 but Estimated vector Check is a 7x1 vector.

4-Vector Check is quite different from my original vector. Even its not approximate version of that.

5-What does this command do? y = filter(b,sum(b),x);

6-Why do we flip our signal i.e., b = flipud(temp(:));

[Star Strider]

The code I quoted is not mine. I just referred you to it, so it is not obvious to me what choices the original author made, or why.

‘5-What does this command do? y = filter(b,sum(b),x);’

That normalises the filter, so that ideally it does not amplify the signal being filtered.

[Sadiq Akbar]

Thanks a lot dear Star Strider for your kind responses. Yes, you are right that you referred me to that code. But still my issue is not resolved. How can I resolve my issue i.e., how to recove my vector approximately from noisy data (AWGN added)?

[Star Strider]

My pleasure!

It is likely not possible to completely recover a signal that has had Gaussian white noise added to it.

Experiment with the sgolayfilt function and wdenoise (wavelet denoising) function. Those are likely the best options.

[Sadiq Akbar]

Thanks a lot for your kind response dear Star Strider. Can you demonstrate it here as I am not too much expert in Matlab?