use double precision

This commit is contained in:
Patrick Brosi 2018-08-10 15:21:27 +02:00
parent 967963e180
commit 486a8136fd
24 changed files with 186 additions and 177 deletions

View file

@ -7,10 +7,10 @@
using pfaedle::osm::BBoxIdx;
// _____________________________________________________________________________
BBoxIdx::BBoxIdx(float padding) : _padding(padding), _size(0) {}
BBoxIdx::BBoxIdx(double padding) : _padding(padding), _size(0) {}
// _____________________________________________________________________________
void BBoxIdx::add(Box<float> box) {
void BBoxIdx::add(Box<double> box) {
// division by 83.000m is only correct here around a latitude deg of 25,
// but should be a good heuristic. 1 deg is around 63km at latitude deg of 44,
// and 110 at deg=0, since we usually dont do map matching in the arctic,
@ -24,21 +24,21 @@ void BBoxIdx::add(Box<float> box) {
size_t BBoxIdx::size() const { return _size; }
// _____________________________________________________________________________
bool BBoxIdx::contains(const Point<float>& p) const {
bool BBoxIdx::contains(const Point<double>& p) const {
return treeHas(p, _root);
}
// _____________________________________________________________________________
util::geo::Box<float> BBoxIdx::getFullWebMercBox() const {
return util::geo::FBox(
util::geo::latLngToWebMerc<float>(_root.box.getLowerLeft().getY(),
util::geo::Box<double> BBoxIdx::getFullWebMercBox() const {
return util::geo::DBox(
util::geo::latLngToWebMerc<double>(_root.box.getLowerLeft().getY(),
_root.box.getLowerLeft().getX()),
util::geo::latLngToWebMerc<float>(_root.box.getUpperRight().getY(),
util::geo::latLngToWebMerc<double>(_root.box.getUpperRight().getY(),
_root.box.getUpperRight().getX()));
}
// _____________________________________________________________________________
bool BBoxIdx::treeHas(const Point<float>& p, const BBoxIdxNd& nd) const {
bool BBoxIdx::treeHas(const Point<double>& p, const BBoxIdxNd& nd) const {
if (!nd.childs.size()) return util::geo::contains(p, nd.box);
for (const auto& child : nd.childs) {
if (util::geo::contains(p, child.box)) return treeHas(p, child);
@ -48,7 +48,7 @@ bool BBoxIdx::treeHas(const Point<float>& p, const BBoxIdxNd& nd) const {
}
// _____________________________________________________________________________
void BBoxIdx::addToTree(const Box<float>& box, BBoxIdxNd* nd, size_t lvl) {
void BBoxIdx::addToTree(const Box<double>& box, BBoxIdxNd* nd, size_t lvl) {
double bestCommonArea = 0;
ssize_t bestChild = -1;

View file

@ -15,9 +15,9 @@ using util::geo::Box;
using util::geo::Point;
struct BBoxIdxNd {
BBoxIdxNd() : box(util::geo::minbox<float>()) {}
explicit BBoxIdxNd(const Box<float>& box) : box(box) {}
Box<float> box;
BBoxIdxNd() : box(util::geo::minbox<double>()) {}
explicit BBoxIdxNd(const Box<double>& box) : box(box) {}
Box<double> box;
std::vector<BBoxIdxNd> childs;
};
@ -26,16 +26,16 @@ struct BBoxIdxNd {
*/
class BBoxIdx {
public:
explicit BBoxIdx(float padding);
explicit BBoxIdx(double padding);
// Add a bounding box to this index
void add(Box<float> box);
void add(Box<double> box);
// Check if a point is contained in this index
bool contains(const Point<float>& box) const;
bool contains(const Point<double>& box) const;
// Return the full total bounding box of this index
util::geo::Box<float> getFullWebMercBox() const;
util::geo::Box<double> getFullWebMercBox() const;
// Return the size of this index
size_t size() const;
@ -46,8 +46,8 @@ class BBoxIdx {
BBoxIdxNd _root;
void addToTree(const Box<float>& box, BBoxIdxNd* nd, size_t lvl);
bool treeHas(const Point<float>& p, const BBoxIdxNd& nd) const;
void addToTree(const Box<double>& box, BBoxIdxNd* nd, size_t lvl);
bool treeHas(const Point<double>& p, const BBoxIdxNd& nd) const;
static const size_t MAX_LVL = 5;
static constexpr double MIN_COM_AREA = 0.0;

View file

@ -136,7 +136,7 @@ void OsmBuilder::read(const std::string& path, const OsmReadOpts& opts,
for (double d : opts.maxSnapDistances) {
for (auto s : orphanStations) {
FPoint geom = *s->pl().getGeom();
DPoint geom = *s->pl().getGeom();
NodePL pl = s->pl();
pl.getSI()->setIsFromOsm(false);
const auto& r = snapStation(g, &pl, &eg, &sng, opts, res, false, d);
@ -374,7 +374,7 @@ void OsmBuilder::readWriteWays(xml::File* i, util::xml::XmlWriter* o,
// _____________________________________________________________________________
NodePL OsmBuilder::plFromGtfs(const Stop* s, const OsmReadOpts& ops) const {
NodePL ret(util::geo::latLngToWebMerc<float>(s->getLat(), s->getLng()),
NodePL ret(util::geo::latLngToWebMerc<double>(s->getLat(), s->getLng()),
StatInfo(ops.statNormzer(s->getName()),
ops.trackNormzer(s->getPlatformCode()), false));
@ -413,7 +413,7 @@ xml::ParserState OsmBuilder::readBBoxNds(xml::File* xml, OsmIdSet* nodes,
double y = util::atof(cur.attrs.find("lat")->second, 7);
double x = util::atof(cur.attrs.find("lon")->second, 7);
if (bbox.contains(Point<float>(x, y))) {
if (bbox.contains(Point<double>(x, y))) {
curId = util::atoul(cur.attrs.find("id")->second);
nodes->add(curId);
}
@ -705,7 +705,7 @@ void OsmBuilder::readNodes(xml::File* xml, Graph* g, const RelLst& rels,
keepAttrs, fl))
.id) {
Node* n = 0;
auto pos = util::geo::latLngToWebMerc<float>(nd.lat, nd.lng);
auto pos = util::geo::latLngToWebMerc<double>(nd.lat, nd.lng);
if (nodes->count(nd.id)) {
n = (*nodes)[nd.id];
n->pl().setGeom(pos);
@ -930,7 +930,7 @@ std::string OsmBuilder::getAttr(const DeepAttrRule& s, osmid id,
// _____________________________________________________________________________
Nullable<StatInfo> OsmBuilder::getStatInfo(Node* node, osmid nid,
const FPoint& pos, const AttrMap& m,
const DPoint& pos, const AttrMap& m,
StAttrGroups* groups,
const RelMap& nodeRels,
const RelLst& rels,
@ -994,7 +994,7 @@ double OsmBuilder::dist(const Node* a, const Node* b) const {
}
// _____________________________________________________________________________
double OsmBuilder::webMercDistFactor(const util::geo::FPoint& a) const {
double OsmBuilder::webMercDistFactor(const util::geo::DPoint& a) const {
// euclidean distance on web mercator is in meters on equator,
// and proportional to cos(lat) in both y directions
@ -1035,7 +1035,7 @@ void OsmBuilder::fixGaps(Graph* g, NodeGrid* ng) const {
otherN = (*nb->getAdjListOut().begin())->getOtherNd(nb);
else
otherN = (*nb->getAdjListIn().begin())->getOtherNd(nb);
FLine l;
DLine l;
l.push_back(*otherN->pl().getGeom());
l.push_back(*n->pl().getGeom());
@ -1057,7 +1057,7 @@ void OsmBuilder::fixGaps(Graph* g, NodeGrid* ng) const {
// _____________________________________________________________________________
EdgeGrid OsmBuilder::buildEdgeIdx(Graph* g, size_t size,
const Box<float>& webMercBox) const {
const Box<double>& webMercBox) const {
EdgeGrid ret(size, size, webMercBox, false);
for (auto* n : *g->getNds()) {
for (auto* e : n->getAdjListOut()) {
@ -1070,7 +1070,7 @@ EdgeGrid OsmBuilder::buildEdgeIdx(Graph* g, size_t size,
// _____________________________________________________________________________
NodeGrid OsmBuilder::buildNodeIdx(Graph* g, size_t size,
const Box<float>& webMercBox,
const Box<double>& webMercBox,
bool which) const {
NodeGrid ret(size, size, webMercBox, false);
for (auto* n : *g->getNds()) {
@ -1084,7 +1084,7 @@ NodeGrid OsmBuilder::buildNodeIdx(Graph* g, size_t size,
// _____________________________________________________________________________
Node* OsmBuilder::depthSearch(const Edge* e, const StatInfo* si,
const util::geo::FPoint& p, double maxD,
const util::geo::DPoint& p, double maxD,
int maxFullTurns, double minAngle,
const SearchFunc& sfunc) const {
// shortcuts
@ -1125,10 +1125,10 @@ Node* OsmBuilder::depthSearch(const Edge* e, const StatInfo* si,
if (cur.fromEdge &&
cur.node->getInDeg() + cur.node->getOutDeg() >
2) { // only intersection angles
const FPoint& toP = *cand->pl().getGeom();
const FPoint& fromP =
const DPoint& toP = *cand->pl().getGeom();
const DPoint& fromP =
*cur.fromEdge->getOtherNd(cur.node)->pl().getGeom();
const FPoint& nodeP = *cur.node->pl().getGeom();
const DPoint& nodeP = *cur.node->pl().getGeom();
if (util::geo::innerProd(nodeP, fromP, toP) < minAngle) fullTurn = 1;
}
@ -1150,24 +1150,24 @@ Node* OsmBuilder::depthSearch(const Edge* e, const StatInfo* si,
// _____________________________________________________________________________
bool OsmBuilder::isBlocked(const Edge* e, const StatInfo* si,
const util::geo::FPoint& p, double maxD,
const util::geo::DPoint& p, double maxD,
int maxFullTurns, double minAngle) const {
return depthSearch(e, si, p, maxD, maxFullTurns, minAngle, BlockSearch());
}
// _____________________________________________________________________________
Node* OsmBuilder::eqStatReach(const Edge* e, const StatInfo* si,
const util::geo::FPoint& p, double maxD,
const util::geo::DPoint& p, double maxD,
int maxFullTurns, double minAngle) const {
return depthSearch(e, si, p, maxD, maxFullTurns, minAngle, EqSearch());
}
// _____________________________________________________________________________
void OsmBuilder::getEdgCands(const FPoint& geom, EdgeCandPQ* ret, EdgeGrid* eg,
void OsmBuilder::getEdgCands(const DPoint& geom, EdgeCandPQ* ret, EdgeGrid* eg,
double d) const {
double distor = webMercDistFactor(geom);
std::set<Edge*> neighs;
Box<float> box = util::geo::pad(util::geo::getBoundingBox(geom), d / distor);
Box<double> box = util::geo::pad(util::geo::getBoundingBox(geom), d / distor);
eg->get(box, &neighs);
for (auto* e : neighs) {
@ -1186,7 +1186,7 @@ std::set<Node*> OsmBuilder::getMatchingNds(const NodePL& s, NodeGrid* ng,
std::set<Node*> ret;
double distor = webMercDistFactor(*s.getGeom());
std::set<Node*> neighs;
Box<float> box =
Box<double> box =
util::geo::pad(util::geo::getBoundingBox(*s.getGeom()), d / distor);
ng->get(box, &neighs);
@ -1204,7 +1204,7 @@ std::set<Node*> OsmBuilder::getMatchingNds(const NodePL& s, NodeGrid* ng,
Node* OsmBuilder::getMatchingNd(const NodePL& s, NodeGrid* ng, double d) const {
double distor = webMercDistFactor(*s.getGeom());
std::set<Node*> neighs;
Box<float> box =
Box<double> box =
util::geo::pad(util::geo::getBoundingBox(*s.getGeom()), d / distor);
ng->get(box, &neighs);
@ -1282,7 +1282,7 @@ std::set<Node*> OsmBuilder::snapStation(Graph* g, NodePL* s, EdgeGrid* eg,
auto ne = g->addEdg(e->getFrom(), n, e->pl());
ne->pl().setLength(webMercDist(n, e->getFrom()));
FLine l;
DLine l;
l.push_back(*e->getFrom()->pl().getGeom());
l.push_back(*n->pl().getGeom());
*ne->pl().getGeom() = l;
@ -1290,7 +1290,7 @@ std::set<Node*> OsmBuilder::snapStation(Graph* g, NodePL* s, EdgeGrid* eg,
auto nf = g->addEdg(n, e->getTo(), e->pl());
nf->pl().setLength(webMercDist(n, e->getTo()));
FLine ll;
DLine ll;
ll.push_back(*n->pl().getGeom());
ll.push_back(*e->getTo()->pl().getGeom());
*nf->pl().getGeom() = ll;

View file

@ -162,7 +162,7 @@ class OsmBuilder {
OsmRel nextRel(xml::File* xml, const OsmFilter& filter,
const AttrKeySet& keepAttrs) const;
Nullable<StatInfo> getStatInfo(Node* node, osmid nid, const FPoint& pos,
Nullable<StatInfo> getStatInfo(Node* node, osmid nid, const DPoint& pos,
const AttrMap& m, StAttrGroups* groups,
const RelMap& nodeRels, const RelLst& rels,
const OsmReadOpts& ops) const;
@ -172,14 +172,14 @@ class OsmBuilder {
void deleteOrphEdgs(Graph* g) const;
double dist(const Node* a, const Node* b) const;
double webMercDist(const Node* a, const Node* b) const;
double webMercDistFactor(const FPoint& a) const;
double webMercDistFactor(const DPoint& a) const;
NodeGrid buildNodeIdx(Graph* g, size_t size,
const util::geo::Box<float>& webMercBox,
const util::geo::Box<double>& webMercBox,
bool which) const;
EdgeGrid buildEdgeIdx(Graph* g, size_t size,
const util::geo::Box<float>& webMercBox) const;
const util::geo::Box<double>& webMercBox) const;
void fixGaps(Graph* g, NodeGrid* ng) const;
void collapseEdges(Graph* g) const;
@ -190,7 +190,7 @@ class OsmBuilder {
uint32_t writeComps(Graph* g) const;
bool edgesSim(const Edge* a, const Edge* b) const;
const EdgePL& mergeEdgePL(Edge* a, Edge* b) const;
void getEdgCands(const FPoint& s, EdgeCandPQ* ret, EdgeGrid* eg,
void getEdgCands(const DPoint& s, EdgeCandPQ* ret, EdgeGrid* eg,
double d) const;
std::set<Node*> getMatchingNds(const NodePL& s, NodeGrid* ng, double d) const;
@ -204,14 +204,14 @@ class OsmBuilder {
// Checks if from the edge e, a station similar to si can be reach with less
// than maxD distance and less or equal to "maxFullTurns" full turns. If
// such a station exists, it is returned. If not, 0 is returned.
Node* eqStatReach(const Edge* e, const StatInfo* si, const FPoint& p,
Node* eqStatReach(const Edge* e, const StatInfo* si, const DPoint& p,
double maxD, int maxFullTurns, double maxAng) const;
Node* depthSearch(const Edge* e, const StatInfo* si,
const util::geo::FPoint& p, double maxD, int maxFullTurns,
const DPoint& p, double maxD, int maxFullTurns,
double minAngle, const SearchFunc& sfunc) const;
bool isBlocked(const Edge* e, const StatInfo* si, const FPoint& p,
bool isBlocked(const Edge* e, const StatInfo* si, const DPoint& p,
double maxD, int maxFullTurns, double minAngle) const;
StatGroup* groupStats(const NodeSet& s) const;