dsi_bitstream/codes/rice.rs
1/*
2 * SPDX-FileCopyrightText: 2023 Sebastiano Vigna
3 *
4 * SPDX-License-Identifier: Apache-2.0 OR LGPL-2.1-or-later
5 */
6
7//! Rice codes.
8//!
9//! Rice codes (AKA Golomb−Rice codes) are a form of approximated [Golomb
10//! codes](crate::codes::golomb) in which the parameter *b* is a power of
11//! two. This restriction makes the code less precise in modeling data with a
12//! geometric distribution, but encoding and decoding can be performed without
13//! any integer arithmetic, and thus much more quickly.
14//!
15//! The implied distribution of a Rice code is [the same as that of a Golomb
16//! code with the same parameter](crate::codes::golomb).
17//!
18//! For natural numbers distributed with a geometric distribution with base *p*,
19//! the base-2 logarithm of the optimal *b* is [⌈log₂(ln((√5 + 1)/2) / ln(1 -
20//! *p*))⌉](log2_b).
21//!
22//! The supported range is [0 . . 2⁶⁴ – 1) for log₂(*b*) in [0 . . 64).
23//!
24//! # References
25//!
26//! Robert F. Rice, “[Some practical universal noiseless coding
27//! techniques](https://ntrs.nasa.gov/api/citations/19790014634/downloads/19790014634.pdf)”.
28//! Jet Propulsion Laboratory, Pasadena, CA, Tech. Rep. JPL-79-22, JPL-83-17,
29//! and JPL-91-3, March 1979.
30//!
31//! Aaron Kiely. “[Selecting the Golomb parameter in Rice
32//! coding](https://tda.jpl.nasa.gov/progress_report/42-159/159E.pdf)”.
33//! Interplanetary Network Progress report 42-159, Jet Propulsion Laboratory,
34//! 2004.
35
36use crate::traits::*;
37
38/// Returns the length of the Rice code for `n` with parameter `log2_b`.
39#[must_use]
40#[inline(always)]
41pub fn len_rice(n: u64, log2_b: usize) -> usize {
42 debug_assert!(log2_b < 64);
43 (n >> log2_b) as usize + 1 + log2_b
44}
45
46/// Returns the optimal value of log₂*b* for a geometric distribution of base
47/// *p*, that is, ⌈log₂(ln((√5 + 1)/2) / ln(1 - *p*))⌉
48pub fn log2_b(p: f64) -> usize {
49 ((-((5f64.sqrt() + 1.0) / 2.0).ln() / (-p).ln_1p()).log2()).ceil() as usize
50}
51
52/// Trait for reading Rice codes.
53pub trait RiceRead<E: Endianness>: BitRead<E> {
54 #[inline(always)]
55 fn read_rice(&mut self, log2_b: usize) -> Result<u64, Self::Error> {
56 debug_assert!(log2_b < 64);
57 Ok((self.read_unary()? << log2_b) + self.read_bits(log2_b)?)
58 }
59}
60
61/// Trait for writing Rice codes.
62pub trait RiceWrite<E: Endianness>: BitWrite<E> {
63 #[inline(always)]
64 fn write_rice(&mut self, n: u64, log2_b: usize) -> Result<usize, Self::Error> {
65 debug_assert!(log2_b < 64);
66
67 let mut written_bits = self.write_unary(n >> log2_b)?;
68 #[cfg(feature = "checks")]
69 {
70 // Clean up n in case checks are enabled
71 let n = n & (1_u128 << log2_b).wrapping_sub(1) as u64;
72 written_bits += self.write_bits(n, log2_b)?;
73 }
74 #[cfg(not(feature = "checks"))]
75 {
76 written_bits += self.write_bits(n, log2_b)?;
77 }
78 Ok(written_bits)
79 }
80}
81
82impl<E: Endianness, B: BitRead<E>> RiceRead<E> for B {}
83impl<E: Endianness, B: BitWrite<E>> RiceWrite<E> for B {}
84
85#[cfg(test)]
86#[test]
87fn test_log2_b() {
88 use crate::prelude::golomb::b;
89
90 let mut p = 1.0;
91 for _ in 0..100 {
92 p *= 0.9;
93 let golomb = b(p);
94 if golomb & -(golomb as i64) as u64 == golomb {
95 assert_eq!(golomb, 1 << log2_b(p));
96 }
97 }
98}