Programmable Devices
CPLDs, FPGAs, SoC FPGAs, Configuration, and Transceivers
20725 Discussions

Division Approximation In Fm Demodulator

Altera_Forum
Honored Contributor II
2,677 Views

i am a post-graduate student at the university of ioannina greece "physics department". my thesis subject has to do about an implement of sdr receiver and demodulation of fm signal. i already have implement the ddc part in an fpga with my i and q signals in two complement's integer format (the same format used by the adc) but i have problem with the demodulation and especially with division.i use a digital phase discriminator algorithm witch produce an output following the next equation (i*dq – qdi)/(i^2 + q^2).i can compute the numerator and the denominator but the problem is with division. i can use a divisor but the result from this division is just a quotient and a remainder, this wont give the desired output in a fixed scale capable to be processed properly by a dac.is there any possible way to implement a divisor with a fixed range quotient in regard with the remainder.in other words is it possible to produce an approximated quotient from the division of two complement format signals. with only the knowledge of the first quotient and the remainder, and if that is possible what is the form of the output?any answer guidance proposed bibliography or links will be very helpful. (http://www.fpga4fun.com/forum/viewtopic.php?t=1669&start=0&postdays=0&postorder=asc&highlight=)

0 Kudos
23 Replies
Altera_Forum
Honored Contributor II
1,443 Views

Hi, 

 

If you mean how to round up to nearest integer value then: 

With 2's complement of any computation result we can do rounding by looking at the remainder's MSB. If this MSB bit is 1 then add 1 to quotient else leave the quotient unchanged. This applies to both positive and negative numbers. 

Be careful about overflowing at this addition of 1 

 

If rounding is not enough then you have to increase the bit resolution. 

 

edit: 

I am assuming your remainder is not signed which is the case in any truncation of signed values so check the format of your remainder.
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

My remainder is a singed two complement number.I use the lpm_divide megafunction and a cyclone III fpga.Do you have an idea how i round up with a remainder in that form. 

 

edit: 

To give a hint a based much more on this....http://rfdesign.com/mag/604rfdf4.pdf but i don't really understand what they mean by the term automatic scaling function if you have an idea....
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Hi, 

 

The easiest way is to have the quotient with more bits(then ignore the remainder), discard the excess bits(this is now your remainder) and check its MSB as above post.
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

The automatic scaling - as explained in your document - is meant to keep the data values at full representation of 16 bits. i.e. there is variable gain control (similar to AGC)depending on input values. Obviously the document doesn't tell how they have implemented that but there might be a comparison and a mult/div stage inside. 

 

edit: 

the diagram obviously indicates *2^n and divide by 2^m by the automatic scaling unit. . the n/m values can be arrived at by checking the values of every bit and its index for num/denum then it is possible to target a result that equals +2^15-1 by shifting num/denum right or left. 

 

It seems you have 16 samples per symbol at divider output(since there is a decimation factor of 16 afterwards). If they are looking at values within one symbol window of time to adjust the representation to 16 bits then it gets more difficult to scale. , if not then all you need is to scale every positive value up and every negative value down. 

This scaling is really vague to me, so is the modulation, is it FM or FSK?
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Another demod option I did 20 years ago is to take the Q/I as sin/cos and compute an arctan(sin/cos) to determine phase. I'd use a cordic for the arctan. You only need the quotient for this value. This gives you a PM demod by the way. 

The next step is to take the difference between phase samples x[n] - x[n-1] which gives you a dθ/dt which is proportional to frequency offset. The frequency value is then equal to ± half the sample rate. For instance, with a 100MHz sample rate your demod range would cover ±50MHz. If this is wider than necessary you can insert pipeline stages between the subtraction terms such then your FM term is now x[n] - x[n-1-m] and now your demod range is ±50MHz/m.
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Hi all, 

 

I have written the following Matlab code to model FM modulation/demodualtion (without upconversion) according to original method. I am sure the second method in the last post should work and be compared. 

The code should throw some light on the issue of scaling. It is clear that symbols are recovered but their range needs scaling. Hope it is going to help. 

 

******************************* 

clear all; 

%generate random zeros and ones 

symb = round(rand(1,50)); 

 

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

%convert symbols to two frequencies(f1,f2) 

for i = 1:50 

if symb(i) == 0 incr(i) = 10; else incr(i) = 20; end 

end 

 

%one cycle sine data 

lut = exp(2*pi*j*[0:255]*1/256); 

 

%modulate using f1,f2 generated with lut pointer from modulo adder 

k = 1; 

n = 0; 

for t = 1:length(incr)  

for i = 1:64 %upsample symbols by 64  

n = n + 1; 

x(n) = lut(k); 

k = k + incr(t);  

if k>256 k=k-256; end;%modulo 256 addition 

end 

end 

 

%check FM signal spectrum 

f = linspace(0,1,length(x)); 

plot(f,20*log10(abs(fft(x)))); 

 

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

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

%demodulate 

 

I = real(x); I_diff = [diff(I) 0]; 

Q = imag(x); Q_diff = [diff(Q) 0]; 

m = (I.*Q_diff - Q.*I_diff)./(I.^2 + Q.^2) ; 

 

%m = decimate(m,64); ?? 

m = m(1:64:end); %decimate by 64 

 

%compare Rx symbols with Tx symbols  

str = sprintf('%d ',symb); 

figure; plot(m) 

title(str) 

***********************************
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Thanks for the answers!!!!Another interesting approach that eliminates the fact that I and Q are similar in amplitude is this http://ieeexplore.ieee.org/ielx5/8452/26620/01186832.pdf?arnumber=1186832 witch gives in a more detailed representation the need and the implementation of scaling.

0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Hi, 

 

Thanks for the link, however I can't get through this link. May be you can email it to me at: 

kadhiem_ayob@yahoo.co.uk 

 

I have added the upconversion/downconversion as well in the code below. 

The recovered I/Q fluctuate between .2 ~ .4 even though the inputs are unitary, they don't go minus/plus as I thought first. 

 

 

///////////////////////////////////////// 

 

clear all; 

%******************* Tx *********************** 

%generate random zeros and ones 

symb = round(rand(1,50)); 

 

%convert symbols to two frequencies(f1,f2) 

%f1 = 10/256,f2 = 20/256 

for i = 1:50 

if symb(i)==0 incr(i)=10; else incr(i)=20; end 

end 

 

%one cycle sine/cosine data 

lut = exp(2*pi*j*[0:255]*1/256); 

 

%modulate using f1,f2  

%f1,f2 generated with lut pointer from modulo adder 

k = 1; 

n = 0; 

for t = 1:length(incr)  

for i = 1:64 %upsample symbols by 64  

n = n + 1; 

x(n) = lut(k); 

k = k + incr(t); 

if k>256 k=k-256; end; %modulo 256 addition 

end 

end 

 

%check FM signal spectrum 

figure;hold 

f = linspace(0,1,length(x)); 

plot(f,20*log10(abs(fft(x)))); 

 

%upconvert to carrier = 100/256 

carrier = exp(j*2*pi*[0:length(x)-1]*100/256); 

x_up = real(carrier .* x); 

 

%check FM signal spectrum 

f = linspace(0,1,length(x_up)); 

plot(f,20*log10(abs(fft(x_up))),'r'); 

 

%******************* Rx *********************** 

%downconvert 

xr = x_up .* real(carrier); 

xi = x_up .* -imag(carrier); 

 

%filter off the high sideband 

h = fir1(40,2*20/256+.05); 

xr = filtfilt(h,1,xr); 

xi = filtfilt(h,1,xi); 

x = complex(xr,xi); 

 

%check FM signal spectrum 

f = linspace(0,1,length(x)); 

plot(f,20*log10(abs(fft(x))),'g--'); 

legend('baseband','upconverted','downconverted'); 

 

%demodulate 

I = real(x); I_diff = [diff(I) 0]; 

Q = imag(x); Q_diff = [diff(Q) 0]; 

m = (I.*Q_diff - Q.*I_diff)./(I.^2 + Q.^2) ; 

 

% decimate by 64 

m = m(32:64:end);  

 

%compare Rx symbols with Tx symbols  

figure; plot(m) 

str = sprintf('%d',symb);legend(str) 

 

////////////////////////////////////////////// 

 

Thanks
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Hi all, 

 

Seems we are one step left. 

I did some playing on my Matlab code(posted earlier) and arrived at this conclusion: 

 

If Reciever I & Q are unitary(swinging +1) then the recovered symbols swing in the range .2 ~ .4 

 

If Reciever I & Q are scaled as in hardware(swinging +32767) then the recovered symbols swing in the range .15 ~ .56(sort of). 

 

The conclusion is this: The information is in the remainder of division only and not in the quotient, the quotient is always zero. Thus this post should be about "How to recover information from remainder?" and not about rounding per se. 

 

I have two suggestions, I hope somebody will help by any contribution: 

 

1) either read the remainder only and scale it to 16 bits. Ignore the quotient except watching for zero-crossing(due to noise ...etc) then saturate the remainder. 

 

or 

2) scale the numerator by shifting n bits towards MSB inserting zeros. Then take the quotient and ignore the remainder and the rounding(adding 1 is trivial for 16 bit full representation). 

 

There is still the issue of "automatic scaling" as mentioned in the first document. Do they really mean it "automatic and variable" or do they just mean "scaling"
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

It's true i think that the amplitude of I^2+Q^2 is always a constant due to the thing that I and Q have a difference in phase of π, so like the common expression cos(f)^2+sin(f)^2=1 the amplitude of I^2+Q^2 must be a constant.So we don"t need to worry about the division part of the algorithm m = (I.dQ - Q.dI)/(I^2 + Q^2)we simply can compute m = (I.dQ - Q.dI).

0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

True, an ideal normalised sin^2+cos^2 = 1. 

In your case you recover I,Q with values scaled and far from ideal so I wonder why the FM document people have added this term to the denum. 

 

If you don't divide you are likely to get very large value yet the symbol will only swing a fraction near the top so you will still need to divide and you will need to read the remainder of division to extract the symbols then decode them.
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

In the general case, I and Q aren't normalized and cover a certain amplitude range. Some kind of normalization, e. g. dividing by I^2+Q^2, is required, if the decoder output has to be correctly scaled. Cause the carrier magnitude can be expected to vary only slowly, other means as an AGC can be used as well. This is the case with usual FM receivers.

0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Thanks for the replays.You are probably right..FvM do you have any idea on how i can implement the AGC you mentioned.there is a same idea its called Automatic scaling function in a work according to this link http://rfdesign.com/mag/604rfdf4.pdf if you have time take a look!

0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

 

--- Quote Start ---  

Thanks for the replays.You are probably right..FvM do you have any idea on how i can implement the AGC you mentioned.there is a same idea its called Automatic scaling function in a work according to this link http://rfdesign.com/mag/604rfdf4.pdf if you have time take a look! 

--- Quote End ---  

 

 

Hi Xristos, 

 

I actually indicated in an above post a bit about automatic scaling mentioned in your document. It multiplies the num by 2^n and the denum by 2^-m without defining n,m. 

 

First you don't need AGC as this is meant for amplitude scaling against a variable input signal. In your case we are looking for phase measurement. All you need is scale and get the division right.(imagine the input signal as just different frequencies with same amplitude, though some advice clipping this swing) 

 

I believe yur document people mean the following for n,m of automatic scaling: 

(remember all values are positive) look for the location of MSB that equals '1' in the 32 bit wide num and denum. then shift num value until you get this '1' in bit(31). Similarly shift the denum towrds LSB until this '1' is at bit(16). 

Now divide and ignore the remainder. 

You should get a signal that swings across 16 bits. 

Personally I don't see any relevance of using the denum. After all we are forcing it to 1 on the 16 bits though its fraction is left over on the other 16 bits infront, either way it works(ignoring denum or not).
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

I was referring to the above discussion about omitting the division, if the signal magnitude is constant. Generally it isn't constant. Automatic scaling still leaves an 1:2 amplitude variation. So if you require the demodulated signal to be scaled correctly, you either should use the division or somethink like ADC. I don't want to suggest a particular implementation, cause a suitable design depends on the overall system structure.

0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Thanks kaz!!!But my I and Q singals have two complement form.Am i able to follow the procedure you mention above, or i have to change them into an other form.Something else by the term Msb you mean the last 16 bits of the number?So i scale the first one of the last 16 bits of the numerator until it goes from it's place to bit 31?All the arithmitic in my system is in two complement form...The singal being digitized by a an adc and has a two complement form.Then i multiply it with a sin and a cosine which have also two complement form.The coefficient of Fir filter are also in two complement form so I and Q have two complement form.In all these multiplications i truncate the result to Msb. 

 

....by the way in similar discussion in another forum someone claim that i don't need the last division by I^2+Q^2 because as they said i measure angles not amplitude. 

 

In fact, in a typical analog FM receiver, amplitude noise can mess up 

reception of the signal so they clip the amplitude to remove the 

noise.
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Hi Xristos, 

 

I am aware that you are doing 2's complement. But since the Matlab model shows positive values only then I think you can safely consider all values to be so (unless proved otherwise). You may therefore either assume all 16 bits as magnitude bits(even though they are meant to cover negative values. Then you will need to offset values down if you want + representation. 

If in doubt you may as well use 15 bits and move the '1' to bit (30) rather than bit(31) for the num. For the denum move the bit to bit(16). In this case you will cover 15 bits only(even after any offset). For offset you can subtract each value by a the mean(you don't need to compute the mean in real time but you can precompute it)
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Hi kaz! 

 

I believe too and i wish that recovered I and Q are positive numbers!!!Another option that sounds interesting is to build a cordic algorithm so i can compute the P[n]=atan(Q[n]/I[n]) and then take the message by derivation m[n]=P[n]-P[n-1]!!!
0 Kudos
Altera_Forum
Honored Contributor II
1,443 Views

Hi Xristos, 

 

the atan method leads to the same phenomenon of small swing just over zero. 

So I don't see any gain using that method. Moreover the cordic method is more difficult to implement compared to other methods. 

Attached a fully functioning Matlab code that uses 3 methods and compares the swing of signal.
0 Kudos
Altera_Forum
Honored Contributor II
1,337 Views

thanks Kaz!!! 

 

I agree with you cordic is not such a good idea!!!Another question... 

In cyclone II data sheet they that the embedded multipliers can achieve performance up to 250 Mhz. That means that i can multiply for example samples from an ADC up to 250 Msamples/s with each sample to be 12 bit in width without any problem? 

http://www.alteraforum.com/proweb/misc/progress.gif (http://www.alteraforum.com/forum/editpost.php?do=editpost&p=19329)
0 Kudos
Reply