usrpL0变采样滤波器的性能优化

    项目是一个从物理层到协议层到应用层的LTE系统，射频端在使用ettus公司的USRP软件无线电平台，嗯，这些就不要说太多了，说这篇文章中的重点。

    FPGA端最高只支持25M的采样速率，而LTE的基带需要的是30.72M的样本速率，所以在L0到PHY之间必须有一个变采样的操作。老师思忖再三，决定还是把这个工作放到L0的PC上处理，而不是在FPGA上就变采样，貌似是在FPGA上做会对FPGA的代码改动比较大。放到L0的PC程序上处理当然没有问题，不过速度就成了必须要解决的问题。

    做过变采样的同学都知道，变采样必须要有一个滤波操作，何况我们这样一个插值抽取倍数都还不算小的系统，运算量是很大的。师兄之前写了一个版本，用10000阶的滤波器进行变采样滤波，1ms的数据平均执行时间在2000多us，明显不够用。后来他使用查表、减少除法运算等等方法，把这个时间压缩到了800us，再加上我用SSE指令集再做一下并行，这个时间就压缩到220us左右了，足够满足1ms的实时操作需求。

    这个版本我们就用了很久，但是后来，由于系统的需求，整个系统的延时要变的更小，而之前的一些设计就给整个系统带来了太大的延时。这个延时出现在哪儿，我们一次又一次的定位：

1、同步与数据的排队(定位了延时的产生原因，找到了解决延时的办法)

1.1 问题描述

    在帧号机制完成之前，由于usrpL0和物理层的程序是两个进程，因此我们无法在两个进程中定位同一子帧。帧号机制设计好之后，我们测量并计算了一子帧数据从usrpL0接收到它，到物理层接收到它所消耗的时间，大概在2500us左右。而这一过程合理的延时最多不应该超过1200us。(其中包含数据本身的时常1000us)

1.2 问题分析

    usrpL0最开始设计的时候，组帧的功能是在usrpL0实现的。后来我们发现，同步的结果(时偏)是由物理层给出的，而将同步结果传递给usrpL0进程，不仅是一种异步的处理，而且还会额外消耗系统资源，于是就把组帧功能放在了物理层。而usrpL0的机制没有改变，依然会等待一个子帧的数据量，再向上提交。

    数据好像排了两次队，那么其中就有一次是多余的。usrpL0提供的一子帧长的数据是非同步的，物理层为了获取一个同步的子帧，必须要等待usrpL0向上提交包含该子帧的一整个长度为30720的数据包。如图1所示，时间轴上方是USRP提交的数据，下方左边两个箭头之间表示的是一个同步后的子帧。但是，为了获得这个子帧，物理层必须要等到USRP提交的第二个包上交，也就是下方最右边箭头指示的时间才能收到。那么下方的右边两个箭头之间的时间就是第二次排队多余的耗时。

图 1 第二次排队详解

   

    至此，问题已经基本分析清楚了，2500us的耗时如图2所示。usrpL0等待一子帧数据之后才进行滤波，这里耗时有1ms；滤波因为是密集运算，本身受系统负载的影响，最差耗时有0.4ms；物理层的第二次排队，如上面分析的，最差情况会耗时近1ms；再加上进程间通信和其他开销0.2ms，那么总的耗时就有近2600ms了。

图2  问题分析

1.3初步方案

    问题分析出来以后，大的方向就很明确了，就是减少一次排队。既然组帧功能放在了物理层，那usrpL0层的排队就可以去掉。由于usrpL0使用的UHD是通过以太网实现射频前段和PC的通信，为了与以太网协议相适应，一个最小的包长是362个sample。为了方便的处理，我们在usrpL0和物理层的接口一侧(采样速率30.72M)取一子帧长度的1/40，也就是768个sample作为一个packet，也就是新的usrpL0向物理层传递的最小粒度，取代原来以一子帧作为最小粒度的方案。这样，第一次排队虽然还存在，但是等待的最长时间就从1ms降到了362个362/25000(以usrpL0的采样速率25M为准)=0.0145ms，大大降低了排队时延。

    我们前面说过，尽管我们优化了很多，对一个子帧的数据进行变采样滤波也要花费300us的时间。而我们的优化目标是整体降到200us以内，那么，排队时延降下来了，变采样滤波有没有可能也优化呢？答案是肯定的。如果我们滤波的基本单位也是packet而不是子帧，那么usrpL0的底层每收到一个packet，就对它执行滤波，紧接着就交付物理层等待组帧，那么从usrpL0的底层从空口收到包含一个同步子帧的最后一个packet，到物理层组帧模块收到这个子帧的全部数据，这段延时就可以从原来对一个子帧滤波的时间，降到对一个packet滤波的时间。大胆假设一个packet的滤波时间是一个子帧的1/40，那么滤波造成的延时就只有7.5us了。虽然1/40是达不到的，但是我们认为降到1/10，也就是30us还是可能的。

    也就是说，如果我们把usrpL0的滤波和向物理层传递的粒度都改成小的packet，那么根据上面的分析，排队时延降至15us，滤波时延降至30us，再加上包括进程间通信等其他流程的固定开销，保守估计100us(上文分析这部分的时延是200us，是因为当时传输的是整个子帧，粒度变小以后就会降低)，总的时延就降至145us，就符合我们之前所提的要求了。

图3 理想情况

2、对小packet滤波的设计和优化(尽可能地优化了运行时间，使得小packet滤波成为可能，并且验证了1中定位的问题)

2.1初步设计

    初步设计的思路很简单。师兄之前写的滤波程序很合理，接口没有写死，只要把输入数据的长度改为原来的1/40就行了。

    我们担心的是，滤波的函数有与处理数据量无关的固定开销，而新的方案会把这部分开销放大40倍，也就是每对一个packet进行滤波，都要消耗这部分的开销。而在串行的处理过程中，如果一个packet的处理时间超过了25us(1ms的1/40)，那么这个系统就无法满足最基本的实时性要求。

2.2内存分配优化

    改完一试，果然不行，一个子帧的滤波时间达到了平均40us。实时性满足不了，整个系统就根本转不起来。那就想办法优化吧。

    滤波的函数乍一看，好像并没有什么可以优化的地方。代码使用查表去掉了大部分的除法和取模运算，而且运算部分也用SSE指令集优化过，SSE指令的使用虽然离极限还有距离，但是也算是精心布置过的…...这时候，程序初始化时候的一条语句引起了注意。

short *in_buffer = new short[75000]();

    这句代码意思也很清楚，就是申请了一个临时变量存储我们的输入数据，用完了会delete掉它。我们处理的信号是short类型的复信号，30720个sample需要61440个short类型的空间，再加上处理的时候额外需要的空间，再加上一点余量(后来看书才发现，其实分配内存的时候是不应该留余量的，这样有助于发现越界错误)，就写到了75000。

    除了设置内存余量这样的不良习惯以外，这句语句还有两个严重的缺陷：第一，没有根据输入数据的长度来决定分配内存空间；第二，密集运算使用的临时内存应该一次性分配而不是动态分配。第一个问题导致了，即使输入的数据长度已经降到了1/40，程序在申请内存的时候还是申请了75000个，而第二个问题导致了，每次执行滤波都要动态分配内存，而这个操作是很耗时的。

    为了解决这个问题，我们在滤波处理所在的类的内部增加了一个数组作为数据成员。滤波程序执行的时候，把要处理的数据直接放到这块内存里面，省去了动态分配内存的步骤。每个packet的滤波时间就降到平均15us了。

2.3多线程并行优化

    到了这里，系统已经可以跑起来了，但是有时能正常运行半个小时，有时候几秒钟以后就会报错。分析日志一看，好像是实时性还是不够，15us已经比25us低很多了啊，这是怎么回事呢？只好继续分析。

    usrpL0对一个packet的处理基本上就分成三步：通过UHD驱动从前段获取一个packet的数据量、变采样滤波、通过socket将数据发给物理层。分块测试时间，我们发现，第一步的平均时间是9us，第二步的平均时间是15us，第三步的平均时间是1~2us。这三步我们是串行处理的，三步加在一起，开销正好是实时性需求的临界25us。在这种情况下，当系统的运行状态不是很稳定，比如说其他负载很大，或者因为夏天散热不畅导致性能下降的时候，开销就会超过临界。

    问题分析清楚以后，优化的思路很快也出来了，就是用多线程把这三步并行起来。并行的效果是，原来处理一个packet的耗时是三步的耗时相加，而并行之后就是三步的耗时中最大的一个。由于底层的数据是25us才到一个packet，所以当我们的处理时间快于25us时，第一步就会等待底层的数据。为了更明显地表述并行优化带来的好处，我们假设底层随时有数据，也就是说，第一步不包含等待数据的时间。我们以3个packet的处理时间为例，如图4所示。

图4 并行流程说明

    在串行的处理模式下，处理一个packet需要9+15+2=26us，3个packet的处理时间为26*3=84us。而在多线程并行处理的情况下，假设第一步不需要等待，第三步需要等待第二步的处理结果，那么3个packet的处理时间为9+15*3+2=56us。假设第一步永远不需要等待，那么当packet的数量趋于无穷时，串行处理一个packet的时间仍然是26us，而并行处理一个packet的时间就仅为15us了。

    在实现的时候，我们将三个步骤分别细分为一个线程，线程之间使用zeromq的inproc机制传递数据。由于第一步和第三步的运算量不大，可以绑定到同一个CPU核，第二步单独绑定到一个CPU核，多线程的并行处理就算完成了。经过测试，这样的处理在我们的服务器上能很好地满足系统的实时要求。

3、进一步提升性能

    小包滤波的功能完成之后，我们的系统已经可以比较稳定地运行了，但是延时性能如何，我们仍然需要测试一下。我们依旧像原来那样，测量了一下接收端一个子帧从usrpL0接收到，到物理层完成对其的组帧的过程的耗时。这个耗时按理来说，应该是1000us加上几十us，其中1000us是子帧本身的持续时间，几十us是系统的处理时间。我们看到大部分(90%)的子帧能在1200us之内处理完，而剩下的小部分子帧的处理时间则超过了1200us，有的甚至超过了2000us。这依然是系统所不能忍受的，因为处理时间过长就意味着这个帧会失效，进而影响到若干个帧的调度和处理。那这一小部分子帧的处理时间过长又是怎么回事呢？我们继续往下看。

3.1日志输出到内存

    首先，我们先定位一下极少部分处理时间超过2000us的子帧。干想是想不出什么线索的，最好的办法还是分析日志；单个数据看不出什么线索，那就把一堆数据连在一起分析，差分，绘图，总会有蛛丝马迹的线索。(这是我这阶段最大的收获之一)

    通过仔细分析日志我们发现，这样超时的子帧通常是连续出现的，而且一般是每5秒出现一次，非常规律。于是我们开始在操作系统的后台处理上下功夫。WHY找到的一条线索很有用，他说linux将内存中的数据(脏页)写入硬盘好像默认是5秒一次，而我们在运行我们的系统时，正好有大量的数据需要从内存移入硬盘，那就是我们的日志。我们最开始想把这个默认值改一下，但是因为日志的数据量确实很大，不管是增大写入频率，还是降低写入频率，好像都不能从根本上解决问题。WHY又提出来一个更直接的办法，把日志直接写到/dev/shm里面去。/dev/shm这个目录不在硬盘上，而是在内存里，这就省去了从内存写入硬盘这一IO操作。而且日志存在这个目录下，只要不重启服务器，存在里面的日志不会丢失，如果我们需要长期保存，只要额外把需要的东西手动移到硬盘就行了。  

    操作起来也很简单。我们的日志都是用重定向的方法生成的，比如a.out>test.txt。现在只要把重定向的目标改成a.out>/dev/shm/test.txt就行了。尝试了一下，再看日志，5秒一次的超时已经不见了。

3.2网卡网卡，奇怪的五指山

    在此之前，我们测试的接口时间(不包括子帧本身的持续时间1ms)已经稳定在了200us以内，这已经可以大概满足系统的实时性要求了，但是距离理想状态(几十us)好像还是有一定距离。而且我们在分析日志的时候，还发现了这样的一件怪事。

 

图 5 奇怪的延时概率分布

    图5反应了延时的概率分布。概率在某些值上出现了峰值，而在某些值上则出现了低谷。我们实际接收到的数据，按理来说只有3种不同的数据：系统消息、数据负载和空数据，按道理说不会有这么多个峰值。如果能定位这个奇怪的现象，说不定就能把延时性能进一步地优化。通过进一步分析日志我们发现，UHD的底层接收连续两个packet之间的耗时是不同的，呈现出一定规律性的变化，即大概每8~9个packet，就会出现一次较大的延时(100多us)，而正常的延时通常是在几个us。如图6所示，第2个和第1个packet之间隔了101us，而第3个和第2个packet之间只隔了9us；第6个和第7个packet之间隔了185us，而第7个和第8个packet之间只隔了9us。

图 6 接收packet的延时

    我们对这个现象的原因做了很多个假设，又一次一次地推翻。最后，师兄发现了Intel系列网卡的InterruptThrottleRate参数可能是影响这个问题的关键。

    InterruptThrottleRate

    该参数用来控制网卡每秒钟能够产生的最大中断数目。参数的有效值为：0,1,3,100-100000 (0=off, 1=dynamic, 3=dynamic conservative) 默认值: 3 。增大中断的意义是显而易见的，降低CPU的使用率。但是同时也增加了数据包处理的延迟，这里的延迟主要指的是在缓冲区的排队延迟。

    这个参数决定了网卡每秒钟最多能产生多少次中断。这个值越小，网卡的两次中断之间就会隔更长的时间。而我们的以太网卡，至少每15us就要发一个数据包，如果网卡两次中断之间的间隔大于这个值，就意味着这个数据包要在网卡缓存，直到网卡产生下一次中断，把这个数据包交给CPU为止。加入网卡两次中断之间的最大间隔是100us，那就意味着，有可能会有6-7个数据包被网卡“暂扣”了，最大的延时就会达到100us。

    我们起初尝试直接对这个数据进行了改动，发现并没有什么改善。但我们仍然觉得这是一条有价值的线索，直到我们找到了这样一句说明

    user can now modifyInterruptThrottleRate settings using ethtool-C rx-usecs N

    有可能是操作系统使用ethtool把InterruptThrottleRate参数封装起来了，因此不能直接操作它，但是操作系统提供了ethtool间接地操作这个参数。我们使用ethtool –c命令查看了一下与usrp连接的以太网卡的配置，发现rx-usecs参数的值是3，和InterruptThrottleRate的默认值是一致的，也就是说，网卡每秒产生的最大中断数目是由系统根据负载情况动态分配的。但是我们的应用需要高的实时性，不应该设定最大的中断数，于是我们执行ethtool –C rx-usecs 0，把rx-usecs的值改成了0，之前那个奇怪的图就不见了，新的延时性能图如图7所示，接口延时被基本控制在了80us以内，比较符合之前的理想预期了。

图 7更改配置之后的概率分布

PS

最后再附上变采样滤波函数的部分代码吧，这个是还没有优化的原始代码，其实是可以照着教材编出来的。如果看官有更佳的方法，请联系我。

// table_len should be decied by the length of the filter output// but a default value of 50000 is adequate for our appvoid buildPolyphaseCache(const int P, const int Q, int *&branch_table, int *&offset_table, const int table_len) {    branch_table = new int[table_len];    offset_table = new int[table_len];    for (int i = 0; i < table_len; ++i) {        branch_table[i] = (i * Q) % P;        offset_table[i] = (i * Q - branch_table[i]) / P;    }}

void polyphaseResample(short *in_data, unsigned int in_sample_count,    short *&out_data, unsigned int &out_sample_count,    short *filter_coeff, unsigned int filter_len,    const int P, const int Q, const unsigned int fix_shift,    int *branch_table, int *offset_table) {    // calculate the sample count of output    // out_sample_count = (int)ceil(in_sample_count * (float)P / (float)Q);    out_sample_count = (int)(in_sample_count * (float)P / (float)Q);    // allocate space for output data, both real and imag    // default initilization make them 0    out_data = new short[out_sample_count * 2]();    short *output_ptr = out_data;    // int output_ix = (filter_len - 1) / 2 / Q;    int output_ix = 0;    short const *input_end = in_data + in_sample_count * 2;    short const *output_end = out_data + out_sample_count * 2;    short const *filter_end = filter_coeff + filter_len;    while (output_ptr < output_end) {        int output_branch = branch_table[output_ix];        int input_offset = offset_table[output_ix];        short *input_ptr = in_data + 2 * input_offset;        short *filter_ptr = filter_coeff + output_branch;                while (input_ptr >= input_end) {            input_ptr -= 2;            filter_ptr += P;        }                int sum_real = 0, sum_imag = 0;        while ((input_ptr >= in_data) && (filter_ptr < filter_end)) {            sum_real += (*input_ptr) * (*filter_ptr);            sum_imag += (*(input_ptr + 1)) * (*filter_ptr);            input_ptr -= 2;            filter_ptr += P;        }        *output_ptr = sum_real >> fix_shift;        *(output_ptr + 1) = sum_imag >> fix_shift;        output_ptr += 2;        output_ix++;    }}