Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with HTTPS or Subversion.

Download ZIP
Browse files

Preserve significant figures of the original GPS coordinates

The coordinates are stored as EXIF rational numbers which hold
the number of significant figures of the data implicitly.  Using
a fixed rational denominator of 10000 for altitude (as was
previously done, for example) implies a precision of 0.1
millimeters, which is more what GPS technology can provide (and
is well beyond what consumer units can do).  The denominator is
now scaled to the number of decimal places in the original XML
file.  Since at least some GPS units will specify the number of
decimal points in relation to the confidence of the location
reading, this preserves the confidence level in the EXIF tags.
  • Loading branch information...
commit e42132fd55d2736451a47fcf1976f00cf7b770ed 1 parent da1540a
@dfandrich dfandrich authored
View
4 RELEASES
@@ -117,8 +117,10 @@ v1.X: TBD
- Dan Fandrich took over maintenance of gpscorrelate
- Fixed to write GPSTimeStamp tag as unsigned rational, as per spec
- Added automatic time zone offset detection by default (the previous
- behaviour can be selected with "-z 0")
+ behaviour can be selected with "-z 0")
- GUI now displays the version number in the title bar
- File loading dialogs now have appropriate filters on file extensions
- Stop documenting the -p option, which never worked (--degmins hasn't
been touched)
+ - Store GPS coordinates in such a way as to preserve the number of
+ significant figures of the original data
View
11 correlate.c
@@ -35,6 +35,8 @@
#include "correlate.h"
#include "unixtime.h"
+#define MIN(a,b) (((a)<(b))?(a):(b))
+
/* Internal functions used to make it work. */
static void Round(const struct GPSPoint* First, struct GPSPoint* Result,
time_t PhotoTime);
@@ -190,8 +192,11 @@ struct GPSPoint* CorrelatePhoto(const char* Filename,
/* This is the point, exactly.
* Copy out the data and return that. */
Actual->Lat = Search->Lat;
+ Actual->LatDecimals = Search->LatDecimals;
Actual->Long = Search->Long;
+ Actual->LongDecimals = Search->LongDecimals;
Actual->Elev = Search->Elev;
+ Actual->ElevDecimals = Search->ElevDecimals;
Actual->Time = Search->Time;
Options->Result = CORR_OK;
@@ -280,8 +285,11 @@ void Round(const struct GPSPoint* First, struct GPSPoint* Result,
/* Copy the numbers over... */
Result->Lat = CopyFrom->Lat;
+ Result->LatDecimals = CopyFrom->LatDecimals;
Result->Long = CopyFrom->Long;
+ Result->LongDecimals = CopyFrom->LongDecimals;
Result->Elev = CopyFrom->Elev;
+ Result->ElevDecimals = CopyFrom->ElevDecimals;
Result->Time = CopyFrom->Time;
/* Done! */
@@ -303,13 +311,16 @@ void Interpolate(const struct GPSPoint* First, struct GPSPoint* Result,
/* Now calculate the Latitude. */
Result->Lat = First->Lat + ((First->Next->Lat - First->Lat) * Scale);
+ Result->LatDecimals = MIN(First->LatDecimals, First->Next->LatDecimals);
/* And the longitude. */
Result->Long = First->Long + ((First->Next->Long - First->Long) * Scale);
+ Result->LongDecimals = MIN(First->LongDecimals, First->Next->LongDecimals);
/* And the elevation. If elevation wasn't set, it should be zero.
* Which works quite fine for us. */
Result->Elev = First->Elev + ((First->Next->Elev - First->Elev) * Scale);
+ Result->ElevDecimals = MIN(First->ElevDecimals, First->Next->ElevDecimals);
/* The time is not interpolated, but matches photo. */
Result->Time = PhotoTime;
View
107 exif-gps.cpp
@@ -47,6 +47,9 @@
#include "gpsstructure.h"
#include "exif-gps.h"
+#define MAX(a,b) (((a)>(b))?(a):(b))
+#define MIN(a,b) (((a)<(b))?(a):(b))
+
/* Debug
int main(int argc, char* argv[])
{
@@ -369,14 +372,12 @@ char* ReadGPSTimestamp(const char* File, char* DateStamp, char* TimeStamp, int*
return Copy;
};
-void ConvertToRational(double Number,
- long int* Numerator, long int* Denominator,
- int Rounding)
+static void ConvertToRational(double Number, int Decimals, char *Buf, int BufSize)
{
// This function converts the given decimal number
// to a rational (fractional) number.
//
- // Examples in comments use Number as 25.12345, Rounding as 4.
+ // Examples in comments use Number as 25.12345, Decimals as 4.
// Split up the number.
double Whole = trunc(Number);
@@ -384,7 +385,7 @@ void ConvertToRational(double Number,
// Calculate the "number" used for rounding.
// This is 10^Digits - ie, 4 places gives us 10000.
- double Rounder = pow(10, Rounding);
+ double Rounder = pow(10, Decimals);
// Round the fractional part, and leave the number
// as greater than 1.
@@ -426,13 +427,47 @@ void ConvertToRational(double Number,
}
// Copy out the numbers.
- *Numerator = (int)NumTemp;
- *Denominator = (int)DenTemp;
+ snprintf(Buf, BufSize, "%ld/%ld", (long)NumTemp, (long)DenTemp);
// And finished...
}
+/* Converts a floating point number with known significant decimal places
+ * into a string representation of a set of latitude or longitude rational
+ * numbers.
+ */
+static void ConvertToLatLongRational(double Number, int Decimals, char *Buf, int BufSize)
+{
+ int Deg, Min, Sec;
+ Deg = (int)floor(fabs(Number)); // Slice off after decimal.
+ Min = (int)floor((fabs(Number) - floor(fabs(Number))) * 60); // Now grab just the minutes.
+ double FracPart = ((fabs(Number) - floor(fabs(Number))) * 60) - (double)Min; // Grab the fractional minute.
+ // Calculate the appropriate denominator based on the number of
+ // significant figures in the original data point. Splitting off the
+ // minutes and integer seconds reduces the number of significant
+ // figures by 3.6 (log10(60*60)), so round it down to 3 in order to
+ // preserve the maximum precision. Cap it at 9 to avoid overflow
+ // in the EXIF rational data type.
+ int Multiplier = powl(10, MAX(0, MIN(Decimals - 3, 9)));
+ Sec = (int)floor(FracPart * 60 * Multiplier); // Convert to seconds.
+ snprintf(Buf, BufSize, "%d/1 %d/1 %d/%d", Deg, Min, Sec, Multiplier);
+ //printf("New style lat/long: %f -> %d/%d/ %d/%d\n", Number, Deg, Min, Sec, Multiplier);
+}
+
+/* Converts a floating point number into a string representation of a set of
+ * latitude or longitude rational numbers, using the older, not as accurate
+ * style, which nobody should really be using any more.
+ */
+static void ConvertToOldLatLongRational(double Number, char *Buf, int BufSize)
+{
+ int Deg, Min;
+ Deg = (int)floor(fabs(Number)); // Slice off after decimal.
+ Min = (int)floor((fabs(Number) - floor(fabs(Number))) * 6000);
+ snprintf(Buf, BufSize, "%d/1 %d/100 0/1", Deg, Min);
+ //printf("Old style lat/long: %f -> %s\n", Number, Buf);
+}
+
int WriteGPSData(const char* File, const struct GPSPoint* Point,
const char* Datum, int NoChangeMtime, int DegMinSecs)
{
@@ -464,9 +499,6 @@ int WriteGPSData(const char* File, const struct GPSPoint* Point,
Exiv2::ExifData &ExifToWrite = Image->exifData();
char ScratchBuf[100];
- long int Nom, Denom;
- long int Deg, Min, Sec;
- double FracPart;
// Do all the easy constant ones first.
// GPSVersionID tag: standard says it should be four bytes: 02 00 00 00
@@ -491,8 +523,9 @@ int WriteGPSData(const char* File, const struct GPSPoint* Point,
ExifToWrite.add(Exiv2::ExifKey("Exif.GPSInfo.GPSAltitudeRef"), Value.get());
// And the actual altitude.
Value = Exiv2::Value::create(Exiv2::unsignedRational);
- ConvertToRational(fabs(Point->Elev), &Nom, &Denom, 4);
- snprintf(ScratchBuf, 100, "%ld/%ld", Nom, Denom);
+ // 3 decimal points is beyond the limit of current GPS technology
+ int Decimals = MIN(Point->ElevDecimals, 3);
+ ConvertToRational(fabs(Point->Elev), Decimals, ScratchBuf, sizeof(ScratchBuf));
/* printf("Altitude: %f -> %s\n", Point->Elev, ScratchBuf); */
Value->read(ScratchBuf);
@@ -525,7 +558,9 @@ int WriteGPSData(const char* File, const struct GPSPoint* Point,
// as the sign is encoded in LatRef.
// Further note: original code did not translate between
// dd.dddddd to dd mm.mm - that's why we now multiply
- // by 6000 - x60 to get minutes, x100 to get to mmmm/100.
+ // by 60*N - x60 to get minutes, xN to get to mmmm/N.
+ // N is 10^S where S is the number of significant decimal
+ // places.
//
// Rereading the EXIF standard, it's quite ok to do DD MM SS.SS
// Which is much more accurate. This is the new default, unless otherwise
@@ -534,18 +569,9 @@ int WriteGPSData(const char* File, const struct GPSPoint* Point,
if (DegMinSecs)
{
- Deg = (int)floor(fabs(Point->Lat)); // Slice off after decimal.
- Min = (int)floor((fabs(Point->Lat) - floor(fabs(Point->Lat))) * 60); // Now grab just the minutes.
- FracPart = ((fabs(Point->Lat) - floor(fabs(Point->Lat))) * 60) - (double)Min; // Grab the fractional minute.
- Sec = (int)floor(FracPart * 6000); // Convert to seconds.
-
- /* printf("New style latitude: %f -> %ld/%ld/ %ld/100\n", Point->Lat, Deg, Min, Sec); */
-
- snprintf(ScratchBuf, 100, "%ld/1 %ld/1 %ld/100", Deg, Min, Sec);
+ ConvertToLatLongRational(Point->Lat, Point->LatDecimals, ScratchBuf, sizeof(ScratchBuf));
} else {
- Deg = (int)floor(fabs(Point->Lat)); // Slice off after decimal.
- Min = (int)floor((fabs(Point->Lat) - floor(fabs(Point->Lat))) * 6000);
- snprintf(ScratchBuf, 100, "%ld/1 %ld/100 0/1", Deg, Min);
+ ConvertToOldLatLongRational(Point->Lat, ScratchBuf, sizeof(ScratchBuf));
}
Value->read(ScratchBuf);
ExifToWrite.add(Exiv2::ExifKey("Exif.GPSInfo.GPSLatitude"), Value.get());
@@ -562,37 +588,14 @@ int WriteGPSData(const char* File, const struct GPSPoint* Point,
// Eastern hemisphere. Where I live.
ExifToWrite["Exif.GPSInfo.GPSLongitudeRef"] = "E";
}
- // Now the actual longitude itself.
- // This is done as three rationals.
- // I choose to do it as:
- // dd/1 - degrees.
- // mmmm/100 - minutes
- // 0/1 - seconds
- // Exif standard says you can do it with minutes
- // as mm/1 and then seconds as ss/1, but its
- // (slightly) more accurate to do it as
- // mmmm/100 than to split it.
- // We also absolute the value (with fabs())
- // as the sign is encoded in LongRef.
- // Further note: original code did not translate between
- // dd.dddddd to dd mm.mm - that's why we now multiply
- // by 6000 - x60 to get minutes, x100 to get to mmmm/100.
+ // Now the actual longitude itself, in the same way as latitude
Value = Exiv2::Value::create(Exiv2::unsignedRational);
if (DegMinSecs)
{
- Deg = (int)floor(fabs(Point->Long)); // Slice off after decimal.
- Min = (int)floor((fabs(Point->Long) - floor(fabs(Point->Long))) * 60); // Now grab just the minutes.
- FracPart = ((fabs(Point->Long) - floor(fabs(Point->Long))) * 60) - (double)Min; // Grab the fractional minute.
- Sec = (int)floor(FracPart * 6000); // Convert to seconds.
-
- /* printf("New style longitude: %f -> %ld/%ld/ %ld/100\n", Point->Long, Deg, Min, Sec); */
-
- snprintf(ScratchBuf, 100, "%ld/1 %ld/1 %ld/100", Deg, Min, Sec);
+ ConvertToLatLongRational(Point->Long, Point->LongDecimals, ScratchBuf, sizeof(ScratchBuf));
} else {
- Deg = (int)floor(fabs(Point->Long)); // Slice off after decimal.
- Min = (int)floor((fabs(Point->Long) - floor(fabs(Point->Long))) * 6000);
- snprintf(ScratchBuf, 100, "%ld/1 %ld/100 0/1", Deg, Min);
+ ConvertToOldLatLongRational(Point->Long, ScratchBuf, sizeof(ScratchBuf));
}
Value->read(ScratchBuf);
ExifToWrite.add(Exiv2::ExifKey("Exif.GPSInfo.GPSLongitude"), Value.get());
@@ -619,7 +622,7 @@ int WriteGPSData(const char* File, const struct GPSPoint* Point,
}
Value = Exiv2::Value::create(Exiv2::unsignedRational);
- snprintf(ScratchBuf, 100, "%d/1 %d/1 %d/1",
+ snprintf(ScratchBuf, sizeof(ScratchBuf), "%d/1 %d/1 %d/1",
TimeStamp.tm_hour, TimeStamp.tm_min,
TimeStamp.tm_sec);
Value->read(ScratchBuf);
@@ -627,7 +630,7 @@ int WriteGPSData(const char* File, const struct GPSPoint* Point,
// And we should also do a datestamp.
Value = Exiv2::Value::create(Exiv2::unsignedRational);
- snprintf(ScratchBuf, 100, "%d/1 %d/1 %d/1",
+ snprintf(ScratchBuf, sizeof(ScratchBuf), "%d/1 %d/1 %d/1",
TimeStamp.tm_year + 1900,
TimeStamp.tm_mon + 1,
TimeStamp.tm_mday);
View
3  gpsstructure.h
@@ -30,8 +30,11 @@
struct GPSPoint {
double Lat;
+ int LatDecimals;
double Long;
+ int LongDecimals;
double Elev;
+ int ElevDecimals;
time_t Time;
int EndOfSegment;
struct GPSPoint* Next;
View
17 gpx-read.c
@@ -39,6 +39,16 @@
static struct GPSPoint* FirstPoint;
static struct GPSPoint* LastPoint;
+/* Returns the number of decimal places in the given decimal number string */
+static int NumDecimals(const char *Decimal)
+{
+ const char *Dec = strchr(Decimal, '.');
+ if (Dec) {
+ return strspn(Dec+1,"0123456789");
+ }
+ return 0;
+}
+
static void ExtractTrackPoints(xmlNodePtr Start)
{
/* The pointer passed to us should be the start
@@ -128,16 +138,22 @@ static void ExtractTrackPoints(xmlNodePtr Start)
/* Clear the structure first... */
LastPoint->Lat = 0;
+ LastPoint->LatDecimals = 0;
LastPoint->Long = 0;
+ LastPoint->LongDecimals = 0;
LastPoint->Elev = 0;
+ LastPoint->ElevDecimals = 0;
LastPoint->Time = 0;
LastPoint->EndOfSegment = 0;
/* Write the data into LastPoint, which should be a new point. */
LastPoint->Lat = atof(Lat);
+ LastPoint->LatDecimals = NumDecimals(Lat);
LastPoint->Long = atof(Long);
+ LastPoint->LongDecimals = NumDecimals(Long);
if (Elev) {
LastPoint->Elev = atof(Elev);
+ LastPoint->ElevDecimals = NumDecimals(Elev);
}
LastPoint->Time = ConvertToUnixTime(Time, GPX_DATE_FORMAT, 0, 0);
@@ -145,6 +161,7 @@ static void ExtractTrackPoints(xmlNodePtr Start)
printf("TrackPoint. Lat %s (%f), Long %s (%f). Elev %s (%f), Time %d.\n",
Lat, atof(Lat), Long, atof(Long), Elev, atof(Elev),
ConvertToUnixTime(Time, GPX_DATE_FORMAT, 0, 0));
+ printf("Decimals %d %d %d\n", LastPoint->LatDecimals, LastPoint->LongDecimals, LastPoint->ElevDecimals);
*/
View
2  gui.c
@@ -898,7 +898,7 @@ void SetListItem(GtkTreeIter* Iter, const char* Filename, char* Time, double Lat
/* Radius of earth ~6000km */
if (Elev > -7000000)
{
- snprintf(ElevScratch, sizeof(ElevScratch), "%fm", Elev);
+ snprintf(ElevScratch, sizeof(ElevScratch), "%.2fm", Elev);
} else {
snprintf(ElevScratch, sizeof(ElevScratch), " ");
}
View
4 main-command.c
@@ -110,10 +110,10 @@ static void ShowFileDetails(const char* File, int MachineReadable)
* it machine readable. */
if (MachineReadable)
{
- printf("\"%s\",\"%s\",%f,%f,%f\n",
+ printf("\"%s\",\"%s\",%f,%f,%.3f\n",
File, Time, Lat, Long, Elev);
} else {
- printf("%s: %s, Lat %f, Long %f, Elevation %f.\n",
+ printf("%s: %s, Lat %f, Long %f, Elevation %.3f.\n",
File, Time, Lat, Long, Elev);
}
} else {
Please sign in to comment.
Something went wrong with that request. Please try again.