In this post I am going to look at some of the raw satellite data directly without using RTKLIB to calculate a solution. Because of the complexity of the full solution, it is easy to start treating RTKLIB as a black box, adjusting input parameters to see how they affect the output without really understanding what is going on inside. It is often possible, however, to simplify the problem down to something we can solve, or at least analyze, with some very simple calculations. Doing this can provide some real insight into what is going on inside the more complex full solution.
In this case, I have a data set from my two M8N receivers, both stationary and within a couple meters of each other with skies partially obscured by trees. There are enough cycle slips in the data to cause RTKLIB some trouble with the solution. The eventual goal here is to improve RTKLIB’s response to cycle slips, but to start with I just want to understand more about what the signals look like during the cycle slips.
Since the two receivers are very close together, double differencing should be very effective at cancelling errors from the ionosphere and troposphere as well as the receiver and satellite clock biases. Since the receivers are stationary, the position component of the double differences will be constant and can be ignored. Double differencing is a very powerful error cancellation technique and is the basis for all differential GPS solutions. I’ll skip a detailed tutorial here, but if you aren’t familiar with how and why it works, a Google search will quickly bring up many explanations better than I could give.
To get the double differences, we first calculate the single difference value for each satellite by subtracting the raw carrier phase measurement from one receiver from the raw carrier phase measurement from the other receiver. Since the satellite clock bias and atmospheric error terms are in both measurements, they will be cancelled or at least minimized by this subtraction.
Next, we subtract the single difference value from one satellite from the single difference value from a second reference satellite. The result is the double difference. This second subtraction will cancel the receiver clock biases, further reducing the errors. Since the receivers are stationary (relative to each other) the position component of the double differences will be constant and so the dominant remaining component will be the phase bias. The phase bias should also be constant unless there is a cycle slip. So, what we expect to see in the double differences is a constant value if there is no cycle slip, and a jump equal to the size of the slip if there is a slip. For this exercise we don’t care about the actual values of the phase biases, only any changes in the values.
The cycle slips are reported by the receiver in the raw data. One thing to be aware of it that the receiver is reporting possible cycle slips, not confirmed cycle slips, so not all the reported slips are real. How many of the reported slips are real is something I hope to be able to tell by looking directly at the double differences.
I wrote a simple matlab script to parse the raw observation files (*.obs). This is quite a bit easier if the observation files are in Rinex 3 format rather than the older Rinex 2 format. This is a configuration option in RTKCONV and CONVBIN to select this. Unfortunately it does not seem to work in CONVBIN so I had to use RTKCONV to create the observation files. Here is a snippet of the Rinex 3 observation file for a typical cycle slip on a single satellite.
The top line is the time stamp, each line below is for one satellite. In this example, the first column is the satellites number, the third column is the carrier phase. A cycle slip shows up as a fourth digit after the decimal point in the carrier phase measurement. In this example there is a cycle slip for satellite R8.
Here is another example, but in this case, the cycle slip is reported for all satellites.
This is very difficult for RTKLIB to deal with since the cycle slips, whether they are real or not, immediately cause it to reset the phase bias estimates. In this case, RTKLIB resets all the bias bias estimates and has to start over completely. We’ll use this example to look more closely at what is going on inside the solution.
Let’s look at the G13 satellite first. We’ll use G12 as the reference satellite since that is what RTKLIB used. Here are the carrier phase measurements for both satellites from the two receivers, pulled from the observation file.
Note the y-axis has a scale of 10^6. Even though the receivers are not moving relative to the earth, the distance between the satellites and the receivers is changing quite quickly. This is because the satellites are moving and because the earth is rotating.
Below is the single difference measurement for G13, calculated by subtracting the rover measurement (yellow line above) from the base measurement (purple line above). The mean is also subtracted to center the plot around zero. The red circles are where cycle slips occurred.
The cancellations that occur from the subtraction have reduced the y-axis scale down from 10^6 to just tens of cycles but the remaining variation is still too large to see a cycle slip. Again, this calculation is done directly on the raw data from the observation file, RTKLIB is not used at all to process the data.
Next we get the double difference by subtracting the single difference measurement for the reference satellite (G12 in this case) from the single difference for G13 plotted above. Again the mean of the result is removed to center the result around zero. This results in the plot on the left below. Again, the red circles are where cycle slips occurred.
We have now cancelled some additional errors and the variation on the y-axis is very small, in fact small enough to see a cycle slip. It is quite encouraging how this very simple mathematical operation on the raw data can give us such meaningful results.
The above plot on the right comes from data in the RTKLIB state output file and shows the kalman filter estimate of the phase-bias for satellite G13. Again the mean is removed so only the variation is relevant. The large excursions occur when the kalman filter states are reset because of the reported cycle slips. The filter states are re-initialized using the pseudorange measurements which are much less accurate than the carrier-phase measurements but give an absolute value rather than just a relative number. As you can see, in this case it takes several minutes for RTKLIB to fully recover the phase bias estimate but when it does, it converges back to the same value as before the reported slip.
Looking back at the plot on the left, we can see for the first cycle slip there is no jump in the measurement. If there was a cycle slip we would expect a jump in this measurement equal to the size of the cycle slip. For the second cycle slip there is a brief half cycle discontinuity which I don’t understand, but since it quickly returns to zero we’ll ignore it, at least for the time being.
What our very simple calculation is telling us is in this case is that neither potential cycle slip reported by the receiver was actually real and RTKLIB would have done much better to not have reset the phase bias estimates. This is also confirmed by RTKLIB since it eventually converges back to the same phase-bias estimate as it had previous to the potential cycle slip.
Let’s look at another satellite from the same reported cycle slip. On the left below is the double difference for G18 calculated directly from the raw data just as we did for G13 above. On the right again is the kalman filter estimate of the phase bias from RTKLIB and in both plots, the reported cycle slips occur at the red circles.
This time the first cycle slip is real and you can see the two cycle jump in the double difference value on the left. On the right you can see that RTKLIB eventually settles on the same two cycle slip, but takes several minutes to figure it out. The second cycle slip again causes a brief half cycle deviation in the double difference but then returns to zero, and again we will ignore it.
Based on the plots on the left which are generated with only three basic subtractions using the raw carrier phase data, it appears that RTKLIB could be doing a significantly better job of handling cycle slips. Not only in determining if potential reported cycle slips are real or not, but also in correcting the cycle slip more quickly when it is real. It is true that the baseline in this example is very short and with longer baselines the errors in the double differences will be greater, making things more difficult, but I still have to believe there is opportunity for improvement here.
In the example above, of the thirteen satellites processed by RTKLIB, all of which were reported in the raw data as having cycle slips, only four of them actually did. Even if RTKLIB reset the phase-bias estimates for the four satellites with slips, it would still handle the slip much better using the other nine satellites than it does now where it resets all thirteen phase-bias estimates. If it could correct the slip as well as correctly identify it, the robustness to slips would be even greater.
I am actually quite surprised how clear cut the cycle slip examples I looked at are. I expect this is not always the case or cycle slip correction would probably have been added to RTKLIB a long time ago. Still, I hope to use this approach to try and add some amount of cycle slip detection and correction to RTKLIB, even if it can only be applied in some subset of cases.
Anyone else have experience trying something like this?