Full Stack GPS Receiver

I performed some quick experiments of search time vs channels found. Here's the results:

Search Bands | Channels found | Time
-------------+----------------+--------
    21       |      6         |  48.202
    31       |      8         |  65.194
    41       |      9         |  92.639
    51       |      9         | 103.231
    71       |      7         |  99.600

It looks like 31 or 41 is a good balance between runtime and performance.

Over the weekend I wrote some simple acquisition code. The method is rather simple, but effective.

The +/-5kHz band is 31 evenly spaced local oscillators - each 330 Hz apart. At any given time, one Space Vehicle is being searched for, over all 31 frequencies. Each of 1023 phase alignments is tested for 2ms, and each test is offset by 1us from the last one. When a 'sniff' of a signal is received, a channel is assigned to it and it tries to achieve lock. If a better 'sniff' is seen, the channel is reset with the new values, and once again tries to achieve lock.

The scheduling of which SV is being acquired changes depending on if lock is achieved, if any 'sniffs' were found, and an initial setting of priorities from a text file. At the moment, this uses about 15% of a CPU Core (out of a current total of around 30%).

One problem was that when a channel starts to lock, it can fall one of two ways - either onto the correct frequency (yay!) or fall to 500Hz to one side (wah!). If it falls to one side the channel's in-phase value alternates between + and - every millisecond, rather than staying steady for 20ms. My fix for this was when you have half a second of bad NAV data, but still have a signal lock then adjust the NCO by 500 Hz, towards (and past) the initial frequency. This jumps the channel pretty much directly into lock.

All in all the results have been pretty good: - It has been able to tease out 8 channels from a test data set:

Update at  195.000    Acquiring 29
Channel status:
SV, WeekNum, FrameOfWeek,     msOfFrame,  earlyPwr, promptPwr,   latePwr, frame, bitErrors
23,    1942,       38729,  2001.1655120,     89461,   1110140,    125161,  12345       0
03,    1942,       92342,  2007.5771226,    271183,   4084861,    152405,  12345       0
14,    1942,       92342,  2006.9831913,    169201,   3817052,    223213,  12345       0
26,    1942,       92342,  2005.1505016,    131473,   2327256,    141497,  12345       0
31,    1942,       92342,  2010.4623747,    303547,   6850993,    251145,  12345       0
04,    1942,       92342,  2006.0406834,    134309,   3619075,    157964,  12345       0
22,    1942,       92342,  2010.4915444,    110075,   6369297,    146286,  12345       0
16,    1942,       92342,  1999.1110117,    123204,    988464,     62537,  12345       0



Space Vehicle Positions:   BAD TIME DETECTED - SV position dropped
sv,            x,            y,            z,            t
23, -10332457.50,  13443819.33, -20441422.71, 554054.00768279
 3, -10332457.50,  13443819.33, -20441422.71, 554054.00768279
14, -16698154.35, -12949614.61, -15924944.22, 554054.00704834
26, -25841408.28,  -4860111.89,  -4072673.60, 554054.00566661
31, -15305881.89,  -4074165.08, -21307244.57, 554054.01025238
 4, -24522414.32,  -7782155.47,  -7385698.32, 554054.00603940
22, -18729989.63,  10711501.63, -15288354.25, 554054.01043638


Solution ECEF:  -4582619.86,    609345.20,  -4379474.77,     0.00394
Solution LLA:     -43.64293,    172.42588,        10.66

I thought up a away to improve the speed - now uses about 3% of a CPU core per channel. The speedup also allows the early/late tracking that is used to follow the Space Vehicle's Gold code to be change from +/- 1 .0 chip (1/1023th of a ms) to +/-0.75 chip.

This allows more accurate tracking, giving more accurate pseudoranges and more precise fixes. Adding an averaging filter removes quite a bit of the remaining noise too.

Now I am in the ballpark of a commercial GPS module!

Here is the console output as it runs, showing it as it tracks 7 channels:

I've got the fast version now generating positioning solutions. Here's a screenshot after processing 400MB of raw samples.

196 seconds of data was processed for seven channels in 88.9 seconds, - 45% CPU usage. Looking Great!

Now to work on to a speedy satellite acquistion code.

GREAT NEWS - A REAL-TIME VERSION WILL BE POSSIBLE

The initial version of the software was slow - it took over 6 hours to process 200 seconds of 16Mb/s data. This was expected, as the system was processing one bit at a time, one channel at a time, and was made that way to allow me to get the techniques correct before focusing on performance.

I've completed a new version of the channel tracking code (in the fast_fsgps directory) that now processes data in 'gulps' of 32 bits at a time.

This greatly speeds things up, and I can now process the same dataset in 110 seconds.

Debugging this code was one of the harder things I've done. I hit odd compiler behavior (you can't shift an 'int' 32 bits) and just my own stupidity a few times, but it now works cleanly and very quick.

Now to incorporate the NAV message code, to get position fixes in real time.

It will also make control loop tuning a lot quicker to test.

PS. I Just notice that I left some debugging code in place, where math was replaced with a loop to help find a fault. It is now 300x faster!

I've been fine-tuning the the process by which the timing and phase of the Gold codes is followed. This has approximately halved the spread of the fixes.

The points are 20ms apart, so an average of a second's worth of points would be less spread - absolute error is about 15m.

There is still a bit of room for refinement of the process. I played around with a PID loop for a few nights, but it doesn't seem to be that suited to the process. It tends to bounce around due to the noise in the system. However the (a-b)/(a+b) discriminator (where 'a' is the late code power and b is the early code power) is great for assessing the relative performance of new methods.

Sometimes it is easier to work backwards. If you have positions for multiple Space Vehicles, and a position in Earth Centered, Earth Fixed (ECEF) coordinates, how can you verify it?

An example position could be (3851787, -78307, 5066308)

At http://www.oc.nps.edu/oc2902w/coord/llhxyz.htm there is a tool to convert from ECEF (in km, not meters) to Lat/Long/Alt.

Converted to Lat/Long/Alt the output you get is:

Latitude,Longitude, Height (ellipsoidal) from ECEF

Latitude : 52.93473 deg N

Longitude : 1.16467 deg E

Height : 174.1 m

This looks intuitively OK (somewhere in England, at a reasonable altitude). So the next step is to check this on a map. For Latitude, N is positive and S is negative. For longitude W is positive and E is negative.

A quick search on Google Maps for "52.93473, -1.16467" will show where the solution is. For these numbers you should be somewhere just south of Nottingham in England.

But what do you do if the final location is wrong?

For debugging, you can also run this process backwards, but you will also need to find your altitude. For me the best tool for this is google Earth.

Once you have latitude, longitude and altitude you can convert your known position to ECEF position using the above website.

You can then run the calculations backwards to see what errors exist in your Space Vehicle positioning.

For example, the five Space Vehicles used had the following X,Y,Z and transmit time:

SV A ( 10796463, -11059512, 21468778) @ 466786.70444515
SV B ( 25624907, 5493917, 6376860) @ 466786.69914821
SV C ( 17562614, 17600801, 9306610) @ 466786.69828776
SV D ( 12436654, 8851501, 21813026) @ 466786.70476175
SV E ( 18466455, 6615577, 17512965) @ 466786.70642941

The distance in each axis can be calculate by subtracting the known position (3851787, -78307, 5066308):

SV A ( 6944676, -10981205, 16402470)
SV B ( 21773120, 5572224, 1310552)
SV C ( 13710827, 17679108, 4240302)
SV D ( 8584867, 8929808, 16746718)
SV E ( 14614668, 6693884, 12446657)

The straight line distance can now be calculated with d = sqrt(x*x+y+y+z*z)

SV A 20,925,018.76
SV B 22,513,018.04
SV C 22,770,985.88
SV D 20,830,121.85
SV E 20,330,171.52

The time of flight for the signal can be calculated by dividing the distance by the speed of light (299,792,458 m/s):

SV A 0.069,798,349
SV B 0.075,095,344
SV C 0.075,955,833
SV D 0.069,481,807
SV E 0.067,814,152

If you add this time of flight to the transmit time they should all be equal:

SV A 466786.774,243,500
SV B 466786.774,243,555
SV C 466786.774,243,593
SV D 466786.774,243,558
SV E 466786.774,243,553

So in this case, the range of times of about +/- 50ns or so - equating to +/- 15m.

It is very hard to find raw RF data for analysis. The only source I was able to find on the internet was http://www.jks.com/gps/gps.html

It contains a 59 MB file, that has the following properties:

• Sample rate - 5.456 MHz
• Sample size - 1 bit per sample
• Sample order - within a byte bit 0 is the earliest sample, and bit 7 is the latest sample
• Intermediate frequency (IF) - 4.092 MHz.

It is important to note that the intermediate frequency is above half the sample rate - this will cause signals to be aliased down to 4.092 - (5.456/2) = 1.364 MHz, and the spectrum to be inverted.

If you download this file, and combine it with the project's C source, you will be able to experiment with this project and obtain a valid position fix. To do this, use the following command line

fsgps -s 5456000 -i 4092000  gps.samples.1bit.I.fs5456.if4092.bin