《一种基于GPU的高阶数字FIR滤波器频域并行处理实现方法.pdf》由会员分享，可在线阅读，更多相关《一种基于GPU的高阶数字FIR滤波器频域并行处理实现方法.pdf（11页珍藏版）》请在专利查询网上搜索。
1、10申请公布号CN102340296A43申请公布日20120201CN102340296ACN102340296A21申请号201110204946622申请日20110721H03H17/02200601G06F9/3820060171申请人东北大学秦皇岛分校地址066004河北省秦皇岛市经济技术开发区泰山路143号72发明人汪晋宽张春宏韩英华宋昕高静74专利代理机构北京理工大学专利中心11120代理人李爱英郭德忠54发明名称一种基于GPU的高阶数字FIR滤波器频域并行处理实现方法57摘要为解决高阶数字FIR滤波器并行处理效率的问题，提出了一种高效的、适合GPU体系结构的高阶数字FIR滤波。
2、器并行处理算法，该方法采用重叠保留方法结合GPU自身结构特点优化实现高阶数字FIR的频域并行化处理。通过计算FIR频率响应系数，将待处理的输入数据传送给GPU；数据重叠搬移；滤波计算处理；数据合并搬移；将合并搬移结果YY0，Y1，YK1传送到主机内存等步骤完成高阶数字FIR滤波器频域并行处理。对比在CPU上单线程所实现的FIR频域重叠保留方法，其吞吐率，即每秒处理样点的数量有着极大地提高，典型的加速比在100倍以上。51INTCL19中华人民共和国国家知识产权局12发明专利申请权利要求书2页说明书5页附图3页CN102340310A1/2页21一种基于GPU的高阶数字FIR滤波器频域并行处理实。
3、现方法，其特征在于根据通用图形处理器GPU众核体系结构的特点，将重叠保留方法的处理过程划分为六个核函数频响计算核函数、重叠搬移核函数、傅里叶变换核函数、乘法核函数、傅里叶逆变换核函数、合并搬移核函数；实现该方法的具体步骤如下步骤一、确定有限冲击响应FIR频率响应系数；当给定的响应系数为FIR的冲击响应系数HH0，H1，HM1时，将FIR的冲击响应系数HH0，H1，HM1，经过尾部填0扩展至长度为N，将扩展后的系数传送给GPU，启动频响计算核函数对FIR的冲击响应系数H进行N点的傅里叶变换，结果为有限冲击响应FIR频率响应系数HH0，H1，HN1，保存在GPU内存中；当给定的响应系数为FIR频率。
4、响应系数HH0，H1，HN1时，将频率响应系数H传给GPU，并保存在GPU内存中，保存时，第一个字节的地址为存储器位宽的整数倍；步骤二、将待处理的输入数据传送给GPU；将一块长度为NBLK的待滤波样点数据XB0，B1，BK1从主机内存中传入GPU的内存中，其中NBLKKL，K为整数，L为每个数据块的长度；其中BICI，DI，0IK1，CI表示在BI中起点为0，而长度为LM1的连续样点数据块，DI表示在BI中起点为LM1，而长度为M1的连续样点数据块；步骤三、数据的重叠搬移；启动重叠搬移核函数完成数据重叠搬移操作，即将待滤波样点数据XB0，B1，BK1以及输入的长度为M1的初始状态数据S0，重叠。
5、搬移为EE0，E1，EK1，其中E0S0，B0，其中，S0为上一次处理待滤波样点数据X时其中的DK1，对于I0的数据块EIDI1，BI，同时将DK1搬移到中S0，作为下一次处理过程的初始状态数据使用；重叠数据的搬移分为三个步骤，第一步，启动M1个线程同时工作，线程TI0IM2完成将S0中第I个数据搬移到重叠搬移结果数据E0的S0中的第I的位置，然后再完成将DK1中第I个数据搬移到S0中的第I的位置；第二步，启动L个线程同时工作，线程TI0IL1完成将待滤波样点数据X的BJ0JK1块中第I个数据搬移到重叠搬移结果数据E的BJ0JK1块中的第I的位置；第三步，启动M1个线程同时工作，线程TI0IM。
6、2完成将待滤波样点数据X的DJ0JK1块中第I个数据搬移到重叠搬移结果数据E的DJ0JK1块中的第I的位置；步骤四、滤波计算处理；启动傅里叶变换核函数，完成对每个EI做N点的傅里叶变换运算，得到运算结果FF0，F1，FK1；然后启动乘法核函数完成FI与系数H相乘操作，即IIFIH；然后启动傅里叶逆变换核函数，完成对每个II做N点的傅里叶逆变换运算，并得结果RR0，R1，RK1；在启动乘法核函数核函数时，启动N个线程同时工作，线程TI0IN1完成将FJ0JK1中第I个数据与频响系数H中的第I个系数相乘，并将结果保存在RJ0JK1块中的第I的位置；权利要求书CN102340296ACN102340。
7、310A2/2页3步骤五、数据的合并搬移；启动合并搬移核函数完成数据的合并搬移操作，即将结果数据RR0，R1，RK1合并搬移为YY0，Y1，YK1，其中RIZI，YI，其中ZI表示在RI中起点为0，长度为M1的连续数据点所组成的数据块，YI表示在RI中起点为M1，长度为L的连续数据点所组成的数据块；在启动合并搬移核函数时，启动L个线程同时工作，线程TI0IL1完成将乘积RJ0JK1中第M1I个数据搬移到结果数据YJ0JK1块中的第I的位置，其中，启动的L个线程在同一时刻搬运长度为L的连续数据块；步骤六、将合并搬移结果YY0，Y1，YK1传送到主机内存，若还有剩余数据需要处理重复步骤二，否则结束。
8、处理过程，完成高阶数字FIR滤波器频域并行处理。权利要求书CN102340296ACN102340310A1/5页4一种基于GPU的高阶数字FIR滤波器频域并行处理实现方法技术领域0001本发明涉及一种高阶数字FIR滤波器频域并行处理实现方法，特别是一种基于GPU的高阶数字FIR滤波器频域并行处理实现方法，属于数字信号处理领域。技术背景0002在数字信号处理系统中，数字有限冲击响应FIR滤波器是最为核心和基础的数字信号处理算法之一。由于数字FIR滤波器，特别是高阶的数字FIR滤波器的计算复杂度是相当高的，因而高效的、并行化的FIR实现方法对于加速FIR的处理是极其至关重要的。0003目前，随着。
9、通用图形处理器GPU技术的迅猛发展，GPU正在被广泛地应用于众多的应用领域之中，由于GPU所具有的众核体系结构，使其能够提供十分强大的计算能力。0004因此，对于FIR滤波器这种计算复杂度较高的算法，通过设计合理高效且适合GPU结构的并行算法，就可以使得FIR处理在GPU上得到很高的加速比，从而极大地缩短其处理时间。因此，研发适合于GPU体系结构的数字FIR滤波器并行处理算法，具有很高的实际价值和现实意义。发明内容0005本发明所要解决的技术问题是高阶数字FIR滤波器并行处理效率的问题，提出了一种高效的、适合GPU体系结构的高阶数字FIR滤波器并行处理算法，该方法采用重叠保留方法结合GPU自身。
10、结构特点优化实现高阶数字FIR的频域并行化处理。0006为了发挥GPU众核的体系结构，需要将重叠保留方法的计算过程加以分解，并合理地将计算负载分配到在GPU上执行的每个线程中，同时还要优化各个线程对于GPU内存的访问，以最大限度地利用GPU所提供的存储带宽。在GPU的编程模型中，核函数KERNEL是由编程人员定义，并可以在GPU上由众多线程加以并行执行的功能单元。作为CPU加速器的GPU，是典型的FORKJION并行模式。在主机端启动核函数时，通过指定调用配置参数，编程人员可以控制启动执行核函数的线程数量以及线程的组织结构。如何将重叠保留方法中的计算过程分解成为不同的核函数，并确定各个核函数所。
11、要完成的处理，对于整个算法的并发执行效率起着至关重要的影响。如果核函数内部控制流比较复杂，则会大大地降低GPU的并行效率；然而如果核函数功能过于简单，将导致核函数的数量增加，从而增大启动核函数的总体时间，因此同样会降低GPU并行效率。0007因此，根据GPU的自身结构特点，以及重叠保留方法所需的数据处理过程特性，提出了优化的核函数划分方法，即将重叠保留方法的处理划分为六个核函数频响计算核函数、重叠搬移核函数、傅里叶变换核函数、乘法核函数、傅里叶逆变换核函数、合并搬移核函数，算法的总体处理过程如下0008根据通用图形处理器GPU众核体系结构的特点，将重叠保留方法的处理过程划分为六个核函数频响计算。
12、核函数、重叠搬移核函数、傅里叶变换核函数、乘法核函数、傅里叶说明书CN102340296ACN102340310A2/5页5逆变换核函数、合并搬移核函数；0009实现该方法的具体步骤如下0010步骤一、确定有限冲击响应FIR频率响应系数；0011当给定的响应系数为FIR的冲击响应系数HH0，H1，HM1时，将FIR的冲击响应系数HH0，H1，HM1，经过尾部填0扩展至长度为N，将扩展后的系数传送给GPU，启动频响计算核函数对FIR的冲击响应系数H进行N点的傅里叶变换，结果为有限冲击响应FIR频率响应系数HH0，H1，HN1，保存在GPU内存中；0012当给定的响应系数为FIR频率响应系数HH0。
13、，H1，HN1时，将频率响应系数H传给GPU，并保存在GPU内存中，保存时，第一个字节的地址为存储器位宽的整数倍；0013步骤二、将待处理的输入数据传送给GPU；0014将一块长度为NBLK的待滤波样点数据XB0，B1，BK1从主机内存中传入GPU的内存中，其中NBLKKL，K为整数，L为每个数据块的长度；其中BICI，DI，0IK1，CI表示在BI中起点为0，而长度为LM1的连续样点数据块，DI表示在BI中起点为LM1，而长度为M1的连续样点数据块；0015步骤三、数据的重叠搬移；0016启动重叠搬移核函数完成数据重叠搬移操作，即将待滤波样点数据XB0，B1，BK1以及输入的长度为M1的初始。
14、状态数据S0，重叠搬移为EE0，E1，EK1，其中E0S0，B0，其中，S0为上一次处理待滤波样点数据X时其中的DK1，对于I0的数据块EIDI1，BI，同时将DK1搬移到中S0，作为下一次处理过程的初始状态数据使用；0017重叠数据的搬移分为三个步骤，第一步，启动M1个线程同时工作，线程TI0IM2完成将S0中第I个数据搬移到重叠搬移结果数据E0的S0中的第I的位置，然后再完成将DK1中第I个数据搬移到S0中的第I的位置；第二步，启动L个线程同时工作，线程TI0IL1完成将待滤波样点数据X的BJ0JK1块中第I个数据搬移到重叠搬移结果数据E的BJ0JK1块中的第I的位置；第三步，启动M1个线。
15、程同时工作，线程TI0IM2完成将待滤波样点数据X的DJ0JK1块中第I个数据搬移到重叠搬移结果数据E的DJ0JK1块中的第I的位置；0018步骤四、滤波计算处理；0019启动傅里叶变换核函数，完成对每个EI做N点的傅里叶变换运算，得到运算结果FF0，F1，FK1；然后启动乘法核函数完成FI与系数H相乘操作，即IIFIH；然后启动傅里叶逆变换核函数，完成对每个II做N点的傅里叶逆变换运算，并得结果RR0，R1，RK1；0020在启动乘法核函数核函数时，启动N个线程同时工作，线程TI0IN1完成将FJ0JK1中第I个数据与频响系数H中的第I个系数相乘，并将结果保存在RJ0JK1块中的第I的位置；。
16、0021步骤五、数据的合并搬移；0022启动合并搬移核函数完成数据的合并搬移操作，即将结果数据RR0，R1，RK1合并搬移为YY0，Y1，YK1，其中RIZI，YI，其中ZI表示在RI中起点为0，说明书CN102340296ACN102340310A3/5页6长度为M1的连续数据点所组成的数据块，YI表示在RI中起点为M1，长度为L的连续数据点所组成的数据块；0023在启动合并搬移核函数时，启动L个线程同时工作，线程TI0IL1完成将乘积RJ0JK1中第M1I个数据搬移到结果数据YJ0JK1块中的第I的位置，其中，启动的L个线程在同一时刻搬运长度为L的连续数据块；0024步骤六、将合并搬移结果。
17、YY0，Y1，YK1传送到主机内存，若还有剩余数据需要处理重复步骤二，否则结束处理过程，完成高阶数字FIR滤波器频域并行处理。0025有益效果0026通过对于GPU内存访问的优化设计，即所设计的重叠搬移以及合并搬移操作方法，较其他非连续的存储访问方法，能够最大限度地避免内存访问竞争的发生频率，从而可以有效地提高内存访问效率，并提高整体算法的执行效率。另外，通过优化线程组织结构以及启动线程的数量，能够达到最大限度地均衡各个流多处理器上的工作负载，从而提高GPU整体使用率。通过这些优化设计手段，较一般非优化设计，可平均提高执行效率10左右。同时保存数据时，第一个字节的地址为存储器位宽的整数倍，以确。
18、保每次访问内存时能读取到最多的数据。0027本发明方法对比在CPU上单线程所实现的FIR频域重叠保留方法，其吞吐率，即每秒处理样点的数量有着极大地提高，典型的加速比在100倍以上。如在GEFORCE580对于复数数据的处理可达到400M样点/秒。附图说明0028图1是通用数据重叠搬移操作图；0029图2是实例数据重叠搬移操作图；0030图3是通用数据合并搬移操作图；0031图4是实例数据合并搬移操作图；0032图5是并行FIR滤波处理流程框图；具体实施方式0033本实施例以1025点FIR为例，说明一个具体的实施方法。0034我们选择傅里叶变换的长度为21120482K，L1024，每次处理的。
19、数据长度NBLK2201M，故KNBLK/L2101024。0035根据前面所述方法，我们创建下列在GPU上执行的核函数00361频响计算核函数0037其输入为1025点的FIR时域冲击响应，所进行的处理是将1024个点的时域冲击响应在其尾部填充1023个0点，并进行2048点的傅里叶变换，从而得到相应的频率响应系数，并将其存储到GPU的共享内存中。00382重叠搬移核函数0039其处理是将输入的220点样本数据，进行重叠搬移。图1为通用的重叠搬移示意图，XB0，B1，BK1是输入的待滤波的数据块，其长度为NBLKKL；X是由K个子块BI组成，0IK1，每个子块的长度为L。每个子块BI则是由长。
20、度为LM1的数据块CI说明书CN102340296ACN102340310A4/5页7和长度为M1的数据块DI组成，即BICI，DI。S0是重叠搬移状态数据块，其长度为M1，初始值为M1个0。重叠搬移的结果是EE0，E1，EK1和S10，其中E0S0，B0，对于I0的数据块EIDI1，BI，S10DK1。而在本具体实例中由于LM11024，因此CI数据块长度为0，即DIBI，具体搬移操作参见图2。重叠搬移过程分为三个步骤0040步骤201启动1024个线程完成将状态数据块S0搬移到E0中，即将S0M搬移到EM，其中0M1023，以及将D1023搬移到S10中。在此步骤中，线程根据下列搬移公式进。
21、行数据搬移0041线程TI需要搬移EIS0I，然后需要搬移S0IX220210I，其中0I1023。例如线程T0将S00搬移到X0，将X220210搬移到S00；0042步骤202启动1024个线程完成将待滤波数据X搬移到E的奇数块中。在此步骤中，线程根据下列搬移公式进行数据搬移0043线程TI需要搬移E2M1210IXM210I，其中0I1023，M0，1，2101。例如线程T0在M0时，将数据X0搬移到E1024，在M1时，数据X1024搬移到E3072；0044步骤203启动1024个线程完成待滤波数据X搬移到E的偶数块中。在此步骤中，线程根据下列搬移公式进行数据搬移0045线程TI需要。
22、搬移E2M1210IXM210I，其中0I1023，M0，1，2101。例如线程T0在M0时，将数据X0搬移到E2048，在M1时，数据X1024搬移到E4096；00463傅里叶变换核函数0047其处理是将经过重叠搬移后的数据，进行傅里叶变换计算，其中傅里叶变换的长度为2112048。00484乘法核函数0049其处理是将傅里叶变换的结果数据，与FIR频响系数相乘。在此过程中，启动2048个线程，并根据下列公式进行计算0050线程TI计算PM211IEM211IHI/211，其中，0I2047，M0，1，2101，E为重叠搬移操作的结果数据。例如线程T0在M0时，计算P0E0H0/211，在。
23、M1时，计算P211E211H0/211；00515傅里叶逆变换核函数0052其处理是将与频响系数相乘所得的结果数据，进行傅里叶逆变换计算，其中傅里叶逆变换的长度为2112048。00536合并搬移核函数0054其处理是将傅里叶逆变换的结果数据，进行合并搬移。图3为通用的合并搬移示意图，将傅里叶逆变换的结果数据RR0，R1，RK1合并搬移为YY0，Y1，YK1，其中RIZI，YI，其中ZI长度为M1的连续数据点，YI则是长度为L的连续数据点。而在本具体实例中，LM12101024，具体搬移操作参见图4。在此过程中，启动1024个线程，并根据下列公式进行数据搬移操作0055线程TI需要搬移YM2。
24、10IR2M1210I，其中，0I1023，M0，1，2101，Y为搬移后的结果数据。例如线程T0在M0时，将数据R210搬移说明书CN102340296ACN102340310A5/5页8到Y0，在M1时，将数据R3072搬移到Y1024；0056上面六个核函数是条件，有了上述的六个可在GPU上并行执行的核函数之后，我们就可以按照下述步骤来完成FIR的并行处理，其中具体的滤波计算处理操作参见图50057步骤501将1025点的FIR冲击响应数据从主机内存拷贝到GPU内存中，启动频响计算核函数，完成FIR频响系数的计算。0058步骤502将一块2201M点的数据从主机内存拷贝到GPU内存中，在。
25、GPU上启动210个线程，用以执行重叠搬移核函数。0059步骤503执行傅里叶变换核函数。0060步骤504在GPU上启动211个线程，用以执行乘法核函数。0061步骤505执行傅里叶逆变换核函数。0062步骤506在GPU上启动210个线程，用以执行合并搬移核函数。0063步骤507将合并搬移的结果从GPU内存中拷贝到主机内存中，如果有待处理的数据，则转到步骤二执行。说明书CN102340296ACN102340310A1/3页9图1说明书附图CN102340296ACN102340310A2/3页10图2图3图4说明书附图CN102340296ACN102340310A3/3页11图5说明书附图CN102340296A。