Geekly articles weekly

Flightradar24 - bagaimana cara kerjanya? Bagian 2, protokol ADS-B

Saya akan menebak dan mengatakan bahwa setiap orang yang teman atau keluarganya pernah terbang dengan pesawat, telah menggunakan Flightradar24 - layanan gratis dan nyaman untuk melacak penerbangan secara real time.

gambar

Pada bagian pertama ide-ide dasar operasi dijelaskan. Sekarang mari kita melangkah lebih jauh dan mencari tahu, data apa yang sebenarnya mengirim dan menerima antara pesawat dan stasiun darat. Kami juga akan mendekode data ini menggunakan Python.

Sejarah


Seharusnya jelas bahwa data pesawat tidak dimaksudkan hanya untuk pengguna, yang ingin melihatnya di smartphone mereka. Sistem ini disebut ADS - B (Siaran pengawasan otomatis - siaran), dan dirancang untuk mentransmisikan data informasi pesawat secara otomatis ke pusat kendali - parameter yang berbeda, seperti koordinat, kecepatan, heading, ketinggian, dan lainnya, sedang dikirim. Sebelumnya, operator hanya bisa melihat titik di layar radar. Dan itu pasti tidak cukup ketika jumlah pesawat meningkat secara drastis.

Secara teknis, ADS-B terdiri dari pemancar di dalam pesawat, yang secara berkala mengirimkan bingkai data informasi pada frekuensi yang relatif tinggi 1090 MHz (ada beberapa mode lain, tetapi mereka tidak begitu menarik bagi kami, karena koordinatnya hanya dikirimkan di sini ) Tentu saja, ada juga penerima di suatu tempat di bandara, tetapi bagi kami, seperti bagi pengguna, penerima kami sendiri lebih menarik.

Sebagai perbandingan, sistem pertama yang dirancang untuk pengguna biasa, Airnav Radarbox, muncul pada 2007, dan berharga sekitar $ 900, dan sekitar $ 250 per tahun dikenakan biaya berlangganan layanan jaringan mereka (ide utama sistem ini adalah untuk mengumpulkan dan berbagi data dari banyak penerima , penerima yang berdiri sendiri relatif tidak berguna).



Sekarang, ketika penerima RTL-SDR menjadi lebih banyak tersedia, perangkat serupa dapat dibuat seharga $ 30. Itu dapat ditemukan di bagian pertama artikel , dan kita akan melangkah lebih jauh dan menggambarkan protokol itu sendiri - mari kita lihat cara kerjanya.

Menerima sinyal


Pertama, kita perlu merekam sampel sinyal. Seluruh sinyal hanya memiliki panjang 120 mikrodetik, dan untuk melihat rinciannya dalam "resolusi" yang baik, lebih baik untuk memiliki radio SDR dengan sampel rate minimal 5 MHz.

gambar

Setelah merekam, kami mendapatkan file WAV dengan sampling rate 5.000.000 sampel / detik, 30 detik dari catatan tersebut memiliki ukuran sekitar 500MB. Tentu saja, tidak ada gunanya untuk mendengarkannya dengan pemutar media favorit Anda - file tersebut tidak mengandung suara, tetapi sinyal radio yang langsung didigitalkan sendiri - ini persis cara Radio Didefinisikan Perangkat Lunak bekerja.

Kami dapat membuka dan memproses file ini dengan Python. Mereka yang ingin melakukan percobaan sendiri, dapat mengunduh sampel dari tautan ini .

Mari kita buka file dan lihat, apa yang ada di dalamnya.

from scipy.io import wavfile import matplotlib.pyplot as plt import numpy as np fs, data = wavfile.read("adsb_20190311_191728Z_1090000kHz_RF.wav") data = data.astype(float) I, Q = data[:, 0], data[:, 1] A = np.sqrt(I*I + Q*Q) plt.plot(A) plt.show()

Hasil: kita melihat beberapa 'impuls' atas kebisingan.



Setiap "impuls" adalah sinyal, yang strukturnya terlihat jelas jika kita meningkatkan resolusi pada grafik.



Seperti yang dapat kita lihat, gambar sepenuhnya konsisten dengan uraian di atas. Kami sekarang dapat memproses data ini.

Decoding


Pertama, kita perlu memiliki sedikit aliran. Sinyal itu sendiri dikodekan dengan pengkodean manchester:



Dari perbedaan setengah gigitan kita dapat dengan mudah mendapatkan "0" dan "1" yang nyata.

  bits_str = "" for p in range(8): pos = start_data + bit_len*p p1, p2 = A[pos: pos + bit_len/2], A[pos + bit_len/2: pos + bit_len] avg1, avg2 = np.average(p1), np.average(p2) if avg1 < avg2: bits_str += "0" elif avg1 > avg2: bits_str += "1"

Struktur sinyal itu sendiri terlihat seperti ini:



Mari kita lihat isinya lebih detail.

DF (Downlink Format, 5 bit) - mendefinisikan jenis pesan. Ada beberapa jenis:


( sumber halaman )

Kami hanya tertarik pada tipe DF17, karena hanya yang ini berisi koordinat pesawat.

ICAO (24 bit) - adalah kode pesawat internasional yang unik. Kami dapat memeriksa pesawat dengan kodenya di situs web ini (sayangnya, penulis telah berhenti memperbarui database, tetapi masih relevan). Misalnya, untuk kode 3c5ee2 kita dapat memiliki informasi berikut:



DATA (56 atau 112 bit) - adalah data itu sendiri, yang akan kami dekode. 5 bit pertama data adalah bidang Kode Jenis , yang berisi subtipe dari data yang disimpan (tidak dicampur dengan bidang DF ). Ada banyak jenis seperti itu:


( sumber tabel )

Mari kita lihat beberapa contoh.

Identifikasi pesawat

Contoh dalam bentuk biner:

00100 011 000101 010111 000111 110111 110001 111000

Bidang data:

 +------+------+------+------+------+------+------+------+------+------+ | TC,5 | EC,3 | C1,6 | C2,6 | C3,6 | C4,6 | C5,6 | C6,6 | C7,6 | C8,6 | +------+------+------+------+------+------+------+------+------+------+

TC = 00100b = 4, dan setiap simbol C1-C8 berisi kode, yang harus cocok dengan indeks dalam string ini:
#ABCDEFGHIJKLMNOPQRSTUVWXYZ ##### _ ################ 0123456789 ######

Setelah decoding, mudah untuk mendapatkan nama pesawat: EWG7184

 symbols = "#ABCDEFGHIJKLMNOPQRSTUVWXYZ#####_###############0123456789######" code_str = "" for p in range(8): c = int(bits_str[8 + 6*p:8 + 6*(p + 1)], 2) code_str += symbols[c] print("Aircraft Identification:", code_str.replace('#', ''))

Posisi di udara

Nama decoding sederhana, tetapi koordinatnya lebih rumit. Mereka ditransmisikan dalam bentuk 2 frame, genap dan ganjil. Kode bidang TC = 01011b = 11.



Contoh frame data genap dan ganjil:

 01011 000 000101110110 00 10111000111001000 10000110101111001 01011 000 000110010000 01 10010011110000110 10000011110001000

Perhitungan koordinat itu sendiri menggunakan rumus agak rumit:


( sumber )

Saya bukan ahli GIS, jadi saya tidak tahu dari mana asalnya. Siapa tahu lebih baik, silakan tulis komentar.

Perhitungan ketinggian lebih sederhana - tergantung pada bit tertentu, dapat direpresentasikan sebagai kelipatan 25 atau 100 kaki.

Kecepatan udara

Dataframe dengan TC = 19. Yang menarik di sini adalah bahwa kecepatan dapat relatif ke tanah (lebih akurat, disebut Kecepatan Tanah), dan Kecepatan Udara, diukur oleh sensor udara pesawat (bisa kurang akurat karena angin). Banyak bidang lain yang berbeda juga ditransmisikan:


( sumber )

Kesimpulan


Seperti yang dapat kita lihat, teknologi ADS-B telah menjadi simbiosis yang menarik, ketika standar dapat digunakan secara luas tidak hanya untuk para profesional, tetapi juga untuk pengguna biasa. Tapi yang pasti, peran utama dalam hal ini dilakukan oleh semakin murahnya teknologi penerima SDR digital, yang memungkinkan untuk menerima sinyal dengan frekuensi lebih tinggi daripada gigahertz pada perangkat yang sangat murah.

Standar itu sendiri, tentu saja, memiliki lebih banyak data. Mereka yang tertarik dapat melihat PDF di halaman ICAO atau mengunjungi situs web mode-s.org yang telah disebutkan di atas. Sepertinya artikel ini tidak akan benar-benar digunakan oleh sebagian besar pembaca, tetapi saya harap, setidaknya gagasan umum tentang cara kerjanya, sekarang lebih jelas.

Omong-omong, dekoder ADS-B Python sudah ada, dapat diselidiki di github . Pemilik penerima SDR juga dapat membangun dan menjalankan dekoder ADS-B yang siap digunakan dari halaman ini , dan (saya akan ulangi lagi) beberapa detail yang kami juga di bagian pertama artikel ini.

Kode sumber parser, dijelaskan di atas, diberikan di bawah spoiler. Ini hanya contoh uji yang tidak berpura-pura dengan kualitas produksi, tetapi secara umum berfungsi, dan dapat mengurai file WAV , yang dicatat di atas.

Kode Sumber (Python)
 from __future__ import print_function from scipy.io import wavfile from scipy import signal import matplotlib.pyplot as plt import numpy as np import math import sys def parse_message(data, start, bit_len): max_len = bit_len*128 A = data[start:start + max_len] A = signal.resample(A, 10*max_len) bits = np.zeros(10*max_len) bit_len *= 10 start_data = bit_len*8 # Parse first 8 bits bits_str = "" for p in range(8): pos = start_data + bit_len*p p1, p2 = A[pos: pos + bit_len/2], A[pos + bit_len/2: pos + bit_len] avg1, avg2 = np.average(p1), np.average(p2) if avg1 < avg2: bits_str += "0" elif avg1 > avg2: bits_str += "1" df = int(bits_str[0:5], 2) # Aircraft address (db - https://junzis.com/adb/?q=3b1c5c ) bits_str = "" for p in range(8, 32): pos = start_data + bit_len * p p1, p2 = A[pos: pos + bit_len / 2], A[pos + bit_len / 2: pos + bit_len] avg1, avg2 = np.average(p1), np.average(p2) if avg1 < avg2: bits_str += "0" elif avg1 > avg2: bits_str += "1" # print "Aircraft address:", bits_str, hex(int(bits_str, 2)) address = hex(int(bits_str, 2)) # Filter specific aircraft (optional) # if address != "0x3c5ee2": # return if df == 16 or df == 17 or df == 18 or df == 19 or df == 20 or df == 21: # print "Pos:", start, "DF:", msg_type # Data (56bit) bits_str = "" for p in range(32, 88): pos = start_data + bit_len*p p1, p2 = A[pos: pos + bit_len/2], A[pos + bit_len/2: pos + bit_len] avg1, avg2 = np.average(p1), np.average(p2) if avg1 < avg2: bits_str += "0" # bits[pos + bit_len / 2] = 50 elif avg1 > avg2: bits_str += "1" # http://www.lll.lu/~edward/edward/adsb/DecodingADSBposition.html # print "Data:" # print bits_str[:8], bits_str[8:20], bits_str[20:22], bits_str[22:22+17], bits_str[39:39+17] # Type Code: tc, ec = int(bits_str[:5], 2), int(bits_str[5:8], 2) # print("DF:", df, "TC:", tc) # 1 - 4 Aircraft identification # 5 - 8 Surface position # 9 - 18 Airborne position (w/ Baro Altitude) # 19 Airborne velocities if tc >= 1 and tc <= 4: # and (df == 17 or df == 18): print("Aircraft address:", address) print("Data:") print(bits_str[:8], bits_str[8:14], bits_str[14:20], bits_str[20:26], bits_str[26:32], bits_str[32:38], bits_str[38:44]) symbols = "#ABCDEFGHIJKLMNOPQRSTUVWXYZ#####_###############0123456789######" code_str = "" for p in range(8): c = int(bits_str[8 + 6*p:8 + 6*(p + 1)], 2) code_str += symbols[c] print("Aircraft Identification:", code_str.replace('#', '')) print() if tc == 11: print("Aircraft address:", address) print("Data: (11)") print(bits_str[:8], bits_str[8:20], bits_str[20:22], bits_str[22:22+17], bits_str[39:39+17]) # Bit 22 contains the F flag which indicates which CPR format is used (odd or even) # First frame has F flag = 0 so is even and the second frame has F flag = 1 so odd # f = bits_str[21:22] # print("F:", int(f, 2)) # Altitude alt1b = bits_str[8:20] if alt1b[-5] == '1': bits = alt1b[:-5] + alt1b[-4:] n = int(bits, 2) alt_ft = n*25 - 1000 print("Alt (ft)", alt_ft) # lat_dec = int(bits_str[22:22+17], 2) # lon_dec = int(bits_str[39:39+17], 2) # print("Lat/Lon:", lat_dec, lon_dec) # http://airmetar.main.jp/radio/ADS-B%20Decoding%20Guide.pdf print() if tc == 19: print("Aircraft address:", address) print("Data:") # print(bits_str) print(bits_str[:5], bits_str[5:8], bits_str[8:10], bits_str[10:13], bits_str[13] ,bits_str[14:24], bits_str[24], bits_str[25:35], bits_str[35:36], bits_str[36:65]) subtype = int(bits_str[5:8], 2) # https://mode-s.org/decode/adsb/airborne-velocity.html spd, hdg, rocd = -1, -1, -1 if subtype == 1 or subtype == 2: print("Velocity Subtype 1: Ground speed") v_ew_sign = int(bits_str[13], 2) v_ew = int(bits_str[14:24], 2) - 1 # east-west velocity v_ns_sign = int(bits_str[24], 2) v_ns = int(bits_str[25:35], 2) - 1 # north-south velocity v_we = -1*v_ew if v_ew_sign else v_ew v_sn = -1*v_ns if v_ns_sign else v_ns spd = math.sqrt(v_sn*v_sn + v_we*v_we) # unit in kts hdg = math.atan2(v_we, v_sn) hdg = math.degrees(hdg) # convert to degrees hdg = hdg if hdg >= 0 else hdg + 360 # no negative val if subtype == 3: print("Subtype Subtype 3: Airspeed") hdg = int(bits_str[14:24], 2)/1024.0*360.0 spd = int(bits_str[25:35], 2) vr_sign = int(bits_str[36], 2) vr = int(bits_str[36:45], 2) rocd = -1*vr if vr_sign else vr # rate of climb/descend print("Speed (kts):", spd, "Rate:", rocd, "Heading:", hdg) print() # print() def calc_coordinates(): def _cprN(lat, is_odd): nl = _cprNL(lat) - is_odd return nl if nl > 1 else 1 def _cprNL(lat): try: nz = 15 a = 1 - math.cos(math.pi / (2 * nz)) b = math.cos(math.pi / 180.0 * abs(lat)) ** 2 nl = 2 * math.pi / (math.acos(1 - a/b)) return int(math.floor(nl)) except: # happens when latitude is +/-90 degree return 1 def floor_(x): return int(math.floor(x)) lat1b, lon1b, alt1b = "10111000111010011", "10000110111111000", "000101111001" lat2b, lon2b, alt2b = "10010011101011100", "10000011000011011", "000101110111" lat1, lon1, alt1 = int(lat1b, 2), int(lon1b, 2), int(alt1b, 2) lat2, lon2, alt2 = int(lat2b, 2), int(lon2b, 2), int(alt2b, 2) # 131072 is 2^17, since CPR lat and lon are 17 bits each cprlat_even, cprlon_even = lat1/131072.0, lon1/131072.0 cprlat_odd, cprlon_odd = lat2/131072.0, lon2/131072.0 print(cprlat_even, cprlon_even) j = floor_(59*cprlat_even - 60*cprlat_odd) print(j) air_d_lat_even = 360.0 / 60 air_d_lat_odd = 360.0 / 59 # Lat lat_even = float(air_d_lat_even * (j % 60 + cprlat_even)) lat_odd = float(air_d_lat_odd * (j % 59 + cprlat_odd)) if lat_even >= 270: lat_even = lat_even - 360 if lat_odd >= 270: lat_odd = lat_odd - 360 # Lon ni = _cprN(lat_even, 0) m = floor_(cprlon_even * (_cprNL(lat_even)-1) - cprlon_odd * _cprNL(lat_even) + 0.5) lon = (360.0 / ni) * (m % ni + cprlon_even) print("Lat", lat_even, "Lon", lon) # Altitude # Q-bit (bit 48) indicates whether the altitude is encoded in multiples of 25 or 100 ft (0: 100 ft, 1: 25 ft) # The value can represent altitudes from -1000 to +50175 ft. if alt1b[-5] == '1': bits = alt1b[:-5] + alt1b[-4:] n = int(bits, 2) alt_ft = n*25 - 1000 print("Alt (ft)", alt_ft) fs, data = wavfile.read("adsb_20190311_191728Z_1090000kHz_RF.wav") T = 1/fs print("Sample rate %f MS/s" % (fs / 1e6)) print("Cnt samples %d" % len(data)) print("Duration: %fs" % (T * len(data))) data = data.astype(float) cnt = data.shape[0] # Processing only part on file (faster): # cnt = 10000000 # data = data[:cnt] print("Processing I/Q...") I, Q = data[:, 0], data[:, 1] A = np.sqrt(I*I + Q*Q) bits = np.zeros(cnt) # To see scope without any processing, uncomment # plt.plot(A) # plt.show() # sys.exit(0) print("Extracting signals...") pos = 0 avg = 200 msg_start = 0 # Find beginning of each signal while pos < cnt - 16*1024: # P1 - message start while pos < cnt - 16*1024: if A[pos] < avg and A[pos+1] > avg and pos - msg_start > 1000: msg_start = pos bits[pos] = 100 pos += 4 break pos += 1 start1, start2, start3, start4 = msg_start, 0, 0, 0 # P2 while pos < cnt - 16*1024: if A[pos] < avg and A[pos+1] > avg: start2 = pos bits[pos] = 90 pos += 1 break pos += 1 # P3 while pos < cnt - 16*1024: if A[pos] < avg and A[pos+1] > avg: start3 = pos bits[pos] = 80 pos += 1 break pos += 1 # P4 while pos < cnt - 16*1024: if A[pos] < avg and A[pos+1] > avg: start4 = pos bits[pos] = 70 pos += 1 break pos += 1 sig_diff = start4 - start1 if 20 < sig_diff < 25: bits[msg_start] = 500 bit_len = int((start4 - start1) / 4.5) # print(pos, start1, start4, ' - ', bit_len) # start = start1 + 8*bit_len parse_message(A, msg_start, bit_len) pos += 450 # For debugging: check signal start # plt.plot(A) # plt.plot(bits) # plt.show()


Semoga bermanfaat, terima kasih sudah membaca.

Source: https://habr.com/ru/post/id447078/

More articles:


All Articles