Sunday, May 13, 2018

Some interesting FSK signals found on HF

The following FSK modulated signals were picked up recently using KiwiSDRs located at DF0KL and at Newport,OR using kiwirecorder.py

All signals use a format similar, but not identical, to STANAG 5065 where frames are delimited by LFSR-generated pseudo-random sequences. These formats may also be related to the patent EP2220799A2.

KiwiSDR location Frequency* (kHz) Date UTC Shift (Hz) Baud Format
Newport 4905.0 4/2/18 08:49 850 50 (3)
Newport 4985.0 4/2/18 08:47 850 50 (1)
DF0KL 6376.3 4/2/18 08:32 850 50 (2)
DF0KL 13419.8 4/5/18 14:46 850 50 (2)
DF0KL 16123 4/2/18 10:35 850 50 (1)
Newport 7455.0 4/2/18 08:56 850 50 (1)
DF0KL 8518.8 4/2/18 08:35 850 50 (2)
*center frequency between mark and shift

The recorded IQ data streams were demodulated using a combination of octave and C++ code, see https://github.com/hcab14/signal-analysis

For testing if a given bit sequence is generated by a LFSR, a simple method was used which is described, e.g., in this thesis. This method uses Gaussian elimination over GF(2) for finding the the generating polynomial (which is not guaranteed to be minimal, unlike the result of the Berlekamp–Massey algorithm).

All signals shown in the table above consist of 7- or 21-bit long frames which are marked by LFSR-generated pseudo-random sequences:

Format (1): 7-bit frame delimited by F1
X X X X X X F1
1 2 3 4 5 6 7


Format (2): 21-bit frame delimited by F2, F3
X X X X X X X X X X X X X F2 X X X X X 1 F3
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21


Format (3): 7-bit frame delimited by F4
X X X X X X F4
1 2 3 4 5 6 7

The pseudo-random sequences F1, F2, F3, F4 are generated the following polynomials:

sequence polynomial period
F1 P1=1+x28+x31 231-1
F2 P2=1+x+x6 26-1
F3 P3=1+x+x7 27-1
F4=F1 P4=(1+x)P1 231-1


For format (2) the pseudo-random sequences could have been found by autocorrelation as they are short, whereas for the other formats the length of the pseudo-random sequence is 231-1, i.e., at 50/7 baud the sequences repeats after about 3480 days only.

Format (2) which can be easily identified as a peak at position 21 in the bit stream autocorrelation is also described, e.g., here and here.

Wednesday, January 24, 2018

Using an extended Kalman filter for TDoA HF geo-location

In this post the precision of geo-location using TDoA on shortwave is explored using the data shown here before.

The idea is to feed the time series of time differences into an extended Kalman filter (EKF) using two models for the propagation time from the transmitter to the receiver:

1) Simple ray-geometry for oblique propagation in the ionosphere:
  • based on arXiv:1104.2248 which uses the Breit and Tuve + Martyn theorems to convert vertical to oblique incidence
  • using vertical electron density profiles from the IRI2016 model as an input
  • group refractive index: 1/sqrt(1-fN2/f2), i.e., no magentoionic effects are taken into account
  • when integrating over the group refractive index in the vertical direction, analytic integration of 1/sqrt(f(h)) is performed in each height interval in order to avoid the spurious numerical singularity at the reflection point
  • spherical coordinates are used (rE=6371km)
2) Alternatively, great-circle distances on the WGS84 ellipsoid are used:

All receiver coordinates are assumed to be in WGS84, however in practice it is suspected that there are deviations from the true locations of up to 20 km, as the KiwiSDR positions are provided by the KiwiSDR owners and are not taken from the GPS position solutions used for time synchronization.


The plot below shows the output of the EKF
  • implementation in octave following ADA285972
  • state = [lat, lon] of the unknown position
  • measurements = time differences between pairs of KiwiSDRs
  • the measurements covariance matrix (constant) is build from the covariances of time differences (and is set to zero for time differences not sharing a common receiver)
  • state covariance matrix = diag(0.01, 0.01)2 deg2 (constant)
  • start covariance = diag(0.1, 0.1)2 deg2
  • there are about 60 measurements which are run 5 times sequentially through the EKF

Extended Kalman Filter (EKF) output; note that the range of the x(y) axis corresponds to approximately 17(67)km.

The simple ionosphere propagation model improves the match to the known position but there is still a remaining bias of about 30km North, 4km West. Clearly, more studies are needed in order to find out if this improvement is by chance or real.


Thursday, January 18, 2018

Cyprus OTHR

Yesterday I found the OTHR from Cyprus operating with 25 Hz pulse repetition rate and μ=500kHz/s sweep rate. The most likely geo-location is just north of the known location:

TDoA plots 17030 kHz 20180113T1126Z

Cross-correlations 173030 kHz 20180113T1126Z
I was thinking about an alternative way of geo-location using an extended Kalman filter, and found this reference: www.dtic.mil/dtic/tr/fulltext/u2/a285972.pdf. It is written very well and I was able to implement it quickly in octave using dfpdp from the optim toolbox for the Jacobian. More on this later.

Wednesday, January 17, 2018

OTHR on 6950 kHz

TDoA maps 6950 KHz 20180112T2138Z 

Cross-correlations 6950 kHz 20180112T2138Z 
The bandwidth of this signal is larger than the 12kHz seen by the KiwiSDR. The auto-correlation for this signal has peaks at multiples of about 0.0462938428 sec but it is interesting that the pulses look like this:
abs(z) vs. time 6950 kHz 20180112T2138Z
so an interesting form of modulation is being used.

Friday, January 5, 2018

De-chirping applied to a HF over the horizon radar signal (2)

This is a follow up to this blog post. Again the signals of the OTHR in Cyprus were analyzed using GPS time-stamped IQ samples from several KiwiSDRs.

In this case the OTHR signal was found between 13965 and 13985 kHz. The slope of the FMCW signal is μ=1MHz/s, and there are 50 sweeps per second. This results in a bandwidth of about 20kHz.

As is well known, the frequencies Δf of the spectrum of the de-chirped signal correspond to time delays τ=-Δf/μ in the case of a stationary target.
  • for this one has to multiply the received signal with the original waveform
  • tuning offsets w.r.t. center frequency of the OTHR signal correspond to a time delay of the original waveform
  • when Δf is positive the delay corresponds to the frequency Δf-fs (fs is the sampling frequency), i.e. the ranges wrap around the x-axis.
Applying a 2nd Fourier transform to the time direction of time vs. Δf data one obtains the Doppler frequency fd vs. Δf where fd = v•λ/2 (for velocity v and wavelength λ). Here one has to be careful to adjust the phases in each bin to correct for the fact that the sampling times are not exactly aligned to multiples of 0.02 seconds.

Plots for 5 KiwiSDR are shown below. The splitting of the ionospheric reflection is most likely due to O- and X-mode propagation in the F layer.

13975 kHz 20180105T1114Z @CS5SEL; the band at fd=0 is caused by an interference

13975 kHz 20180105T1114Z @DF0KL

13975 kHz 20180105T1114Z @G8JNJ

13975 kHz 20180105T1114Z @Izhvesk

13975 kHz 20180105T1114Z @Julisdalen
Thanks to all KiwiSDR owners who make their receivers available and do connect a GPS antenna.

Sunday, December 17, 2017

TDoA measurements using GPS time-stamped IQ samples from KiwiSDRs (8)

Today I picked up a HF radar signal on 13500+-20 kHz. These are the time-differences between four KiwiSDRs:
TDoA 13490 kHz

This is most likely the Calypso CODAR radar system:
TDoA maps for 14390 kHz
Within the TDoA resolution the three sites in Ta'Sopu (Gozo), Ta'Barkat (Malta) and Pozzallo Harbour (Sicily) cannot be distinguished.

It is still puzzling for me that the TDoA maps seem to provide reasonable results, as for now great-circle distances are used, i.e., delays due to skywave propagation are not taken into account.