diff --git a/.gitignore b/.gitignore index 453dfd0..c6d576d 100644 --- a/.gitignore +++ b/.gitignore @@ -5,6 +5,9 @@ __pycache__/ *.py[cod] .pytest_cache/ +#cython +.c + # C extensions *.so diff --git a/Makefile b/Makefile new file mode 100644 index 0000000..70cc54d --- /dev/null +++ b/Makefile @@ -0,0 +1,19 @@ +install: + pip install . --upgrade + +uninstall: + pip uninstall pyModeS -y + +ext: + python setup.py build_ext --inplace + +test: + python -m pytest + +clean: + find pyModeS/decoder -type f -name '*.c' -delete + find pyModeS/decoder -type f -name '*.so' -delete + find . | grep -E "(__pycache__|\.pyc|\.pyo$$)" | xargs rm -rf + rm -rf *.egg-info + rm -rf .pytest_cache + rm -rf build/* diff --git a/README.rst b/README.rst index b6fc7c0..fa0cb09 100644 --- a/README.rst +++ b/README.rst @@ -321,8 +321,19 @@ Here is an example: Unit test --------- -To perform unit tests. First, install ``tox`` through pip. Then, run the following commands: +To perform unit tests, ``pytest`` must be install first. -.. code:: bash +Build Cython extensions +:: - $ tox + $ make ext + +Run unit tests +:: + + $ make test + +Clean build files +:: + + $ make clean diff --git a/pyModeS/__init__.py b/pyModeS/__init__.py index bc491d9..a2014df 100644 --- a/pyModeS/__init__.py +++ b/pyModeS/__init__.py @@ -3,11 +3,16 @@ from __future__ import absolute_import, print_function, division import os import warnings -from .decoder.common import * +try: + from .decoder import c_common as common + from .decoder.c_common import * +except: + from .decoder import common + from .decoder.common import * + from .decoder import tell from .decoder import adsb from .decoder import commb -from .decoder import common from .decoder import bds from .extra import aero from .extra import tcpclient diff --git a/pyModeS/decoder/acas.py b/pyModeS/decoder/acas.py index c275030..e22759b 100644 --- a/pyModeS/decoder/acas.py +++ b/pyModeS/decoder/acas.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - """ Decoding Air-Air Surveillance (ACAS) DF=0/16 diff --git a/pyModeS/decoder/adsb.py b/pyModeS/decoder/adsb.py index a84dc1b..d28d77a 100644 --- a/pyModeS/decoder/adsb.py +++ b/pyModeS/decoder/adsb.py @@ -1,18 +1,3 @@ -# Copyright (C) 2015 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - """ADS-B Wrapper. The ADS-B wrapper also imports functions from the following modules: @@ -38,7 +23,7 @@ from pyModeS.decoder import uncertainty from pyModeS.decoder.bds.bds05 import ( airborne_position, airborne_position_with_ref, - altitude, + altitude as altitude05, ) from pyModeS.decoder.bds.bds06 import ( surface_position, @@ -143,18 +128,13 @@ def altitude(msg): if tc < 5 or tc == 19 or tc > 22: raise RuntimeError("%s: Not a position message" % msg) - if tc >= 5 and tc <= 8: + elif tc >= 5 and tc <= 8: # surface position, altitude 0 return 0 - msgbin = common.hex2bin(msg) - q = msgbin[47] - if q: - n = common.bin2int(msgbin[40:47] + msgbin[48:52]) - alt = n * 25 - 1000 - return alt else: - return None + # airborn position + return altitude05(msg) def velocity(msg, rtn_sources=False): @@ -177,7 +157,7 @@ def velocity(msg, rtn_sources=False): as reference, 'mag_north' for magnetic north as reference), rate of climb/descent source ('Baro' for barometer, 'GNSS' for GNSS constellation). - + In the case of surface messages, None will be put in place for vertical rate and its respective sources. """ diff --git a/pyModeS/decoder/allcall.py b/pyModeS/decoder/allcall.py index 6978e42..1901698 100644 --- a/pyModeS/decoder/allcall.py +++ b/pyModeS/decoder/allcall.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - """ Decoding all call replies DF=11 diff --git a/pyModeS/decoder/bds/bds05.py b/pyModeS/decoder/bds/bds05.py index 15a64fe..dce5742 100644 --- a/pyModeS/decoder/bds/bds05.py +++ b/pyModeS/decoder/bds/bds05.py @@ -1,19 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - - # ------------------------------------------ # BDS 0,5 # ADS-B TC=9-18 diff --git a/pyModeS/decoder/bds/bds06.py b/pyModeS/decoder/bds/bds06.py index 08dccdd..d882569 100644 --- a/pyModeS/decoder/bds/bds06.py +++ b/pyModeS/decoder/bds/bds06.py @@ -1,19 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - - # ------------------------------------------ # BDS 0,6 # ADS-B TC=5-8 diff --git a/pyModeS/decoder/bds/bds08.py b/pyModeS/decoder/bds/bds08.py index cd8068c..bdd7106 100644 --- a/pyModeS/decoder/bds/bds08.py +++ b/pyModeS/decoder/bds/bds08.py @@ -1,19 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - - # ------------------------------------------ # BDS 0,8 # ADS-B TC=1-4 diff --git a/pyModeS/decoder/bds/bds09.py b/pyModeS/decoder/bds/bds09.py index 1b36d64..1a26ca6 100644 --- a/pyModeS/decoder/bds/bds09.py +++ b/pyModeS/decoder/bds/bds09.py @@ -1,19 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - - # ------------------------------------------ # BDS 0,9 # ADS-B TC=19 diff --git a/pyModeS/decoder/bds/bds10.py b/pyModeS/decoder/bds/bds10.py index 1bea8cb..02cb658 100644 --- a/pyModeS/decoder/bds/bds10.py +++ b/pyModeS/decoder/bds/bds10.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 1,0 # Data link capability report diff --git a/pyModeS/decoder/bds/bds17.py b/pyModeS/decoder/bds/bds17.py index 07a58af..5614fce 100644 --- a/pyModeS/decoder/bds/bds17.py +++ b/pyModeS/decoder/bds/bds17.py @@ -1,19 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - - # ------------------------------------------ # BDS 1,7 # Common usage GICB capability report diff --git a/pyModeS/decoder/bds/bds20.py b/pyModeS/decoder/bds/bds20.py index 63f28d6..32ced8c 100644 --- a/pyModeS/decoder/bds/bds20.py +++ b/pyModeS/decoder/bds/bds20.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 2,0 # Aircraft identification diff --git a/pyModeS/decoder/bds/bds30.py b/pyModeS/decoder/bds/bds30.py index cb4d94c..59a7a9a 100644 --- a/pyModeS/decoder/bds/bds30.py +++ b/pyModeS/decoder/bds/bds30.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 3,0 # ACAS active resolution advisory diff --git a/pyModeS/decoder/bds/bds40.py b/pyModeS/decoder/bds/bds40.py index 181149f..a9891ad 100644 --- a/pyModeS/decoder/bds/bds40.py +++ b/pyModeS/decoder/bds/bds40.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 4,0 # Selected vertical intention diff --git a/pyModeS/decoder/bds/bds44.py b/pyModeS/decoder/bds/bds44.py index 50fc0cc..7bd59eb 100644 --- a/pyModeS/decoder/bds/bds44.py +++ b/pyModeS/decoder/bds/bds44.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 4,4 # Meteorological routine air report diff --git a/pyModeS/decoder/bds/bds45.py b/pyModeS/decoder/bds/bds45.py index 3dd5c0a..adf3197 100644 --- a/pyModeS/decoder/bds/bds45.py +++ b/pyModeS/decoder/bds/bds45.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 4,5 # Meteorological hazard report diff --git a/pyModeS/decoder/bds/bds50.py b/pyModeS/decoder/bds/bds50.py index ae65786..91ec924 100644 --- a/pyModeS/decoder/bds/bds50.py +++ b/pyModeS/decoder/bds/bds50.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 5,0 # Track and turn report diff --git a/pyModeS/decoder/bds/bds53.py b/pyModeS/decoder/bds/bds53.py index 962e260..6673489 100644 --- a/pyModeS/decoder/bds/bds53.py +++ b/pyModeS/decoder/bds/bds53.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 5,3 # Air-referenced state vector diff --git a/pyModeS/decoder/bds/bds60.py b/pyModeS/decoder/bds/bds60.py index 8998644..75f2927 100644 --- a/pyModeS/decoder/bds/bds60.py +++ b/pyModeS/decoder/bds/bds60.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - # ------------------------------------------ # BDS 6,0 # Heading and speed report diff --git a/pyModeS/decoder/c_common.pxd b/pyModeS/decoder/c_common.pxd new file mode 100644 index 0000000..e311ddf --- /dev/null +++ b/pyModeS/decoder/c_common.pxd @@ -0,0 +1,23 @@ +# cython: language_level=3 + +cdef int char_to_int(unsigned char binstr) +cdef unsigned char int_to_char(unsigned char i) + +cpdef str hex2bin(str hexstr) +cpdef long bin2int(str binstr) +cpdef long hex2int(str binstr) + +cpdef unsigned char df(str msg) +cpdef long crc(str msg, bint encode=*) + +cpdef long floor(double x) +cpdef str icao(str msg) +cpdef bint is_icao_assigned(str icao) + +cpdef int typecode(str msg) +cpdef int cprNL(double lat) +cpdef str idcode(str msg) +cpdef int altcode(str msg) + +cdef str data(str msg) +cpdef bint allzeros(str msg) diff --git a/pyModeS/decoder/c_common.pyx b/pyModeS/decoder/c_common.pyx new file mode 100644 index 0000000..2c76cf7 --- /dev/null +++ b/pyModeS/decoder/c_common.pyx @@ -0,0 +1,431 @@ +# cython: language_level=3 + +cimport cython +from cpython cimport array +from cpython.bytes cimport PyBytes_GET_SIZE +from cpython.bytearray cimport PyByteArray_GET_SIZE + +from libc.math cimport cos, acos, fabs, M_PI as pi, floor as c_floor + + +cdef int char_to_int(unsigned char binstr): + if 48 <= binstr <= 57: # 0 to 9 + return binstr - 48 + if 97 <= binstr <= 102: # a to f + return binstr - 97 + 10 + if 65 <= binstr <= 70: # A to F + return binstr - 65 + 10 + return 0 + +cdef unsigned char int_to_char(unsigned char i): + if i < 10: + return 48 + i # "0" + i + return 97 - 10 + i # "a" - 10 + i + +@cython.boundscheck(False) +@cython.overflowcheck(False) +cpdef str hex2bin(str hexstr): + """Convert a hexdecimal string to binary string, with zero fillings.""" + # num_of_bits = len(hexstr) * 4 + cdef hexbytes = bytes(hexstr.encode()) + cdef Py_ssize_t len_hexstr = PyBytes_GET_SIZE(hexbytes) + # binstr = bin(int(hexbytes, 16))[2:].zfill(int(num_of_bits)) + cdef bytearray _binstr = bytearray(4 * len_hexstr) + cdef unsigned char[:] binstr = _binstr + cdef unsigned char int_ + cdef Py_ssize_t i + for i in range(len_hexstr): + int_ = char_to_int(hexbytes[i]) + binstr[4*i] = int_to_char((int_ >> 3) & 1) + binstr[4*i+1] = int_to_char((int_ >> 2) & 1) + binstr[4*i+2] = int_to_char((int_ >> 1) & 1) + binstr[4*i+3] = int_to_char((int_) & 1) + return _binstr.decode() + +@cython.boundscheck(False) +cpdef long bin2int(str binstr): + """Convert a binary string to integer.""" + # return int(binstr, 2) + cdef bytearray binbytes = bytearray(binstr.encode()) + cdef Py_ssize_t len_ = PyByteArray_GET_SIZE(binbytes) + cdef long cumul = 0 + cdef unsigned char[:] v_binstr = binbytes + for i in range(len_): + cumul = 2*cumul + char_to_int(v_binstr[i]) + return cumul + +@cython.boundscheck(False) +cpdef long hex2int(str hexstr): + """Convert a binary string to integer.""" + # return int(hexstr, 2) + cdef bytearray binbytes = bytearray(hexstr.encode()) + cdef Py_ssize_t len_ = PyByteArray_GET_SIZE(binbytes) + cdef long cumul = 0 + cdef unsigned char[:] v_hexstr = binbytes + for i in range(len_): + cumul = 16*cumul + char_to_int(v_hexstr[i]) + return cumul + +@cython.boundscheck(False) +cpdef unsigned char df(str msg): + """Decode Downlink Format vaule, bits 1 to 5.""" + cdef str dfbin = hex2bin(msg[:2]) + # return min(bin2int(dfbin[0:5]), 24) + cdef long df = bin2int(dfbin[0:5]) + if df > 24: + return 24 + return df + +# the CRC generator +# G = [int("11111111", 2), int("11111010", 2), int("00000100", 2), int("10000000", 2)] +cdef array.array _G = array.array('l', [0b11111111, 0b11111010, 0b00000100, 0b10000000]) + +@cython.boundscheck(False) +@cython.wraparound(False) +cpdef long crc(str msg, bint encode=False): + """Mode-S Cyclic Redundancy Check. + + Detect if bit error occurs in the Mode-S message. When encode option is on, + the checksum is generated. + + Args: + msg (string): 28 bytes hexadecimal message string + encode (bool): True to encode the date only and return the checksum + Returns: + int: message checksum, or partity bits (encoder) + + """ + # the CRC generator + # G = [int("11111111", 2), int("11111010", 2), int("00000100", 2), int("10000000", 2)] + # cdef array.array _G = array.array('l', [0b11111111, 0b11111010, 0b00000100, 0b10000000]) + cdef long[4] G = _G + + # msgbin_split = wrap(msgbin, 8) + # mbytes = list(map(bin2int, msgbin_split)) + cdef bytearray _msgbin = bytearray(hex2bin(msg).encode()) + cdef unsigned char[:] msgbin = _msgbin + + cdef Py_ssize_t len_msgbin = PyByteArray_GET_SIZE(_msgbin) + cdef Py_ssize_t len_mbytes = len_msgbin // 8 + cdef Py_ssize_t i + + if encode: + for i in range(len_msgbin - 24, len_msgbin): + msgbin[i] = 0 + + cdef array.array _mbytes = array.array( + 'l', [bin2int(_msgbin[8*i:8*i+8].decode()) for i in range(len_mbytes)] + ) + + cdef long[:] mbytes = _mbytes + + cdef long bits, mask + cdef Py_ssize_t ibyte, ibit + + for ibyte in range(len_mbytes - 3): + for ibit in range(8): + mask = 0x80 >> ibit + bits = mbytes[ibyte] & mask + + if bits > 0: + mbytes[ibyte] = mbytes[ibyte] ^ (G[0] >> ibit) + mbytes[ibyte + 1] = mbytes[ibyte + 1] ^ ( + 0xFF & ((G[0] << 8 - ibit) | (G[1] >> ibit)) + ) + mbytes[ibyte + 2] = mbytes[ibyte + 2] ^ ( + 0xFF & ((G[1] << 8 - ibit) | (G[2] >> ibit)) + ) + mbytes[ibyte + 3] = mbytes[ibyte + 3] ^ ( + 0xFF & ((G[2] << 8 - ibit) | (G[3] >> ibit)) + ) + + cdef long result = (mbytes[len_mbytes-3] << 16) | (mbytes[len_mbytes-2] << 8) | mbytes[len_mbytes-1] + + return result + + + +cpdef long floor(double x): + """Mode-S floor function. + + Defined as the greatest integer value k, such that k <= x + For example: floor(3.6) = 3 and floor(-3.6) = -4 + + """ + return c_floor(x) + +cpdef str icao(str msg): + """Calculate the ICAO address from an Mode-S message. + + Applicable only with DF4, DF5, DF20, DF21 messages. + + Args: + msg (String): 28 bytes hexadecimal message string + + Returns: + String: ICAO address in 6 bytes hexadecimal string + + """ + cdef unsigned char DF = df(msg) + cdef long c0, c1 + + if DF in (11, 17, 18): + addr = msg[2:8] + elif DF in (0, 4, 5, 16, 20, 21): + c0 = crc(msg, encode=True) + c1 = hex2int(msg[-6:]) + addr = "%06X" % (c0 ^ c1) + else: + addr = None + + return addr + + +cpdef bint is_icao_assigned(str icao): + """Check whether the ICAO address is assigned (Annex 10, Vol 3).""" + if (icao is None) or (not isinstance(icao, str)) or (len(icao) != 6): + return False + + cdef long icaoint = hex2int(icao) + + if 0x200000 < icaoint < 0x27FFFF: + return False # AFI + if 0x280000 < icaoint < 0x28FFFF: + return False # SAM + if 0x500000 < icaoint < 0x5FFFFF: + return False # EUR, NAT + if 0x600000 < icaoint < 0x67FFFF: + return False # MID + if 0x680000 < icaoint < 0x6F0000: + return False # ASIA + if 0x900000 < icaoint < 0x9FFFFF: + return False # NAM, PAC + if 0xB00000 < icaoint < 0xBFFFFF: + return False # CAR + if 0xD00000 < icaoint < 0xDFFFFF: + return False # future + if 0xF00000 < icaoint < 0xFFFFFF: + return False # future + + return True + +@cython.boundscheck(False) +@cython.wraparound(False) +cpdef int typecode(str msg): + """Type code of ADS-B message + + Args: + msg (string): 28 bytes hexadecimal message string + + Returns: + int: type code number + """ + if df(msg) not in (17, 18): + return -1 + # return None + + cdef str tcbin = hex2bin(msg[8:10]) + return bin2int(tcbin[0:5]) + +@cython.cdivision(True) +cpdef int cprNL(double lat): + """NL() function in CPR decoding.""" + + if lat == 0: + return 59 + + if lat == 87 or lat == -87: + return 2 + + if lat > 87 or lat < -87: + return 1 + + cdef int nz = 15 + cdef double a = 1 - cos(pi / (2 * nz)) + cdef double b = cos(pi / 180.0 * fabs(lat)) ** 2 + cdef double nl = 2 * pi / (acos(1 - a / b)) + NL = floor(nl) + return NL + +@cython.boundscheck(False) +@cython.wraparound(False) +cpdef str idcode(str msg): + """Compute identity (squawk code). + + Applicable only for DF5 or DF21 messages, bit 20-32. + credit: @fbyrkjeland + + Args: + msg (String): 28 bytes hexadecimal message string + + Returns: + string: squawk code + + """ + if df(msg) not in [5, 21]: + raise RuntimeError("Message must be Downlink Format 5 or 21.") + + cdef bytearray _mbin = bytearray(hex2bin(msg).encode()) + cdef unsigned char[:] mbin = _mbin + + cdef bytearray _idcode = bytearray(4) + cdef unsigned char[:] idcode = _idcode + + cdef unsigned char C1 = mbin[19] + cdef unsigned char A1 = mbin[20] + cdef unsigned char C2 = mbin[21] + cdef unsigned char A2 = mbin[22] + cdef unsigned char C4 = mbin[23] + cdef unsigned char A4 = mbin[24] + # _ = mbin[25] + cdef unsigned char B1 = mbin[26] + cdef unsigned char D1 = mbin[27] + cdef unsigned char B2 = mbin[28] + cdef unsigned char D2 = mbin[29] + cdef unsigned char B4 = mbin[30] + cdef unsigned char D4 = mbin[31] + + # byte1 = int(A4 + A2 + A1, 2) + # byte2 = int(B4 + B2 + B1, 2) + # byte3 = int(C4 + C2 + C1, 2) + # byte4 = int(D4 + D2 + D1, 2) + + idcode[0] = int_to_char((char_to_int(A4)*2 + char_to_int(A2))*2 + char_to_int(A1)) + idcode[1] = int_to_char((char_to_int(B4)*2 + char_to_int(B2))*2 + char_to_int(B1)) + idcode[2] = int_to_char((char_to_int(C4)*2 + char_to_int(C2))*2 + char_to_int(C1)) + idcode[3] = int_to_char((char_to_int(D4)*2 + char_to_int(D2))*2 + char_to_int(D1)) + + return _idcode.decode() + + #return str(byte1) + str(byte2) + str(byte3) + str(byte4) + + +@cython.boundscheck(False) +@cython.wraparound(False) +cpdef int altcode(str msg): + """Compute the altitude. + + Applicable only for DF4 or DF20 message, bit 20-32. + credit: @fbyrkjeland + + Args: + msg (String): 28 bytes hexadecimal message string + + Returns: + int: altitude in ft + + """ + if df(msg) not in [0, 4, 16, 20]: + raise RuntimeError("Message must be Downlink Format 0, 4, 16, or 20.") + + # Altitude code, bit 20-32 + cdef bytearray _mbin = bytearray(hex2bin(msg).encode()) + cdef unsigned char[:] mbin = _mbin + + cdef char mbit = mbin[25] # M bit: 26 + cdef char qbit = mbin[27] # Q bit: 28 + cdef int alt = 0 + cdef bytearray vbin + cdef bytearray _graybytes = bytearray(11) + cdef unsigned char[:] graybytes = _graybytes + + if mbit == 48: # unit in ft, "0" -> 48 + if qbit == 49: # 25ft interval, "1" -> 49 + vbin = _mbin[19:25] + _mbin[26:27] + _mbin[28:32] + alt = bin2int(vbin.decode()) * 25 - 1000 + if qbit == 48: # 100ft interval, above 50175ft, "0" -> 48 + graybytes[8] = mbin[19] + graybytes[2] = mbin[20] + graybytes[9] = mbin[21] + graybytes[3] = mbin[22] + graybytes[10] = mbin[23] + graybytes[4] = mbin[24] + # _ = mbin[25] + graybytes[5] = mbin[26] + # cdef char D1 = mbin[27] # always zero + graybytes[6] = mbin[28] + graybytes[0] = mbin[29] + graybytes[7] = mbin[30] + graybytes[1] = mbin[31] + # graybytes = D2 + D4 + A1 + A2 + A4 + B1 + B2 + B4 + C1 + C2 + C4 + + alt = gray2alt(_graybytes.decode()) + + if mbit == 49: # unit in meter, "1" -> 49 + vbin = _mbin[19:25] + _mbin[26:31] + alt = int(bin2int(vbin.decode()) * 3.28084) # convert to ft + + return alt + + + +cpdef int gray2alt(str codestr): + cdef str gc500 = codestr[:8] + cdef int n500 = gray2int(gc500) + + # in 100-ft step must be converted first + cdef str gc100 = codestr[8:] + cdef int n100 = gray2int(gc100) + + if n100 in [0, 5, 6]: + return -1 + #return None + + if n100 == 7: + n100 = 5 + + if n500 % 2: + n100 = 6 - n100 + + alt = (n500 * 500 + n100 * 100) - 1300 + return alt + + +cdef int gray2int(str graystr): + """Convert greycode to binary.""" + cdef int num = bin2int(graystr) + num ^= num >> 8 + num ^= num >> 4 + num ^= num >> 2 + num ^= num >> 1 + return num + + +cdef str data(str msg): + """Return the data frame in the message, bytes 9 to 22.""" + return msg[8:-6] + + +cpdef bint allzeros(str msg): + """Check if the data bits are all zeros. + + Args: + msg (String): 28 bytes hexadecimal message string + + Returns: + bool: True or False + + """ + d = hex2bin(data(msg)) + + if bin2int(d) > 0: + return False + else: + return True + + +def wrongstatus(data, sb, msb, lsb): + """Check if the status bit and field bits are consistency. + + This Function is used for checking BDS code versions. + + """ + # status bit, most significant bit, least significant bit + status = int(data[sb - 1]) + value = bin2int(data[msb - 1 : lsb]) + + if not status: + if value != 0: + return True + + return False diff --git a/pyModeS/decoder/common.py b/pyModeS/decoder/common.py index 17f55c6..b1edf35 100644 --- a/pyModeS/decoder/common.py +++ b/pyModeS/decoder/common.py @@ -15,32 +15,11 @@ def hex2int(hexstr): return int(hexstr, 16) -def int2hex(n): - """Convert a integer to hexadecimal string.""" - # strip 'L' for python 2 - return hex(n)[2:].rjust(6, "0").upper().rstrip("L") - - def bin2int(binstr): """Convert a binary string to integer.""" return int(binstr, 2) -def bin2hex(hexstr): - """Convert a hexdecimal string to integer.""" - return int2hex(bin2int(hexstr)) - - -def bin2np(binstr): - """Convert a binary string to numpy array.""" - return np.array([int(i) for i in binstr]) - - -def np2bin(npbin): - """Convert a binary numpy array to string.""" - return np.array2string(npbin, separator="")[1:-1] - - def df(msg): """Decode Downlink Format value, bits 1 to 5.""" dfbin = hex2bin(msg[:2]) @@ -102,7 +81,7 @@ def crc_legacy(msg, encode=False): ) ng = len(generator) - msgnpbin = bin2np(hex2bin(msg)) + msgnpbin = np.array([int(i) for i in hex2bin(msg)]) if encode: msgnpbin[-24:] = [0] * 24 @@ -116,7 +95,9 @@ def crc_legacy(msg, encode=False): msgnpbin[i : i + ng] = np.bitwise_xor(msgnpbin[i : i + ng], generator) # last 24 bits - reminder = bin2int(np2bin(msgnpbin[-24:])) + msgbin = np.array2string(msgnpbin[-24:], separator="")[1:-1] + reminder = bin2int(msgbin) + return reminder @@ -148,7 +129,7 @@ def icao(msg): addr = msg[2:8] elif DF in (0, 4, 5, 16, 20, 21): c0 = crc(msg, encode=True) - c1 = hex2int(msg[-6:]) + c1 = int(msg[-6:], 16) addr = "%06X" % (c0 ^ c1) else: addr = None @@ -161,7 +142,7 @@ def is_icao_assigned(icao): if (icao is None) or (not isinstance(icao, str)) or (len(icao) != 6): return False - icaoint = hex2int(icao) + icaoint = int(icao, 16) if 0x200000 < icaoint < 0x27FFFF: return False # AFI diff --git a/pyModeS/decoder/surv.py b/pyModeS/decoder/surv.py index 45f219d..91b8c38 100644 --- a/pyModeS/decoder/surv.py +++ b/pyModeS/decoder/surv.py @@ -1,18 +1,3 @@ -# Copyright (C) 2018 Junzi Sun (TU Delft) - -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. - -# This program 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 . - """ Warpper for short roll call surveillance replies DF=4/5 diff --git a/setup.py b/setup.py index e24ce41..9f62029 100644 --- a/setup.py +++ b/setup.py @@ -15,6 +15,13 @@ Steps for deploying a new version: # Always prefer setuptools over distutils from setuptools import setup, find_packages +# Compile some parts +from setuptools.extension import Extension +from Cython.Build import cythonize + +extensions = [Extension("pyModeS.decoder.c_common", ["pyModeS/decoder/c_common.pyx"])] + + # To use a consistent encoding from codecs import open from os import path @@ -57,6 +64,7 @@ setup( "Programming Language :: Python :: 2", "Programming Language :: Python :: 3", ], + ext_modules=cythonize(extensions), # What does your project relate to? keywords="Mode-S ADS-B EHS ELS Comm-B", # You can just specify the packages manually here if your project is diff --git a/tests/benchmark.py b/tests/benchmark.py new file mode 100644 index 0000000..c7e5276 --- /dev/null +++ b/tests/benchmark.py @@ -0,0 +1,130 @@ +import sys +import time +import pandas as pd +from tqdm import tqdm +from pyModeS.decoder import adsb + +fin = sys.argv[1] + +df = pd.read_csv(fin, names=["ts", "df", "icao", "msg"]) +df_adsb = df[df["df"] == 17].copy() + +total = df_adsb.shape[0] + + +def native(): + + from pyModeS.decoder import common + + # airborne position + m_air_0 = None + m_air_1 = None + + # surface position + m_surf_0 = None + m_surf_1 = None + + for i, r in tqdm(df_adsb.iterrows(), total=total): + ts = r.ts + m = r.msg + + downlink_format = common.df(m) + crc = common.crc(m) + icao = adsb.icao(m) + tc = adsb.typecode(m) + + if 1 <= tc <= 4: + category = adsb.category(m) + callsign = adsb.callsign(m) + if tc == 19: + velocity = adsb.velocity(m) + + if 5 <= tc <= 8: + if adsb.oe_flag(m): + m_surf_1 = m + t1 = ts + else: + m_surf_0 = m + t0 = ts + + if m_surf_0 and m_surf_1: + position = adsb.surface_position( + m_surf_0, m_surf_1, t0, t1, 50.01, 4.35 + ) + altitude = adsb.altitude(m) + + if 9 <= tc <= 18: + if adsb.oe_flag(m): + m_air_1 = m + t1 = ts + else: + m_air_0 = m + t0 = ts + + if m_air_0 and m_air_1: + position = adsb.position(m_air_0, m_air_1, t0, t1) + altitude = adsb.altitude(m) + + +def cython(): + + from pyModeS.decoder import c_common as common + + # airborne position + m_air_0 = None + m_air_1 = None + + # surface position + m_surf_0 = None + m_surf_1 = None + + for i, r in tqdm(df_adsb.iterrows(), total=total): + ts = r.ts + m = r.msg + + downlink_format = common.df(m) + crc = common.crc(m) + icao = adsb.icao(m) + tc = adsb.typecode(m) + + if 1 <= tc <= 4: + category = adsb.category(m) + callsign = adsb.callsign(m) + if tc == 19: + velocity = adsb.velocity(m) + + if 5 <= tc <= 8: + if adsb.oe_flag(m): + m_surf_1 = m + t1 = ts + else: + m_surf_0 = m + t0 = ts + + if m_surf_0 and m_surf_1: + position = adsb.surface_position( + m_surf_0, m_surf_1, t0, t1, 50.01, 4.35 + ) + altitude = adsb.altitude(m) + + if 9 <= tc <= 18: + if adsb.oe_flag(m): + m_air_1 = m + t1 = ts + else: + m_air_0 = m + t0 = ts + + if m_air_0 and m_air_1: + position = adsb.position(m_air_0, m_air_1, t0, t1) + altitude = adsb.altitude(m) + + +if __name__ == "__main__": + t1 = time.time() + native() + dt1 = time.time() - t1 + + t2 = time.time() + cython() + dt2 = time.time() - t2 diff --git a/tests/sample_run_adsb.py b/tests/sample_run_adsb.py index c6c7c16..f0d5134 100644 --- a/tests/sample_run_adsb.py +++ b/tests/sample_run_adsb.py @@ -1,44 +1,46 @@ -from __future__ import print_function -from pyModeS import adsb, ehs +import sys +import time +import csv + +if len(sys.argv) > 1 and sys.argv[1] == "cython": + from pyModeS.c_decoder import adsb +else: + from pyModeS.decoder import adsb + +print("===== Decode ADS-B sample data=====") + +f = open("tests/data/sample_data_adsb.csv", "rt") + +msg0 = None +msg1 = None + +tstart = time.time() +for i, r in enumerate(csv.reader(f)): + + ts = int(r[0]) + m = r[1].encode() + + icao = adsb.icao(m) + tc = adsb.typecode(m) + + if 1 <= tc <= 4: + print(ts, m, icao, tc, adsb.category(m), adsb.callsign(m)) + if tc == 19: + print(ts, m, icao, tc, adsb.velocity(m)) + if 5 <= tc <= 18: + if adsb.oe_flag(m): + msg1 = m + t1 = ts + else: + msg0 = m + t0 = ts + + if msg0 and msg1: + pos = adsb.position(msg0, msg1, t0, t1) + alt = adsb.altitude(m) + print(ts, m, icao, tc, pos, alt) -# === Decode sample data file === +dt = time.time() - tstart - -def adsb_decode_all(n=None): - print("===== Decode ADS-B sample data=====") - import csv - - f = open("tests/data/sample_data_adsb.csv", "rt") - - msg0 = None - msg1 = None - - for i, r in enumerate(csv.reader(f)): - if n and i > n: - break - - ts = r[0] - m = r[1] - icao = adsb.icao(m) - tc = adsb.typecode(m) - if 1 <= tc <= 4: - print(ts, m, icao, tc, adsb.category(m), adsb.callsign(m)) - if tc == 19: - print(ts, m, icao, tc, adsb.velocity(m)) - if 5 <= tc <= 18: - if adsb.oe_flag(m): - msg1 = m - t1 = ts - else: - msg0 = m - t0 = ts - - if msg0 and msg1: - pos = adsb.position(msg0, msg1, t0, t1) - alt = adsb.altitude(m) - print(ts, m, icao, tc, pos, alt) - - -if __name__ == "__main__": - adsb_decode_all(n=100) +print("Execution time: {} seconds".format(dt)) diff --git a/tests/test_c_common.py b/tests/test_c_common.py new file mode 100644 index 0000000..71599b9 --- /dev/null +++ b/tests/test_c_common.py @@ -0,0 +1,60 @@ +from pyModeS.decoder import c_common as common + + +def test_conversions(): + assert common.hex2bin("6E406B") == "011011100100000001101011" + + +def test_crc_decode(): + + assert common.crc("8D406B902015A678D4D220AA4BDA") == 0 + assert common.crc("8d8960ed58bf053cf11bc5932b7d") == 0 + assert common.crc("8d45cab390c39509496ca9a32912") == 0 + assert common.crc("8d74802958c904e6ef4ba0184d5c") == 0 + assert common.crc("8d4400cd9b0000b4f87000e71a10") == 0 + assert common.crc("8d4065de58a1054a7ef0218e226a") == 0 + + assert common.crc("c80b2dca34aa21dd821a04cb64d4") == 10719924 + assert common.crc("a800089d8094e33a6004e4b8a522") == 4805588 + assert common.crc("a8000614a50b6d32bed000bbe0ed") == 5659991 + assert common.crc("a0000410bc900010a40000f5f477") == 11727682 + assert common.crc("8d4ca251204994b1c36e60a5343d") == 16 + assert common.crc("b0001718c65632b0a82040715b65") == 353333 + + +def test_crc_encode(): + parity = common.crc("8D406B902015A678D4D220AA4BDA", encode=True) + assert parity == 11160538 + + +def test_icao(): + assert common.icao("8D406B902015A678D4D220AA4BDA") == "406B90" + assert common.icao("A0001839CA3800315800007448D9") == "400940" + assert common.icao("A000139381951536E024D4CCF6B5") == "3C4DD2" + assert common.icao("A000029CFFBAA11E2004727281F1") == "4243D0" + + +def test_modes_altcode(): + assert common.altcode("A02014B400000000000000F9D514") == 32300 + + +def test_modes_idcode(): + assert common.idcode("A800292DFFBBA9383FFCEB903D01") == "1346" + + +def test_graycode_to_altitude(): + assert common.gray2alt("00000000010") == -1000 + assert common.gray2alt("00000001010") == -500 + assert common.gray2alt("00000011011") == -100 + assert common.gray2alt("00000011010") == 0 + assert common.gray2alt("00000011110") == 100 + assert common.gray2alt("00000010011") == 600 + assert common.gray2alt("00000110010") == 1000 + assert common.gray2alt("00001001001") == 5800 + assert common.gray2alt("00011100100") == 10300 + assert common.gray2alt("01100011010") == 32000 + assert common.gray2alt("01110000100") == 46300 + assert common.gray2alt("01010101100") == 50200 + assert common.gray2alt("11011110100") == 73200 + assert common.gray2alt("10000000011") == 126600 + assert common.gray2alt("10000000001") == 126700 diff --git a/tests/test_common.py b/tests/test_common.py index df2ac36..571b0cf 100644 --- a/tests/test_common.py +++ b/tests/test_common.py @@ -1,10 +1,8 @@ -from pyModeS import common +from pyModeS.decoder import common def test_conversions(): assert common.hex2bin("6E406B") == "011011100100000001101011" - assert common.bin2hex("011011100100000001101011") == "6E406B" - assert common.int2hex(11160538) == "AA4BDA" def test_crc_decode(): @@ -28,7 +26,7 @@ def test_crc_decode(): def test_crc_encode(): parity = common.crc("8D406B902015A678D4D220AA4BDA", encode=True) - assert common.int2hex(parity) == "AA4BDA" + assert parity == 11160538 def test_icao(): diff --git a/tox.ini b/tox.ini deleted file mode 100644 index 9a0ca74..0000000 --- a/tox.ini +++ /dev/null @@ -1,11 +0,0 @@ -[tox] -toxworkdir = /tmp/tox -envlist = py2,py3 - -[testenv] -deps = - pytest - numpy - pyzmq - pyrtlsdr -commands = py.test