BPSK Demodulation

[Altera_Forum]

Hello again, I am still struggling with BPSK demodulation :) 

 OK, now I am trying to implement Costas Loop for carrier synchronization, which involves Inphase and Quadrature Phase paths. Inphase is looking good, and works as expected, however, the Q component looks awkward. Take a look at these waveforms:  http://www.alteraforum.com/forum/attachment.php?attachmentid=10363&stc=1   I understand that multiplying sin * cos will produce a high frequency component, which is a sin modulated by half amplitude of the data signal. BPSK(t) * sin(2*pi*f*t) = 0.5 * [DATA(t) * sin(0) + DATA(t) * sin(2*pi*f*t)] => 0.5 * DATA(t) * sin(2*pi*f*t) => LPF{0.5 * DATA(t) * sin(2*pi*f*t)} = 0  I feel that problem is in my LPF. I used the same Matched RRC filter that I used in the in-phase arm. Papers I read didn't mention that these two filters should be different, that's why I used the same filter. (I think those are unpractical DSP writers, KAZ :D )  So, should I use a separate Loop Filter for Q component? How do I determine it's parameters?

[Altera_Forum]

I found out the reason behind that error... NCO_model produces the same sample value for a single sample of phase input, I need to pass a vector instead. 

Code now is as follows:   load LPF.mat N = length(rf_signal); phi_inc = zeros(1, N); bb_sig_I = zeros(1, N); bb_sig_Q = zeros(1, N); I_f = zeros(1, N); Q_f = zeros(1, N); error_sig = zeros(1, N); err = zeros(1, N); % NCO fs = 1.0E8; f = 5.0E6; phi_inc = ((2^32)*f/fs).*ones(1, N); for i = 1:N = NCO_model(phi_inc); % Downconverting to Baseband bb_sig_I(i) = rf_signal(i)*cos_out(i); bb_sig_Q(i) = rf_signal(i)*sin_out(i); % Filtering I_f(i) = filter(LPF, bb_sig_I(i)); I_f(i) = filter(LPF, I_f(i)); Q_f(i) = filter(LPF, bb_sig_Q(i)); Q_f(i) = filter(LPF, Q_f(i)); % Error error_sig(i) = I_f(i).*Q_f(i); % Loop Filter (Leaky Integrator) alpha = 0.1; err(i) = filter(alpha, , error_sig(i)); phi_inc = phi_inc - err(i); end   phi_inc is a vector now, passed to NCO_model on each iteration, and produces new vectors for sin & cos. I only multiply with the index i of both sin_out & cos_out. However, as you stated... filtering is not correct.

[Altera_Forum]

I solved it Kazem! 

Here is the code:   load LPF.mat N = length(rf_signal); cos_out = zeros(1, N); sin_out = zeros(1, N); phi_inc = zeros(1, N); bb_sig_I = zeros(1, N); bb_sig_Q = zeros(1, N); I_f = zeros(1, N); Q_f = zeros(1, N); error_sig = zeros(1, N); alpha = 0.1; t = 1:N; fs = 1.0E8; fc = 5.0E6; for i = 1:N cos_out(i) = cos(2*pi*fc*t(i)/fs + phi_inc(i)); sin_out(i) = sin(2*pi*fc*t(i)/fs + phi_inc(i)); % Downconverting to Baseband bb_sig_I(i) = rf_signal(i)*cos_out(i); bb_sig_Q(i) = rf_signal(i)*sin_out(i); % Filtering tmp = filter(LPF, filter(LPF, bb_sig_I(1:i))); I_f(i) = tmp(i); tmp = filter(LPF, filter(LPF, bb_sig_Q(1:i))); Q_f(i) = tmp(i); % Error error_sig(i) = (5*10^-5)*pi*sign(I_f(i)*Q_f(i)); phi_inc(i+1) = phi_inc(i) - error_sig(i); end   I replaced my NCO with a cos & sin functions to speed up the execution. Filtering is done as you suggested. Error signal is calculated in a slightly different way I found online, but this doesn't require IIR or Loop Filter.  https://www.alteraforum.com/forum/attachment.php?attachmentid=10382   https://www.alteraforum.com/forum/attachment.php?attachmentid=10383   Well, I think this is the desired result... but I don't understand error signal calculation. Can you please explain it a little bit and how different is it from the previous one? And I believe this is only theoretical, and can't be implemented in FPGA. Right? How should I modify it to simulate hardware?

[Altera_Forum]

Hi Siraj, 

 I don't want to discourage you bu I think scaling your error till it is dead is no good.  It looks like qpsk/QAM carrier tracking that I did years ago was based on different error circuitry involving slicing as well so we better ignore its details.  I thought of your BPSK design and is practically new to me but it looks like it doesn't specifically apply to BPSK. It looks like a basic form of PLL.  You might have to do matlab modelling to prove the concept in a simpler way away from BPSK as follows:  1)make your input a clean tone say at 0.1 Fs 2) apply NCO at -0.1 and do the rest as usual at the error detector I and Q should be dc (if all ok) and error should be zero (or settle to zero)  repeat above with NCO initial freq set to -0.11 and repeat since we don't expect I and Q to be dc but rather at a negative offset we should have the loop push them back to dc  repeat above with nco at -0.09 I and Q would be at positive offset and the loop should push them back to dc. the error should change orientation.  Thus the whole issue is finding how to decide if a complex frequency(I/Q) is positive,0 or negative. The actual values of error is not that important but its sense on the NCO is  You need to be aware of some other issues in modelling your loop but can do that later such as effect of filter delay (group delay) and actual design latency between error and nco update.  But for now the crucial point is finding a way to have good error detection (don't worry about its scaling)  I believe you need a lot more work to get this loop right and hopefully tell me.

[Altera_Forum]

how do I apply negative frequency to my nco? do you mean in regard to its sin & cos output?

[Altera_Forum]

--- Quote Start ---  how do I apply negative frequency to my nco? do you mean in regard to its sin & cos output?  --- Quote End ---    If you are using the nco model then just swap sin/cos (or use negative tuning word if supported)  Alternatively you can use the following line:  nco = exp(j*2*pi*i*freq);  inside your loop  set freq to any value +/- 0 ~ 0.5 and I suggest using -0.1, -0.11,-0.09 (not -1.1 that was typing error)  you rf input can be just: rf = cos(2*pi*(0:N-1)*.1); fixed for all three tests

[Altera_Forum]

By the way, you can apply all positive frequencies due to mirroring of RF i.e. RF can be +0.1 and nco +0.1, +0.11,+0.09

[Altera_Forum]

Alright... I just want to note that the second plots in the above posts are for the phase, NOT the error.

[Altera_Forum]

Test fc = -0.09 Error signal: 

https://www.alteraforum.com/forum/attachment.php?attachmentid=10391   I feel they are reversed from what you stated they should be. Should I invert the sign of imag(carrier) in code?

[Altera_Forum]

Looks like inverting the sign of imag(carrier) makes everything as you stated! 

First, does that mean the model is a correct Costas Loop? Second, why multiply with -imag(carrier) rather than +imag(carrier)?

[Altera_Forum]

There is no indication of lock.  

You have killed the error by /1e5 then (5*10^-5) , this adds almost nothing to nco tuning word whose range is 0 ~ +/-0.5  for .1 case error should settle zero and filtered signal should be dc (I think you got that because error is forced to nothing) for .09 and .11 case error should settle at .01 (either plus or minus) and filtered signal should be dc (not sine wave)

[Altera_Forum]

--- Quote Start ---  There is no indication of lock.  You have killed the error by /1e5 then (5*10^-5) , this adds almost nothing to nco tuning word whose range is 0 ~ +/-0.5  for .1 case error should settle zero and filtered signal should be dc (I think you got that because error is forced to nothing) for .09 and .11 case error should settle at .01 (either plus or minus) and filtered signal should be dc (not sine wave)  --- Quote End ---    correction:  can you set RF signal to 0.1, 0.11,0.09 and keep nco at 0.1  also why you have nested your filter function twice?

[Altera_Forum]

I rerun the code without killing the error and bringing it back to life (as Celine Dion says :P ) 

But nothing changed... exactly same result. So, it's not behaving as Costas Loop. What do you suggest? I surrender! (Also as Celine Dion says :D )

[Altera_Forum]

--- Quote Start ---  correction:  can you set RF signal to 0.1, 0.11,0.09 and keep nco at 0.1  --- Quote End ---    I will try now.    --- Quote Start ---   also why you have nested your filter function twice?  --- Quote End ---    Because you said to filter 200 taps x 2... remember? By the way I tried to filter once... I just got a high frequency component.

[Altera_Forum]

--- Quote Start ---  correction:  can you set RF signal to 0.1, 0.11,0.09 and keep nco at 0.1  --- Quote End ---    Nothing changed also!

[Altera_Forum]

--- Quote Start ---  Nothing changed also!  --- Quote End ---    Now what I said should apply but as follows:  phy_inc should settle to 0.1, 0.11 and 0.09 for the three cases. I & Q should be dc only in all three cases especially towards end of vector if lock occurs. error may oscillate initially.  Do not kill the error by that scaling, you need a target of .01 be added or subtracted from nco value of 0.1 initially. Once the nco is set to this correct value then the error I think will tend to drop towards zero.  200 x 2 does not mean nesting. I simply meant resource of two flters each 200 taps so it 200 x 2 (one per branch) use y = filter(LPF,1,x(1:i));  if it doesn't cut the high band then check response.

[Altera_Forum]

well I finally tried my muscles. There is some life in it but it requires a lot of tweaks. It lloks it is tending to lock but requires very long stime so I tried 6000 samples. 

  clear all; close all; clc; N = 6000; h = fir1(100,0.3); rf_signal = cos(2*pi*(0:N-1)*0.1); %.1,.11,.09 bb_sig = zeros(1, N); y_f = zeros(1, N); error_sig = zeros(1, N); error_f = zeros(1, N); %PLL fs = 1; f = .1; phi_inc = f/fs; for i = 1:N %Downconverting to Baseband bb_sig(i) = rf_signal(i)*exp(j*2*pi*i*phi_inc); %Filtering y_f = filter(h,1,bb_sig(1:i)); %Error error(i) = real(y_f(i))*imag(y_f(i))/2^16; phi_inc = phi_inc - error(i); temp(i) = phi_inc; end plot(temp,'.-'); figure; plot(error); figure; plot(real(y_f(end-1000:end)),'.-');

[Altera_Forum]

1- I ran your code for 12000 samples, output signal got distorted, why? 

2- You are changing the frequency rather than shifting phase. Is it what you intend to do? My code controls phase not frequency. 3- Why are you dividing by 2^16? 4- Phase, or more precisely, frequency is settling at 0.1, 0.11, 0.09 as you mentioned, but I don't see I&Q as DC... they are sinusoidal waves.

[Altera_Forum]

--- Quote Start ---  1- I ran your code for 12000 samples, output signal got distorted, why? 2- You are changing the frequency rather than shifting phase. Is it what you intend to do? My code controls phase not frequency.  --- Quote End ---    when you run nco you control the tuning word (phi_inc). That is what I am doing. phase and frequency are related. In design the nco will keep running nonstop from phase 1 and after some clocks (latency) its tuning word is updated.  In the code the index i is the phase and it is incremented steadily.  If you can run modelsim it gets much easier to check long sims and do any tweaks.  The important thing initial is to watch the phi_in value tending towards the rf centre frequency. once it settles the signal will settle at a constant (dc).  /2^16 is the tweek that turned out to control the loop better for our setup. you may require modifying it in final sign off of your radio !  output may look sinusoidal but check dynamic range (is it around zero or just osscilation around a constant level)

[Altera_Forum]

the filter needs to cut off the high sidband. 

rf signal has two frequencies at +0.1 and -0.1 when we mix it with 0.1 (or you can choose =-0.1) we get the sum and the difference: 0.1 + 0.1 = 0.2 (high sideband) 0.1 - 0.1 = 0 (our dream dc)  the filter must pass dc but cut highband at 0.2 (of Fs) the function fir1 that I used cuts off at 0.3 of Nyquist (i.e. 0.3/2 of Fs = 0.15 so is just right for single tone).  bwt this costas loop seems horrible to control compared to my previous work on qpsk/qam which was predictable. in this loop even if reverse the sense of error I get same result.  so surely this is not how it should be done. If you make research and find out I will be interested to know.  Should we use filter for the error (I tried IIR and IIR + a second accum and lost control).

[Altera_Forum]

Thank you so much for bearing with me Mr Kazem. 

I really got tired today, I'll see what I can do tomorrow and let you know.

[Altera_Forum]

Hello Mr Kazem. Back to work :) 

   --- Quote Start ---   the filter must pass dc but cut highband at 0.2 (of Fs) the function fir1 that I used cuts off at 0.3 of Nyquist (i.e. 0.3/2 of Fs = 0.15 so is just right for single tone).  --- Quote End ---    1- I checked the frequency response of this filter using freqz. But it shows the cutoff frequency at 0.3 not 0.15.  2- Regarding the control of NCO, I just have this thought and please explain. If we control frequency as you suggest, suppose that the RF carrier is 0.2 and I tuned my radio to 0.1, the loop will push until it reaches 0.2, right? Will the loop span all the range from 0.1 to 0.2 to lock? Suppose there is another RF signal with carrier 0.15, don't you think it might lock on it while it's looking for 0.2? I just have this thought nagging on me for a while, I hope you explain it.  3- You told me that variable 'i' in the 'for' loop represents phase, why not time?  4- I attached a simulation of Costas Loop I found on MATLAB FileExchange. I hope you check it.