Calculating polynomial division result as well as remainder (CRC)

Question

I'm trying to write a table-based CRC routine for receiving Mode S uplink interrogator messages. On the downlink side, the CRC is just the 24-bit CRC based on polynomial P=0x1FFF409. So far, so good -- I wrote a table-based implementation that follows the usual byte-at-a-time convention, and it's working fine.

On the uplink side, though, things get weird. The protocol specification says that calculating the target uplink address is by finding:

U' = x^24 * U / G(x)

...where U is the received message and G(x) is the encoding polynomial 0x1FFF409, resulting in:

U' = x^24 * m(x) + A(x) + r(x) / G(x)

...where m(x) is the original message, A(x) is the address, and r(x) is the remainder. I want the low-order quotient A(x); e.g., the result of the GF(2) polynomial division operation instead of the remainder. The remainder is effectively discarded. The target address is encoded with the transmitted checksum such that the receiving aircraft can validate the checksum by comparing it with its address.

This is great and all, and I have a bitwise implementation which follows from the above. Please ignore the weird shifting of the polynomial and checksum, this has been cribbed from this Pascal implementation (on page 15) which assumes 32-bit registers and makes optimizations based on that assumption. In reality the message and checksum come as a single, 56-bit transmission.

#This is the reference bit-shifting implementation. It is slow.
def uplink_bitshift_crc():
    p = 0xfffa0480 #polynomial (0x1FFF409 shifted left 7 bits)
    a = 0x00000000 #rx'ed uplink data (32 bits)
    adr = 0xcc5ee900 #rx'ed checksum (24 bits, shifted left 8 bits)
    ad = 0 #will hold division result low-order bits

    for j in range(56):
        #if MSBit is 1, xor w/poly
        if a & 0x80000000:
            a = a ^ p
        #shift off the top bit of A (we're done with it),
        #and shift in the top bit of adr
        a = ((a << 1) & 0xFFFFFFFF) + ((adr >> 31) & 1)
        #shift off the top bit of adr
        adr = (adr << 1) & 0xFFFFFFFF

        if j > 30:
            #shift ad left 1 bit and shift in the msbit of a
            #this extracts the LS 24bits of the division operation
            #and ignores the remainder at the end
            ad = ad + ((a >> 31) & 1)
            ad = ((ad << 1) & 0xFFFFFFFF)

    #correct the ad
    ad = ad >> 2
    return ad

The above is of course slower than molasses in software and I'd really like to be able to construct a lookup table that would allow similar byte-at-a-time calculation of the received address, or massage the remainder (which is quickly calculated) into a quotient.

TL;DR: Given a message, the encoding polynomial, and the remainder (calculated by the normal CRC method), is there a faster way to obtain the quotient of the polynomial division operation than by using shift registers to do polynomial division "longhand"?

Answer 1

You might take a look at the PyCRC library, I guess this may answer your questions.

Answer 2

Too late for the OP, but I'm posting this for others that might see this question. You can generate two tables to operate a byte at a time. The first 256 by 8 bit table is indexed by the current leading 8 bits of the dividend (message), and the 8 bit values are the quotients. The second 256 by 32 bit table is indexed by the 8 bit quotient and the 32 bit values are the 32 bit product of the 8 bit quotient times the 25 bit polynomial (since this is a carryless multiply, the product is 32 bits, (x^7 * x^24 = x^31)), which you xor to the upper 32 bits of the dividend, which will zero out the upper 8 bits of the dividend. Then loop back for the next 8 bits of the dividend.

A modern X86 cpu has the carryless multiply instruction, PCLMULQDQ that operates on 128 bit xmm registers, performing a 64 bit by 64 bit multiply to produce a 128 bit product (since it's a carryless multiply bit 127 is always 0, so it's really a 127 bit product). A multiply of the 56 bit message by the 41 bit constant 2^64/G(x) will produce a 96 bit product, of which the upper 32 bits will be the quotient (lower 64 bits are not used).