跳转至

16 fft-快速傅里叶变换

孟伟

一、概述

fft 模块提供高性能的快速傅里叶变换(FFT)和逆快速傅里叶变换(IFFT)功能,支持 float32 和 q15 定点两种内核,适用于信号处理、频谱分析等场景。

二、核心示例

1、核心示例是指:使用本库文件提供的核心 API,开发的基础业务逻辑的演示代码;

2、核心示例的作用是:帮助开发者快速理解如何使用本库,所以核心示例的逻辑都比较简单;

3、更加完整和详细的 demo,请参考 LuatOS 仓库 中各个产品目录下的 demo/fft

function test_fft_fun()
    if not fft then
        -- 如果不支持FFT库,进入死循环并提示错误
        while 1 do
            sys.wait(1000) -- 等待1秒
            log.info("bsp", "此BSP不支持fft库,请检查") -- 输出错误信息
        end
    end

    -- FFT 参数配置段

    -- 设置FFT基本参数
    local N = 2048 -- FFT 点数:决定频率分辨率和计算复杂度,必须是2的幂次方
    local fs = 2000 -- 采样频率 (Hz):根据奈奎斯特定理,可分析的最高频率为fs/2=1000Hz
    local freq = 200 -- 测试信号频率 (Hz):生成一个200Hz的正弦波作为测试信号

    -- 输出测试开始信息和参数
    log.info("fft", "q15 测试开始", "N=" .. N, "fs=" .. fs, "freq=" .. freq)


    -- 内存分配段:为FFT计算准备缓冲区
    -- 分配 zbuff 缓冲区:
    -- 使用16位整数缓冲区来存储Q15格式的数据
    -- 每个复数点需要2个16位整数(实部和虚部),所以总大小为 N * 2 * 2 字节
    local real_i16 = zbuff.create(N * 2) -- 实部缓冲区:存储时域信号的实部
    local imag_i16 = zbuff.create(N * 2) -- 虚部缓冲区:存储时域信号的虚部(初始为0)


    -- 测试信号生成段:生成一个200Hz的正弦波测试信号
    -- 生成 U12 整数正弦波并写入缓冲区(避免浮点预处理干扰)
    -- U12格式:12位无符号整数,范围0-4095,2048对应0电平
    for i = 0, N - 1 do
        -- 计算时间点:i/fs 表示第i个样本的时间(秒)
        local t = i / fs

        -- 生成200Hz正弦波:sin(2π * 频率 * 时间)
        local x = math.sin(2 * math.pi * freq * t)

        -- 将浮点数转换为U12格式:
        -- 2048为直流偏置(中间值),2047为幅度范围
        -- +0.5是为了四舍五入到最接近的整数
        local val = math.floor(2048 + 2047 * x + 0.5)

        -- 数值范围限制:确保在U12的有效范围内(0-4095)
        if val < 0 then val = 0 end
        if val > 4095 then val = 4095 end

        -- 将数据写入缓冲区:
        -- seek定位到正确的位置,每个样本占2字节
        -- writeU16写入16位无符号整数(低12位有效)
        real_i16:seek(i * 2, zbuff.SEEK_SET); real_i16:writeU16(val) -- 实部写入正弦波数据
        imag_i16:seek(i * 2, zbuff.SEEK_SET); imag_i16:writeU16(0) -- 虚部写入0(实数信号)
    end


    -- 旋转因子生成段:为Q15 FFT计算准备旋转因子
    -- 生成 Q15 旋转因子到 zbuff(避免任何浮点旋转因子)
    -- 旋转因子用于FFT计算中的复数乘法,长度是FFT点数的一半
    local Wc_q15 = zbuff.create((N // 2) * 2) -- 余弦旋转因子缓冲区
    local Ws_q15 = zbuff.create((N // 2) * 2) -- 正弦旋转因子缓冲区(实际存储-sin值)

    -- 生成Q15格式的旋转因子表
    fft.generate_twiddles_q15_to_zbuff(N, Wc_q15, Ws_q15)


    -- 使用定点Q15内核执行FFT
    -- 执行 Q15 FFT(输入为 U12),显式传入 Q15 旋转因子
    local t0 = mcu.ticks() -- 记录开始时间

    -- 执行FFT计算:
    -- real_i16, imag_i16: 输入输出缓冲区(原地计算)
    -- N: FFT点数
    -- Wc_q15, Ws_q15: Q15格式旋转因子
    -- {core = "q15", input_format = "u12"}: 使用Q15定点内核,输入格式为U12
    fft.run(real_i16, imag_i16, N, Wc_q15, Ws_q15, { core = "q15", input_format = "u12" })

    local t1 = mcu.ticks() -- 记录结束时间
    -- 输出Q15 FFT计算耗时
    log.info("fft", "q15 FFT 完成", "耗时:" .. (t1 - t0) .. "ms")


    -- 分析FFT结果,找到主要频率成分
    -- 扫描前半部分频谱(0 ~ fs/2),寻找主峰
    -- 由于频谱的对称性,只需要分析前N/2个点
    local peak_k, peak_pow = 1, -1 -- peak_k: 峰值位置, peak_pow: 峰值功率

    -- 遍历所有频率bin(跳过直流分量k=0)
    for k = 1, (N // 2) - 1 do
        -- 定位到第k个频率点的实部和虚部
        real_i16:seek(k * 2, zbuff.SEEK_SET)
        imag_i16:seek(k * 2, zbuff.SEEK_SET)

        -- 读取实部和虚部值(Q15格式)
        local rr = real_i16:readI16() -- 实部
        local ii = imag_i16:readI16() -- 虚部

        -- 计算功率谱:|X[k]|² = real² + imag²
        -- 功率谱表示该频率成分的能量大小
        local p = rr * rr + ii * ii

        -- 更新峰值信息
        if p > peak_pow then
            peak_pow = p -- 更新最大功率值
            peak_k = k -- 更新峰值位置
        end
    end

    -- 计算峰值对应的实际频率:频率 = bin索引 * 频率分辨率
    -- 频率分辨率 = 采样频率 / FFT点数
    local peak_freq = (peak_k) * fs / N

    -- 输出主峰频率信息
    -- peak_k: 峰值所在的bin索引
    -- peak_freq: 计算出的实际频率(应该接近200Hz)
    log.info("fft", "主峰(Hz/bin)", string.format("%.2f", peak_freq), peak_k)


    -- 使用浮点F32内核进行相同计算,对比性能
    -- 比较:使用 f32 内核处理相同输入(验证 q15 与 f32 结果一致性)
    -- 复制相同的 U12 输入到新的 zbuff(浮点格式)

    -- 为浮点FFT分配缓冲区:每个浮点数占4字节
    local real_f32 = zbuff.create(N * 4) -- 实部缓冲区(浮点)
    local imag_f32 = zbuff.create(N * 4) -- 虚部缓冲区(浮点)

    -- 将Q15数据复制到浮点缓冲区
    for i = 0, N - 1 do
        -- 从Q15缓冲区读取U12数据
        real_i16:seek(i * 2, zbuff.SEEK_SET)
        local val = real_i16:readU16()

        -- 写入浮点缓冲区(直接写入U12原始值,不进行格式转换)
        real_f32:seek(i * 4, zbuff.SEEK_SET)
        imag_f32:seek(i * 4, zbuff.SEEK_SET)
        real_f32:writeF32(val) -- 实部写入原始U12值
        imag_f32:writeF32(0.0) -- 虚部写入0.0
    end

    -- 生成 f32 旋转因子(浮点格式)
    -- 返回Lua table格式的旋转因子,而不是zbuff
    local Wc, Ws = fft.generate_twiddles(N)

    -- 执行 f32 FFT(输入为 U12)
    local t2 = mcu.ticks() -- 记录开始时间

    -- 执行浮点FFT计算:
    -- 使用默认的f32内核,输入格式为u12
    fft.run(real_f32, imag_f32, N, Wc, Ws, { input_format = "u12" })

    local t3 = mcu.ticks() -- 记录结束时间
    local dt_f32 = (t3 - t2) -- 计算浮点FFT耗时

    -- 输出浮点FFT计算耗时
    log.info("fft", "f32 FFT 完成", "耗时:" .. dt_f32 .. "ms")


    -- 结果总结段:对比两种实现的性能


    -- 总结耗时对比:Q15定点 vs F32浮点
    -- 在没有硬件浮点加速的嵌入式设备上,Q15会比F32快很多,
    -- 而780和8000系列都没有硬件浮点加速,建议使用q15计算提高速度,但如果追求计算精度,仍然可以用浮点计算
    log.info("fft", "对比(q15 vs f32, ms)", string.format("%d / %d", (t1 - t0), dt_f32))
end

sys.taskInit(test_fft_fun)

三、常量详解

核心库常量,顾名思义是由合宙 LuatOS 内核固件中定义的、不可重新赋值或修改的固定值,在脚本代码中不需要声明,可直接调用;

每个常量对应的常量取值仅做日志打印时查询使用,不要将这个常量取值用做具体的业务逻辑判断,因为LuatOS内核固件可能会变更每个常量对应的常量取值;

如果用做具体的业务逻辑判断,一旦常量取值发生改变,业务逻辑就会出错;

fft 模块没有常量。

四、函数详解

fft.generate_twiddles(N)

功能

生成 float32 旋转因子表。

参数

N

参数含义:FFT 点数;
数据类型:number
取值范围:必须为 2 的幂次方,且 N < 65536
是否必选:必须传入此参数;
注意事项:推荐范围为 N  16384(已验证稳定运行);

返回值

local Wc, Ws = fft.generate_twiddles(N)

Wc

含义说明:cos 旋转因子表;
数据类型:table
取值范围:长度为 N/2  Lua 数组;
注意事项:返回的是 float32 精度的旋转因子;

Ws

含义说明:-sin 旋转因子表;
数据类型:table
取值范围:长度为 N/2  Lua 数组;
注意事项:返回的是 float32 精度的旋转因子;

示例

local N = 2048
local Wc, Ws = fft.generate_twiddles(N)

fft.generate_twiddles_q15_to_zbuff(N, Wc_zb, Ws_zb)

功能

生成 q15 定点旋转因子到 zbuff(零浮点)。

参数

N

参数含义:FFT 点数;
数据类型:number
取值范围:必须为 2 的幂次方,且 N < 65536
是否必选:必须传入此参数;

Wc_zb

参数含义:输出缓冲(cos旋转因子),用于存放 int16 Q15 格式的 cos 旋转因子;
数据类型:zbuff
取值范围:长度至少为 (N/2)*2 字节;
是否必选:必须传入此参数;

Ws_zb

参数含义:输出缓冲(-sin旋转因子),用于存放 int16 Q15 格式的 -sin 旋转因子(前向FFT用);
数据类型:zbuff
取值范围:长度至少为 (N/2)*2 字节;
是否必选:必须传入此参数;

返回值

无返回值

示例

local N = 2048
local Wc_q15 = zbuff.create((N//2)*2)
local Ws_q15 = zbuff.create((N//2)*2)
fft.generate_twiddles_q15_to_zbuff(N, Wc_q15, Ws_q15)

fft.run(real, imag, N, Wc, Ws[, opts])

功能

原地 FFT 计算。

参数

real

参数含义:实部容器;
数据类型:table  zbuff
取值范围:float32 路径:Lua 数组或 zbuff(float32)
         q15 路径:zbuff(int16)。当 opts.core="q15"  opts.input_format  "u12"/"u16"/"s16" 时生效;
是否必选:必须传入此参数;

imag

参数含义:虚部容器;
数据类型:table  zbuff
取值范围:float32 路径:Lua 数组或 zbuff(float32)
         q15 路径:zbuff(int16)。当 opts.core="q15"  opts.input_format  "u12"/"u16"/"s16" 时生效;
是否必选:可选;

N

参数含义:FFT 点数;
数据类型:number
取值范围:必须为 2 的幂次方;
是否必选:必须传入此参数;

Wc

参数含义:旋转因子 cos
数据类型:table  zbuff
取值范围:float32 路径:Lua 数组或 zbuff(float32)
         q15 路径:zbuff(int16),推荐配合 fft.generate_twiddles_q15_to_zbuff 生成;
是否必选:必须传入此参数;

Ws

参数含义:旋转因子 -sin
数据类型:table  zbuff
取值范围:float32 路径:Lua 数组或 zbuff(float32)
         q15 路径:zbuff(int16),推荐配合 fft.generate_twiddles_q15_to_zbuff 生成;
是否必选:必须传入此参数;

opts

参数含义:可选配置表;
数据类型:table;
是否必选:可选;
参数格式:
{
参数含义:核心计算路径;
数据类型:string;
是否必选:可选;
取值范围:"f32"(默认)、"q15";
f32浮点内核,精度高(32位),计算稳定,适合精密分析;
q15定点内核,速度快(16位整数),内存省,适合实时处理但精度略低;
参数示例:"f32"
参数名称:opts.core

参数含义:输入数据格式;
数据类型:string;
是否必选:q15 路径时必填;
取值范围:"f32"、"u12"、"u16"、"s16";f32 路径时无需配置;
参数示例:"f32"
注意事项:--f32: 标准浮点输入,适用于已处理的信号数据;
--u12: 12位无符号整数(0~4095),常见于ADC采样,自动去直流分量;
--u16: 16位无符号整数(0~65535),适用于高精度ADC或预处理数据;
--s16: 16位有符号整数(-32768~32767),适用于已去直流的差分信号;
参数名称:opts.input_format
}

返回值

无返回值

示例

_-- f32 路径示例(zbuff float32)_
local N=2048
local real=zbuff.create(N*4);
local imag=zbuff.create(N*4)
local Wc,Ws=fft.generate_twiddles(N)
fft.run(real, imag, N, Wc, Ws)

_-- q15 路径示例(U12 整数输入)_
local N=2048
local real_i16=zbuff.create(N*2);
local imag_i16=zbuff.create(N*2)
local Wc_q15=zbuff.create((N//2)*2);
local Ws_q15=zbuff.create((N//2)*2)
fft.generate_twiddles_q15_to_zbuff(N, Wc_q15, Ws_q15)
_-- 写入 U12 数据到 real_i16 后:_
fft.run(real_i16, imag_i16, N, Wc_q15, Ws_q15, {core="q15", input_format="u12"})

fft.ifft(real, imag, N, Wc, Ws[, opts])

功能

原地 IFFT 计算。

注意事项

就地修改 real/imag,并在 f32 路径下包含 1/N 归一化

参数

real

参数含义:实部容器;
数据类型:table  zbuff
取值范围:float32 路径:Lua 数组或 zbuff(float32)
         q15 路径:zbuff(int16)。当 opts.core="q15"  opts.input_format  "u12"/"u16"/"s16" 时生效;
是否必选:必须传入此参数;

imag

参数含义:虚部容器;
数据类型:table  zbuff
取值范围:float32 路径:Lua 数组或 zbuff(float32)
         q15 路径:zbuff(int16)。当 opts.core="q15"  opts.input_format  "u12"/"u16"/"s16" 时生效;
是否必选:可选;

N

参数含义:FFT 点数;
数据类型:number
取值范围:必须为 2 的幂次方;
是否必选:必须传入此参数;

Wc

参数含义:旋转因子 cos
数据类型:table  zbuff
取值范围:float32 路径:Lua 数组或 zbuff(float32)
         q15 路径:zbuff(int16),推荐配合 fft.generate_twiddles_q15_to_zbuff 生成;
是否必选:必须传入此参数;

Ws

参数含义:旋转因子 -sin
数据类型:table  zbuff
取值范围:float32 路径:Lua 数组或 zbuff(float32)
         q15 路径:zbuff(int16),推荐配合 fft.generate_twiddles_q15_to_zbuff 生成;
是否必选:必须传入此参数;

opts

参数含义:可选配置表;
数据类型:table;
是否必选:可选;
参数格式:
{
参数含义:核心计算路径;
数据类型:string;
是否必选:可选;
取值范围:"f32"(默认)、"q15";
--f32浮点内核,精度高(32位),计算稳定,适合精密分析;
--q15定点内核,速度快(16位整数),内存省,适合实时处理但精度略低;
参数示例:"f32"
参数名称:opts.core

参数含义:输入数据格式;
数据类型:string;
是否必选:q15 路径时必填;
取值范围:"f32"、"u12"、"u16"、"s16";f32 路径时无需配置;
参数示例:"f32"
注意事项:--f32: 标准浮点输入,适用于已处理的信号数据;
--u12: 12位无符号整数(0~4095),常见于ADC采样,自动去直流分量;
--u16: 16位无符号整数(0~65535),适用于高精度ADC或预处理数据;
--s16: 16位有符号整数(-32768~32767),适用于已去直流的差分信号;
参数名称:opts.input_format
}

返回值

无返回值

示例

_-- 示例:使用 IFFT 重构信号_
local N=2048
local real=zbuff.create(N*4);
local imag=zbuff.create(N*4)
local Wc,Ws=fft.generate_twiddles(N)
_-- 先执行 FFT 处理..._
fft.run(real, imag, N, Wc, Ws)
_-- 处理频域数据...-- 然后执行 IFFT 重构_
fft.ifft(real, imag, N, Wc, Ws)

fft.fft_integral(real, imag, n, df)

功能

频域积分(1/(jω))。

参数

real

参数含义:real实部
数据类型:table  zbuff
取值范围:float32 精度;
是否必选:必须传入此参数;

imag

参数含义:imag虚部
数据类型:table  zbuff
取值范围:float32 精度;
是否必选:必须传入此参数;

n

参数含义:点数;
数据类型:number
取值范围:必须为 2 的幂次方;
是否必选:必须传入此参数;

df

参数含义:频率分辨率;
数据类型:number
取值范围:fs/nfs为采样率);
是否必选:必须传入此参数;

返回值

无返回值

示例

_-- 先完成 FFT 得到谱 (real, imag),再调用积分:_
fft.fft_integral(real, imag, N, fs/N)

五、产品支持说明

支持 LuatOS 开发的所有产品都支持 fft 核心库。