dump1090/crc.c

556 lines
17 KiB
C

// Part of dump1090, a Mode S message decoder for RTLSDR devices.
//
// crc.h: Mode S CRC calculation and error correction.
//
// Copyright (c) 2014,2015 Oliver Jowett <oliver@mutability.co.uk>
//
// This file is free software: you may copy, redistribute and/or modify it
// under the terms of the GNU General Public License as published by the
// Free Software Foundation, either version 2 of the License, or (at your
// option) any later version.
//
// This file is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program. If not, see <http://www.gnu.org/licenses/>.
#include "dump1090.h"
#include <assert.h>
// Errorinfo for "no errors"
static struct errorinfo NO_ERRORS;
// Generator polynomial for the Mode S CRC:
#define MODES_GENERATOR_POLY 0xfff409U
// CRC values for all single-byte messages;
// used to speed up CRC calculation.
static uint32_t crc_table[256];
// Syndrome values for all single-bit errors;
// used to speed up construction of error-
// correction tables.
static uint32_t single_bit_syndrome[112];
static void initLookupTables()
{
int i;
uint8_t msg[112/8];
for (i = 0; i < 256; ++i) {
uint32_t c = i << 16;
int j;
for (j = 0; j < 8; ++j) {
if (c & 0x800000)
c = (c<<1) ^ MODES_GENERATOR_POLY;
else
c = (c<<1);
}
crc_table[i] = c & 0x00ffffff;
}
memset(msg, 0, sizeof(msg));
for (i = 0; i < 112; ++i) {
msg[i/8] ^= 1 << (7 - (i & 7));
single_bit_syndrome[i] = modesChecksum(msg, 112);
msg[i/8] ^= 1 << (7 - (i & 7));
}
}
uint32_t modesChecksum(uint8_t *message, int bits)
{
uint32_t rem = 0;
int i;
int n = bits/8;
assert(bits % 8 == 0);
assert(n >= 3);
for (i = 0; i < n-3; ++i) {
rem = (rem << 8) ^ crc_table[message[i] ^ ((rem & 0xff0000) >> 16)];
rem = rem & 0xffffff;
}
rem = rem ^ (message[n-3] << 16) ^ (message[n-2] << 8) ^ (message[n-1]);
return rem;
}
static struct errorinfo *bitErrorTable_short;
static int bitErrorTableSize_short;
static struct errorinfo *bitErrorTable_long;
static int bitErrorTableSize_long;
// compare two errorinfo structures
static int syndrome_compare(const void *x, const void *y) {
struct errorinfo *ex = (struct errorinfo*)x;
struct errorinfo *ey = (struct errorinfo*)y;
return (int)ex->syndrome - (int)ey->syndrome;
}
// (n k), the number of ways of selecting k distinct items from a set of n items
static int combinations(int n, int k)
{
int result = 1, i;
if (k == 0 || k == n)
return 1;
if (k > n)
return 0;
for (i = 1; i <= k; ++i) {
result = result * n / i;
n = n - 1;
}
return result;
}
// Recursively populates an errorinfo table with error syndromes
//
// in:
// table: the table to fill
// n: first entry to fill
// maxSize: max size of table
// offset: start bit offset for checksum calculation
// startbit: first bit to introduce errors into
// endbit: (one past) last bit to introduce errors info
// base_entry: template entry to start from
// error_bit: how many error bits have already been set
// max_errors: maximum total error bits to set
// out:
// returns: the next free entry in the table
// table: has been populated between [n, return value)
static int prepareSubtable(struct errorinfo *table, int n, int maxsize, int offset, int startbit, int endbit, struct errorinfo *base_entry, int error_bit, int max_errors)
{
int i = 0;
if (error_bit >= max_errors)
return n;
for (i = startbit; i < endbit; ++i) {
assert(n < maxsize);
table[n] = *base_entry;
table[n].syndrome ^= single_bit_syndrome[i + offset];
table[n].errors = error_bit+1;
table[n].bit[error_bit] = i;
++n;
n = prepareSubtable(table, n, maxsize, offset, i + 1, endbit, &table[n-1], error_bit + 1, max_errors);
}
return n;
}
static int flagCollisions(struct errorinfo *table, int tablesize, int offset, int startbit, int endbit, uint32_t base_syndrome, int error_bit, int first_error, int last_error)
{
int i = 0;
int count = 0;
if (error_bit > last_error)
return 0;
for (i = startbit; i < endbit; ++i) {
struct errorinfo ei;
ei.syndrome = base_syndrome ^ single_bit_syndrome[i + offset];
if (error_bit >= first_error) {
struct errorinfo *collision = bsearch(&ei, table, tablesize, sizeof(struct errorinfo), syndrome_compare);
if (collision != NULL && collision->errors != -1) {
++count;
collision->errors = -1;
}
}
count += flagCollisions(table, tablesize, offset, i+1, endbit, ei.syndrome, error_bit + 1, first_error, last_error);
}
return count;
}
// Allocate and build an error table for messages of length "bits" (max 112)
// returns a pointer to the new table and sets *size_out to the table length
static struct errorinfo *prepareErrorTable(int bits, int max_correct, int max_detect, int *size_out)
{
int maxsize, usedsize;
struct errorinfo *table;
struct errorinfo base_entry;
int i, j;
assert (bits >= 0 && bits <= 112);
assert (max_correct >=0 && max_correct <= MODES_MAX_BITERRORS);
assert (max_detect >= max_correct);
if (!max_correct) {
*size_out = 0;
return NULL;
}
maxsize = 0;
for (i = 1; i <= max_correct; ++i) {
maxsize += combinations(bits - 5, i); // space needed for all i-bit errors
}
#ifdef CRCDEBUG
fprintf(stderr, "Preparing syndrome table to correct up to %d-bit errors (detecting %d-bit errors) in a %d-bit message (max %d entries)\n", max_correct, max_detect, bits, maxsize);
#endif
table = malloc(maxsize * sizeof(struct errorinfo));
base_entry.syndrome = 0;
base_entry.errors = 0;
for (i = 0; i < MODES_MAX_BITERRORS; ++i)
base_entry.bit[i] = -1;
// ignore the first 5 bits (DF type)
usedsize = prepareSubtable(table, 0, maxsize, 112 - bits, 5, bits, &base_entry, 0, max_correct);
#ifdef CRCDEBUG
fprintf(stderr, "%d syndromes (expected %d).\n", usedsize, maxsize);
fprintf(stderr, "Sorting syndromes..\n");
#endif
qsort(table, usedsize, sizeof(struct errorinfo), syndrome_compare);
#ifdef CRCDEBUG
{
// Show the table stats
fprintf(stderr, "Undetectable errors:\n");
for (i = 1; i <= max_correct; ++i) {
int j, count;
count = 0;
for (j = 0; j < usedsize; ++j)
if (table[j].errors == i && table[j].syndrome == 0)
++count;
fprintf(stderr, " %d undetectable %d-bit errors\n", count, i);
}
}
#endif
// Handle ambiguous cases, where there is more than one possible error pattern
// that produces a given syndrome (this happens with >2 bit errors).
#ifdef CRCDEBUG
fprintf(stderr, "Finding collisions..\n");
#endif
for (i = 0, j = 0; i < usedsize; ++i) {
if (i < usedsize-1 && table[i+1].syndrome == table[i].syndrome) {
// skip over this entry and all collisions
while (i < usedsize && table[i+1].syndrome == table[i].syndrome)
++i;
// now table[i] is the last duplicate
continue;
}
if (i != j)
table[j] = table[i];
++j;
}
if (j < usedsize) {
#ifdef CRCDEBUG
fprintf(stderr, "Discarded %d collisions.\n", usedsize - j);
#endif
usedsize = j;
}
// Flag collisions we want to detect but not correct
if (max_detect > max_correct) {
int flagged;
#ifdef CRCDEBUG
fprintf(stderr, "Flagging collisions between %d - %d bits..\n", max_correct+1, max_detect);
#endif
flagged = flagCollisions(table, usedsize, 112 - bits, 5, bits, 0, 1, max_correct+1, max_detect);
#ifdef CRCDEBUG
fprintf(stderr, "Flagged %d collisions for removal.\n", flagged);
#else
#endif
if (flagged > 0) {
for (i = 0, j = 0; i < usedsize; ++i) {
if (table[i].errors != -1) {
if (i != j)
table[j] = table[i];
++j;
}
}
#ifdef CRCDEBUG
fprintf(stderr, "Discarded %d flagged collisions.\n", usedsize - j);
#endif
usedsize = j;
}
}
if (usedsize < maxsize) {
#ifdef CRCDEBUG
fprintf(stderr, "Shrinking table from %d to %d..\n", maxsize, usedsize);
table = realloc(table, usedsize * sizeof(struct errorinfo));
#endif
}
*size_out = usedsize;
#ifdef CRCDEBUG
{
// Check the table.
unsigned char *msg = malloc(bits/8);
for (i = 0; i < usedsize; ++i) {
int j;
struct errorinfo *ei;
uint32_t result;
memset(msg, 0, bits/8);
ei = &table[i];
for (j = 0; j < ei->errors; ++j) {
msg[ei->bit[j] >> 3] ^= 1 << (7 - (ei->bit[j]&7));
}
result = modesChecksum(msg, bits);
if (result != ei->syndrome) {
fprintf(stderr, "PROBLEM: entry %6d/%6d syndrome %06x errors %d bits ", i, usedsize, ei->syndrome, ei->errors);
for (j = 0; j < ei->errors; ++j)
fprintf(stderr, "%3d ", ei->bit[j]);
fprintf(stderr, " checksum %06x\n", result);
}
}
free(msg);
// Show the table stats
fprintf(stderr, "Syndrome table summary:\n");
for (i = 1; i <= max_correct; ++i) {
int j, count, possible;
count = 0;
for (j = 0; j < usedsize; ++j)
if (table[j].errors == i)
++count;
possible = combinations(bits-5, i);
fprintf(stderr, " %d entries for %d-bit errors (%d possible, %d%% coverage)\n", count, i, possible, 100 * count / possible);
}
fprintf(stderr, " %d entries total\n", usedsize);
}
#endif
return table;
}
// Precompute syndrome tables for 56- and 112-bit messages.
void modesChecksumInit(int fixBits)
{
initLookupTables();
switch (fixBits) {
case 0:
bitErrorTable_short = bitErrorTable_long = NULL;
bitErrorTableSize_short = bitErrorTableSize_long = 0;
break;
case 1:
// For 1 bit correction, we have 100% coverage up to 4 bit detection, so don't bother
// with flagging collisions there.
bitErrorTable_short = prepareErrorTable(MODES_SHORT_MSG_BITS, 1, 1, &bitErrorTableSize_short);
bitErrorTable_long = prepareErrorTable(MODES_LONG_MSG_BITS, 1, 1, &bitErrorTableSize_long);
break;
default:
// Detect out to 4 bit errors; this reduces our 2-bit coverage to about 65%.
// This can take a little while - tell the user.
fprintf(stderr, "Preparing error correction tables.. ");
bitErrorTable_short = prepareErrorTable(MODES_SHORT_MSG_BITS, 2, 4, &bitErrorTableSize_short);
bitErrorTable_long = prepareErrorTable(MODES_LONG_MSG_BITS, 2, 4, &bitErrorTableSize_long);
fprintf(stderr, "done.\n");
break;
}
}
// Given an error syndrome and message length, return
// an error-correction descriptor, or NULL if the
// syndrome is uncorrectable
struct errorinfo *modesChecksumDiagnose(uint32_t syndrome, int bitlen)
{
struct errorinfo *table;
int tablesize;
struct errorinfo ei;
if (syndrome == 0)
return &NO_ERRORS;
assert (bitlen == 56 || bitlen == 112);
if (bitlen == 56) { table = bitErrorTable_short; tablesize = bitErrorTableSize_short; }
else { table = bitErrorTable_long; tablesize = bitErrorTableSize_long; }
if (!table)
return NULL;
ei.syndrome = syndrome;
return bsearch(&ei, table, tablesize, sizeof(struct errorinfo), syndrome_compare);
}
// Given a message and an error-correction descriptor,
// apply the error correction to the given message.
void modesChecksumFix(uint8_t *msg, struct errorinfo *info)
{
int i;
if (!info)
return;
for (i = 0; i < info->errors; ++i)
msg[info->bit[i] >> 3] ^= 1 << (7 - (info->bit[i] & 7));
}
#ifdef CRCDEBUG
int main(int argc, char **argv)
{
int shortlen, longlen;
int i;
struct errorinfo *shorttable, *longtable;
if (argc < 3) {
fprintf(stderr, "syntax: crctests <ncorrect> <ndetect>\n");
return 1;
}
initLookupTables();
shorttable = prepareErrorTable(MODES_SHORT_MSG_BITS, atoi(argv[1]), atoi(argv[2]), &shortlen);
longtable = prepareErrorTable(MODES_LONG_MSG_BITS, atoi(argv[1]), atoi(argv[2]), &longlen);
// check for DF11 correction syndromes where there is a syndrome with lower 7 bits all zero
// (which would be used for DF11 error correction), but there's also a syndrome which has
// the same upper 17 bits but nonzero lower 7 bits.
// empirically, with ncorrect=1 ndetect=2 we get no ambiguous syndromes;
// for ncorrect=2 ndetect=4 we get 11 ambiguous syndromes:
/*
syndrome 1 = 000C00 bits=[ 44 45 ]
syndrome 2 = 000C1B bits=[ 30 43 ]
syndrome 1 = 001400 bits=[ 43 45 ]
syndrome 2 = 00141B bits=[ 30 44 ]
syndrome 1 = 001800 bits=[ 43 44 ]
syndrome 2 = 00181B bits=[ 30 45 ]
syndrome 1 = 001800 bits=[ 43 44 ]
syndrome 2 = 001836 bits=[ 29 42 ]
syndrome 1 = 002400 bits=[ 42 45 ]
syndrome 2 = 00242D bits=[ 29 30 ]
syndrome 1 = 002800 bits=[ 42 44 ]
syndrome 2 = 002836 bits=[ 29 43 ]
syndrome 1 = 003000 bits=[ 42 43 ]
syndrome 2 = 003036 bits=[ 29 44 ]
syndrome 1 = 003000 bits=[ 42 43 ]
syndrome 2 = 00306C bits=[ 28 41 ]
syndrome 1 = 004800 bits=[ 41 44 ]
syndrome 2 = 00485A bits=[ 28 29 ]
syndrome 1 = 005000 bits=[ 41 43 ]
syndrome 2 = 00506C bits=[ 28 42 ]
syndrome 1 = 006000 bits=[ 41 42 ]
syndrome 2 = 00606C bits=[ 28 43 ]
*/
// So in the DF11 correction logic, we just discard messages that require more than a 1 bit fix.
fprintf(stderr, "checking %d syndromes for DF11 collisions..\n", shortlen);
for (i = 0; i < shortlen; ++i) {
if ((shorttable[i].syndrome & 0xFF) == 0) {
int j;
// all syndromes with the same first 17 bits should sort immediately after entry i,
// so this is fairly easy
for (j = i+1; j < shortlen; ++j) {
if ((shorttable[i].syndrome & 0xFFFF80) == (shorttable[j].syndrome & 0xFFFF80)) {
int k;
int mismatch = 0;
// we don't care if the only differences are in bits that lie in the checksum
for (k = 0; k < shorttable[i].errors; ++k) {
int l, matched = 0;
if (shorttable[i].bit[k] >= 49)
continue; // bit is in the final 7 bits, we don't care
for (l = 0; l < shorttable[j].errors; ++l) {
if (shorttable[i].bit[k] == shorttable[j].bit[l]) {
matched = 1;
break;
}
}
if (!matched)
mismatch = 1;
}
for (k = 0; k < shorttable[j].errors; ++k) {
int l, matched = 0;
if (shorttable[j].bit[k] >= 49)
continue; // bit is in the final 7 bits, we don't care
for (l = 0; l < shorttable[i].errors; ++l) {
if (shorttable[j].bit[k] == shorttable[i].bit[l]) {
matched = 1;
break;
}
}
if (!matched)
mismatch = 1;
}
if (mismatch) {
fprintf(stderr,
"DF11 correction collision: \n"
" syndrome 1 = %06X bits=[",
shorttable[i].syndrome);
for (k = 0; k < shorttable[i].errors; ++k)
fprintf(stderr, " %d", shorttable[i].bit[k]);
fprintf(stderr, " ]\n");
fprintf(stderr,
" syndrome 2 = %06X bits=[",
shorttable[j].syndrome);
for (k = 0; k < shorttable[j].errors; ++k)
fprintf(stderr, " %d", shorttable[j].bit[k]);
fprintf(stderr, " ]\n");
}
} else {
break;
}
}
}
}
free(shorttable);
free(longtable);
return 0;
}
#endif