5

My objective is to show the results of:

  • How do I go about intentionally aliasing a signal?
  • undersampling
  • Nyquist rate sampling and oversampling

I'm first going to get the MATLAB code working well, then I'll rewrite it for the MSP430.

First I need to downsample the .wav file to get an incomplete or impartial data stream that I can then reconstuct. Here the flow chart:
analog -> sampling analog filter -> ADC -> resample down -> resample up -> DAC -> reconstruction analog filter

This is my first signal processing project using MATLAB (v2010), so be nice. Original code below; updated code here.

    %Play decimated file ( soundsc(y,fs) ) 
%Play Original file ( soundsc(play,fs ) )
%Play reconstucted File ( soundsc(final,fs) )

[piano,fs]=wavread('piano.wav'); % loads piano
play=piano(:,1); % Renames the file as "play"

t = linspace(0,time,length(play));          % Time vector
x = play;
y = decimate(x,25);

stem(x(1:30)), axis([0 30 -2 2])   % Original signal
title('Original Signal')
figure
stem(y(1:30))                        % Decimated signal
title('Decimated Signal')

%changes the sampling rate

fs1 = fs/2;
fs2 = fs/3;
fs3 = fs/4;
fs4 = fs*2;
fs5 = fs*3;
fs6 = fs*4;

wavwrite(y,fs/25,'PianoDecimation');


%------------------------------------------------------------------

%Downsampled version of piano is now upsampled to the original
[PianoDecimation,fs]=wavread('PianoDecimation.wav'); % loads piano
play2=PianoDecimation(:,1); % Renames the file as "play

%upsampling
UpSampleRatio = 2;  % 2*fs = nyquist rate sampling
play2Up=zeros(length(PianoDecimation)*UpSampleRatio, 1);
play2Up(1:UpSampleRatio:end) = play2; % fill in every N'th sample

%low pass filter

ResampFilt = firpm(44, [0 0.39625 0.60938 1], [1 1 0 0]);


fsUp = (fs*UpSampleRatio)*1;
wavwrite(play2Up,fsUp,'PianoUpsampled');

%Plot2
%data vs time plot
time=(1/44100)*length(play2);
t=linspace(0,time,length(play2));
stem(t,play2)
title('Upsampled graph of piano')
xlabel('time(sec)');
ylabel('relative signal strength')



[PianoUpsampled,fs]=wavread('PianoUpsampled.wav'); % loads piano
final=PianoUpsampled(:,1); % Renames the file as "play"


%-------------------------------------------------------------
%resampleing
[piano,fs]=wavread('piano.wav'); % loads piano
x=piano(:,1); % Renames the file as "play"
m = resample(x,3,2);

**%this is were i would need to sample it at the different frequecys (both above and below and at) nyquist frequency.*I think.***

soundsc(left,fs) % shows the resaultant audio file , which is the same as original ( only at or above nyquist frequency however)

How can I improve this? How do I sample at different frequencies?

Note that I've crossposted this on SO, here.

Community
  • 1
Martin
  • 109
  • 1
  • 7
  • 2
    I don't think this question fits E&R. You'd be better off at Stack Overflow. – tyblu Dec 28 '10 at 00:31
  • i posted it in all the forums hoping to get a wide spread of answers, as everyones help is wanted. Its funny i posted in stack overflow first and yet to get an answer were as both the other sites i have people have atleast looked at it. – Martin Dec 28 '10 at 00:36
  • 2
    This is a pure EE question. Better to present it here than to the mathematically illiterate programmers at SO. – markrages Dec 28 '10 at 01:31
  • Im getting yelled at were ever i post this. mark is this what you were saying earlier about doing it in matlab first then reconstructing it for msp430? – Martin Dec 28 '10 at 02:58
  • 3
    @Martin Cross posting is not viewed nicely among the stack exchange communities. It especially seems fishy when you create different account (and names) for each site rather then using the same name and linking all sites together. – Kellenjb Dec 28 '10 at 05:28
  • @markrages E&R is not a EE site, nor is SO a computer science site. It may so happen that EEs tend more toward E&R, but that doesn't necessarily mean that all topics that EEs go over are appropriate for this site. However, with that said, I personally am fine with DSP questions being here. – Kellenjb Dec 28 '10 at 05:31
  • 2
    @Martin "using nyquist" is not really a proper use of the word. It almost sounds like nyquist could be a function that you call. Nyquist is the name of a guy when used by itself. It is more commonly said as Nyquist Frequency, or Nyquist–Shannon sampling theorem, or folding frequency. Either way, Nyquist isn't a process. Just some tips on how to word it better. – Kellenjb Dec 28 '10 at 05:42
  • Here @Martin: this may be of use to get better answers. – tyblu Dec 28 '10 at 13:03
  • Discussion about relevance goes here: http://meta.electronics.stackexchange.com/questions/126/ – endolith Dec 29 '10 at 16:54
  • I am getting the error in your updated code " time is undeclared" –  Nov 28 '12 at 12:44

2 Answers2

3

If you want to downsample in a way that produces aliasing, just throw away some fraction of the samples, ie, only keep every 2nd sample, etc.

If you want to downsample in a way that doesn't alias, precede the decimation with a lowpass filter having a cutoff below the nyquist frequency of the sample rate you'll have after you decimate. Since decimating involves throwing away samples and an FIR lowpass filter has no future dependence on its results, you don't have to bother doing a convolution dot product for the outputs you'd only discard in the decimation, instead you can just skip to the next one.

To upsample, insert zeroes between the input samples (use zeroes, do not repeate the input values). This produces alias frequencies, so follow with a lowpass filter below the nyquist rate of the original sample rate. Again, you can design the convolution computation to skip over the zero inputs.

For small integer ratios of sample rates, discarding or zero inserting samples and FIR filters work very well. For larger ratios, CIC filters may be preferred as they can be efficiently computed in hardware using wrapping binary computation (you intentionally let the calculations overflow!), though they have a sinc function-shaped rather than rectangular frequency response and so need a clean up filter to select only the rectangular area near the tip of the sinc response and/or compensate for the droop of the sinc.

Chris Stratton
  • 33,491
  • 3
  • 44
  • 90
  • It might be helpful if you expand on why you can't just repeat the input values when he up-samples. – Kellenjb Dec 28 '10 at 05:35
  • Repeating samples is comparable to doing the right thing - alternating input samples and zeroes, and then feeding that through an FIR kernel consisting of a few ones. That is to say, convolving it with a rectangular pulse in time, which is a sinc function in the frequency domain. This gives the frequency response a sinc-function shape rather than flat passband up to the desired response frequency. More intuitively, the repetition boosts the low frequencies but starts to destructively interfere with the higher ones. Besides, it's cheaper to multiply by zero... – Chris Stratton Dec 28 '10 at 05:57
  • standard practice is just to edit the improved information into your answer. It is great to edit your answer, make it better, add section headings and such. – Kortuk Dec 28 '10 at 18:17
  • @Kortuk, I intentionally left it out of the original answer, because I think it's a divergence that does not help with the base explanation. I'm happy to have it as supplementary material though. – Chris Stratton Dec 28 '10 at 18:23
  • Awesome, just trying to spread knowledge. – Kortuk Dec 28 '10 at 18:48
3

I am going to start simple, and build up.

F= Frequency

F(Hz=1/s) E.x. 100Hz = 1000 (Cyc/sec) F(s)= 1/(2f)


Example problem: 1000 hz = Highest frequency 1/2(1000hz) = 1/2000 = 
5x10(-3) sec/cyc or a sampling rate of 5ms

If your highest frequency is 1000Hz and you sample with a period of 5mS(500Hz) your signal has already lost a large amount of its information. You are sampling at 1/4th of the nyquist rate. For a 1000Hz wave, you must sample at at least 2000Hz, and that is in a perfect world. You really need to go a bit faster than this due to technology limitations, but this is beyond what you need to know.

What is my nyquist frequency?

Wait, do you mean Nyquist Frequency or Nyquist rate?

This is probably seeming a pretty trivial distinction, but many people are referencing different things when they switch terms, some people can control the "Symbol Rate" some can control the "Sampling Rate," some can control both(lucky engineers).

You need to sample at a frequency that is greater than double your highest frequency.

If you go lower than this your signal will alias, this is irreversible once it has happened, the best you can do is filter off the corrupted part of your signal, but this does not recover information, just removes the corrupted frequencies.

Downsampling in Matlab

Lets say that you have your signal, you have avoided aliasing, and you want to downsample to half the points. This is very easy in matlab. Take a look at the decimate function. This will handle the mess of making sure your down-sampling does not introduce aliasing, at a cost. When it filters off the higher frequencies, they are lost. This is what you ask for when you half your datapoints. You can just use the downsample function. If you had the higher frequencies, it will alias. If you know there should not be those frequencies there, you can always use decimate to make sure that noise does not sneak into your system.

I would like to go into more detail, but I had to cut myself short as I am out of time for today. Leave me a comment if there is more that would help you, and I will come back and add more information. I am sure I did not cover everything you wanted.

Kortuk
  • 13,412
  • 8
  • 62
  • 85
  • I have been talking with someone in stack flow that has a similar point to yours. So first i would decimate the signal using matlab, then reconstruct the .wav file then take the new .wav file and sample it at nyquist rate? – Martin Dec 28 '10 at 19:00
  • 1
    Martin, it is still unclear what your are trying to do. Are you trying to determine the highest frequency component? Look at the FFT of the original file. Are you trying to subsample below Nyquist? Are you trying to do everything "by the book", or are you trying to undersample on purpose to see aliasing? – markrages Dec 28 '10 at 21:20
  • i want to sample at nyquist frequency to test the principle of it, and yes i want to undersample it on purpose. if you want heres my code as of now. http://pastebin.com/D0QzMkzy . and my audio clip http://www.4shared.com/audio/11xvNmkd/piano.html – Martin Dec 28 '10 at 21:30
  • @Martin, this was not clear in your question, I would have said, "How do I go about intentionally aliasing a signal?" I will add more later. – Kortuk Dec 29 '10 at 02:51
  • Ok ill go update it. But its kinda of both, just showing the effects of undersampling which inturn will intentionally aliasing the signal. – Martin Dec 29 '10 at 03:11
  • @Martin, down-sampling and under-sampling are very different terms, they often happen together, but it is not implied that one happened when the other did. – Kortuk Dec 29 '10 at 16:11