SDR-related utility to compute channels for PFB channelizers.

channelize_example.py

"""
Computes channelization scheme for a set of input frequencies.

Example output:

$ python channelize.py
Input frequency spacing:
        freq[00] 854037500
        freq[01] 854212500      dfreq=175000.000000
        freq[02] 854312500      dfreq=100000.000000
        freq[03] 854362500      dfreq=50000.000000
        freq[04] 854512500      dfreq=150000.000000
        freq[05] 854962500      dfreq=450000.000000
        freq[06] 855487500      dfreq=525000.000000
        freq[07] 855737500      dfreq=250000.000000
        freq[08] 856212500      dfreq=475000.000000
        freq[09] 856237500      dfreq=25000.000000
        freq[10] 856587500      dfreq=350000.000000
        freq[11] 856687500      dfreq=100000.000000
        freq[12] 856962500      dfreq=275000.000000
        freq[13] 857212500      dfreq=250000.000000
        freq[14] 857237500      dfreq=25000.000000
        freq[15] 858212500      dfreq=975000.000000
        freq[16] 858237500      dfreq=25000.000000
        freq[17] 858687500      dfreq=450000.000000
        freq[18] 859212500      dfreq=525000.000000
        freq[19] 859237500      dfreq=25000.000000
*** Forcing number of channels to be odd by adding a channel.
*** NOTE: Center channel (ch[19], pfb_ch[00])  may have issues if there is I/Q DC offset leakage.
group[0]: low: 854025000.000000, high: 855000000.000000, center: 854512500.000000, samp_rate: 975000.000000, num_chans: 39, chan_width: 25000
        ch[00], pfb_ch[20]: 854037500.000000
        ch[07], pfb_ch[27]: 854212500.000000
        ch[11], pfb_ch[31]: 854312500.000000
        ch[13], pfb_ch[33]: 854362500.000000
        ch[19], pfb_ch[00]: 854512500.000000
        ch[37], pfb_ch[18]: 854962500.000000
group[1]: low: 855475000.000000, high: 857250000.000000, center: 856362500.000000, samp_rate: 1775000.000000, num_chans: 71, chan_width: 25000
        ch[00], pfb_ch[36]: 855487500.000000
        ch[10], pfb_ch[46]: 855737500.000000
        ch[29], pfb_ch[65]: 856212500.000000
        ch[30], pfb_ch[66]: 856237500.000000
        ch[44], pfb_ch[09]: 856587500.000000
        ch[48], pfb_ch[13]: 856687500.000000
        ch[59], pfb_ch[24]: 856962500.000000
        ch[69], pfb_ch[34]: 857212500.000000
        ch[70], pfb_ch[35]: 857237500.000000
*** Forcing number of channels to be odd by adding a channel.
group[2]: low: 858200000.000000, high: 859275000.000000, center: 858737500.000000, samp_rate: 1075000.000000, num_chans: 43, chan_width: 25000
        ch[00], pfb_ch[22]: 858212500.000000
        ch[01], pfb_ch[23]: 858237500.000000
        ch[19], pfb_ch[41]: 858687500.000000
        ch[40], pfb_ch[19]: 859212500.000000
        ch[41], pfb_ch[20]: 859237500.000000

"""

M = 1.0e6

# Example P25 system:
# http://www.radioreference.com/apps/db/?sid=8285
# Not included: control channels.
INPUT_FREQS = [854.0375 * M,
               854.2125 * M,
               854.3125 * M,
               854.3625 * M,
               854.5125 * M,
               854.9625 * M,
               855.4875 * M,
               855.7375 * M,
               856.2125 * M,
               856.2375 * M,
               856.5875 * M,
               856.6875 * M,
               856.9625 * M,
               857.2125 * M,
               857.2375 * M,
               858.2125 * M,
               858.2375 * M,
               858.6875 * M,
               859.2125 * M,
               859.2375 * M]

# See note below about defining groups by hand.  I did this because I liked my grouping
# better than what define_groups produced.
GROUPS = [[], [], []]
# group0: 37 channels w/25kHz spacing
GROUPS[0] = [854.0375 * M,
             854.2125 * M,
             854.3125 * M,
             854.3625 * M,
             854.5125 * M,
             854.9625 * M]

# group1: 71 channels w/25kHz spacing
GROUPS[1] = [855.4875 * M,
             855.7375 * M, 
             856.2125 * M, 
             856.2375 * M,
             856.5875 * M, 
             856.6875 * M, 
             856.9625 * M, 
             857.2125 * M,
             857.2375 * M]
 
# group2: 42 channels w/25kHz spacing
GROUPS[2] = [858.2125 * M, 
             858.2375 * M, 
             858.6875 * M,
             859.2125 * M, 
             859.2375 * M]



MAX_BW_HZ = 2.54 * M # RTL max reliable sample rate = 2.56 MS/s

def test_chan_width(freqs, chan_width):
    """
    Given a sorted list of frequencies and possible channel width, returns
    True if some multiple of chan_width lands exactly on each frequency in the
    list.
    """
    freq_idx = 0
    current_freq = freqs[0]
    while(True):
        if(freqs[freq_idx] == current_freq):
            freq_idx += 1
        else:
            # Try next channel
            current_freq += chan_width

        # Check for done
        # Matched all freqs
        if(freq_idx == len(freqs)):
            return True

        # Stop searching after maximum frequency
        if(current_freq > freqs[freq_idx]):
            return False


def compute_chan_width(freqs):
    """
    Given a list of frequencies, returns the largest channel size that would
    result in all frequencies being at the center of a channel.
    """
    freqs.sort()
    last_freq = None
    dfreqs = []
    for freq in freqs:
        if(last_freq != None):
            dfreq = (freq - last_freq)
            dfreqs.append(dfreq)
            # print freq, dfreq
        # else:
            # print freq, 0
        last_freq = freq

    # Largest possible channel width that could work is the minimum spacing
    # between frequencies.  This is the starting point.
    chan_width = min(dfreqs)

    # Iteratively search for a channel spacing that works.  Might end up at 1Hz!
    while((test_chan_width(freqs, chan_width) is False) and (chan_width > 0)):
        chan_width -= 1
    num_chans = int((freqs[-1] - freqs[0]) / chan_width) + 1
    # print "chan_width: %f, num_chans: %d" % (chan_width, num_chans)

    return chan_width, num_chans


def define_groups(freqs, max_bw_hz):
    """
    Greedily divide list of input frequencies into groups of at most max_bw_hz.
    """

    freqs.sort()
    groups = []

    low_freq = freqs[0]
    group = []
    for freq in freqs:
    # for i, freq in enumerate(freqs):
        if((freq - low_freq) <= max_bw_hz):
            group.append(freq)
        else:
            # Append previous group
            groups.append(group)
            group = [freq]
            low_freq = freq

    if(len(group) != 0):
        groups.append(group)

    msg = "Partitioned %d input frequencies into %d groups "
    msg += "based on your requirement for maximum group b/w of %f Hz:"
    msg %= (len(freqs), len(groups), MAX_BW_HZ)
    print msg

    for i, group in enumerate(groups):
        msg = "group[%d]: %d frequencies spanning %f Hz: %s"
        msg %= (i, len(group), group[-1] - group[0], repr(group))
        print msg

    return groups

def show_spacing(freqs):
    """
    Prints freqs list and delta frequency between each item in the list.
    Assumes list is sorted.
    """
    print "Input frequency spacing:"
    last_freq = None
    for i, freq in enumerate(freqs):
        if(last_freq != None):
            dfreq = freq - last_freq
            msg = "\tfreq[%02d] %d\tdfreq=%f"
            msg %= (i, freq, dfreq)
            print msg
        else:
            msg = "\tfreq[%02d] %d"
            msg %= (i, freq)
            print msg
        last_freq = freq


def main():
    show_spacing(INPUT_FREQS)

    # Originally the program divided the frequency list into groups based on
    # the max b/w of each RTL dongle, but it's not as good as doing it
    # manually.  The define_groups routine should be re-written such that it
    # divides the groups at the largest gaps, to reduce overall system b/w.
    # But for now they are manually allocated at the top of the file.

    # groups = define_groups(INPUT_FREQS, MAX_BW_HZ)
    groups = GROUPS

    for i, group in enumerate(groups):
        chan_width, num_chans = compute_chan_width(group)

        # Get channel numbers - depending on whether the center frequency is
        # used these may change.  Also these numbers are zero-based but the
        # PFB channelizer in gnuradio starts with zero in the middle of the
        # tuned band, so they always need to be remapped.

        # Preliminary channel list
        base = group[0]
        chans = []
        for freq in group:
            chan = int((freq - base)/chan_width)
            chans.append(chan)

        if(num_chans & 1 == 0):
            # Even number of channels: make odd by adding another channel so
            # that center of channels lines up with desired signal in PFB
            # output.
            num_chans += 1
            print "*** Forcing number of channels to be odd by adding a channel."

        # Compute shift for channel remapping for PFB
        shift = (num_chans+1)/2
        pfb_chan_map = range(num_chans)[shift:] + range(num_chans)[:shift]

        center_chan = (num_chans - 1) / 2
        if(center_chan in chans):
            # FIXME: Perhaps do something to avoid using the enter channel
            # because of the I/Q DC leakage artifact?
            msg = "*** NOTE: Center channel (ch[%02d], pfb_ch[%02d]) "
            msg += " may have issues if there is I/Q DC offset leakage."
            msg %= (center_chan, pfb_chan_map[center_chan])
            print msg

        tune_low = group[0] - (chan_width / 2.0)
        tune_high = tune_low + (num_chans * chan_width)
        tune_freq = (tune_high + tune_low) / 2.0
        samp_rate = tune_high - tune_low
        msg = 'group[%d]: low: %f, high: %f, center: %f, samp_rate: %f, num_chans: %d, chan_width: %d'
        msg %= (i, tune_low, tune_high, tune_freq, samp_rate, num_chans, chan_width)
        print msg

        # remapped_chans = []
        # for chan in chans:
        for freq in group:
            chan = int((freq - base)/chan_width)
            # chans.append(chan)
            # remapped_chans.append(pfb_chan_map[chan])
            print '\tch[%02d], pfb_ch[%02d]: %f' % (chan, pfb_chan_map[chan], freq)


if __name__ == '__main__':
    main()