2
$\begingroup$

I currently have a .wav signal that was recorded at a 48 kHz sample rate, with a central frequency of 5.260 MHz and bandwidth of 4 kHz. I'm trying to get some information from this signal, specifically the location of the frequency peaks using matlab, which I would expect to be around 5.260 MHz. The code I'm using to generate the power spectral density is as follows:

path = '5.260MHz.wav';
f0 = 5.260e6; % frequency 5.260 MHz

% sample properties
[x, Fs] = wavread(path); % Fs = 48 kHz
alias = 48e3*floor(f0/48e3);      % starting frequency of target alias

%estimate spectrum
[psd, f] = pwelch(x, 512, [], [], Fs);

% If the frequency lies in one of the mirror frequency bands,
% we have to rotate about half the sample rate.
if (mod(f0, 48e3) > 24e3)  
    f = Fs - f;
end

plot(f + alias, psd)

This script produces the following plot

enter image description here

The bandwidth does seem to be 4 kHz as expected, but the frequency peaks appear to be in the wrong place. I would expect them to occur between 5.256 and 5.264 MHz, or maybe 5.258 and 5.262 MHz, not between 5.276 and 5.280 MHz.

However, when I generated my own files with the same central frequency and bandwidth, I got the following when I ran the above script:

% signal properties
f0 = 5.260e6;           % frequency 5.260 MHz
sf = 4e3;               % bandwidth 4 kHz

% sample properties
fs = 24e6;              % sampling frequency 24 MHz
N = 10e6;               % number of samples

% generate random frequency modulated sinusoidal signal
fi = smooth(randn(N, 1), 11);
fi = fi / std(fi) * sf + f0;
x = cos(2 * pi * cumsum(fi) / fs);

% downsample to 48 kHz (factor 500)
x = x(1 :500: end);
fs = fs / 500;

wavwrite(x, fs, 'test_5.260MHz.wav');

enter image description here

This plot is exactly as expected, using the same script as the first plot. I'm wondering if I should be treating the first input different in some way to have it graphed correctly.

So far I've tried testing some other generated cosine signals using the same code as above, and they all graph correctly. However, every data signal I've tried has been off to some degree (I can post different sample plots if that would be useful). It could be an error with the data collection instrument, but it's far more likely to be a bug in my code.

$\endgroup$
  • $\begingroup$ I'm new to the site, so if I can make the question clearer or add additional details let me know. $\endgroup$ – Jake Dec 20 '13 at 21:44
  • $\begingroup$ Some minor things: It's "kHz" and "MHz", the unit is "Hz" and the prefixes are "k" and "M". What you are estimating is not the "fft" but the power spectral density. The Fourier transform is only one of the elements needed to estimate the psd, and "fft" is just a fast implementation of the Fourier transform. $\endgroup$ – A. Donda Jan 6 '14 at 15:15
  • $\begingroup$ I think there's an error, it should be alias = 48e3*floor(f0/48e3); $\endgroup$ – A. Donda Jan 6 '14 at 15:16
  • $\begingroup$ @A.Donda Thanks, I'm still new to the world of signals. Also good find, that was definitely a bug I overlooked. I updated my code and question with all of the suggestions. $\endgroup$ – Jake Jan 9 '14 at 0:32
  • $\begingroup$ If you could post your data file, I could have a look at it. $\endgroup$ – A. Donda Jan 9 '14 at 20:27
1
$\begingroup$

You cannot capture a 5 MHz signal by sampling at 48 kHz. That doesn't mean there won't be anything in your sample, it means it will represent the signal incorrectly.

What happens if you use a too small sampling frequency is called aliasing: the original frequency shows up as another frequency altogether. The translation from the real to the sampled frequency is like a process of wrapping or folding the infinite frequency axis into your limited frequency interval $[0 ; \frac12 f_s]$, or $[-\frac12 f_s ; \frac12 f_s]$ for complex-valued signals. $\frac12 f_s$ is also called Nyquist frequency.

In your example, not just the main frequency but the total bandwidth of the signal of 800 kHz is much larger than the sampling frequency. Because of this, I wouldn't expect any clearly defined peak to be present. Since you seem to find peaks, my guess is that that is due to your very primitive approach to spectral estimation. If you use a better method, e.g. Welch's modified periodogram, which is implemented in Matlab as pwelch, I'd expect that you see a very flat (white) spectrum.


Addition based on nimrodm's comment: My initial statement has to be relativized if the bandwidth of the signal is smaller than $\frac12 f_s$. In that case, if there is additional knowledge about the true signal frequency, this knowledge can be used to resolve the ambiguity introduced by aliasing, and perfect reconstruction of the true signal is possible. See undersampling. However, it is not possible to determine the true frequency from the undersampled signal and nothing else.


A simulated example: Matlab code to generate a frequency-modulated signal with peak at 5 MHz and bandwidth 400 kHz, initially sampled at 24 MHz.

% signal properties
f0 = 5e6;               % frequency 5 MHz
sf = 0.4e6;             % bandwidth 400 kHz

% sample properties
fs = 24e6;              % sampling frequency 24 MHz
N = 10e6;               % number of samples

% generate random frequency modulated sinusoidal signal
fi = smooth(randn(N, 1), 11);
fi = fi / std(fi) * sf + f0;
x = cos(2 * pi * cumsum(fi) / fs);

Estimating the spectrum using pwelch

% estimate spectrum
[psd, f] = pwelch(x, 512, [], [], fs);
figure
plot(f, psd)
xlabel('f [Hz]')
ylabel('normalized power spectral density [1/Hz]')

gives this result:

If the signal is downsampled to 48 kHz

% downsample to 48 kHz (factor 500)
fs = fs / 500;
x = x(1 :500: end);

the estimated spectrum looks like this:

This is, as expected, basically a white (flat) spectrum, and the many small "peaks" are random fluctuations due to the fact that this spectrum was estimated based on only a 500th of the original 10 million samples. The estimate could be improved by repeating the procedure and averaging over the different results.

If the bandwidth is changed to 4 kHz

sf = 4e3;               % bandwidth 4 kHz

the results of spectal estimation on original and downsampled signal look like this:

The narrower peak fits nicely into the range 0 to 24 kHz of the downsampled spectral estimate, it just necessarily appears at the wrong place, in this case 8 kHz.

8 kHz = 5 MHz – 104 · 48 KHz is the smallest alias of 5 MHz with respect to a sampling rate of 48 kHz, the next largest is 48 kHz – 8 kHz = 40 kHz (mirroring around the Nyquist frequency of 24 kHz). All the larger aliases derive from these two by adding integer multiples of 48 kHz: 56, 88, 104, 136 kHz up to 4952, 4984, 5000, 5032 kHz and so on.

Since the spectral estimation cannot distinguish these aliases, the frequency axis of the spectrum plot becomes ambiguous: instead of representing frequency components of the signal from 0 to 24 kHz, it could as well represent frequencies from 48 to 24 kHz (backwards because of the mirroring!),

or from 48 to 72 kHz, and so on.

If we know the frequency content of the signal lies around 5 MHz, we can select the corresponding aliased frequency interval, 4992 to 5016 kHz,

and recover the correct location of the peak at 5 GHz.

This trick does not work, however, if the bandwidth of the signal is larger than $\frac12 f_s$, because the true frequency content of the signal has been aliased across several of these frequency intervals.

$\endgroup$
  • $\begingroup$ Correct but ... The signal bandwidth is 800 KHz (+/- 0.4 MHz) and not 5 MHz. Had the bandwidth been less than 24 KHz, sampling the 5 MHz carrier would not have created aliasing (this is sometimes called "IF sampling" as in "Intermediate Frequency" and its effectiveness depends on the imperfections of the ADC) $\endgroup$ – nimrodm Dec 21 '13 at 18:43
  • 2
    $\begingroup$ @nimrodm, thanks for your comment. I think it is not correct to say that in that case aliasing does not occur. However, you are right that the ambiguity introduced by aliasing can be resolved, using prior knowledge about the signal frequency. I have updated my answer to include information about "undersampling". $\endgroup$ – A. Donda Dec 21 '13 at 23:20
  • $\begingroup$ @A. Donna, Thanks for the response. I think I understand the concept of aliasing better now. So if I were to have the same signal of 5MHz, but a bandwidth of say 4Khz, I would be able to perfectly reconstruct the signal using a process called "IF sampling"? How would I have to alter the above matlab script above to do this and have the correct units on the x-axis? $\endgroup$ – Jake Jan 2 '14 at 22:11
  • $\begingroup$ @Jake, are you interested in reconstructing the true signal itself, or just in a correct estimate of its spectrum? $\endgroup$ – A. Donda Jan 3 '14 at 17:02
  • $\begingroup$ I extended the answer explaining what aliasing does and how to compensate for it in spectral estimation. $\endgroup$ – A. Donda Jan 3 '14 at 18:56
1
$\begingroup$

You can capture a higher frequency signal by a using a sample rate much lower than twice the signal/spectral content frequency. It's called undersampling. To be useful, one has to know approximately how many "folds" the signal frequency is above the sample rate, there has to be no baseband (or other folded spectrum) signal to alias against, thus the bandwidth of the signal has to be narrower than Fs/2, and the spectrum of the signal can't cross any multiple of Fs/2. The jitter of the sampling clock also has to be quite low (as low as the max absolute time jitter required for a sampler that is clocked at well above the signal frequency.)

This is no different from baseband sampling, as one needs a-priori knowledge (that the signal was, in fact, correctly low-pass-filtered baseband below Fs/2) to reconstruct the original signal correctly,

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.