STM32 FFT之后计算模值比较慢的问题

xfdr0805 · 发表于 2017-3-8 22:28:01

本帖最后由 xfdr0805 于 2017-3-8 22:29 编辑

STM32F103c8t6用DSP库进行 256点 FFT运算，仿真发现FFT运算倒挺快的，跟官方速度差不多，但是计算幅值的时间，时间特别长，有差不多40多ms，造成3264点阵刷新比较慢，速度差不多20帧吧，反正也过的去，这里想请教一下大家，这种情况求模运算比FFT运算慢那么多是怎么回事？有没有其它办法提高速度呢，请做过这方面的朋友给指点一下，谢谢了

void dsp_asm_powerMag(void)
{
  s16 lX,lY;
  u32 i;
  for(i=0;i<NPT/2;i++)
  {
lX  = (lBUFOUT << 16) >> 16;
lY  = (lBUFOUT >> 16);
{
float X = NPT * ((float)lX) /32768;
float Y = NPT * ((float)lY) /32768;
float Mag = sqrt(X*X + Y*Y)/NPT;
            if(i==0)
            {
                     lBUFMAG = (u32)(Mag* 32768*31/256);
            }
            else
            {
                     lBUFMAG = (u32)(Mag* 65536*31/150);
            }
}
  }
}

yijinxiaoyou · 发表于 2017-3-8 22:28:02

xfdr0805 发表于 2017-3-9 19:27
等会下班了验证一下,最后一个计算是正确的，但还是慢，还没原来的快呢

/*****************************************************************************************
平方根
******************************************************************************************/
float mySqrt(float x)
{
float a = x;
unsigned int i = *(unsigned int *)&x;
i = (i + 0x3f76cf62) >> 1;
x = *(float *)&i;
x = (x + a / x) * 0.5;
return x;
}
不会啊快的不只一点

xfdr0805 · 发表于 2017-3-9 07:23:58

这段函数能不能查表替换掉？

wmaxwell · 发表于 2017-3-9 08:41:12

少用除法，比如除以2，改为乘以0.5

yijinxiaoyou · 发表于 2017-3-9 08:55:34

本帖最后由 yijinxiaoyou 于 2017-3-9 09:00 编辑

1、F1没有浮点核用q15整点算吧
2、除法2的幂用移位
3、求根网上一大堆快速算法
4、dsp库的频谱泄露真的很严重用MATLAB对比下就知道了

yijinxiaoyou · 发表于 2017-3-9 09:02:38

1024点360us 103跑起来没问题的

ieeloleei · 发表于 2017-3-9 13:56:29

yijinxiaoyou 发表于 2017-3-9 08:55
1、F1没有浮点核用q15整点算吧
2、除法2的幂用移位
3、求根网上一大堆快速算法

第四点不赞同，按照你的说法，“频谱泄露”是ARM DSP库算法造成的，让人觉得这是ARM DSP库的缺陷，然而并不是这样。请仔细理解“频谱泄露”的含义，仔细查看FFT函数的说明。
---
无论是用labview还是matlab验证，计算结果都是一样的，这是我的经验。除非你使用的不对。

yijinxiaoyou · 发表于 2017-3-9 16:18:34

本帖最后由 yijinxiaoyou 于 2017-3-9 16:33 编辑

ieeloleei 发表于 2017-3-9 13:56
第四点不赞同，按照你的说法，“频谱泄露”是ARM DSP库算法造成的，让人觉得这是ARM DSP库的缺陷，然而并 ...

不知道你测的是什么信号我测发时候同样的数组经过dsp库与MATLAB比较高次有泄露后来就没深究后来用自己的算法没有什么问题

yuzeyuan1 · 发表于 2017-3-9 17:18:26

我有一大堆矩阵运算，前阵子算了下计算速度，有些太慢了，就分别试用了有单精度浮点运算单元的片子和有有单精度浮点运算单元的片子F767。发现用单精度片子时，库里自带的函数确实比手敲的快太多了。
感觉原因是这样的：float型乘float型，理论上需要用double型进行存储，所以可能是用float做乘除法时它自动转成了double型进行计算，就没法用FPU了，所以速度会慢，比没有FPU的片子略快一些，但跟库函数这种直接用FPU进行运算的东西比还是慢了很多。

xfdr0805 · 发表于 2017-3-9 19:22:54

根据大家回复，找到了一个快速开根号算法

xfdr0805 · 发表于 2017-3-9 19:23:48

xfdr0805 发表于 2017-3-9 19:22
根据大家回复，找到了一个快速开根号算法

这里给出实现32位无符号整数开方得到16位无符号整数的C语言代码。

//------------------------------------------------------------------------------- -
/****************************************/
/*Function: 开根号处理 */
/*入口参数：被开方数，长整型 */ /*出口参数：开方结果，整型 */
/****************************************/
unsigned int sqrt_16(unsigned long M)
{
unsigned int N, i;
unsigned long tmp, ttp; // 结果、循环计数
if (M == 0) // 被开方数，开方结果也为0
return 0;
N = 0;
tmp = (M >> 30); // 获取最高位：B[m-1]
M <<= 2;
if (tmp > 1) // 最高位为1
{
N ++; // 结果当前位为1，否则为默认的0
tmp -= N;
}
for (i=15; i>0; i--) // 求剩余的15位
{
N <<= 1; // 左移一位
tmp <<= 2;
tmp += (M >> 30); // 假设
ttp = N;
ttp = (ttp<<1)+1;
M <<= 2;
if (tmp >= ttp) // 假设成立
{
tmp -= ttp; N ++;
}
}
return N;
}

xfdr0805 · 发表于 2017-3-9 19:25:12

xfdr0805 发表于 2017-3-9 19:23
这里给出实现32位无符号整数开方得到16位无符号整数的C语言代码。

/*还有一个*/
/* 来至 Quake 3 的源码 */
float CarmSqrt(float x){
union
{
      int intPart;
      float floatPart;
} convertor;
union
{
      int intPart;
      float floatPart;
} convertor2;
convertor.floatPart = x;
convertor2.floatPart = x;
convertor.intPart = 0x1FBCF800 + (convertor.intPart >> 1);
convertor2.intPart = 0x5f3759df - (convertor2.intPart >> 1);
return 0.5f*(convertor.floatPart + (x * convertor2.floatPart));
}

xfdr0805 · 发表于 2017-3-9 19:26:45

本帖最后由 xfdr0805 于 2017-3-9 21:09 编辑

xfdr0805 发表于 2017-3-9 19:22
根据大家回复，找到了一个快速开根号算法

  还有一个带神秘常数的
float Q_rsqrt( float number )
  {
   long i;
   float x2, y;
   const float threehalfs = 1.5F;

   x2 = number * 0.5F;
   y    = number;
   i    = * ( long * ) &y;    // evil floating point bit level hacking
   i    = 0x5f3759df - ( i >> 1 ); // what the func?
   y    = * ( float * ) &i;
   y    = y * ( threehalfs - ( x2 * y * y ) ); // 1st iteration
   // y    = y * ( threehalfs - ( x2 * y * y ) ); // 2nd iteration, this can be removed

   return y;
  }