Frequency detection

What I have:

A digital signal sampled at Fs composed of a QPSK modulation at Fmod plus some additive white gaussian noise.

What I want:

$$modPower/signalPower$$

The idea behind is that I want to detect if there is a communication at Fmod on the channel or if there's just noise.

What I do:

Over the observation window n I compute: $$signalPower = \sum\limits_{i=1}^n s(i)^2$$ $$modPowerCos= \sum\limits_{i=1}^n s(i)*cos(2\pi*Fmod*i/Fs)$$ $$modPowerSin= \sum\limits_{i=1}^n s(i)*sin(2\pi*Fmod*i/Fs)$$ $$modPower = modPowerCos^2 + modPowerSin^2$$ $$ratio = modPower / signalPower$$

My code:

% Ns = Number of samples per symbol
% Nsyn = Number of symbols of the window
% Nt = Length of the window
% Fmod = Modulation Frequency
% Fs = Sampling Frequency
Nt = Nsym*Ns;
signalPower = zeros(1, length(obj.signal));     % Signal Powe
modPowerCos = zeros(1, length(obj.signal));     % Fmod cos correlation
modPowerSin = zeros(1, length(obj.signal));     % Fmod sin correlation
for ii = obj.Ns:length(obj.signal)
for jj = 0:Nt-1
if ii-jj > 0 && ii-jj <= length(obj.IQrx)
signalPower(ii) = signalPower(ii) + ...
obj.signal(ii-jj)^2;
modPowerCos(ii) = modPowerCos(ii) + ...
obj.signal(ii-jj) * ...
cos(2*pi*Fmod*(ii-jj)/Fs);
modPowerSin(ii) = modPowerSin(ii) + ...
obj.signal(ii-jj) * ...
sin(2*pi*Fmod*(ii-jj)/Fs);
end
end
modPowerCos(ii) = abs(modPowerCos(ii));
modPowerSin(ii) = abs(modPowerSin(ii));
end


My problem:

I find a ratio higher than 1 and I think it comes from a wrong formula but I don't know which one.

Example:

In the picture bellow in blue is the computed signal power signalPower while in magenta is the computed Fmod power modPowerCos + modPowerSin. You can see that modPower is higher than signalPower, which bugs me.

• The main issue here seems to be program debugging and not signal processing, and as such, this is off-topic for dsp.SE. – Dilip Sarwate Jun 27 '13 at 19:43
• @DilipSarwate: I disagree. I see this specific question as algorithm debugging, which is definitely on-topic here. Sure, it's expressed as a script, but it's at the algorithmic level (i.e. it's not calling any library functions). – Peter K. Jul 1 '13 at 11:53

There's a problem with your method of computing powers. If you compute the signal power as

$$P_s=\sum_{i=1}^ns_i^2$$

then the quantities

$$\left(\sum_{i=1}^ns_i\cos(2\pi if_m/f_s)\right)^2\text{ and } \left(\sum_{i=1}^ns_i\sin(2\pi if_m/f_s)\right)^2$$

cannot be related to $P_s$ in any way that would be meaningful for your problem of detecting the presence of the QPSK signal (even though in your program you compute the magnitudes instead of the squares shown in your question).

What I think you were trying to do is demodulate the QPSK signal and measure the power of the demodulated signal. However, there are two problems in the way you did it. First, you must compute the power by taking the sum of squared values, and not by squaring the sum. And second, it is necessary to lowpass filter the signals after multiplication by sine and cosine, because the multiplication by sine and cosine does not only shift the spectrum down to DC but also up to twice the modulation frequency, and you need to get rid of that latter part.

What I've said so far is only the explanation of what I think you were trying to do and why it didn't work. A much simpler way of achieving the same would be to apply a bandpass filter with center frequency $f_m$ to the signal, measure the power of the output, and relate it to the total signal power $P_s$.