DDC and Polyphase Channelizer

GNU Radio
1

Hi,:

I’m learning about polyphase channelizer, and I have been reading fred
harris’ paper “Digital Receivers and Transmitters Using Polyphase Filter
Banks for Wireless Communications”. I have also played with the Matlab
programs at the end of the paper. What I’m trying to do next is to use
GNURadio to do the channelization. I made two attempts, one uses
freq_xlating_fir_filter_ccf to implement a DDC channelizer, another uses
pfb_channelizer_ccf, running against the data in appendix III of harris’
paper. I’m hoping the channelizer result from GNURadio code would match
the result from Matlab, but I failed on both, I wonder if someone could
help me spot the problem.

  1. DDC Channelizer:
    1.1 Code: The following is the core part of the code (input parsing and
    various other stuff are omitted), the value for the variables can be
    found in 1.2
    self._filter_taps = gr.firdes.low_pass(1, self._input_rate,
    self._cutoff_freq, self._trans_width,

window=gr.firdes.WIN_BLACKMAN_hARRIS)
print "Number of taps: _filter_taps = ",
len(self._filter_taps)

     # Blocks
     self.ddc_filter = list()
     # i = 0 to M-1
     for i in xrange(self._M):
         offset = self._channel_bandwidth * (float(self._M)/2 - i)
         print "Channel ", i, "'s offset (Hz): ", offset

self.ddc_filter.append(gr.freq_xlating_fir_filter_ccf(self._resampling_ratio,

self._filter_taps, offset, self._input_rate))
self.connect(self.head, self.ddc_filter[i],
self.output_files[i])

1.2 Variable values:
Input sampling freq (Hz): _input_rate = 1120000000.0
Resampling ratio: _resampling_ratio = 28
Output sampling freq (Hz): _output_rate = 40000000.0
Transition width (Hz): _trans_width = 2000000.0
Number of channels: _M = 40
Number of samples: _N = 11200
Oversample rate: _oversample_rate = 1.42857142857
Channel bandwidth (HZ): _channel_bandwidth = 28000000.0
Cutoff freq (HZ): _cutoff_freq = 13000000.0
Number of taps: _filter_taps = 561
Channel 0 ‘s offset (Hz): 560000000.0
Channel 1 ‘s offset (Hz): 532000000.0
…
1.3 Result: The result from GNURadio code is saved to output files, one
for each channel. Then the results time-domain are plotted in Matlab
(-10 to +10, time 1 to 200), and compared to the result plot from
harris’ code, they are similar, but not the same. Channel 0 to 5, 17,
18, 34 to 39 do not have any strong signal, this is reflected in both
plots, channel 20’s signal is exactly the same, but all other channels’
plot are different.

  1. Polyphase Channelizer:
    2.1 Code:
    self._filter_taps = gr.firdes.low_pass(1, self._input_rate,
    self._cutoff_freq, self._trans_width,

window=gr.firdes.WIN_BLACKMAN_hARRIS)
# Calculate the number of taps per channel for our own
information
self._taps_per_channel =
scipy.ceil(float(len(self._filter_taps)) / float(self._M))
print "Number of taps: _filter_taps = ",
len(self._filter_taps)
print "Taps per channel: _taps_per_channel = ",
self._taps_per_channel

     # Construct the channelizer filter
     self.pfb = blks2.pfb_channelizer_ccf(self._M,

self._filter_taps, self._oversample_rate)

     # Connect the blocks
     self.connect(self.head, self.pfb)
     # i = 0 to M-1
     for i in xrange(self._M):
         self.connect((self.pfb, i), self.output_files[i])

2.2 Variable values:
Input sampling freq (Hz): _input_rate = 1120000000.0
Resampling ratio: _resampling_ratio = 28
Output sampling freq (Hz): _output_rate = 40000000.0
Transition width (Hz): _trans_width = 2000000.0
Number of channels: _M = 40
Number of samples: _N = 11200
Oversample rate: _oversample_rate = 1.42857142857
Channel bandwidth (HZ): _channel_bandwidth = 28000000.0
Cutoff freq (HZ): _cutoff_freq = 13000000.0
Number of taps: _filter_taps = 561
Taps per channel: _taps_per_channel = 15.0

2.3 Result: The GNURadio program failed with:
"terminate called after throwing an instance of ‘std::invalid_argument’
what(): gr_pfb_channelizer: oversample rate must be N/i for i in [1,
N]

This application has requested the Runtime to terminate it in an unusual
way.
Please contact the application’s support team for more information."

This doesn’t seem to make sense, since in this case N = 40, oversample
rate =1.42857142857, i = 28

Thanks

2

On Mon, Dec 27, 2010 at 3:48 AM, Jimmy R. [email protected]
wrote:

result from Matlab, but I failed on both, I wonder if someone could help me
spot the problem.

  1. DDC Channelizer:

Using the frequency xlating filters here is not recommended. You’ll
need to use much longer filters than when using the channelizer. That
said, I’m not sure what would cause the problem you identified below.
In my opinion, it’s not worth doing it this way.

self.ddc_filter = list()
1.2 Variable values:
Channel 0 's offset (Hz): 560000000.0
2.1 Code:

2.2 Variable values:
Taps per channel: _taps_per_channel = 15.0
=1.42857142857, i = 28

Thanks

The “oversample_rate” is not exactly what you are specifying here. The
default behavior of this PFB channelizer implementation is to be
critically sampled; that is, the sample rate equals the bandwidth. For
many applications, you want to have the output sample rate be some
other rate, most notably, 2x the bandwidth. The oversample_rate allows
you to specify a number real number, but this cannot be an arbitrary
number in this case. Because of the way the samples are fed to the
filters, the oversample rate, O, has to be set such that the number of
filters, N, can be evenly divided by an integer to to get O; that is,
N/i = O for i in [1,N]. Take the example in in fred’s paper where he
has 64 channels and a 1/48 resampling rate. So he’s using a 48-to-1
downsampling in 64 channels, so his rate is 4/3, which means that the
output sampling rate is (4/3)*bw, where bw is the channel bandwidth.

If you specified that the rate was 2.0, you would feed 32 of the 64
channels and have a sample rate that is twice the bandwidth. I notice
in the remarks in fred’s paper that he’s outputting the channels at 2
samples/symbol, so the oversample_rate here should be 2.

I hope that makes sense and you can figure out from there what you
need to do to get it to work. I’ve worked with fred for a while now,
but I find his Matlab to be almost unreadable (I’ve told him this
before), so I can’t quite tell you from a quick look over the code he
has in this paper what’s really going on.

Also, keep in mind that fred always develops he own filters straight
from the remez algorithm. If you are using gr.firdes.low_pass, you’re
not going to get the same output filters. Closer would be to use the
blcks2.optfir.low_pass, which uses the PM algorithm (which in turn
uses Remez), but I couldn’t give you the specs you need to build the
same exact filter as in fred’s paper. But you’ll get close enough
that the basic shapes of the signals will be the same, if not every
point.

Tom

3

Hi, Tom:

Please see my comments below.

Thanks

On 12/29/2010 12:58 AM, Tom R. wrote:

paper. I’m hoping the channelizer result from GNURadio code would match the
result from Matlab, but I failed on both, I wonder if someone could help me
spot the problem.

  1. DDC Channelizer:

Using the frequency xlating filters here is not recommended. You’ll
need to use much longer filters than when using the channelizer. That

Could you clarify? Does this mean I need to change the parameter to
gr.firdes.low_pass, maybe reduce _cutoff_freq?

said, I’m not sure what would cause the problem you identified below.
In my opinion, it’s not worth doing it this way.

Ok, the reason I did this is because I want a way to verify the result
from pfb channelizer is correct, since I’m not sure if I’m using it
correctly.

    # Blocks

Number of taps: _filter_taps = 561
2. Polyphase Channelizer:
self._taps_per_channel

Number of taps: _filter_taps = 561
This doesn’t seem to make sense, since in this case N = 40, oversample rate
filters, the oversample rate, O, has to be set such that the number of
filters, N, can be evenly divided by an integer to to get O; that is,
N/i = O for i in [1,N]. Take the example in in fred’s paper where he
has 64 channels and a 1/48 resampling rate. So he’s using a 48-to-1
downsampling in 64 channels, so his rate is 4/3, which means that the
output sampling rate is (4/3)*bw, where bw is the channel bandwidth.

Sorry, false alarm here, it looks like a problem with MinGW’s compiler,
it didn’t compile double comparison with 0.0 correctly or something. I
tried the program under Ubuntu and the error disappeared. I think the
number I used is the correct number since in my case N=40, i = 28, and O
= N/i = 1.42857142857

If you specified that the rate was 2.0, you would feed 32 of the 64
channels and have a sample rate that is twice the bandwidth. I notice
in the remarks in fred’s paper that he’s outputting the channels at 2
samples/symbol, so the oversample_rate here should be 2.

I’m using fred’s example in appendix III, not the 64 channels example.
The comment summaries the example as below:
% receiver_40z is a demo of a 40 channel receiver, demodulating 30
channels,
% of nominal symbol rate 20 MHz, separated by 28 MHz centers (1.4 times
symbol rate)
% input sample rate is 40*28 = 1120 MHz.
% receiver performs a 40 point transform on the output of a 40-stage
polyphase filter
% the polyphase filter operates at input rate but outputs at 2
samples/symbol or
% 40 MHz. The resampling rate is 1120/40 = 28-to-1, thus output from
40-channels are
% computed once for every 28 input samples. channelizer is not matched
filter,
% prototype filter is 10% wider than two sided bandwidth of input signal
to accommodate
% frequency uncertainty of separate channel centers.

I hope that makes sense and you can figure out from there what you
need to do to get it to work. I’ve worked with fred for a while now,
but I find his Matlab to be almost unreadable (I’ve told him this
before), so I can’t quite tell you from a quick look over the code he
has in this paper what’s really going on.

I got passed the error, thanks, but still have problems with
channelizer, see below…

Also, keep in mind that fred always develops he own filters straight
from the remez algorithm. If you are using gr.firdes.low_pass, you’re
not going to get the same output filters. Closer would be to use the
blcks2.optfir.low_pass, which uses the PM algorithm (which in turn
uses Remez), but I couldn’t give you the specs you need to build the
same exact filter as in fred’s paper. But you’ll get close enough
that the basic shapes of the signals will be the same, if not every
point.

Thanks for the suggestion, I modified filter statement to:

     self._filter_taps = blks2.optfir.low_pass(1, self._input_rate,
                                               self._cutoff_freq -

self._trans_width/2, self._cutoff_freq + self._trans_width/2,
0.9, 60, 400)

Note I have to use 400 as extra taps otherwise remez won’t work (too
many extremals – cannot continue)

I got the result, but still not the same as the one from Matlab:

  1. The channels seem to be shifted by 20 then wrapped around, i.e.
    assuming channel # starts with 1, the pfb result’s channel 1 seems to be
    matlab’s channel 21, and channel 20 in pfb result seems to be matlab’s
    channel 40, channel 21 in pfb result seems to be matlab’s channel 1.
  2. Assuming the channel shift/wrap is done, the time-domain plot still
    are not the same, pfb channel 1 (matlab’s channel 21)‘s plot seems to be
    the same, and no signal channels do fit, but other channels’ plot are
    not the same (i.e. the situation is similar to my result of DDC).

Is this because the filter taps in GNURadio is different from Matlab? Is
there a way I can create a filter in GNURadio from a list of numbers
exported from Matlab?

4

On Wed, Dec 29, 2010 at 4:32 AM, Jimmy R. [email protected]
wrote:

the
gr.firdes.low_pass, maybe reduce _cutoff_freq?
I was just saying that, in general, you’ll be doing much more work
this way. But since you say later that you’re only using this for
comparisons, then you already know that.

in 1.2
offset = self._channel_bandwidth * (float(self._M)/2 - i)
Input sampling freq (Hz): _input_rate = 1120000000.0
Channel 1 ‘s offset (Hz): 532000000.0
exactly the same, but all other channels’ plot are different.
float(self._M))

i = 0 to M-1

Oversample rate: _oversample_rate = 1.42857142857
way.
critically sampled; that is, the sample rate equals the bandwidth. For

Sorry, false alarm here, it looks like a problem with MinGW’s compiler, it
didn’t compile double comparison with 0.0 correctly or something. I tried
the program under Ubuntu and the error disappeared. I think the number I
used is the correct number since in my case N=40, i = 28, and O = N/i =
1.42857142857

Ok. I didn’t calculate the number myself, just assumed it was wrong if
you were being given that error message :slight_smile:

On a side note, I don’t like that problem happening on MinGW.
Honestly, doing an == on floating point numbers isn’t a good idea, so
I may need to rethink that check.

If you specified that the rate was 2.0, you would feed 32 of the 64
channels and have a sample rate that is twice the bandwidth. I notice
in the remarks in fred’s paper that he’s outputting the channels at 2
samples/symbol, so the oversample_rate here should be 2.

I’m using fred’s example in appendix III, not the 64 channels example. The
comment summaries the example as below:

Yeah, I know. I was using fred’s 64-channel example from the paper
because the numbers are easier to handle/explain.

% computed once for every 28 input samples. channelizer is not matched

point.

Thanks for the suggestion, I modified filter statement to:

self._filter_taps = blks2.optfir.low_pass(1, self._input_rate,
self._cutoff_freq -
self._trans_width/2, self._cutoff_freq + self._trans_width/2,
0.9, 60, 400)

Note I have to use 400 as extra taps otherwise remez won’t work (too many
extremals – cannot continue)

Hm, that looks wrong. You shouldn’t need to use the extra taps input
just to get the algorithm to converge. It’s really meant to be able to
give you minor control over the number of taps you use, specifically
if you need a filter with an even or odd number of taps. Using 400 is
way too many.

The real problem I think comes from your in-band ripple or 0.9. That’s
a value in dB and specifies the magnitude of the peak-to-peak
fluctuations in the passband. I was able to use your settings above
but specified a ripple of 0.1, which produces a filter of 764 taps.

I got the result, but still not the same as the one from Matlab:

  1. The channels seem to be shifted by 20 then wrapped around, i.e. assuming
    channel # starts with 1, the pfb result’s channel 1 seems to be matlab’s
    channel 21, and channel 20 in pfb result seems to be matlab’s channel 40,
    channel 21 in pfb result seems to be matlab’s channel 1.
  2. Assuming the channel shift/wrap is done, the time-domain plot still are
    not the same, pfb channel 1 (matlab’s channel 21)‘s plot seems to be the
    same, and no signal channels do fit, but other channels’ plot are not the
    same (i.e. the situation is similar to my result of DDC).

This strikes me as a simple naming issue; which channel you decide is
0 versus N. The way that I have the output of the channelizer working,
channel 0 is the DC channel, channel 1 is the next one to the “right”
(the next one in moving in the positive frequency direction). Channel
floor(N/2) will either span the last part of the positive spectrum and
the first part of the negative spectrum (if N is even) or be the last
channel in the positive spectrum (if N is odd). Channel N-1 is
directly “left” of channel 0.

It’s hard to know what you mean by the signals being “the same.” Are
they completely different?

Is this because the filter taps in GNURadio is different from Matlab? Is
there a way I can create a filter in GNURadio from a list of numbers
exported from Matlab?

It could be a filter issue. If you want to, generate you prototype
filter in Matlab, then you can just copy-and-paste those into your GNU
Radio program as a list and pass that to the pfb_channelizer.

Tom

5

Hi, Tom:

Please see my comments below.

Thanks

On 1/4/2011 12:36 AM, Tom R. wrote:

wrote:

    self._filter_taps = blks2.optfir.low_pass(1, self._input_rate,

way too many.
0.1, 60)
it, but because it’s doing a curve fit, thinking about it in terms of
one filter being more “relaxed” than another doesn’t quite work. Sorry
I can’t give you more details or math about why this is true. Maybe
That is alright, I think I can take it as a quirk of the algorithm, it’s
nice to know :slight_smile:

are
directly “left” of channel 0.
appendix_III_pfb.matlabfilter.result_time.jpg: This is the result from
probably the shift you’re seeing.
Yes, the lengths are different, optfir is 1524, matlab is 600

  1. Comparing the result from matlab and pfb_channelizer_ccf: Channel 0, 5,
    10, 30, 35 are the same, the others are different, very strange.
    Indeed. Very strange. My guess is that there is a misconception
    somewhere in the code about sample rate. I can’t quite see it in my
    head, but I’m guessing that the channel spacing isn’t exactly the
    channel spacing you think it is.

You mean the oversample rate O = N/i is wrong? As far as I can see,
that’s the only number that could go wrong, since pfb_channelizer_ccf
only takes 3 parameters, # of channels, taps and oversample rate, and it
doesn’t look like the other two can be wrong.

I checked the calculation of O again, but couldn’t find the problem. In
the Matlab code, harris uses loops “for nn=1:28:5600-28” to do the
processing, so it does seems a 28-to-1 downsampling in 40 channels.

Have you looked at the examples in gnuradio-examples/python/pfb?
Specifically, channelize.py, chirp_channelize.py, and synth_to_chan.py.

I checked the first two samples when I wrote the channelizer, although
it seems they didn’t use the oversample rate parameter. I couldn’t find
the last one (synth_to_chan.py) though.

6

On Sun, Jan 2, 2011 at 11:57 PM, Jimmy R. [email protected]
wrote:

wrote:

Thanks for the suggestion, I modified filter statement to:
just to get the algorithm to converge. It’s really meant to be able to
to:
self._filter_taps = blks2.optfir.low_pass(1, self._input_rate,
self._cutoff_freq -
self._trans_width/2, self._cutoff_freq + self._trans_width/2,
0.1, 60)
and no extra taps are needed, although the new filter has 1524 taps.

But I’m a bit confused why this works, I thought smaller ripple is better,
and thus harder to get, right? And the ripple parameter is the maximum
ripple allowed for the filter? If so, 0.9 should be more relaxed than 0.1?

It’s not quite that straight-forward. I’ve found issues when designing
filters with an out-of-band attenuation of 80 that work fine but an
attenuation of 60 fail to converge. The PM algorithm is essentially an
optimization procedure that tries to fit a curve under the provided
conditions. It’s been a while since I’ve looked at the specifics of
it, but because it’s doing a curve fit, thinking about it in terms of
one filter being more “relaxed” than another doesn’t quite work. Sorry
I can’t give you more details or math about why this is true. Maybe
someone else here has gone through this before and has more to offer?

same (i.e. the situation is similar to my result of DDC).
Ok, thanks for the explanation, this cleared things up, I added a shift to
appendix_III_pfb.matlabfilter.result_time.jpg: This is the result from
pfb_channelizer_ccf, using filter generated from matlab code

Thanks, that helps.

Here’s my observations:

  1. The two results from pfb_channelizer_ccf are the same, minus some time
    shift, so I think filter is no the problem here.

Due to group delay in the filters, almost certainly. What are the
filter lengths for the optfir and matlab filters? The difference is
probably the shift you’re seeing.

  1. Comparing the result from matlab and pfb_channelizer_ccf: Channel 0, 5,
    10, 30, 35 are the same, the others are different, very strange.

Indeed. Very strange. My guess is that there is a misconception
somewhere in the code about sample rate. I can’t quite see it in my
head, but I’m guessing that the channel spacing isn’t exactly the
channel spacing you think it is.

Have you looked at the examples in gnuradio-examples/python/pfb?
Specifically, channelize.py, chirp_channelize.py, and synth_to_chan.py.

Tom

7

On Mon, Jan 3, 2011 at 11:25 PM, Jimmy R. [email protected]
wrote:

I checked the calculation of O again, but couldn’t find the problem. In the
Matlab code, harris uses loops “for nn=1:28:5600-28” to do the processing,
so it does seems a 28-to-1 downsampling in 40 channels.

No, not the oversample rate, just the sample rate. I’m thinking that
there is some misunderstanding somewhere about the actual sample rate
and the therefore the channel bandwidths such that the channels you
are pulling out are not covering the same frequencies that you think
they are.

On the other hand, it could be a miscalculation in the
pfb_channelizer, although all of my tests with it came out fine.

Have you looked at the examples in gnuradio-examples/python/pfb?
Specifically, channelize.py, chirp_channelize.py, and synth_to_chan.py.

I checked the first two samples when I wrote the channelizer, although it
seems they didn’t use the oversample rate parameter. I couldn’t find the
last one (synth_to_chan.py) though.

Tom

The synth_to_chan.py is not in the tarball release. It was introduced
a little while ago when I added the synthesis filterbank. You’d have
to get it from the git repo.

Tom

8

Hi, Tom:

Please see my comments below.

Thanks

On 1/4/2011 11:44 PM, Tom R. wrote:

I checked the calculation of O again, but couldn’t find the problem. In the
Matlab code, harris uses loops “for nn=1:28:5600-28” to do the processing,
so it does seems a 28-to-1 downsampling in 40 channels.
No, not the oversample rate, just the sample rate. I’m thinking that
there is some misunderstanding somewhere about the actual sample rate
and the therefore the channel bandwidths such that the channels you
are pulling out are not covering the same frequencies that you think
they are.

Sorry, maybe I didn’t get how the algorithm works, but it looks to me
the actual sample rate does not enter into the equation, except in the
filter generation part? If we use the same filter on both methods, then
the actual sample rate should not cause the difference?

On the other hand, it could be a miscalculation in the
pfb_channelizer, although all of my tests with it came out fine.

Is the test cases included in the GNURadio source code? Also if I want
to understand the logic inside pfb_channelizer implementation, are there
any book or paper I should review besides harris’ paper and book?

Have you looked at the examples in gnuradio-examples/python/pfb?
Specifically, channelize.py, chirp_channelize.py, and synth_to_chan.py.
I checked the first two samples when I wrote the channelizer, although it
seems they didn’t use the oversample rate parameter. I couldn’t find the
last one (synth_to_chan.py) though.

Tom
The synth_to_chan.py is not in the tarball release. It was introduced
a little while ago when I added the synthesis filterbank. You’d have
to get it from the git repo.

Ok, I’ll take a look.

9

On Wed, Jan 5, 2011 at 4:14 AM, Jimmy R. [email protected]
wrote:

like the other two can be wrong.
are pulling out are not covering the same frequencies that you think
they are.

Sorry, maybe I didn’t get how the algorithm works, but it looks to me the
actual sample rate does not enter into the equation, except in the filter
generation part? If we use the same filter on both methods, then the actual
sample rate should not cause the difference?

No, the sample rate as a real number doesn’t enter into the code, but
it does abstractly in where the signals exist, what the channel
bandwidths mean, etc. I’m not sure what else to say here, really, just
to think closely about the signals and the channelizer. Again, it
could be that there is a bug in the code that only worked for my
cases.

On the other hand, it could be a miscalculation in the
pfb_channelizer, although all of my tests with it came out fine.

Is the test cases included in the GNURadio source code? Also if I want to
understand the logic inside pfb_channelizer implementation, are there any
book or paper I should review besides harris’ paper and book?

No, they aren’t. It shouldn’t be difficult to modify the example codes
that are there, though, to see what happens and make sure it makes
sense.

Tom

10

Hi, Tom:

Thank you for the suggestion, I’ll review the code and examples again to
see what I missed.

Thanks