整数形式算术编码参数调整及频率转换问题求解

问题描述 #

基于指定PDF实现了整数形式的算术编码模型以实现无限精度，参考GitHub上的代码实现时遇到以下问题：

• 计算得到的概率是double类型，但API要求传入uint64_t类型的频率计数，目前采用将概率乘以接近1e9的基数转换为整数的方式
• 基数过小会丢失小概率值的精度（比如示例中的pMin=0.000000000913767），基数过大则会导致溢出
• 将频率类型从uint32_t改为uint64_t后，频率总和超过了模型的maximumTotal上限，调整numBits参数也无法解决
• 当前临时方案是使用1e9基数并分配最小概率区间，寻求更优雅的解决方案

相关代码片段：

概率转频率代码 #

double p = calc_prob(i, lap);
uint64_t f = (uint64_t)( p * base);

模型中total和maximumTotal定义 #

private: std::uint64_t total;
protected: std::uint64_t maximumTotal;

numBits与maximumTotal关联代码 #

ArithmeticCoderBase::ArithmeticCoderBase(int numBits) {
    if (numBits < 1 || numBits > 63)
        throw std::domain_error("State size out of range");
    numStateBits = numBits;
    fullRange = static_cast<decltype(fullRange)>(1) << numStateBits;
    halfRange = fullRange >> 1;  // Non-zero
    quarterRange = halfRange >> 1;  // Can be zero
    minimumRange = quarterRange + 2;  // At least 2
    maximumTotal = std::min(std::numeric_limits<decltype(fullRange)>::max() / fullRange, minimumRange);
    stateMask = fullRange - 1;
    low = 0;
    high = stateMask;
}

解决方案 #

1. 动态基数+归一化缩放 #

放弃固定基数，根据所有概率的分布动态计算合适的基数，确保总和不超过maximumTotal：

// 先收集所有概率
std::vector<double> probs;
for (int i = 0; i < symbolCount; ++i) {
    probs.push_back(calc_prob(i, lap));
}

// 计算最大可能的基数，使得所有概率*基数的总和 <= maximumTotal
double sumProb = std::accumulate(probs.begin(), probs.end(), 0.0);
uint64_t maxBase = static_cast<uint64_t>(maximumTotal / sumProb);
maxBase = std::max(maxBase, 1ULL); // 避免基数为0

// 转换为频率并处理小概率
std::vector<uint64_t> freqs;
uint64_t totalFreq = 0;
for (double p : probs) {
    uint64_t f = static_cast<uint64_t>(p * maxBase);
    f = std::max(f, 1ULL); // 保证小概率至少有1的频率
    freqs.push_back(f);
    totalFreq += f;
}

// 修正总和溢出问题，从高频符号中扣减溢出部分
if (totalFreq > maximumTotal) {
    uint64_t overflow = totalFreq - maximumTotal;
    while (overflow > 0) {
        int maxIdx = std::distance(probs.begin(), std::max_element(probs.begin(), probs.end()));
        if (freqs[maxIdx] > 1) {
            freqs[maxIdx]--;
            overflow--;
        }
        probs[maxIdx] = -1; // 标记为已调整，避免重复选择
    }
}

2. 最小概率映射优化 #

针对极小概率直接映射为最小频率1，同时将剩余额度分配给高频符号，保证总和不超限：

const uint64_t minFreq = 1;
uint64_t remainingTotal = maximumTotal - symbolCount * minFreq;
std::vector<uint64_t> freqs;
double sumAdjustedProb = 0.0;

// 计算扣除最小概率后的剩余概率总和
for (double p : probs) {
    sumAdjustedProb += std::max(p - (minFreq * 1.0 / maximumTotal), 0.0);
}

// 分配频率
for (double p : probs) {
    if (p < (minFreq * 1.0 / maximumTotal)) {
        freqs.push_back(minFreq);
    } else {
        double adjustedP = p - (minFreq * 1.0 / maximumTotal);
        uint64_t f = static_cast<uint64_t>((adjustedP / sumAdjustedProb) * remainingTotal);
        freqs.push_back(f + minFreq);
    }
}

// 修正总和不足的情况，补到最大概率符号上
uint64_t totalFreq = std::accumulate(freqs.begin(), freqs.end(), 0ULL);
if (totalFreq < maximumTotal) {
    int maxIdx = std::distance(probs.begin(), std::max_element(probs.begin(), probs.end()));
    freqs[maxIdx] += (maximumTotal - totalFreq);
}

3. 调整maximumTotal计算逻辑 #

当前maximumTotal被minimumRange限制了上限，若编码逻辑允许更大的total值，可修改这部分逻辑：

// 原代码：
// maximumTotal = std::min(std::numeric_limits<decltype(fullRange)>::max() / fullRange, minimumRange);
// 修改为：
uint64_t maxSafeTotal = std::numeric_limits<decltype(fullRange)>::max() / fullRange;
// 去掉minimumRange限制，或保留少量余量避免极端溢出
maximumTotal = maxSafeTotal - 1000;

注意：修改后需确保编码/解码的区间更新逻辑能适配更大的total值，避免区间计算溢出。

4. 高精度概率转换 #

避免double的精度损失，直接用分数形式表示概率再转换为频率：

// 假设calc_prob返回分数形式（分子n，分母d）
std::pair<uint64_t, uint64_t> get_prob_fraction(int i, lap);

uint64_t maxScale = maximumTotal / get_prob_fraction(i, lap).second;
uint64_t f = get_prob_fraction(i, lap).first * maxScale;
f = std::max(f, 1ULL); // 保证最小频率

若无法直接获取分数形式，可使用高精度计算库（如GMP）处理概率的精确值，再转换为整数频率。

内容的提问来源于stack exchange，提问作者stackdon'tknow