Function cbm_rsq_encode

src/semantic/rotsq.c:44–96

/* ── Encode ─────────────────────────────────────────────────────────── */

void cbm_rsq_encode(const float *v, cbm_rsq_code_t *out) {
    rsq_init_diag();

    float rot[CBM_RSQ_DIM];
    for (int d = 0; d < CBM_RSQ_IN_DIM; d++) {
        rot[d] = v[d] * g_rsq_diag[d];
    }
    for (int d = CBM_RSQ_IN_DIM; d < CBM_RSQ_DIM; d++) {
        rot[d] = 0.0F;
    }
    rsq_fwht(rot);
    /* Normalize the transform so the rotation is orthonormal (H/√D): keeps
     * the coordinates in a data-independent range and makes the estimated IP
     * directly comparable to the pre-rotation IP. */
    const float inv_sqrt_d = 1.0F / 32.0F; /* 1/√1024 */
    float lo = rot[0] * inv_sqrt_d;
    float hi = lo;
    for (int d = 0; d < CBM_RSQ_DIM; d++) {
        rot[d] *= inv_sqrt_d;
        if (rot[d] < lo) {
            lo = rot[d];
        }
        if (rot[d] > hi) {
            hi = rot[d];
        }
    }

    /* Per-vector scalar quantization over [lo, hi] at CBM_RSQ_BITS. */
    float range = hi - lo;
    float step = range > 0.0F ? range / (float)CBM_RSQ_LEVELS : 1.0F;
    out->offset = lo;
    out->scale = step;

    int32_t sum = 0;
    memset(out->codes, 0, sizeof(out->codes));
    for (int d = 0; d < CBM_RSQ_DIM; d++) {
        float q = (rot[d] - lo) / step;
        int32_t c = (int32_t)(q + 0.5F);
        if (c < 0) {
            c = 0;
        }
        if (c > CBM_RSQ_LEVELS) {
            c = CBM_RSQ_LEVELS;
        }
        sum += c;
        if (d & 1) {
            out->codes[d >> 1] |= (uint8_t)(c << 4);
        } else {
            out->codes[d >> 1] |= (uint8_t)c;
        }
    }
    out->code_sum = sum;
}

/* ── Estimated inner product from two codes ─────────────────────────── */