Investigation of quarter-wave folded balun operation

For a about a month, I've been trying to understand how the quarter-wave folded balun works. There are some explanations online, such as this one at Stack Exchange or this one at antenna-theory.com. There are many others as well.

I kinda sorta understand what those explanations are getting at, but I still felt like I didn't really understand what is going on. The explanations seem to rely on experience and intuition of people who have worked with these things. I don't have that that experience or intuition, so I tried modelling the situation to gain insight. All computations in this article are made with scikit-fem.

So what is a quarter-wave folded balun and what is it used for? In short, a balun is a device for converting single-ended signals to differential signals, and back again. Though, what exactly a single ended signal or a differential signal means differs from situation to situation. Talking about signal as voltage, then typically a single ended signal has two conductors, one of which carries the signal, while the other one is a reference level (like ground). In a differential signal one has three conductors (or just two, because of a caveat), in which one conductor carries a positive signal, one carries a negative signal and the third one is the reference. The caveat in this case is, that the reference level can be inferred as the average of the two signal levels, and the reference conductor is not strictly needed.

Coaxial cables are often used as feedlines for antennas. They have the nice feature, that when the currents in the coax are differential, i.e. the current in the shield is equal but opposite to the current in the inner conductor, and when the shield is at ground potential, then the eletromagnetic fields carrying the signal are entirely contained within the coax. This means that passing the cable near metals, or touching the cable by hand doesn't affect the signal it carries. If the currents are not balanced, though, then external fields are produced. These fields will interact with the surroundings, changing the signal carried. The fields also radiate, causing the coax itself to act as an antenna, changing the overall behavior of the system.

So how is it possible that the currents in the coax are not equal? I mean, where else could the current from the inner conductor go other than to the shield? Well, the currents can return through some other conductor, such as the ground. For instance, in the case of an antenna the currents will vanish into thin air and pass to ground through the radiation resistance. If the impedance to ground from the terminals of the antenna are equal, then driving the antenna differentially (with respect to ground) will result in an equal but opposite current in the terminals, and thus also in the feed line. However, since the coax shield needs to remain at ground potential to avoid an external electric field, it isn't possible to have the shield directly drive one side of a differential signal. This is where the balun steps in. The balun will take in a single-ended signal from the coax, and transform it into a differential signal, which is suitable for the load.

To begin the investigation, let's look at the bare coax setup below (figure modified from the work of Omegatron on Wikimedia Commons).

Baseline setup without the balun

The source is single ended, which feeds the inner conductor on the left side of the coaxial cable. The frequency from the source is 868 MHz. The source also has an output impedance of 50 ohms, which is the same as the coax characteristic impedance. The coaxial cable shield is connected to ground at the left side. The load is on the right side of the coaxial cable. It consists of two 25 ohm resistors: one is connected between the coaxial inner conductor and ground and the other is connected between the coaxial shield and ground. This load looks like a 50 ohm resistor for differential signals and like a 25 ohm resistor for single-ended signals. Further, to control stray capacitances and inductances, and to act as a reference, the cable is placed over a ground plane.

The situation is modeled as a multiconductor transmission line. You should check out my previous post about the subject, which is way more involved regarding the mathematics and modelling. The transmission line parameters were computed using FEM. Dimensions and materials used were those of RG-316 coax suspended 20 mm above the ground plane.

Simulating this configuration we notice that there is a significant dependence on the length of the coaxial cable. I don't think this comes as a huge surprise for anyone with experience with transmission lines and coaxial cables. It turns out, that by choosing the coax length properly, the load can be matched to the coax to get a very good power transfer. However, it doesn't seem possible to both match the load to the coax and prevent external fields.

As a first example at 860 mm length there is very poor matching of the load. No current is passing in the shield-side load resistor. This results in a large current imbalance between the conductors, and leads to large induced voltages in the shield. There are thus both large magnetic and electric fields external to the coax. Some quantitative parameters to characterize this situation:

• Maximum common-mode current amplitude in coax: 14.47 mA
• Maximum voltage amplitude in coax shield: 6.92 V
• Common-mode current to load: 14.47 mA
• SWR: 2.00

Voltage and current distribution in coax of length 860mm

The situation is quite different at 946mm length. Here current flows to the load from both the inner conductor and the shield, and the currents are approximately equal and opposite, which results in a mostly balanced current in the coaxial cable. However, the differential signal at the load is achieved through a voltage oscillation in the shield, which leads to a significant electric field external to the coax. The observed parameters at this situation are

• Maximum common-mode current amplitude in coax: 1.05 mA
• Maximum voltage amplitude in coax shield: 0.50 V
• Common-mode current to load: 0.053 mA
• SWR: 1.00

Voltage and current distribution in coax of length 946mm

Computing the monitored parameters over a range of coax lengths reveals a pattern, which repeats every half wavelength. We also observe, that there is no length which would get rid of the significant voltage in the shield. The best one can do is what is represented by the 946mm coax length case.

Observed parameters against coax length

Next, we'll look at the same situation with the quarter-wave folded balun. See the illustration below (figure again modified from the work of Omegatron on Wikimedia Commons).

Quarter-wave folded balun setup

The feed and load remain exactly as in the first case. What has changed, however, is the addition of a shorting stub, which connects the inner conductor of the coax to the shield of the coax. This shorting piece is made of the same type of coax, but doesn't have its inner conductor connected anywhere. This situation can be modelled as two multiconductor transmission lines connected at the shorting wire interface. We have a 2-conductor (plus ground) line on the left side and a 3-conductor (plus ground) line on the right side. At the interface, we require compatibility conditions of the voltages and currents:

• Inner conductor voltage is continuous
• Inner conductor current is continuous
• Shield voltage is continuous
• Left shield current = right shield current + shorting current
• Shorting voltage equals shield voltage

The suitable shorting wire piece is approximately a quarter of a wavelength at the desired operating frequency (here 868 MHz), thus the reason for calling it a quarter-wave folded balun. However, the velocity in the inner of the coax and in the shorting wire are not the same. The shorting wire is surrounded in large by air, thus the velocity is close to the speed of light, whereas the inner has a velocity factor of about 0.7. A quarter wavelength at the speed of light is about 86 millimeters, while at velocity factor 0.7 it is about 60 millimeters. It is thus not immediately clear what the actual length for the shorting wire should be. To understand how the shorting length affects the performance, we can compute the same observed parameters for a sweep over the length.

Sweep over shorting length for a coax with overall length 860mm

Sweep over shorting length for a coax with overall length 946mm

The common mode currents and the shield voltage all decrease as the length of the shorting wire is decreased. However, at very short shorting lengths, the SWR becomes very bad. Though it seems that the design is quite forgiving when it comes to the exact shorting length. That is, we get a decent SWR over a large range of lengths. There is a trade-off, however, between achieving the best SWR and getting a lower common-mode current.

The optimal SWR of very close to unity is achieved with a shorting length of 61.5mm, which is slightly more than quarter wavelength at the coax propagation speed. At this shorting length, we observe the following parameters for the 860mm coax length:

• Maximum common-mode current amplitude in coax: 0.17 mA
• Maximum voltage amplitude in coax shield: 0.080 V
• Common-mode current to load: 0.22 mA
• SWR: 1.00

For 946mm, we observe the following:

• Maximum common-mode current amplitude in coax: 0.016 mA
• Maximum voltage amplitude in coax shield: 0.0077 V
• Common-mode current to load: 0.040 mA
• SWR: 1.00

Allowing for a slightly worse SWR and choosing the shorting length as 40mm gives considerably better common-mode behavior. At 860mm:

• Maximum common-mode current amplitude in coax: 0.058 mA
• Maximum voltage amplitude in coax shield: 0.028 V
• Common-mode current to load: 0.080 mA
• SWR: 1.32

And at 946mm:

• Maximum common-mode current amplitude in coax: 0.0056 mA
• Maximum voltage amplitude in coax shield: 0.0027 V
• Common-mode current to load: 0.024 mA
• SWR: 1.31

The correct shorting length is thus not at all clear. However, as long as the shorting length is sensible, the balun does seem to work as intended, reducing the coax common-mode current and shield voltage down to about one one-hundredth of just the bare coax. Now we just need to understand how it works.

The following are visualizations of system with a shorting length of 40mm.

Sweep of overall coax length with shorting length at 40mm

Voltage and current profiles in the system with 860mm coax length and 40mm shorting length

Voltage and current profiles in the system with 946mm coax length and 40mm shorting length

The following are computed with shorting length 61.5mm.

Sweep of overall coax length with shorting length at 61.5mm

Voltage and current profiles in the system with 860mm coax length and 61.5mm shorting length

Voltage and current profiles in the system with 946mm coax length and 61.5mm shorting length

At 61.5mm shorting length, the shorting wire appears to only interact very weakly with the inner conductor: almost no current is flowing through the connection between the inner and the shorting. This gives some insight to the operating mechanism of the balun. It shows that the inner conductor plays only a minor role in the operation, and the main interaction is between the shield and the shorting wire.

It seems to work by having the shorting wire and coax shield form a resonator through capacitive and inductive coupling. When current oscillates back and forth in the shorting-shield resonator, the voltage at the load end of the system oscillates differentially. This oscillation is sustained with only tiny kicks of current from the inner conductor. The closer the shorting length is to the quarter wavelength, the less of a kick is needed from the inner conductor and the closer to zero the shorting current is on the right. The oscillation voltage is differential with respect to the ground plane due to symmetry of the system.

Although the inner conductor is very loosely coupled with the resonator, all the energy in the resonator needs to still come from it. This is because the resonator is shorted at the left end, and while there is a lot of current flowing, the voltage there is zero. On the right end, however, the short-inner connection has a non-zero voltage and also nothing is preventing current from flowing. The only reason current is not flowing in the 61.5mm shorting length example is because the system is driven at the resonant frequency, any lower or higher and we'd see power transferring in and out during the cycle.

When feeding the system with a wave packet, we can observe the resonator starting and stopping oscillation. The following shows a Gaussian wave packet with a characteristic time of 1ns at 868MHz passing through a 946mm long coax with a 61.5mm shorting wire.

In the video, you can see the resonator picking up speed a bit slower than the wave packet would otherwise move, and also the resonator keeps going for a bit longer than the wave packet would. Though this these aren't as clearly visible in the video as I had hoped for.

To really see how the resonator ramps ups and down, let's take a closer look at the power driven to the resonator. This power is computed as the shorting wire voltage times shorting wire current conjugate at the right end. Plotted below is the power against time in the simulation shown in the video above. We see that when the oscillation is ramping up, there is a positive power transfer to the resonator, and when the oscillation is dying out, power is transferred from the resonator back to the coax.

Power to (from) the resonator as a function of time in the simulation shown earlier

To be honest, I still don't fully understand how the balun works. I do not have enough insight or intuition of the situation that I could come up with the design myself, and although I can show it working, I can't explain it in simple terms. In conclusion, it's all a bit of black magic and voodoo.

Coaxial cables as multiconductor transmission lines

When using coaxial cables, it is often assumed that all of the current passed through the inner conductor returns through the shield. This seems like a good assumption, since where else would the current go? It turns out, however, that if the conductors have a finite impedance to some third conductor (for example the ground), this will cause currents to pass to that third conductor and not necessarily return symmetrically through the coax. Such a situation may occur for instance when the coaxial cable is used to feed an antenna. The unbalanced currents (common mode currents) in the coax in turn cause the feedline itself to act as an antenna, and to radiate power unpredictably. Antennas are thus typically connected to the feedline through some form of balun, which presents itself as a large series impedance for any unbalanced currents.

There is a lot of discussion regarding unbalanced currents in coaxial cables on amateur radio forums, and Modern Antenna Design by Thomas A. Milligan quickly mentions three-wire transmission lines in their discussion of baluns. However, I haven't found any computational examples of what is going on in these situations. The following is my attempt at that.

Assume a coaxial cable running on top of a perfectly conducting ground plane. Assume also, that the cable isn't too far away from the plane, so that the fields between the coax and the ground plane can be considered to propagate in the TEM mode. This situation can be thought of as a multiconductor transmission line, with two signal conductors and a reference conductor. We can draw the equivalent circuit of a short line segment of length \(\Delta x\) as shown in the figure below.

Inductances \(L_1\) and \(L_2\) are the self-inductances of conductor 1, i.e. the inner conductor with return through the ground plane, and of conductor 2, i.e. the shield with return through the ground plane. The inductance \(L_{12}\) is the mutual inductance between the two signal-bearing conductors. The capacitances bear a similar nomenclature, and their meaning is more clear from the figure. However, note that the capacitance \(C_1\) from the inner conductor to the ground plane is zero, because the electric field would need to pass through the shield first.

In the limit of \(\Delta x \rightarrow 0\) we get the telegrapher's equations for the situation as

\( \begin{align} \frac{\partial V_1}{\partial x} &= -L_1 \frac{\partial I_1}{\partial t} - L_{12} \frac{\partial I_2}{\partial t} \\ \frac{\partial I_1}{\partial x} &= -C_{12} \frac{\partial V_1}{\partial t} + C_{12} \frac{\partial V_2}{\partial t} \\ \\ \frac{\partial V_2}{\partial x} &= -L_2 \frac{\partial I_2}{\partial t} - L_{12} \frac{\partial I_1}{\partial t} \\ \frac{\partial I_2}{\partial x} &= -(C_2 + C_{12}) \frac{\partial V_2}{\partial t} + C_{12} \frac{\partial V_1}{\partial t} \\ \end{align} \)

Equivalently in matrix form as

\( \begin{align} \frac{\partial \boldsymbol{I}}{\partial x} &= -\boldsymbol{C} \frac{\partial \boldsymbol{V}}{\partial t} \\ \frac{\partial \boldsymbol{V}}{\partial x} &= -\boldsymbol{L} \frac{\partial \boldsymbol{I}}{\partial t} \\ \end{align} \)

where \(\boldsymbol{V} = ( V_1, V_2 )\), \(\boldsymbol{I} = (I_1, I_2)\), while \(\boldsymbol{L}\) is the inductance matrix of the system and \(\boldsymbol{C}\) is the capacitance matrix of the system

\begin{align} \boldsymbol{L} &= \begin{pmatrix} L_1 & L_{12} \\ L_{12} & L_2 \end{pmatrix} \\ \boldsymbol{C} &= \begin{pmatrix} C_{12} & -C_{12} \\ -C_{12} & C_2 + C_{12} \end{pmatrix} \end{align}

With a bit of manipulation, we can derive a vector valued wave equation for the voltage as

\( \begin{align} \boldsymbol{L} \boldsymbol{C} \frac{\partial^2 \boldsymbol{V}}{\partial t^2} &= \frac{\partial^2 \boldsymbol{V}}{\partial x^2} \end{align} \)

This describes two waves propagating at velocities \(1/\sqrt{\lambda_1}\) and \(1/\sqrt{\lambda_2}\), where \(\lambda_1\) and \(\lambda_2\) are the eigenvalues of matrix \(\boldsymbol{L} \boldsymbol{C}\).

Determining the capacitance and inductance parameters requires solving the electrostatic and magnetostatic problems in the appropriate geometries. The approach I took, was to use FEM to compute the electric potential for the electrostatic problem, and the magnetic vector potential for the magnetostatic problem. Computing the capacitance is straightforward from relating the electric field energy of the FEM computation to the energy stored in a capacitor. That is, we have the electric field energy from FEM as

\( \begin{align} W = \int_\Omega \boldsymbol{E} \cdot \boldsymbol{D}\ \mathrm{d}x. \end{align} \)

On the other hand, we have the energy stored in a capacitor as

\( \begin{align} W = \frac{1}{2} CU^2. \end{align} \)

This is then simply computed for the two cases:

1. inner conductor at potential \(V\), shield and ground at potential \(0\)
2. inner conductor and shield at potential \(V\), ground at potential \(0\)

Computing the inductance also follows a similar idea, but there is some additional complication to compute the mutual inductance. The energy in the magnetic field is given by

\( \begin{align} W = \int_\Omega \boldsymbol{H} \cdot \boldsymbol{B}\ \mathrm{d}x. \end{align} \)

On the other hand, we have the energy stored in an inductor as

\( \begin{align} W = \frac{1}{2} LI^2. \end{align} \)

Here we have three different cases that need to be considered.

1. Current \(I\) through inner conductor, current \(0\) through shield and current \(-I\) through ground plane.
2. Current \(0\) through inner conductor, current \(I\) through shield and current \(-I\) through ground plane.
3. Current \(I\) through inner conductor, current \(-I\) through shield and current \(0\) through ground plane.

The inductances \(L_1\) and \(L_2\) in the multiconductor transmission line are the inductance of the inner conductor with return through the ground plane and the inductance of the shield with return through the ground plane, respectively. They can thus be directly computed through the magnetic field energy of cases 1 and 2. The inductance \(L_{12}\), however, is the mutual inductance between the inner and outer conductors. To compute that, we need to consider case 3, where we connect the inner conductor to the shield at the far end. This is equivalent to connecting inductances \(L_1\) and \(L_2\) in series together complete with the mutual inductance. Since the currents are opposing, the mutual inductance acts to decrease the total inductance, and we can write the inductance of the series systems as

\( \begin{align} L_\mathrm{series} = L_1 + L_2 - 2 L_{12} \end{align} \)

Thus we can compute \(L_\mathrm{series}\), \(L_1\) and \(L_2\) from the magnetic field energy, and then solve the mutual inductance \(L_{12}\) from the above equation.

Computation domain for the parameters in the example

Placing a coax of roughly RG316 dimensions and materials over a groundplane at a height of 20mm gives the following values (note: the skin effect on the current distribution is not perfectly accounted for, that's a whole new can of worms).

\( \begin{align} L_1 & = 1023.02 \mathrm{nH}/\mathrm{m} \\ L_2 & = 800.58 \mathrm{nH}/\mathrm{m} \\ L_{12} & = 800.58 \mathrm{nH}/\mathrm{m} \\ C_2 & = 14.05 \mathrm{pF}/\mathrm{m} \\ C_{12} & = 105.05 \mathrm{pF}/\mathrm{m} \end{align} \)

Interestingly, the mutual inductance and the inductance of the shield evaluate as very nearly the same value. This means that the voltage induced in the shield due to changing differential currents is basically zero! This can't simply be a lucky coincidence, but I haven't yet worked out the details.

For a visualization of the results, below is a solution of the telegrapher's equations (again using FEM) for a cable length of 2.1 meters. On the left side, the inner conductor is driven with a time harmonic input of amplitude 1V oscillating at 210 MHz, while the shield is grounded. On the right side, the inner conductor is connected to the shield through a 47 ohm resistance. No connection to ground is made on the right side.

Contrary to intuition, the shield really does remain at ground potential! This is in fact exactly what I've observed every time I've tried to measure the potential on the coax shield.

So, what happens if we drive the shield at the left side and ground the inner conductor? On the right side, inner remains terminated to the shield through a 47 ohm resistance without any other connections.

Even in this case, the signal propagates mostly in the differential mode. Also, what is clearly visible here is that the wavelength in the shield is longer than in the inner conductor. The propagation velocity of the common mode signal is thus greater than the differential mode signal. This can also be determined from the eigendecomposition of matrix \(\boldsymbol{L} \boldsymbol{C}\). The eigenmodes in fact almost exactly correspond with the differential and common mode currents (another happy almost coincidence), and the corresponding propagation velocities are

\( \begin{align} v_1 & = 1/\sqrt{\lambda_1} \approx 0.690c \\ v_2 & = 1/\sqrt{\lambda_2} \approx 0.995c \end{align} \)

where \(c\) is the speed of light.

How about the initial motivation? What if there is a finite impedance from both conductors to ground? Some references say that it is unsymmetric impedance to ground that causes inbalance in the currents. I'm not sure what is meant exactly, but the following happens when both the inner and the shield are connected to ground through 23 ohm resistances, with no other connections on the right side. On the left side, the shield is grounded and the inner conductor is driven.

So even if the load is symmetric with respect to ground, it may cause common mode currents to flow in the coax. Also, since in such a situation the coax is not terminated at its characteristic impedance, the length of the cable will also affect the behavior.

A quick look into GPS, and reception using MAX2769

How does GPS work? I mean, I know about the basic priciple of GNSS positioning. However, in addition to that, how does a GPS receiver actually determine the pseudoranges used in the computation. What do the signals look like and how do you receive those signals?

There are a lot of very good resources available on the Internet to learn all about GPS signals, such as ESA Navipedia as well as the actual Navstar GPS interface control documents. This post is mostly a reference for myself about the things I've learned, which might be useful for other people as well.

Radios in general

Before we can talk about GPS specifically, we need to look at radio reception and sampling a bit. So, say we have a signal s(t) occurring around frequency f.

A(t) is the information, which is transmitted. It is much lower in frequency than f, and we can consider it a constant in the time scale of f.

Most radio receivers use a techique called heterodyne as their operating principle. This trick allows the receiver to reduce the frequency of received transmission through the use of a high frequency local oscillator signal. Consider the radio receiver generates two local oscillator signals at frequency f_0, which have a 90 degree phase shift. Let's denote these as

These are both mixed with the received signal in separate mixers, to obtain

 

To simplify the math, we can think of r_I(t) and r_Q(t) as being the real and imaginary components of a complex signal r(t) as

It also useful to write s(t) in terms of the exponential function as

Since s(t) is real, we can write the mixer outputs as

This allows us to think of m_I(t) and m_Q(t) also as the real and imaginary components of a complex signal m(t)

Consider that f_0 is selected to be close to f. So what we're left with after the mixing is f-f_0, which is a low frequency, and f+f_0, which is a high frequency. We're ultimately not interested in f+f_0, and it just occurs as a side product of the mixing. Since f+f_0 is much higher in frequency than f-f_0, it is straightforward to filter out, and were left with just the low frequency complex signal n(t), which is

A GPS receiver will digitize these signals (the signals denoted as the real part and imaginary part separately) and provide the samples for digital signal processing. It remains particularly useful to consider the signal as complex valued in the DSP. In the digital domain, we can apply the heterodyne concept again, to move the signal to 0 Hz. This allows us to correct for effects caused by Doppler shift for instance.

All this boils down to is taking the actual information A(t) from a high frequency carrier down to very close to 0 Hz. Ideally, we'd like to get A(t) exactly, but without exact timing knowledge, the best we can do is to get A(t) exp(j phi), where phi is some phase error. The information contained in A(t) must be designed in a way, which allows recovering the phase, or it needs to such that the phase error doesn't matter.

GPS signals

With that out of the way, then what and how are the GPS satellites actually transmitting? Well, it turns out that every GPS satellite transmits on the same frequency at the same time. Each of them is basically yelling out their name all the time. Although this results in a cacophony of signals at the receiver side, that name is the key to telling them apart.

The name is actually a pseudo random number (PRN) sequence, or PRN code. Each sequence is 1023 bits long, and is transmitted at a 1.023 MHz bit rate, thus repeating every millisecond. The satellites are synchronized so that bit 0 of their respective PRN codes occurs at the same time. If you observe the code of two satellites to not be synchronized to each other, you know that the time it took for the signal to arrive to you from each satellite must be different. This also tells you something (though not everything) about how much farther away one satellite is compared to the other.

The PRN codes above are related to the so-called C/A code transmitted on 1575.42 MHz with BPSK modulation. There are other codes and other frequencies used also in GPS, let alone other GNSS, but I'm not interested in those at the moment. In the following, we'll just happily ignore the other codes, bands and systems.

BPSK modulation simply means, that for each bit in the transmitted bit stream, the A(t) of the transmission has a value of either -1 or 1. For instance A(t)=-1 when the bit transmitted at time t is a 0 and A(t)=1 when the transmitted bit is a 1. The sign will be ambiguous in the reception anyway, so the choice doesn't make a difference.

PRN codes 1 through 32 encoded as a 1023x32 sized image. Unsurprisingly they mostly look like random noise.

The PRN sequences have a strong autocorrelation at 0 lag, while having a very weak autocorrelation at any other lag. In addition, the PRN sequence of one satellite has a very weak cross-correlation to the PRN sequence of any other satellite at any lag. In other words, the PRN sequences have a kind of near-orthogonality with everything else but a coherent version of themselves.

It is also useful to consider the PRN code as a sequence of -1 and 1, instead of 0 and 1. This choice results in the strong correlation having a value of 1023, and the weak correlations being close to 0.

Autocorrelation of PRN code 20. Correlation has a very strong peak at 0 lag, while it is almost 0 elsewhere. Autocorrelations of other PRN codes look similar.

Cross correlation between PRN codes 20 and 21. Correlation is almost 0 everywhere.

It is this near-orthogonality that allows the receiver to separate all the incoming signals from one another. Also, since the correlation exists only at a very precise lag, finding the signal also results in knowing the code phase, which itself is dependent on the time it took for the signal to reach the receiver.

Satellite acquisition

A GPS receiver generates locally the PRN sequence of the satellite it wants to receive, and computes the correlation between the received signal and the generated sequence. If there is a strong correlation, then it knows that it is receiving a satellite with that particular PRN sequence, and that the generated PRN sequence is coherent with the received signal. If there is a weak correlation, then either the satellite isn't visible, or our PRN isn't coherent with the satellite. When searching for a particular satellite, the GPS receiver needs to actually scan through different lags to hopefully find a strong correlation.

Unfortunately there is also a second parameter, which the GPS receiver needs to take into account, when finding a satellite. This is the Doppler shift of the received signal (as well as receiver clock frequency error). Although each satellite is transmitting on very precisely 1575.42 MHz in orbit, the motion of the satellite relative to the receiver causes up to 5 kHz difference on the frequency received on ground. Down-conversion without taking the signal Doppler shift into account leaves the receiver with the signal

That is, the phase of the information in A(t) rotates with frequency f_Doppler. In order to properly compute the correlation between the receiver generated PRN and the received signal, this frequency error needs to be quite small. As example, say the correlation integration time is 10 milliseconds, which is pretty typical. If the frequency is off by just 100 Hz, then the phase error makes one complete revolution during the integration time. This results in the latter half of the integration cancelling the result of the first half, and thus will produce a low correlation even if the correlation would otherwise be high. Shorter integration times suffer less from frequency error, but on the other hand suffer more from noise. The receiver thus needs to scan over both the PRN code lags and the Doppler shifts, when searching for a satellite.

Note that the Doppler shift (and receiver clock error) also causes an error in the bit rate seen in the received signal. This could become an issue with very long integration times. This, however, is not a big problem, as the integration time would need to be around 150 milliseconds for half a bit of accumulated error with the maximum expected Doppler shift. In practice, the integration time is limited to 20 ms due to the navigation message bit rate (more on this later).

MAX2769 based receiver

Data collection setup. Bus Pirate for SPI communication, logic analyzer for data recording. Blue coax leads to a GPS antenna with integrated amplifier. Notice bodge through hole inductor providing bias to active antenna near the antenna SMA connector.

Now that I understood much more about the GPS signals, I needed to acquire some. There are resources online from which you can get pre-recorded GPS signals to play around with, but where's the fun in that? (Though in all fairness, I did download some samples to play around with before attempting reception of my own).

I set to receive the satellites with a MAX2769 evaluation board, which was donated to me some time ago (thank you Tero!). The MAX2769 is a complete software defined radio frontend meant for receiving GNSS signals. It has an LNA for amplifying the incoming signal. It has a PLL for producing the local oscillator frequency. It has an IQ mixer and the required intermediate frequency filters. Finally, it has the IQ sampling ADC for sampling that intermediate frequency.

I've had a decent amount of problems using the MAX2769, most of which have to do with configuring the IC properly. Many settings are not properly documented in the freely available documentation. Looking around the internet, I found one previous hobby project using the MAX2769. They wanted to go with a direct conversion, having the IF frequency at 0 Hz. I opted not to do that since the documentation is particularly lacking in that part. It's not really clear how anything behaves near DC. In the end, I mostly ended up ignoring the configuration the other project used, and instead made my own conclusions based on the documentation, default values and some testing with an RF generator and a spectrum analyzer. Nothing is set in stone though, and its more than likely I will still change the settings in the future.

The configuration I ended up using in my first reception test is the following:

Register 0: 0xA2919A3

Register 1: 0x8550888

Register 2: 0xEAFF1DC

Register 3: 0x9CC0008

Register 4: 0x0C00080

Register 5: 0x8000070

Register 6: 0x8000000

Register 7: 0x10061B2

This sets the IF to 4.092 MHz with 2.5 MHz IF bandwidth. It enables IQ sampling with 2 bits for both, and sets output mode for unsigned values and signal levels for CMOS. Sample rate is set to 8.184 MHz. To write the configuration, I used a Bus Pirate connected to the SPI port of the IC.

Note that with these settings my IF occurs at the Nyquist frequency. However due to IQ sampling, I can separate positive frequencies from negative frequencies, and thus I'm not actually losing any data (though I might be actually recording the complex conjugate of what I want to record).

The data rate produced by the MAX2769 is somewhat high, producing 4 bits of data at 8.184 MHz. I decided that a logic analyzer would be the most straightforward way to record the data. I used my 8 channel USB connected logic analyzer for the purpose (one of the ubiquitous Cypress FX2LP based logic analyzers).

With everything set, on 2020/03/18 at 22:15 local time in Espoo, Finland, I recorded about 500 milliseconds of data. I had the antenna placed next to a window facing South-South-East. Using a satellite tracking tool I knew there were 12 GPS satellites above the horizon at that time. However some were close to the horizon, and some were not in view from the window.

Satellites visible in the sky from Espoo, Finland on 2020/03/18 22:15 local time.

If you're interested in my data, you can download the data as CSV or binary. Both of these formats are horribly inefficient in storing the data, but are somewhat convenient to use.

The CSV format is simply I,Q per line. Both I and Q get values 0,1,2 or 3. To decode a pair as a complex value, use the formula n_k = (I_k-1.5)  + j (Q_k-1.5). Subtracting 1.5 is to reduce the DC content of the signal.

The binary format is simply a dump of numpy.complex type values. To read it, use
n = numpy.fromfile("data_20200318.bin", dtype=numpy.complex)
The contents are the same as for the CSV file and the described decoding.

If you do use the data in either of the formats, remember that it contains a digitization of the intermediate frequency. You should downconvert the signal from 4.092 MHz down to 0. The Fourier transform of the correctly downconverted signal should show power within an about 2.5 MHz band around 0 Hz and look similar to the figure below.

Spectrum of recorded data after digital domain downconverison.

A thing to note is that IF filter center frequency is not exactly at 4.092 MHz even if that's what it is configured to. There should still be enough width to let the GPS signals pass. Also, I haven't made completely sure that the sign of the frequency axis is correct. This would affect at velocity computations, since Doppler shifts would occur with an incorrect sign. Since this data set is simply too short for such things anyway, I won't worry about it too much.

First look

The satellite transmitting PRN10 happened to be quite high in the sky and in the correct direction to not be blocked by the house. In order to find the signal, we can go through the same process a GPS receiver would do to find the satellite. A typical GPS receiver would look at different Doppler shift and PRN phase parameters successively in time. With the recorder data, however, we can go through the same piece of data over and over again until we find the correct parameters. The following figures are produced by computing the correlation of the PRN code against the first 10 milliseconds of recorded data, while stepping through the Doppler from -5 and 5 kHz with 512 steps and through the PRN phases with 8184 steps. The actual correlation result is a complex number because of the phase error we are not yet correcting. The value shown in the figures is the absolute value of the correlation.

By their very nature, the correlation peaks are very narrow. This makes it difficult to properly visualize them. Nevertheless, a peak can be seen in the figure below (near phase 200, Doppler 0 Hz). I managed to record signals from a GPS satellite!

Results of PRN10 acquisition. Do you see the correlation peak?

Not only that, in fact I managed to receive signals from multiple GPS satellites! Below is the same plot computed for PRN code 27 (near phase 550, Doppler -2 kHz). This satellite was also in a favourable position in the sky.

Results of PRN10 acquisition. This peak is a bit more difficult to find.

Below is a zoomed-in view of the PRN10 plot - just to enjoy the correlation peak to the fullest 😃

PRN10 peak. Now you definitely see it. It's beautiful!

Navigation message

The PRN code phase alone isn't enough to determine the pseudoranges to the satellites. The PRN code simply repeats too often. The distance difference between satellites is easily more than one light-millisecond (which is about 300 km). Additional information is still needed to determine the delay between the satellites as a number of full milliseconds.

Of course it should come as no surprise, that such information is available in the GPS signal. Turns out that, in addition to the carrier BPSK modulation with the PRN sequence, the PRN sequence itself is BPSK modulated with a bit stream called the navigation message.

GPS modulation with PRN code and navigation message (figure by P.F. Lammertsma CC BY-SA 3.0)
https://en.wikipedia.org/wiki/GPS_signals#/media/File:GPS_signal_modulation_scheme.svg

Essentially the navigation message will flip the polarity of the PRN code depending on the bit of the navigation message. In other words, the GPS signal A(t) will still be a sequence of -1 and 1 according to the PRN sequence, but which value denotes a 1 and which denotes a 0 gets flipped depending on the navigation message bit.

The navigation message has a bit period of exactly 20 PRN codes, and the bit transitions only occur at PRN code phase 0. Thus the integration periods used for the correlation should always start at code phase 0. I'm not sure if GPS receivers actually do this, but I've used this approach in my implementation.

The navigation message bit stream contains all sorts of useful information that a GPS receiver might needs to know. Such as the aforementioned sychronization data, to allow the proper evaluation of pseudoranges to the satellites. It also contains time information as well as orbit parameters of the satellites. It also includes a short message service type field (GPS Special Message), which I have no idea what it contains. Specification says it can be a short sequence of ASCII characters.

Navigation message in the recorded data

Assume that the Doppler shift is perfectly corrected, the received signal then looks like

Assuming further that A_NAV(t) is constant within each integration period, then correlating this signal against a coherent locally generated PRN code then gives for each integration period i the correlation c_i as

where the amplitude has also been normalized.

These are of course idealized formulas. In reality there is noise on top of the signal also. This noise comes from all typical noise sources, but also from the other satellites since the PRN codes are not perfectly orthogonal.

The correlation thus evaluates as a complex value on the unit circle, and jumps between a point and its antipode as the navigation message bit changes. A GPS receiver will typically constantly adjust for the phase error to align the correlation on the real axis, and thus the value will jump between 1 and -1.

The plots below show what the correlation looks like in time (after meticulous hand-tuning of the Doppler frequency and residual phase error) as well as the effect of using different integration times on the correlation result.

It thus seems obvious that the longest possible integration time is the way to go. Since the navigation message bit period is 20ms, that limits the practical maximum integration time. For best signal to noise ratio, that is the value to use then. Unfortunately it is not possible to know a priori where exactly are the bit transition boundaries of the navigation message. If the 20ms integration period ends up with 10ms of positive polarity and 10ms of negative polarity, then the output is zero.

In the plots above, I've hand-tuned the integration times to occur exactly during the bit periods. The following plots show how the output degrades, if there is an error in synchronizing the bit periods and integration times.

The integration is thus pretty sensitive to the amount of offset. What I understand GPS receiver usually do, is they use a shorter integration time first, and use that to find the bit periods. Then when they know where the transitions occur, they switch to 20ms integration for maximum signal to noise ratio.

For now, the amount of data I have is not enough to do anything useful with the navigation message. There aren't enough bits to determine the time difference between the satellites. For that I'll need to record some more data.

Tracking

The last thing I wanted to cover here is the concept of signal tracking. It is a huge amount of computational effort to go through the data in the way I did. Acquiring the signal does require sweeping over Doppler frequencies and PRN code phases. However, after the signal has been acquired, you don't need to keep looking for it for each successive integration time, since the parameters are quite stable. There are two things to track, the PRN code phase and the signal phase.

Tracking the code phase is typically done by computing the correlations with three phases. These are called the early, prompt and late correlations. The prompt correlation is the best guess of what the PRN code phase is at any given time. The early and late correlations are computed with the code phase being delayed and advanced by a tiny bit (not too much, or we'll lose the correlation completely). If the receiver. for instance, notices that the early correlation is stronger than the prompt correlation, it would adjust the code phase accordingly. This will keep the code phase in lock, as long as it doesn't suddenly change by a large amount.

The code phase drifts because the Doppler shift (and clock errors) cause the PRN bit frequency to look slightly different than 1.023 MHz. As mentioned earlier, at the maximal 5 kHz Doppler shift, this error is half a PRN bit period over 150 milliseconds. The PRN code phase will thus drift by one bit every 300 milliseconds in the most extreme case.

The signal phase error phi needs to be controlled so that it at least stays approximately constant. However, it is useful to control the phase so that the PRN correlation output aligns on the real axis. That is to have phi = 0 or phi = pi. In fact, without the information content in the navigation message, we can't tell them apart.

Given the correlation as

We want to find a correction to signal phase phi, which brings c to the real axis without getting confused about A_NAV changing sign.

The method I used in my tracking experiments was based on the following expression of the phase correction

,

where it is important to note that we deliberately get rid of quadrant information, hence arctan and not the complex angle.

Thinking about the phase correction as

,

the arctan correction has two important properties. If c is aligned on the real axis (i.e. when the imaginary component of c is zero) then the correction is also zero.

Second, if c isn't on the real axis, the correction turns c back onto the real axis toward either the positive real direction or negative real direction, depending on which one is closer.

This correction is optimal, in the sense that the corrected phi will turn c exactly to the real axis. On the other hand, it is quite expensive to compute. Due to noise in the signal it is not a good idea to correct for the entire observed error. For stability, the correction should in fact be made as slow as possible. It thus isn't very important if the computed correction is exact or approximate. Approximate corrections can be much faster to compute, and be just as good.

Conclusion

This is about as far as I've managed to get so far in my quest to understand GPS. In addition to Python code for playing around with the recorded data, I've worked on implementing the Doppler shift correction, PRN code generation and integration in Verilog.

It is still some way away from working in practice. However, I can feed the recorded data to the model in Icarus Verilog, and get out the correlation values correctly. However, I haven't yet decided how to output the correlation data from an actual piece of hardware. I think a simple serial line is too slow even for the correlation output.

In addition to IO questions, the Verilog implementation doesn't do any tracking (PRN phase, signal phase) of the signals by itself, and instead needs to be controlled from the outside. In my concept, the FPGA would do only the fast stuff (handling data at 8.184+ Msamples/s), and output a 1 kHz correlation output based on the configuration given to it. This correlation output would be fed to software, which processes it and gives control parameters back to the FPGA to keep the signal locked.

So for this to work in practice, I need to first solve the IO questions, get that data fed to a processor and implement signal tracking software. A soft core might be a viable option, but the FPGA kits I have at hand are quite small.