Q15格式全四象限定点正余弦近似实现方案咨询

首先对齐基础参数约定：

• 输入theta为Q15格式，取值范围-1~+1对应角度-π~+π（即-180°~+180°），覆盖完整圆周4个象限
• 输出sin、cos均为Q15格式，取值范围-1~+1，满幅0x7FFF对应数值1，0x8000对应数值-1
• 核心逻辑可通过修改定点位宽、移位参数快速适配其他Q格式

核心实现思路 #

全象限覆盖的核心是利用三角函数的对称性做象限折叠，只需要实现0~π/2（第一象限）的三角函数近似，再通过符号调整、函数映射覆盖其余三个象限，不需要存储或计算完整周期的数值，能大幅降低资源占用。

步骤1：象限折叠与归一化 #

首先对输入theta做预处理，提取象限标志，把任意输入角度折叠到0~π/2的计算区间，映射规则如下：

四个象限映射规则（折叠后的归一化值x范围为0~π/2）：

• 第一象限（theta∈[0, π/2)，对应Q15输入值[0x0000, 0x4000)）：sin符号为正，cos符号为正，直接取x=theta
• 第二象限（theta∈[π/2, π)，对应Q15输入值[0x4000, 0x7FFF]）：sin符号为正，cos符号为负，取x=π - theta（对应Q15计算为0x8000 - theta）
• 第三象限（theta∈[-π, -π/2)，对应Q15输入值[0x8000, 0xC000)）：sin符号为负，cos符号为负，取x=theta + π（对应Q15计算为theta - 0x8000）
• 第四象限（theta∈[-π/2, 0)，对应Q15输入值[0xC000, 0xFFFF]）：sin符号为负，cos符号为正，取x=-theta

步骤2：0~π/2区间的三角函数近似 #

两种工程上验证过的可行方案，按需选择即可：

方案1：查表+线性插值（速度最快，资源占用均衡，工业场景首选） #

• 只需要预存第一象限的sin值查找表，不需要单独存cos表（cos(x)=sin(π/2 -x)，用同一套查表逻辑即可）
• 表长可按需配置：256点长度搭配线性插值时，全象限误差小于1个Q15 LSB（即误差<1/32768），完全满足电机控制、信号生成、音频处理等绝大多数场景需求；如果资源极度紧张，缩到64点加插值也能做到1%以内的误差
• 插值计算全为定点运算：取折叠后x的高位作为查表索引，低位作为插值系数，用相邻两个表项做线性加权，计算过程用32位中间变量缓存避免溢出

方案2：Minimax多项式近似（无查找表，适合存储资源极紧的场景） #

• 不要用泰勒展开，同阶数下Minimax多项式的最大误差比泰勒展开小一个数量级以上，0~π/2区间用5阶Minimax多项式即可做到Q15精度下误差小于1个LSB
• 适配Q15格式的近似公式（x为0_{π/2区间归一化到0}1的Q15值）：
  sin(x) ≈ 1.57031*x - 0.646*x³ + 0.0791*x⁵
• 所有系数提前转为Q15定点数，乘法后做对应移位修正，cos值直接将π/2 -x代入同一公式计算即可，不需要额外的系数集

步骤3：定点异常处理 #

定点实现必须做边界防护，避免出现溢出跳变：

• 所有乘加运算的中间结果统一用32位整型存储，计算完成后做饱和截断：数值超过0x7FFF钳位到0x7FFF，低于-0x8000钳位到-0x8000，禁止出现补码溢出反转
• 初始化时校验三个边界点：theta=0时sin=0、cos=32767；theta=0.5（π/2）时sin=32767、cos=0；theta=-1（-π）时sin=0、cos=-32768，确保跨象限边界无跳变
• 格式适配：如果需要切换到Q31等其他定点格式，只需要修改移位位数、饱和阈值、查找表位宽即可，核心的象限折叠逻辑完全复用

参考核心实现代码（Q15格式，查表插值方案） #

#include <stdint.h>
// Q15饱和截断宏
#define Q15_SAT(v)  ((v) > 0x7FFF ? 0x7FFF : ((v) < (int16_t)0x8000 ? (int16_t)0x8000 : (int16_t)(v)))

// 256点第一象限sin查找表，Q15格式，实际使用时填充0~pi/2的等间隔采样值，最后一个点补0x7FFF
static const int16_t sin_lut[257] = {0, 201, 402, /* 省略中间采样点 */, 32767};

/**
 * @brief Q15格式同步计算正弦余弦值
 * @param theta_q15 输入角度，范围-1~+1对应-pi~+pi
 * @param sin_out 正弦输出指针，范围-1~+1
 * @param cos_out 余弦输出指针，范围-1~+1
 */
void q15_sincos(int16_t theta_q15, int16_t *sin_out, int16_t *cos_out)
{
    uint16_t abs_theta, x, x_cos;
    int sin_sign = 1, cos_sign = 1;
    uint16_t idx, idx_cos, frac, frac_cos;
    int32_t sin_val, cos_val;

    // 象限判断与折叠
    if (theta_q15 < 0) {
        abs_theta = (uint16_t)(-theta_q15);
        if (abs_theta >= 0x4000) {
            // 第三象限
            x = abs_theta - 0x4000;
            sin_sign = -1;
            cos_sign = -1;
        } else {
            // 第四象限
            x = abs_theta;
            sin_sign = -1;
        }
    } else {
        abs_theta = (uint16_t)theta_q15;
        if (abs_theta >= 0x4000) {
            // 第二象限
            x = 0x8000 - abs_theta;
            cos_sign = -1;
        } else {
            // 第一象限
            x = abs_theta;
        }
    }

    // 计算sin的插值索引与系数，256点对应x右移7位取高8位索引，低7位为插值权重
    idx = x >> 7;
    frac = x & 0x7F;
    sin_val = ((int32_t)sin_lut[idx] * (128 - frac) + (int32_t)sin_lut[idx+1] * frac) >> 7;

    // 计算cos值：用pi/2 -x做同样插值，比三角恒等式计算精度更高
    x_cos = 0x4000 - x;
    idx_cos = x_cos >> 7;
    frac_cos = x_cos & 0x7F;
    cos_val = ((int32_t)sin_lut[idx_cos] * (128 - frac_cos) + (int32_t)sin_lut[idx_cos+1] * frac_cos) >> 7;

    // 符号处理+饱和输出
    *sin_out = Q15_SAT(sin_sign * sin_val);
    *cos_out = Q15_SAT(cos_sign * cos_val);
}

内容的提问来源于stack exchange，提问作者gari