fix misnamed turning_circle parameter

pfaedle.cfg

@@ -5,7 +5,7 @@
 [tram, bus, coach, subway, rail, gondola, funicular, ferry]
 
 routing_transition_penalty_fac:     0.0083
-routing_station_move_penalty_fac:	0.00087
+routing_station_move_penalty_fac:	0.002
 
 # Regular expressions and station comparision is
 # always case insensitive!
@@ -717,8 +717,8 @@ osm_filter_station:
 	highway=bus_stop
 	amenity=bus_station
 
-osm_filter_turning_cycle:
+osm_filter_turning_circle:
-	highway=turning_cycle
+	highway=turning_circle
 	highway=turning_loop
 	junction=roundabout
 	highway=mini_roundabout
@@ -791,7 +791,7 @@ routing_line_station_from_unmatched_time_penalty:   1.1
 # routing_no_lines_penalty_fac:  1
 
 # If the station name does not match, add this penalty
-routing_station_unmatched_penalty: 0.1
+routing_station_unmatched_penalty: 0.2
 
 # Punishment (in seconds) to add to the distance
 # function if a vehicle performans a full turn

src/pfaedle/config/MotConfigReader.cpp

@@ -152,9 +152,9 @@ void MotConfigReader::parse(const std::vector<std::string>& paths,
       }
     }
 
-    if (p.hasKey(secStr, "osm_filter_turning_cycle")) {
+    if (p.hasKey(secStr, "osm_filter_turning_circle")) {
       for (const auto& kvs :
-           p.getStrArr(sec.first, "osm_filter_turning_cycle", ' ')) {
+           p.getStrArr(sec.first, "osm_filter_turning_circle", ' ')) {
         auto fRule = getFRule(kvs);
         cfg.osmBuildOpts.turnCycleFilter[fRule.kv.first].insert(
             osm::AttrFlagPair(fRule.kv.second, getFlags(fRule.flags)));

src/shapevl/Collector.cpp

@@ -62,10 +62,10 @@ double Collector::add(const Trip* oldT, const Shape* oldS, const Trip* newT,
   double unmatchedSegmentsLength;  // total _an. length of unmatched segments
 
   std::vector<double> oldDists;
-  LINE oldL = getWebMercLine(oldS, &oldDists);
+  LINE oldL = getLine(oldS, &oldDists);
 
   std::vector<double> newDists;
-  LINE newL = getWebMercLine(newS, &newDists);
+  LINE newL = getLine(newS, &newDists);
 
   std::vector<std::pair<double, double>> newLenDists;
   std::vector<std::pair<double, double>> oldLenDists;
@@ -88,13 +88,8 @@ double Collector::add(const Trip* oldT, const Shape* oldS, const Trip* newT,
     newLCut.insert(newLCut.end(), newL.begin(), newL.end());
   }
 
-  // determine the scale factor between the distance in projected
+  // convert (roughly) to degrees
-  // coordinates and the real-world distance in meters
+  double SEGL = 15.0 / util::geo::M_PER_DEG;
-  auto avgY =
-      (oldSegs.front().front().getY() + oldSegs.back().back().getY()) / 2;
-  double fac = cos(2 * atan(exp(avgY / 6378137.0)) - 1.5707965);
-
-  double SEGL = 10;
 
   auto old = _dCache.find(oldLCut);
   if (old != _dCache.end()) {
@@ -102,11 +97,11 @@ double Collector::add(const Trip* oldT, const Shape* oldS, const Trip* newT,
     if (match != old->second.end()) {
       fd = match->second;
     } else {
-      fd = util::geo::accFrechetDistC(oldLCut, newLCut, SEGL / fac) * fac;
+      fd = util::geo::accFrechetDistCHav(oldLCut, newLCut, SEGL);
       _dCache[oldLCut][newLCut] = fd;
     }
   } else {
-    fd = util::geo::accFrechetDistC(oldLCut, newLCut, SEGL / fac) * fac;
+    fd = util::geo::accFrechetDistCHav(oldLCut, newLCut, SEGL);
     _dCache[oldLCut][newLCut] = fd;
   }
 
@@ -115,7 +110,7 @@ double Collector::add(const Trip* oldT, const Shape* oldS, const Trip* newT,
   unmatchedSegmentsLength = dA.second;
 
   double totL = 0;
-  for (auto l : oldSegs) totL += util::geo::len(l) * fac;
+  for (auto l : oldSegs) totL += util::geo::latLngLen(l);
 
   // filter out shapes with a length of under 5 meters - they are most likely
   // artifacts
@@ -171,8 +166,8 @@ std::vector<LINE> Collector::segmentize(
   size_t i = 0;
   for (auto st : t->getStopTimes()) {
     cuts.push_back(std::pair<POINT, double>(
-        util::geo::latLngToWebMerc<PFDL_PREC>(st.getStop()->getLat(),
+          {st.getStop()->getLng(),
-                                              st.getStop()->getLng()),
+                                              st.getStop()->getLat()},
         st.getShapeDistanceTravelled()));
     i++;
   }
@@ -204,7 +199,7 @@ std::vector<LINE> Collector::segmentize(
                              t->getStopTimes()[i - 1].getStop()->getLng(),
                              t->getStopTimes()[i].getStop()->getLat(),
                              t->getStopTimes()[i].getStop()->getLng());
-    double len = util::geo::webMercLen(curL);
+    double len = util::geo::latLngLen(curL);
     lenDist.push_back({dist, len});
 
     ret.push_back(curL);
@@ -215,12 +210,11 @@ std::vector<LINE> Collector::segmentize(
 }
 
 // _____________________________________________________________________________
-LINE Collector::getWebMercLine(const Shape* s, std::vector<double>* dists) {
+LINE Collector::getLine(const Shape* s, std::vector<double>* dists) {
   LINE ret;
 
   for (size_t i = 0; i < s->getPoints().size(); i++) {
-    ret.push_back(util::geo::latLngToWebMerc<PFDL_PREC>(s->getPoints()[i].lat,
+    ret.push_back({s->getPoints()[i].lng, s->getPoints()[i].lat});
-                                                        s->getPoints()[i].lng));
     (*dists).push_back(s->getPoints()[i].travelDist);
   }
 
@@ -233,19 +227,6 @@ const std::set<Result>& Collector::getResults() const { return _results; }
 // _____________________________________________________________________________
 double Collector::getAvgDist() const { return _fdSum / _results.size(); }
 
-// _____________________________________________________________________________
-std::vector<double> Collector::getBins(double mind, double maxd, size_t steps) {
-  double bin = (maxd - mind) / steps;
-  double curE = mind + bin;
-
-  std::vector<double> ret;
-  while (curE <= maxd) {
-    ret.push_back(curE);
-    curE += bin;
-  }
-  return ret;
-}
-
 // _____________________________________________________________________________
 void Collector::printCsv(std::ostream* os,
                          const std::set<Result>& result) const {
@@ -404,32 +385,28 @@ std::pair<size_t, double> Collector::getDa(const std::vector<LINE>& a,
   assert(a.size() == b.size());
   std::pair<size_t, double> ret{0, 0};
 
-  // euclidean distance on web mercator is in meters on equator,
+  // convert (roughly) to degrees
-  // and proportional to cos(lat) in both y directions
+  double SEGL = 15.0 / util::geo::M_PER_DEG;
-
-  double fac = webMercDistFactor(a.front().front());
-
-  double SEGL = 10;
 
   for (size_t i = 0; i < a.size(); i++) {
-    double fd = 0;
+    double fdMeter = 0;
     auto old = _dACache.find(a[i]);
     if (old != _dACache.end()) {
       auto match = old->second.find(b[i]);
       if (match != old->second.end()) {
-        fd = match->second;
+        fdMeter = match->second;
       } else {
-        fd = util::geo::frechetDist(a[i], b[i], SEGL / fac) * fac;
+        fdMeter = util::geo::frechetDistHav(a[i], b[i], SEGL);
-        _dACache[a[i]][b[i]] = fd;
+        _dACache[a[i]][b[i]] = fdMeter;
       }
     } else {
-        fd = util::geo::frechetDist(a[i], b[i], SEGL / fac) * fac;
+        fdMeter = util::geo::frechetDistHav(a[i], b[i], SEGL);
-        _dACache[a[i]][b[i]] = fd;
+        _dACache[a[i]][b[i]] = fdMeter;
     }
 
-    if (fd >= 50) {
+    if (fdMeter >= 50) {
       ret.first++;
-      ret.second += util::geo::webMercLen(a[i]) * 100;
+      ret.second += util::geo::latLngLen(a[i]);
     }
   }
 

src/shapevl/Collector.h

@@ -31,8 +31,10 @@ struct lineCmp {
     }
 
     for (size_t i = 0; i < a.size(); i++) {
-      if (util::geo::dist(a[i], b[i]) > .1) {
+      if (util::geo::dist(a[i], b[i]) > .000001) {
-        return (a[i].getX() < b[i].getX()) || (a[i].getX() == b[i].getX() && a[i].getY() < b[i].getY());;
+        return (a[i].getX() < b[i].getX()) ||
+               (a[i].getX() == b[i].getX() && a[i].getY() < b[i].getY());
+        ;
       }
     }
 
@@ -83,7 +85,7 @@ class Collector {
   // Return the averaged average frechet distance
   double getAvgDist() const;
 
-  static LINE getWebMercLine(const Shape* s, std::vector<double>* dists);
+  static LINE getLine(const Shape* s, std::vector<double>* dists);
 
   double getAcc() const;
 
@@ -112,13 +114,11 @@ class Collector {
   std::ostream* _reportOut;
 
   std::pair<size_t, double> getDa(const std::vector<LINE>& a,
-                                         const std::vector<LINE>& b);
+                                  const std::vector<LINE>& b);
 
   static std::vector<LINE> segmentize(
       const Trip* t, const LINE& shape, const std::vector<double>& dists,
       std::vector<std::pair<double, double>>& lenDist);
-
-  static std::vector<double> getBins(double mind, double maxd, size_t steps);
 };
 
 }  // namespace eval

src/util/geo/Geo.h

@@ -56,7 +56,6 @@ const static double AVERAGING_STEP = 20;
 
 const static double M_PER_DEG = 111319.4;
 
-
 // _____________________________________________________________________________
 template <typename T>
 inline Box<T> pad(const Box<T>& box, double padding) {
@@ -1778,6 +1777,29 @@ inline RotatedBox<T> getOrientedEnvelopeAvg(MultiLine<T> ml) {
   return rbox;
 }
 
+// _____________________________________________________________________________
+template <typename T>
+inline double haversine(T lat1, T lon1, T lat2, T lon2) {
+  lat1 *= RAD;
+  lat2 *= RAD;
+
+  const double dLat = lat2 - lat1;
+  const double dLon = (lon2 - lon1) * RAD;
+
+  const double sDLat = sin(dLat / 2);
+  const double sDLon = sin(dLon / 2);
+
+  const double a = (sDLat * sDLat) + (sDLon * sDLon) * cos(lat1) * cos(lat2);
+  return 6378137.0 * 2.0 * asin(sqrt(a));
+}
+
+// _____________________________________________________________________________
+template <typename T>
+inline double haversine(const Point<T>& a, const Point<T>& b) {
+  return haversine(a.getY(), a.getX(), b.getY(), b.getX());
+}
+
+
 // _____________________________________________________________________________
 template <typename T>
 inline Line<T> densify(const Line<T>& l, double d) {
@@ -1878,6 +1900,81 @@ inline double accFrechetDistC(const Line<T>& a, const Line<T>& b, double d) {
   return ca[p.size() * q.size() - 1];
 }
 
+// _____________________________________________________________________________
+template <typename T>
+inline double frechetDistCHav(size_t i, size_t j, const Line<T>& p,
+                              const Line<T>& q, std::vector<float>& ca) {
+  // based on Eiter / Mannila
+  // http://www.kr.tuwien.ac.at/staff/eiter/et-archive/cdtr9464.pdf
+
+  if (ca[i * q.size() + j] > -1)
+    return ca[i * q.size() + j];
+  else if (i == 0 && j == 0)
+    ca[i * q.size() + j] = haversine(p[0], q[0]);
+  else if (i > 0 && j == 0)
+    ca[i * q.size() + j] =
+        std::max(frechetDistCHav(i - 1, 0, p, q, ca), haversine(p[i], q[0]));
+  else if (i == 0 && j > 0)
+    ca[i * q.size() + j] =
+        std::max(frechetDistCHav(0, j - 1, p, q, ca), haversine(p[0], q[j]));
+  else if (i > 0 && j > 0)
+    ca[i * q.size() + j] =
+        std::max(std::min(std::min(frechetDistCHav(i - 1, j, p, q, ca),
+                                   frechetDistCHav(i - 1, j - 1, p, q, ca)),
+                          frechetDistCHav(i, j - 1, p, q, ca)),
+                 haversine(p[i], q[j]));
+  else
+    ca[i * q.size() + j] = std::numeric_limits<float>::infinity();
+
+  return ca[i * q.size() + j];
+}
+
+// _____________________________________________________________________________
+template <typename T>
+inline double frechetDistHav(const Line<T>& a, const Line<T>& b, double d) {
+  // based on Eiter / Mannila
+  // http://www.kr.tuwien.ac.at/staff/eiter/et-archive/cdtr9464.pdf
+
+  auto p = densify(a, d);
+  auto q = densify(b, d);
+
+  std::vector<float> ca(p.size() * q.size(), -1.0);
+  double fd = frechetDistCHav(p.size() - 1, q.size() - 1, p, q, ca);
+
+  return fd;
+}
+
+// _____________________________________________________________________________
+template <typename T>
+inline double accFrechetDistCHav(const Line<T>& a, const Line<T>& b, double d) {
+  auto p = densify(a, d);
+  auto q = densify(b, d);
+
+  assert(p.size());
+  assert(q.size());
+
+  std::vector<float> ca(p.size() * q.size(), 0);
+
+  for (size_t i = 0; i < p.size(); i++)
+    ca[i * q.size() + 0] = std::numeric_limits<float>::infinity();
+  for (size_t j = 0; j < q.size(); j++)
+    ca[j] = std::numeric_limits<float>::infinity();
+  ca[0] = 0;
+
+  for (size_t i = 1; i < p.size(); i++) {
+    for (size_t j = 1; j < q.size(); j++) {
+      float d =
+          util::geo::haversine(p[i], q[j]) * util::geo::dist(p[i], p[i - 1]);
+      ca[i * q.size() + j] =
+          d + std::min(ca[(i - 1) * q.size() + j],
+                       std::min(ca[i * q.size() + (j - 1)],
+                                ca[(i - 1) * q.size() + (j - 1)]));
+    }
+  }
+
+  return ca[p.size() * q.size() - 1];
+}
+
 // _____________________________________________________________________________
 template <typename T>
 inline Point<T> latLngToWebMerc(double lat, double lng) {
@@ -1890,7 +1987,7 @@ inline Point<T> latLngToWebMerc(double lat, double lng) {
 
 // _____________________________________________________________________________
 template <typename T>
-//TODO: rename to lngLat
+// TODO: rename to lngLat
 inline Point<T> latLngToWebMerc(Point<T> lngLat) {
   return latLngToWebMerc<T>(lngLat.getY(), lngLat.getX());
 }
@@ -1904,28 +2001,6 @@ inline Point<T> webMercToLatLng(double x, double y) {
   return Point<T>(lon, lat);
 }
 
-// _____________________________________________________________________________
-template <typename T>
-inline double haversine(T lat1, T lon1, T lat2, T lon2) {
-  lat1 *= RAD;
-  lat2 *= RAD;
-
-  const double dLat = lat2 - lat1;
-  const double dLon = (lon2 - lon1) * RAD;
-
-  const double sDLat = sin(dLat / 2);
-  const double sDLon = sin(dLon / 2);
-
-  const double a = (sDLat * sDLat) + (sDLon * sDLon) * cos(lat1) * cos(lat2);
-  return 6378137.0 * 2.0 * asin(sqrt(a));
-}
-
-// _____________________________________________________________________________
-template <typename T>
-inline double haversine(const Point<T>& a, const Point<T>& b) {
-  return haversine(a.getY(), a.getX(), b.getY(), b.getX());
-}
-
 // _____________________________________________________________________________
 template <typename T>
 inline double webMercMeterDist(const Point<T>& a, const Point<T>& b) {
@@ -1959,6 +2034,15 @@ inline double webMercLen(const Line<T>& g) {
   return ret;
 }
 
+// _____________________________________________________________________________
+template <typename T>
+inline double latLngLen(const Line<T>& g) {
+  double ret = 0;
+  for (size_t i = 1; i < g.size(); i++) ret += haversine(g[i - 1], g[i]);
+  return ret;
+}
+
+
 // _____________________________________________________________________________
 template <typename G>
 inline double webMercDistFactor(const G& a) {