ldpc_toolbox/codes/
ccsds.rs

1//! CCSDS TM Synchronization and Channel Coding LDPC codes.
2//!
3//! This module contains the AR4JA LDPC codes and the C2 code described in the
4//! TM Synchronization and Channel Coding Blue Book.
5//!
6//! ## References
7//! \[1\] [CCSDS 131.0-B-5 TM Synchronization and Channel Coding Blue Book](https://public.ccsds.org/Pubs/131x0b5.pdf).
8
9use crate::sparse::SparseMatrix;
10use enum_iterator::Sequence;
11
12/// AR4JA code definition.
13#[derive(Copy, Clone, Debug, Eq, PartialEq, Hash)]
14pub struct AR4JACode {
15    rate: AR4JARate,
16    k: AR4JAInfoSize,
17}
18
19/// AR4JA code rate.
20#[derive(Copy, Clone, Debug, Eq, PartialEq, Hash, Sequence)]
21pub enum AR4JARate {
22    /// Rate 1/2.
23    R1_2,
24    /// Rate 2/3.
25    R2_3,
26    /// Rate 4/5.
27    R4_5,
28}
29
30/// AR4JA information block size `k`.
31#[derive(Copy, Clone, Debug, Eq, PartialEq, Hash, Sequence)]
32pub enum AR4JAInfoSize {
33    /// k = 1024
34    K1024,
35    /// k = 4096,
36    K4096,
37    /// k = 16384,
38    K16384,
39}
40
41impl AR4JACode {
42    /// Creates an AR4JA code definition.
43    pub fn new(rate: AR4JARate, information_block_size: AR4JAInfoSize) -> AR4JACode {
44        AR4JACode {
45            rate,
46            k: information_block_size,
47        }
48    }
49
50    /// Constructs the parity check matrix for the code.
51    pub fn h(&self) -> SparseMatrix {
52        let m = 1 << self.m().log2();
53        let extra_column_blocks = match self.rate {
54            AR4JARate::R1_2 => 0,
55            AR4JARate::R2_3 => 2,
56            AR4JARate::R4_5 => 6,
57        };
58        let extra_columns = m * extra_column_blocks;
59        let mut h = SparseMatrix::new(3 * m, extra_columns + 5 * m);
60
61        // fill common part (H_1/2)
62        for i in 0..m {
63            // block(0,2) = I_M
64            h.insert(i, extra_columns + 2 * m + i);
65            // block(0,4) = I_M + Pi_1
66            h.insert(i, extra_columns + 4 * m + i);
67            h.toggle(i, extra_columns + 4 * m + self.pi(1, i));
68            // block(1,0) = I_M
69            h.insert(m + i, extra_columns + i);
70            // block(1,1) = I_M
71            h.insert(m + i, extra_columns + m + i);
72            // block(1,3) = I_M
73            h.insert(m + i, extra_columns + 3 * m + i);
74            // block(1,4) = Pi_2 + Pi_3 + Pi_4
75            h.insert(m + i, extra_columns + 4 * m + self.pi(2, i));
76            h.toggle(m + i, extra_columns + 4 * m + self.pi(3, i));
77            h.toggle(m + i, extra_columns + 4 * m + self.pi(4, i));
78            // block(2,0) = I_M
79            h.insert(2 * m + i, extra_columns + i);
80            // block(2,1) = Pi_5 + Pi_6
81            h.insert(2 * m + i, extra_columns + m + self.pi(5, i));
82            h.toggle(2 * m + i, extra_columns + m + self.pi(6, i));
83            // block(2,3) = Pi_7 + Pi_8
84            h.insert(2 * m + i, extra_columns + 3 * m + self.pi(7, i));
85            h.toggle(2 * m + i, extra_columns + 3 * m + self.pi(8, i));
86            // block(2,4) = I_M
87            h.insert(2 * m + i, extra_columns + 4 * m + i);
88        }
89
90        if !matches!(self.rate, AR4JARate::R1_2) {
91            // fill specific H_2/3 part
92            let extra_columns = match self.rate {
93                AR4JARate::R1_2 => unreachable!(),
94                AR4JARate::R2_3 => 0,
95                AR4JARate::R4_5 => 4 * m,
96            };
97            for i in 0..m {
98                // block(1,0) = Pi_9 + Pi_10 + Pi_11
99                h.insert(m + i, extra_columns + self.pi(9, i));
100                h.toggle(m + i, extra_columns + self.pi(10, i));
101                h.toggle(m + i, extra_columns + self.pi(11, i));
102                // block(1,1) = I_M
103                h.insert(m + i, extra_columns + m + i);
104                // block(2,0) = I_M
105                h.insert(2 * m + i, extra_columns + i);
106                // block(2,1) = Pi_12 + Pi_13 + Pi_14
107                h.insert(2 * m + i, extra_columns + m + self.pi(12, i));
108                h.toggle(2 * m + i, extra_columns + m + self.pi(13, i));
109                h.toggle(2 * m + i, extra_columns + m + self.pi(14, i));
110            }
111        }
112
113        if matches!(self.rate, AR4JARate::R4_5) {
114            // fill specific H_4/5 part
115            for i in 0..m {
116                // block(1,0) = Pi_21 + Pi_22 + Pi_23
117                h.insert(m + i, self.pi(21, i));
118                h.toggle(m + i, self.pi(22, i));
119                h.toggle(m + i, self.pi(23, i));
120                // block(1,1) = I_M
121                h.insert(m + i, m + i);
122                // block(1,2) = Pi_15 + Pi_16 + Pi_17
123                h.insert(m + i, 2 * m + self.pi(15, i));
124                h.toggle(m + i, 2 * m + self.pi(16, i));
125                h.toggle(m + i, 2 * m + self.pi(17, i));
126                // block(1,3) = I_M
127                h.insert(m + i, 3 * m + i);
128                // block(2,0) = I_M
129                h.insert(2 * m + i, i);
130                // block(2,1) = Pi_24 + Pi_25 + Pi_26
131                h.insert(2 * m + i, m + self.pi(24, i));
132                h.toggle(2 * m + i, m + self.pi(25, i));
133                h.toggle(2 * m + i, m + self.pi(26, i));
134                // block(2,2) = I_M
135                h.insert(2 * m + i, 2 * m + i);
136                // block(2,3) = Pi_18 + Pi_19 + Pi_20
137                h.insert(2 * m + i, 3 * m + self.pi(18, i));
138                h.toggle(2 * m + i, 3 * m + self.pi(19, i));
139                h.toggle(2 * m + i, 3 * m + self.pi(20, i));
140            }
141        }
142
143        h
144    }
145
146    // Table 7.2 in [1]
147    fn m(&self) -> M {
148        match (self.rate, self.k) {
149            (AR4JARate::R1_2, AR4JAInfoSize::K1024) => M::M512,
150            (AR4JARate::R2_3, AR4JAInfoSize::K1024) => M::M256,
151            (AR4JARate::R4_5, AR4JAInfoSize::K1024) => M::M128,
152            (AR4JARate::R1_2, AR4JAInfoSize::K4096) => M::M2048,
153            (AR4JARate::R2_3, AR4JAInfoSize::K4096) => M::M1024,
154            (AR4JARate::R4_5, AR4JAInfoSize::K4096) => M::M512,
155            (AR4JARate::R1_2, AR4JAInfoSize::K16384) => M::M8192,
156            (AR4JARate::R2_3, AR4JAInfoSize::K16384) => M::M4096,
157            (AR4JARate::R4_5, AR4JAInfoSize::K16384) => M::M2048,
158        }
159    }
160
161    // Table 7-3 and 7-4 in [1]
162    fn theta(k: usize) -> usize {
163        assert!((1..=26).contains(&k));
164        THETA_K[k - 1].into()
165    }
166
167    // Table 7-3 and 7-4 in [1]
168    fn phi(&self, k: usize, j: usize) -> usize {
169        assert!((1..=26).contains(&k));
170        assert!((0..4).contains(&j));
171        let m_index = self.m().log2() - M::M128.log2();
172        PHI_K[j][k - 1][m_index]
173    }
174
175    // Section 7.4.2.4 in [1]
176    fn pi(&self, k: usize, i: usize) -> usize {
177        let m_log2 = self.m().log2();
178        let m = 1 << m_log2;
179        let j = 4 * i / m;
180        // & 0x3 gives mod 4
181        let a = (Self::theta(k) + j) & 0x3;
182        let m_div_4 = 1 << (m_log2 - 2);
183        // & (m_div_4 - 1) gives mod M/4
184        let b = (self.phi(k, j) + i) & (m_div_4 - 1);
185        // << (m_log2 - 2) gives * M/4
186        (a << (m_log2 - 2)) + b
187    }
188}
189
190enum M {
191    M128,
192    M256,
193    M512,
194    M1024,
195    M2048,
196    M4096,
197    M8192,
198}
199
200impl M {
201    fn log2(&self) -> usize {
202        match self {
203            M::M128 => 7,
204            M::M256 => 8,
205            M::M512 => 9,
206            M::M1024 => 10,
207            M::M2048 => 11,
208            M::M4096 => 12,
209            M::M8192 => 13,
210        }
211    }
212}
213
214static THETA_K: [u8; 26] = [
215    3, 0, 1, 2, 2, 3, 0, 1, 0, 1, 2, 0, 2, 3, 0, 1, 2, 0, 1, 2, 0, 1, 2, 1, 2, 3,
216];
217
218// Table 7-3 and 7-4 in [1]
219static PHI_K: [[[usize; 7]; 26]; 4] = [
220    // j = 0
221    [
222        [1, 59, 16, 160, 108, 226, 1148],
223        [22, 18, 103, 241, 126, 618, 2032],
224        [0, 52, 105, 185, 238, 404, 249],
225        [26, 23, 0, 251, 481, 32, 1807],
226        [0, 11, 50, 209, 96, 912, 485],
227        [10, 7, 29, 103, 28, 950, 1044],
228        [5, 22, 115, 90, 59, 534, 717],
229        [18, 25, 30, 184, 225, 63, 873],
230        [3, 27, 92, 248, 323, 971, 364],
231        [22, 30, 78, 12, 28, 304, 1926],
232        [3, 43, 70, 111, 386, 409, 1241],
233        [8, 14, 66, 66, 305, 708, 1769],
234        [25, 46, 39, 173, 34, 719, 532],
235        [25, 62, 84, 42, 510, 176, 768],
236        [2, 44, 79, 157, 147, 743, 1138],
237        [27, 12, 70, 174, 199, 759, 965],
238        [7, 38, 29, 104, 347, 674, 141],
239        [7, 47, 32, 144, 391, 958, 1527],
240        [15, 1, 45, 43, 165, 984, 505],
241        [10, 52, 113, 181, 414, 11, 1312],
242        [4, 61, 86, 250, 97, 413, 1840],
243        [19, 10, 1, 202, 158, 925, 709],
244        [7, 55, 42, 68, 86, 687, 1427],
245        [9, 7, 118, 177, 168, 752, 989],
246        [26, 12, 33, 170, 506, 867, 1925],
247        [17, 2, 126, 89, 489, 323, 270],
248    ],
249    // j = 1
250    [
251        [0, 0, 0, 0, 0, 0, 0],
252        [27, 32, 53, 182, 375, 767, 1822],
253        [30, 21, 74, 249, 436, 227, 203],
254        [28, 36, 45, 65, 350, 247, 882],
255        [7, 30, 47, 70, 260, 284, 1989],
256        [1, 29, 0, 141, 84, 370, 957],
257        [8, 44, 59, 237, 318, 482, 1705],
258        [20, 29, 102, 77, 382, 273, 1083],
259        [26, 39, 25, 55, 169, 886, 1072],
260        [24, 14, 3, 12, 213, 634, 354],
261        [4, 22, 88, 227, 67, 762, 1942],
262        [12, 15, 65, 42, 313, 184, 446],
263        [23, 48, 62, 52, 242, 696, 1456],
264        [15, 55, 68, 243, 188, 413, 1940],
265        [15, 39, 91, 179, 1, 854, 1660],
266        [22, 11, 70, 250, 306, 544, 1661],
267        [31, 1, 115, 247, 397, 864, 587],
268        [3, 50, 31, 164, 80, 82, 708],
269        [29, 40, 121, 17, 33, 1009, 1466],
270        [21, 62, 45, 31, 7, 437, 433],
271        [2, 27, 56, 149, 447, 36, 1345],
272        [5, 38, 54, 105, 336, 562, 867],
273        [11, 40, 108, 183, 424, 816, 1551],
274        [26, 15, 14, 153, 134, 452, 2041],
275        [9, 11, 30, 177, 152, 290, 1383],
276        [17, 18, 116, 19, 492, 778, 1790],
277    ],
278    // j = 2
279    [
280        [0, 0, 0, 0, 0, 0, 0],
281        [12, 46, 8, 35, 219, 254, 318],
282        [30, 45, 119, 167, 16, 790, 494],
283        [18, 27, 89, 214, 263, 642, 1467],
284        [10, 48, 31, 84, 415, 248, 757],
285        [16, 37, 122, 206, 403, 899, 1085],
286        [13, 41, 1, 122, 184, 328, 1630],
287        [9, 13, 69, 67, 279, 518, 64],
288        [7, 9, 92, 147, 198, 477, 689],
289        [15, 49, 47, 54, 307, 404, 1300],
290        [16, 36, 11, 23, 432, 698, 148],
291        [18, 10, 31, 93, 240, 160, 777],
292        [4, 11, 19, 20, 454, 497, 1431],
293        [23, 18, 66, 197, 294, 100, 659],
294        [5, 54, 49, 46, 479, 518, 352],
295        [3, 40, 81, 162, 289, 92, 1177],
296        [29, 27, 96, 101, 373, 464, 836],
297        [11, 35, 38, 76, 104, 592, 1572],
298        [4, 25, 83, 78, 141, 198, 348],
299        [8, 46, 42, 253, 270, 856, 1040],
300        [2, 24, 58, 124, 439, 235, 779],
301        [11, 33, 24, 143, 333, 134, 476],
302        [11, 18, 25, 63, 399, 542, 191],
303        [3, 37, 92, 41, 14, 545, 1393],
304        [15, 35, 38, 214, 277, 777, 1752],
305        [13, 21, 120, 70, 412, 483, 1627],
306    ],
307    // j = 3
308    [
309        [0, 0, 0, 0, 0, 0, 0],
310        [13, 44, 35, 162, 312, 285, 1189],
311        [19, 51, 97, 7, 503, 554, 458],
312        [14, 12, 112, 31, 388, 809, 460],
313        [15, 15, 64, 164, 48, 185, 1039],
314        [20, 12, 93, 11, 7, 49, 1000],
315        [17, 4, 99, 237, 185, 101, 1265],
316        [4, 7, 94, 125, 328, 82, 1223],
317        [4, 2, 103, 133, 254, 898, 874],
318        [11, 30, 91, 99, 202, 627, 1292],
319        [17, 53, 3, 105, 285, 154, 1491],
320        [20, 23, 6, 17, 11, 65, 631],
321        [8, 29, 39, 97, 168, 81, 464],
322        [22, 37, 113, 91, 127, 823, 461],
323        [19, 42, 92, 211, 8, 50, 844],
324        [15, 48, 119, 128, 437, 413, 392],
325        [5, 4, 74, 82, 475, 462, 922],
326        [21, 10, 73, 115, 85, 175, 256],
327        [17, 18, 116, 248, 419, 715, 1986],
328        [9, 56, 31, 62, 459, 537, 19],
329        [20, 9, 127, 26, 468, 722, 266],
330        [18, 11, 98, 140, 209, 37, 471],
331        [31, 23, 23, 121, 311, 488, 1166],
332        [13, 8, 38, 12, 211, 179, 1300],
333        [2, 7, 18, 41, 510, 430, 1033],
334        [18, 24, 62, 249, 320, 264, 1606],
335    ],
336];
337
338/// C2 code definition.
339///
340/// This C2 code is the basic (8176, 7156) LDPC code. Expurgation, shortening
341/// and extension used to construct the (8160, 7136) code should be handled
342/// separately.
343#[derive(Copy, Clone, Debug, Eq, PartialEq, Hash, Default)]
344pub struct C2Code {}
345
346impl C2Code {
347    /// Creates a C2 code definition.
348    pub fn new() -> C2Code {
349        C2Code::default()
350    }
351
352    /// Constructs the parity check matrix for the code.
353    pub fn h(&self) -> SparseMatrix {
354        const N: usize = 511;
355        let mut h = SparseMatrix::new(Self::ROW_BLOCKS * N, Self::COL_BLOCKS * N);
356        for (row, circs) in C2_CIRCULANTS.iter().enumerate() {
357            for (col, circs) in circs.iter().enumerate() {
358                for &circ in circs.iter() {
359                    let circ = usize::from(circ);
360                    for j in 0..N {
361                        h.insert(row * N + j, col * N + (j + circ) % N);
362                    }
363                }
364            }
365        }
366        h
367    }
368
369    const ROW_BLOCKS: usize = 2;
370    const COL_BLOCKS: usize = 16;
371    const BLOCK_WEIGHT: usize = 2;
372}
373
374// Table 7-1 in CCSDS 131.0-B-5
375static C2_CIRCULANTS: [[[u16; C2Code::BLOCK_WEIGHT]; C2Code::COL_BLOCKS]; C2Code::ROW_BLOCKS] = [
376    [
377        [0, 176],
378        [12, 239],
379        [0, 352],
380        [24, 431],
381        [0, 392],
382        [151, 409],
383        [0, 351],
384        [9, 359],
385        [0, 307],
386        [53, 329],
387        [0, 207],
388        [18, 281],
389        [0, 399],
390        [202, 457],
391        [0, 247],
392        [36, 261],
393    ],
394    [
395        [99, 471],
396        [130, 473],
397        [198, 435],
398        [260, 478],
399        [215, 420],
400        [282, 481],
401        [48, 396],
402        [193, 445],
403        [273, 430],
404        [302, 451],
405        [96, 379],
406        [191, 386],
407        [244, 467],
408        [364, 470],
409        [51, 382],
410        [192, 414],
411    ],
412];
413
414#[cfg(test)]
415mod test {
416    use super::*;
417
418    fn pi_k_model(code: &AR4JACode, k: usize, i: usize) -> usize {
419        let m = 1 << code.m().log2();
420        let theta_k = AR4JACode::theta(k);
421        let phi_k = code.phi(k, 4 * i / m);
422        m / 4 * ((theta_k + (4 * i / m)) % 4) + (phi_k + i) % (m / 4)
423    }
424
425    // Checks that AR4JACode::pi matches the simpler (but less efficient)
426    // implementation given in pi_k_model.
427    #[test]
428    fn pi_k() {
429        for rate in enum_iterator::all() {
430            for info_k in enum_iterator::all() {
431                let code = AR4JACode::new(rate, info_k);
432                let m = 1 << code.m().log2();
433                for k in 1..=26 {
434                    for i in 0..m {
435                        assert_eq!(code.pi(k, i), pi_k_model(&code, k, i));
436                    }
437                }
438            }
439        }
440    }
441}