GPX file processing

badidea · Post by **badidea** » Jun 28, 2020 22:53

Has anyone some experience in working with GPX files in Freebasic? Maybe some code to share?
The only related topic I see is this one from 2007 for reading the files: libxml usage.

Things I want to do:
- Read / write the files
- Convert latitude & longitude to flat 2D position space
- Display the tracks with color depending on speed
- Detect and filter out suspicious points (I sometimes reach 200 km/h on my bicycle)
- Filter the height variation (almost none here, but sometimes my phone thinks otherwise)
- Combine tracks
- Remove sections, especially start and end (where I walk around in my backyard)
- Etc.

I know that there are tools for (some of) this, but it seems very possible to do this myself, which is more fun.

badidea · Post by **badidea** » Jul 06, 2020 23:20

First step done: Reading a gpx-file into memory with libxml.

Code: Select all

#Include Once "libxml/xmlreader.bi"
#Define NULL 0

Type trkpt_type
	Dim As Double lat, lon 'latitude and longitude
	Dim As Double dateTime 'UTC
	Dim As Double ele 'elevation
	'dim as double ext_dir, ext_g_spd, ext_h_acc, ext_v_acc
End Type

Type trkpt_list
	Dim As trkpt_type pt(Any)
	Declare Function Add() As Integer
	Declare Function size() As Integer
	Declare Function check() As Integer
End Type

Function trkpt_list.add() As Integer
	Dim As Integer ub = UBound(pt)
	ReDim Preserve pt(ub + 1)
	Return ub + 1
End Function

Function trkpt_list.size() As Integer
	Return UBound(pt) + 1
End Function

'simple data validity check 
Function trkpt_list.check() As Integer
	Dim As Integer errorCount = 0
	For i As Integer = 0 To UBound(pt)
		If pt(i).lat = 0 Then errorCount += 1
		If pt(i).lon = 0 Then errorCount += 1
		If pt(i).ele = 0 Then errorCount += 1
		If pt(i).dateTime = 0 Then errorCount += 1
	Next
	Return errorCount
End Function

'===============================================================================

#Include "vbcompat.bi"

Enum E_DATE_FORMAT
	DATE_FORMAT_ISO
	DATE_FORMAT_NL
	DATE_FORMAT_DE
	DATE_FORMAT_US
End Enum

Function dateString(dateTimeVal As Double, dateFormat As E_DATE_FORMAT) As String
	Select Case dateFormat
	Case DATE_FORMAT_ISO
		Return Format(dateTimeVal, "yyyy-mm-dd")
	Case DATE_FORMAT_NL
		Return Format(dateTimeVal, "dd-mm-yyyy")
	Case DATE_FORMAT_DE
		Return Format(dateTimeVal, "dd.mm.yyyy")
	Case DATE_FORMAT_US
		Return Format(dateTimeVal, "mm/dd/yyyy")
	End Select
	Return "<invalid>"
End Function

'day, month, year format must be: dd, mm, yyyy (e.g. not 7/5/20 or 2020-7-5)
Function dateValue2(dateStr As String, dateFormat As E_DATE_FORMAT) As Long
	If Len(dateStr) <> 10 Then Return -1
	Dim As String y, m, d 'years, months, days
	Select Case dateFormat
	Case DATE_FORMAT_ISO
		If Mid(dateStr, 5, 1) <> "-" Then Return -1
		If Mid(dateStr, 8, 1) <> "-" Then Return -1
		y = Mid(dateStr, 1, 4)
		m = Mid(dateStr, 6, 2)
		d = Mid(dateStr, 9, 2)
	Case DATE_FORMAT_NL
		If Mid(dateStr, 3, 1) <> "-" Then Return -1
		If Mid(dateStr, 6, 1) <> "-" Then Return -1
		y = Mid(dateStr, 7, 4)
		m = Mid(dateStr, 4, 2)
		d = Mid(dateStr, 1, 2)
	Case DATE_FORMAT_DE
		If Mid(dateStr, 3, 1) <> "." Then Return -1
		If Mid(dateStr, 6, 1) <> "." Then Return -1
		y = Mid(dateStr, 7, 4)
		m = Mid(dateStr, 4, 2)
		d = Mid(dateStr, 1, 2)
	Case DATE_FORMAT_US
		If Mid(dateStr, 3, 1) <> "/" Then Return -1
		If Mid(dateStr, 6, 1) <> "/" Then Return -1
		y = Mid(dateStr, 7, 4)
		m = Mid(dateStr, 1, 2)
		d = Mid(dateStr, 4, 2)
	Case Else
		Return -1
	End Select
	Return DateSerial(Val(y), Val(m), Val(d))
End Function

'Convert 2020-05-03T16:28:56Z to 43954.68675925926
Function gpxDateTimeValue(gpxDateTimeStr As String) As Double
	Return dateValue2(Mid(gpxDateTimeStr, 1, 10), DATE_FORMAT_ISO) _
		+ TimeValue(Mid(gpxDateTimeStr, 12, 8))
End Function

Function DiffDateTimeSec(dt1 As Double, dt2 As Double) As Double
	Return (dt2 - dt1) * 24 * 3600
End Function 

'===============================================================================

'returns < 0 or error
Function readGpxFile(fileName As String, track As trkpt_list) As Integer
	Dim As xmlTextReaderPtr pReader = xmlReaderForFile(filename, NULL, 0)
	If (pReader = NULL) Then
		Print "readGpxFile: Unable to open: "; fileName
		Return -1
	End If
	
	Dim As Const ZString Ptr pConstName, pConstValue
	Dim As Const ZString Ptr pAttrName, pAttrValue
	Dim As Integer nodetype, iTrkpt = -1

	Dim As Integer ret = xmlTextReaderRead(pReader)
	While (ret = 1)
		pConstName = xmlTextReaderConstName(pReader)
		pConstValue = xmlTextReaderConstValue(pReader)
		nodetype = xmlTextReaderNodeType(pReader)
		If (nodetype = 1) Then
			If (*pConstName = "trkpt") Then
				iTrkpt = track.add()
				If xmlTextReaderHasAttributes(pReader) = 1 Then
					If xmlTextReaderMoveToAttributeNo(pReader, 0) Then
						pAttrName = xmlTextReaderConstName(pReader)
						pAttrvalue = xmlTextReaderConstValue(pReader)
						If *pAttrName = "lat" Then
							track.pt(iTrkpt).lat = CDbl(*pAttrvalue)
							'print "lat: " & *pAttrvalue
						End If
					End If
					If xmlTextReaderMoveToAttributeNo(pReader, 1) Then
						pAttrName = xmlTextReaderConstName(pReader)
						pAttrvalue = xmlTextReaderConstValue(pReader)
						If *pAttrName = "lon" Then
							track.pt(iTrkpt).lon = CDbl(*pAttrvalue)
							'print "lon: " & *pAttrvalue
						End If
					End If
				End If
			End If
			If (*pConstName = "time") Then
				'actual data is in next element
				ret = xmlTextReaderRead(pReader) 'read next
				pConstName = xmlTextReaderConstName(pReader)
				pConstValue = xmlTextReaderConstValue(pReader)
				If (ret = 1) Then
					track.pt(iTrkpt).dateTime = gpxDateTimeValue(*Cast(ZString Ptr, pConstValue))
					'print "datetime: " & *pConstValue
				End If
			End If
			If (*pConstName = "ele") Then
				'actual data is in next element
				ret = xmlTextReaderRead(pReader) 'read next
				pConstName = xmlTextReaderConstName(pReader)
				pConstValue = xmlTextReaderConstValue(pReader)
				If (ret = 1) Then
					track.pt(iTrkpt).ele = CDbl(*pConstValue)
					'print "ele: " & *pConstValue
					'print
					'getkey()
				End If
			End If
		End If
		'getkey()
		ret = xmlTextReaderRead(pReader)
	Wend

	xmlFreeTextReader(pReader)

	If (ret <> 0) Then Print "readGpxFile: Failed to parse: "; filename

	xmlCleanupParser()
	xmlMemoryDump()
End Function

'-------------------------------------------------------------------------------

Dim As trkpt_list track
Dim As String fileName = "Skiing-downhill.gpx"
Dim As Integer result = readGpxFile(fileName, track)

Print "fileName: " & fileName
Print "readGpxFile: " & IIf(result >= 0, "Ok", "error")
Print "Num track points: " & track.size()
Print "Track error count: " & track.check()
Print

'show first and last 5 data points:
Print " Index", "yyyy-mm-dd - hh:mm:ss", "Latitude", "Longitude", "Elevation"
If track.size() > 10 Then
	For i As Integer = 0 To 4
		With track.pt(i)
			Print i, Format(.dateTime, "yyyy-mm-dd - hh:mm:ss"), Format(.lat, "#.000000"), Format(.lon, "#.000000"), Format(.ele, "#.000")
		End With
	Next
	Print " ..."
	For i As Integer = track.size() - 5 To track.size() - 1
		With track.pt(i)
			Print i, Format(.dateTime, "yyyy-mm-dd - hh:mm:ss"), Format(.lat, "#.000000"), Format(.lon, "#.000000"), Format(.ele, "#.000")
		End With
	Next
End If

GetKey()

Example gpx-file (recorded by me): Skiing-downhill.gpx

badidea · Post by **badidea** » Jul 09, 2020 20:26

Some progress, with the same data-file (as above), with simple spherical approximation of the earth, some statistics generated.
Code (concatenated):

Code: Select all

#Include Once "libxml/xmlreader.bi"
#Define NULL 0

'===============================================================================

Type dbl2d
	As Double x, y
	Declare Constructor
	Declare Constructor(x As Double, y As Double)
End Type

Constructor dbl2d
End Constructor

Constructor dbl2d(x As Double, y As Double)
	This.x = x : This.y = y
End Constructor

'===============================================================================

Type trkpt_type
	'DATA FROM FILE:
	Dim As Double lat, lon 'latitude and longitude
	Dim As Double dateTime 'UTC
	Dim As Double ele 'elevation
	'dim as double ext_dir, ext_g_spd, ext_h_acc, ext_v_acc
	'DERIVED DATA:
	Dim As Double speed 'scalar
	Dim As dbl2d p 'position x,y relative to pt(0)
End Type

Type trkpt_list
	Dim As trkpt_type pt(Any)
	Declare Function Add() As Integer
	Declare Function size() As Integer
	Declare Function check() As Integer
End Type

Function trkpt_list.add() As Integer
	Dim As Integer ub = UBound(pt)
	ReDim Preserve pt(ub + 1)
	Return ub + 1
End Function

Function trkpt_list.size() As Integer
	Return UBound(pt) + 1
End Function

'simple data validity check 
Function trkpt_list.check() As Integer
	Dim As Integer errorCount = 0
	For i As Integer = 0 To UBound(pt)
		If pt(i).lat = 0 Then errorCount += 1
		If pt(i).lon = 0 Then errorCount += 1
		If pt(i).ele = 0 Then errorCount += 1
		If pt(i).dateTime = 0 Then errorCount += 1
	Next
	Return errorCount
End Function

'===============================================================================

#Include "vbcompat.bi"

Enum E_DATE_FORMAT
	DATE_FORMAT_ISO
	DATE_FORMAT_NL
	DATE_FORMAT_DE
	DATE_FORMAT_US
End Enum

Function dateString(dateTimeVal As Double, dateFormat As E_DATE_FORMAT) As String
	Select Case dateFormat
	Case DATE_FORMAT_ISO
		Return Format(dateTimeVal, "yyyy-mm-dd")
	Case DATE_FORMAT_NL
		Return Format(dateTimeVal, "dd-mm-yyyy")
	Case DATE_FORMAT_DE
		Return Format(dateTimeVal, "dd.mm.yyyy")
	Case DATE_FORMAT_US
		Return Format(dateTimeVal, "mm/dd/yyyy")
	End Select
	Return "<invalid>"
End Function

'day, month, year format must be: dd, mm, yyyy (e.g. not 7/5/20 or 2020-7-5)
Function dateValue2(dateStr As String, dateFormat As E_DATE_FORMAT) As Long
	If Len(dateStr) <> 10 Then Return -1
	Dim As String y, m, d 'years, months, days
	Select Case dateFormat
	Case DATE_FORMAT_ISO
		If Mid(dateStr, 5, 1) <> "-" Then Return -1
		If Mid(dateStr, 8, 1) <> "-" Then Return -1
		y = Mid(dateStr, 1, 4)
		m = Mid(dateStr, 6, 2)
		d = Mid(dateStr, 9, 2)
	Case DATE_FORMAT_NL
		If Mid(dateStr, 3, 1) <> "-" Then Return -1
		If Mid(dateStr, 6, 1) <> "-" Then Return -1
		y = Mid(dateStr, 7, 4)
		m = Mid(dateStr, 4, 2)
		d = Mid(dateStr, 1, 2)
	Case DATE_FORMAT_DE
		If Mid(dateStr, 3, 1) <> "." Then Return -1
		If Mid(dateStr, 6, 1) <> "." Then Return -1
		y = Mid(dateStr, 7, 4)
		m = Mid(dateStr, 4, 2)
		d = Mid(dateStr, 1, 2)
	Case DATE_FORMAT_US
		If Mid(dateStr, 3, 1) <> "/" Then Return -1
		If Mid(dateStr, 6, 1) <> "/" Then Return -1
		y = Mid(dateStr, 7, 4)
		m = Mid(dateStr, 1, 2)
		d = Mid(dateStr, 4, 2)
	Case Else
		Return -1
	End Select
	Return DateSerial(Val(y), Val(m), Val(d))
End Function

'Convert 2020-05-03T16:28:56Z to 43954.68675925926
Function gpxDateTimeValue(gpxDateTimeStr As String) As Double
	Return dateValue2(Mid(gpxDateTimeStr, 1, 10), DATE_FORMAT_ISO) _
		+ TimeValue(Mid(gpxDateTimeStr, 12, 8))
End Function

Function DiffDateTimeSec(dt1 As Double, dt2 As Double) As Double
	Return (dt2 - dt1) * 24 * 3600
End Function 

'===============================================================================

Const As Double PI = 4 * Atn(1)

'returns < 0 or error
Function readGpxFile(fileName As String, track As trkpt_list) As Integer
	Dim As xmlTextReaderPtr pReader = xmlReaderForFile(filename, NULL, 0)
	If (pReader = NULL) Then
		'print "readGpxFile: Unable to open: "; fileName
		Return -1
	End If
	
	Dim As Const ZString Ptr pConstName, pConstValue
	Dim As Const ZString Ptr pAttrName, pAttrValue
	Dim As Integer nodetype, iTrkpt = -1

	Dim As Integer ret = xmlTextReaderRead(pReader)
	While (ret = 1)
		pConstName = xmlTextReaderConstName(pReader)
		pConstValue = xmlTextReaderConstValue(pReader)
		nodetype = xmlTextReaderNodeType(pReader)
		If (nodetype = 1) Then
			If (*pConstName = "trkpt") Then
				iTrkpt = track.add()
				If xmlTextReaderHasAttributes(pReader) = 1 Then
					If xmlTextReaderMoveToAttributeNo(pReader, 0) Then
						pAttrName = xmlTextReaderConstName(pReader)
						pAttrvalue = xmlTextReaderConstValue(pReader)
						If *pAttrName = "lat" Then
							track.pt(iTrkpt).lat = CDbl(*pAttrvalue)
							'print "lat: " & *pAttrvalue
						End If
					End If
					If xmlTextReaderMoveToAttributeNo(pReader, 1) Then
						pAttrName = xmlTextReaderConstName(pReader)
						pAttrvalue = xmlTextReaderConstValue(pReader)
						If *pAttrName = "lon" Then
							track.pt(iTrkpt).lon = CDbl(*pAttrvalue)
							'print "lon: " & *pAttrvalue
						End If
					End If
				End If
			End If
			If (*pConstName = "time") Then
				'actual data is in next element
				ret = xmlTextReaderRead(pReader) 'read next
				pConstName = xmlTextReaderConstName(pReader)
				pConstValue = xmlTextReaderConstValue(pReader)
				If (ret = 1) Then
					track.pt(iTrkpt).dateTime = gpxDateTimeValue(*Cast(ZString Ptr, pConstValue))
					'print "datetime: " & *pConstValue
				End If
			End If
			If (*pConstName = "ele") Then
				'actual data is in next element
				ret = xmlTextReaderRead(pReader) 'read next
				pConstName = xmlTextReaderConstName(pReader)
				pConstValue = xmlTextReaderConstValue(pReader)
				If (ret = 1) Then
					track.pt(iTrkpt).ele = CDbl(*pConstValue)
					'print "ele: " & *pConstValue
					'print
					'getkey()
				End If
			End If
		End If
		'getkey()
		ret = xmlTextReaderRead(pReader)
	Wend

	xmlFreeTextReader(pReader)

	If (ret <> 0) Then Print "readGpxFile: Failed to parse: "; filename

	xmlCleanupParser()
	xmlMemoryDump()
End Function

'-------------------------------------------------------------------------------

'simple spherical earth approximation
Function deltaPos(pt1 As trkpt_type, pt2 As trkpt_type) As dbl2d
	Dim As dbl2d dPos
	dPos.x = (pt1.lon - pt2.lon) * 4e7 * Cos((pt1.lat + pt2.lat) * PI / 360) / 360
	dPos.y = (pt1.lat - pt2.lat) * 4e7 / 360
	Return dPos
End Function

'-------------------------------------------------------------------------------

Type track_stats
	Dim As Double totalDist 'm
	Dim As Double totalTime 's
	Dim As Double maxSpeed, avgSpeed 'm/s
	Dim As Double minXpos, maxXpos 'm
	Dim As Double minYpos, maxYpos 'm
	Dim As Double minZpos, maxZpos 'm
	Declare Sub reset_()
	Declare Sub print_()
	Declare Sub update(track As trkpt_list)
End Type

Sub track_stats.reset_()
	totalDist = 0 : totalTime = 0
	maxSpeed = 0 : avgSpeed = 0
	minXpos = 0 : maxXpos = 0
	minYpos = 0 : maxYpos = 0
	minZpos = 0 : maxZpos = 0
End Sub

Sub track_stats.print_()
	Print !"\nTrack stats:"
	Print " totalDist [km]: " & Format(totalDist / 1000, "0.000")
	Print " totalTime [s]: " & Format(totalTime, "0.0")
	Print " avgSpeed [km/h]: " & Format(avgSpeed * 3.6, "0.000") 
	Print " maxSpeed [km/h]: " & Format(maxSpeed * 3.6, "0.000")
	Print " X-range [km]: " & Format(minXpos / 1000, "0.000") & " ... " & Format(maxXpos / 1000, "0.000")
	Print " Y-range [km]: " & Format(minYpos / 1000, "0.000") & " ... " & Format(maxYpos / 1000, "0.000")
	Print " Z-range [m]: " & Format(minZpos, "0.0") & " ... " & Format(maxZpos, "0.0")
End Sub

Sub track_stats.update(track As trkpt_list)
	Dim As dbl2d dPos
	Dim As Double dist, timeDiff
	Dim As Integer last = track.size() - 1
	reset_()
	totalTime = DiffDateTimeSec(track.pt(0).dateTime, track.pt(last).dateTime)
	'compare sequential track points
	track.pt(0).speed = 0
	For i As Integer = 1 To last
		dPos = deltaPos(track.pt(i - 1), track.pt(i))
		dist = Sqr(dPos.x * dPos.x + dPos.y * dPos.y)
		totalDist += dist
		timeDiff = DiffDateTimeSec(track.pt(i - 1).dateTime, track.pt(i).dateTime)
		track.pt(i).speed = dist / timeDiff
		If track.pt(i).speed > maxSpeed Then maxSpeed = track.pt(i).speed
	Next
	avgSpeed = totalDist / totalTime
	'compare track points against first point
	track.pt(0).p = dbl2d(0, 0)
	minZpos = track.pt(0).ele
	maxZpos = track.pt(0).ele
	For i As Integer = 1 To last
		Dim As dbl2d dPos = deltaPos(track.pt(0), track.pt(i))
		track.pt(i).p = dPos
		If dPos.x < minXpos Then minXpos = dPos.x 
		If dPos.x > maxXpos Then maxXpos = dPos.x 
		If dPos.y < minYpos Then minYpos = dPos.y 
		If dPos.y > maxYpos Then maxYpos = dPos.y
		Dim As Double ele = track.pt(i).ele
		If ele < minZpos Then minZpos = ele 
		If ele > maxZpos Then maxZpos = ele
	Next
End Sub

'-------------------------------------------------------------------------------

Dim As trkpt_list track
Dim As track_stats stats
Dim As String fileName = "gpx_files/Skiing-downhill.gpx"
Dim As Integer result = readGpxFile(fileName, track)

Print "fileName: " & fileName
Print "readGpxFile: " & IIf(result >= 0, "Ok", "error")
Print "Num track points: " & track.size()
Print "Track error count: " & track.check()

If (result >= 0) And (track.size() > 0) And (track.check() = 0) Then

	stats.update(track)
	stats.print_()

	track.pt(0).p = dbl2d(0, 0)
	For i As Integer = 1 To track.size() - 1
		Dim As dbl2d dPos = deltaPos(track.pt(i), track.pt(0))
		Dim As Double dist = Sqr(dPos.x * dPos.x + dPos.y * dPos.y)
		track.pt(i).p = dPos
		'track.pt(i).speed = dPos / time
	Next
	
	Print !"\nFirst and last 10 data points:"
	'print " Index", "yyyy-mm-dd - hh:mm:ss", "Latitude", "Longitude", "Elevation", "x-pos", "y-pos"
	Print " Index", "yyyy-mm-dd - hh:mm:ss", "Elevation [m]", "X-pos [m]", "Y-pos [m]", "Speed [m/s]"
	If track.size() > 20 Then
		For i As Integer = 0 To track.size() - 1
			If i > 9 And i < track.size() - 10 Then Continue For
			With track.pt(i)
				Print i, Format(.dateTime, "yyyy-mm-dd - hh:mm:ss"), _
					Format(.ele, "#.000"), _
					Format(.p.x, "#.000"), Format(.p.y, "#.000"), _
					Format(.speed, "#.000")
					'format(.lat, "#.000000"), format(.lon, "#.000000"), _
			End With
		Next
	End If
End If

GetKey()

Outputs:

Code: Select all

fileName: gpx_files/Skiing-downhill.gpx
readGpxFile: Ok
Num track points: 680
Track error count: 0

Track stats:
 totalDist [km]: 6.771
 totalTime [s]: 2843.0
 avgSpeed [km/h]: 8.573
 maxSpeed [km/h]: 45.632
 X-range [km]: -0.890 ... 1.234
 Y-range [km]: -1.894 ... 0.387
 Z-range [m]: 1544.5 ... 2440.0

First and last 10 data points:
 Index        yyyy-mm-dd - hh:mm:ss       Elevation [m] X-pos [m]     Y-pos [m]     Speed [m/s]
 0            2015-02-24 - 13:29:48       2089.000      0.000         0.000         0.000
 1            2015-02-24 - 13:29:50       2089.000      .815          -2.404        1.269
 2            2015-02-24 - 13:29:58       2089.400      -.871         -3.429        .247
 3            2015-02-24 - 13:30:15       2093.300      -3.442        -6.540        .237
 4            2015-02-24 - 13:31:10       2097.300      -5.312        -9.252        .060
 5            2015-02-24 - 13:31:56       2099.700      -.940         -8.967        .095
 6            2015-02-24 - 13:32:07       2095.900      -.909         -7.247        .156
 7            2015-02-24 - 13:32:18       2096.200      3.403         -1.152        .679
 8            2015-02-24 - 13:32:30       2096.400      2.116         -1.393        .109
 9            2015-02-24 - 13:32:37       2096.600      .329          -1.622        .257
 670          2015-02-24 - 14:15:56       1544.500      837.960       1154.363      4.662
 671          2015-02-24 - 14:15:58       1544.800      846.892       1153.298      4.498
 672          2015-02-24 - 14:16:00       1545.700      855.286       1152.424      4.219
 673          2015-02-24 - 14:16:02       1546.000      863.162       1152.122      3.941
 674          2015-02-24 - 14:16:04       1548.000      870.214       1151.724      3.531
 675          2015-02-24 - 14:16:06       1550.000      874.688       1151.069      2.260
 676          2015-02-24 - 14:16:09       1552.000      874.474       1151.451      .146
 677          2015-02-24 - 14:16:16       1552.800      868.437       1153.863      .929
 678          2015-02-24 - 14:16:26       1549.700      861.666       1156.332      .721
 679          2015-02-24 - 14:17:11       1549.700      861.666       1156.332      0.000

Preview of the graphical display (not yet in code above):

Red = slow, Green = fast
There are clearly a few gaps in the logging and outlier data point(s) halfway.

Same file with a different viewer:

badidea · Post by **badidea** » Jul 11, 2020 22:38

Still working on the viewer part, code follows later. But already I see some interesting details. My current phone sometimes records speeds of 75 km/h on a cycling trip, which is not very realistic. I expected a position error at certain moments, but this is not the case. The time seems to be the problem. Possibly the phone clock runs slower then the GPS clock and is corrected at certain intervals. On the image below, the green sections are high speed anomalies, most points are 1 second apart.

TJF · Post by **TJF** » Jul 12, 2020 7:51

badidea wrote:My current phone sometimes records speeds of 75 km/h on a cycling trip, which is not very realistic. I expected a position error at certain moments, but this is not the case. The time seems to be the problem. Possibly the phone clock runs slower then the GPS clock and is corrected at certain intervals. On the image below, the green sections are high speed anomalies, most points are 1 second apart.

Which software do you use to record the track points? Does it trigger on time or on distance? Anyhow, the software should use the GPS-time, not the phone time.

For me it looks like the tracking software is blocked by other (background) processes. There's not enough CPU-time left, and the tracker doesn't compensate that situation sufficiently.

For track recording I use MapFactor Navigator. I never saw such gaps in my tracks.

Regards

badidea · Post by **badidea** » Jul 12, 2020 9:37

TJF wrote:Which software do you use to record the track points? Does it trigger on time or on distance? Anyhow, the software should use the GPS-time, not the phone time.

For me it looks like the tracking software is blocked by other (background) processes. There's not enough CPU-time left, and the tracker doesn't compensate that situation sufficiently.

I use Laufhelden for Sailfish OS. I see that the developer stopped and that someone forked it.
But it seems that my phone was not using GPS at all. Quite amazing to get these results at all, I think.
I'll make another cycling trip today, hopefully better logging now with GPS on.

Edit 13-07-2020:
That was not the problem. I reached 245 km/h according to my program and when I upload the track sports-tracker.com
The weird thing is that the application on my phone reports a maximum speed of 33.2 km/h which is much more plausible.
I'll have to look a the exact data-points to see what is going on, maybe some auto-pause feature that goes wrong.

2nd edit:
I see nothing special in the logging:

Code: Select all

<trkpt lat="51.55826859" lon="5.60881288">
	<time>2020-07-12T15:31:33Z</time>
	<ele>77.4937744140625</ele>
	<extensions>
		<dir>138.199996948242</dir>
		<g_spd>3.4499990940094</g_spd>
		<h_acc>10.7200002670288</h_acc>
		<v_acc>8</v_acc>
	</extensions>
</trkpt>
<trkpt lat="51.55776789" lon="5.60937846">
	<time>2020-07-12T15:31:34Z</time>
	<ele>75.4686279296875</ele>
	<extensions>
		<dir>139</dir>
		<g_spd>5.02999877929688</g_spd>
		<h_acc>11.7919998168945</h_acc>
		<v_acc>8</v_acc>
	</extensions>
</trkpt>

These 2 data-points are 1 second apart, but the distance corresponds to 68 meter (together give 245 km/h).
Must be something weird in the application.
The <extensions> I ignore for now.

TJF · Post by **TJF** » Jul 14, 2020 14:41

It's a bad idea to differentiate the speed as distance vs. time. The position may have a deviation of more then ten meters each. In your case it seems to be even more.

I don't understand your function deltaPos(), where the elevation is unconsidered. I get

Code: Select all

FUNCTION distance( _
  BYVAL lat1 AS DOUBLE, BYVAL lon1 AS DOUBLE _
, BYVAL lat2 AS DOUBLE, BYVAL lon2 AS DOUBLE _
, BYVAL ele1 AS DOUBLE = 0, BYVAL ele2 AS DOUBLE = 0) AS DOUBLE
    STATIC AS DOUBLE fac = ATN(1) / 45 ' = PI/180
    VAR flat1 = lat1 * fac _
      , flat2 = lat2 * fac _
      , cosX = SIN(flat1) * SIN(flat2) _
             + COS(flat1) * COS(flat2) * COS((lon1 - lon2) * fac)
    RETURN ACOS(cosX) * (3443.89849 _ ' earth radius (nautical miles)
                      + .5 * (ele1 + ele2))
END FUNCTION

VAR miles = distance(51.55826859, 5.60881288, _
                     51.55776789, 5.60937846, 77.4937744140625,75.4686279296875) _
  , meter = miles * 1852 _
  , speed = meter * 3.6

?"Distance: " & miles & " nm = " & meter & " m (" & speed & " km/h)"

resulting in

Code: Select all

Distance: 0.03759270038931566 nm = 68.4187147085545 m (246.3073729507962 km/h)

Anyhow, it's much more accurate (and more easy) to use the speed coming from the sensor (in <extensions>)

...
<g_spd>3.4499990940094</g_spd>
...
<g_spd>5.02999877929688</g_spd>
...

g_spd is the speed over ground in nautic miles, so

3.4499990940094 = 6.389398322 km/h
5.02999877929688 = 9,315557739 km/h

Regards

badidea · Post by **badidea** » Aug 02, 2021 22:41

While looking at this project again, I found this website which might be useful for others: https://www.thunderforest.com/docs/static-maps-api/
You can generate an earth map image of a given position, zoom, size, type, etc. (see link).
To use it, you can create your own API-key for free and without sharing personal details other than an e-mail address.

Examples (with my API-key removed), centered view on the Dom tower in Utrecht (NL):
https://tile.thunderforest.com/static/a ... y=xxxxxxxx with gives:

And https://tile.thunderforest.com/static/l ... y=xxxxxxxx gives:

Now to combine this with my GPX track reader...

badidea · Post by **badidea** » Aug 17, 2021 23:30

I succeeded in plotting the GPX tracks on top of a map image from this Thunderforrest website.
See: https://github.com/verybadidea/gpx_track_plotter
Libraries used:
- libxml2 (http://xmlsoft.org/)
- SNC (viewtopic.php?f=7&t=23421) for web get.
- FBImage (viewtopic.php?t=24105) for png support.
Added GPX example files are recorded 5 years ago. You will not find me anymore at the start or end of those tracks :-)
The track color depends on speed, green is faster (cycling) and red is slower (walking or running).
I still have to implement the filtering (which was the actual reason I started this project).

badidea · Post by **badidea** » Aug 22, 2021 14:20

I added a height smoothing filter. The large variation in height caused an overestimation of the distance and speed.
In addition, one can zoom in on an area and get the corresponding map image after key press.
(there is a bug when zooming out + panning, due to conversion error form longitude/latitude to x/y [Cartesian] position and back)
More info: https://github.com/verybadidea/gpx_track_plotter

GPX file processing

GPX file processing

Re: GPX file processing

Re: GPX file processing

Re: GPX file processing

Re: GPX file processing

Re: GPX file processing

Re: GPX file processing

Re: GPX file processing

Re: GPX file processing

Re: GPX file processing