07_a_graph_route_example.cpp 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387
  1. // Boost.Geometry (aka GGL, Generic Geometry Library)
  2. // Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
  3. // Copyright (c) 2008-2012 Bruno Lalande, Paris, France.
  4. // Copyright (c) 2009-2012 Mateusz Loskot, London, UK.
  5. // Use, modification and distribution is subject to the Boost Software License,
  6. // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  7. // http://www.boost.org/LICENSE_1_0.txt)
  8. //
  9. // Example showing Boost.Geometry combined with Boost.Graph, calculating shortest routes
  10. // input: two WKT's, provided in subfolder data
  11. // output: text, + an SVG, displayable in e.g. Firefox)
  12. #include <iostream>
  13. #include <fstream>
  14. #include <iomanip>
  15. #include <limits>
  16. #include <boost/tuple/tuple.hpp>
  17. #include <boost/foreach.hpp>
  18. #include <boost/graph/adjacency_list.hpp>
  19. #include <boost/graph/dijkstra_shortest_paths.hpp>
  20. #include <boost/geometry/geometry.hpp>
  21. #include <boost/geometry/geometries/linestring.hpp>
  22. #include <boost/geometry/io/wkt/read.hpp>
  23. // For output:
  24. #include <boost/geometry/io/svg/svg_mapper.hpp>
  25. // For distance-calculations over the Earth:
  26. //#include <boost/geometry/extensions/gis/geographic/strategies/andoyer.hpp>
  27. // Read an ASCII file containing WKT's, fill a vector of tuples
  28. // The tuples consist of at least <0> a geometry and <1> an identifying string
  29. template <typename Geometry, typename Tuple, typename Box>
  30. void read_wkt(std::string const& filename, std::vector<Tuple>& tuples, Box& box)
  31. {
  32. std::ifstream cpp_file(filename.c_str());
  33. if (cpp_file.is_open())
  34. {
  35. while (! cpp_file.eof() )
  36. {
  37. std::string line;
  38. std::getline(cpp_file, line);
  39. Geometry geometry;
  40. boost::trim(line);
  41. if (! line.empty() && ! boost::starts_with(line, "#"))
  42. {
  43. std::string name;
  44. // Split at ';', if any
  45. std::string::size_type pos = line.find(";");
  46. if (pos != std::string::npos)
  47. {
  48. name = line.substr(pos + 1);
  49. line.erase(pos);
  50. boost::trim(line);
  51. boost::trim(name);
  52. }
  53. Geometry geometry;
  54. boost::geometry::read_wkt(line, geometry);
  55. Tuple tuple(geometry, name);
  56. tuples.push_back(tuple);
  57. boost::geometry::expand(box, boost::geometry::return_envelope<Box>(geometry));
  58. }
  59. }
  60. }
  61. }
  62. // Code to define properties for Boost Graph's
  63. enum vertex_bg_property_t { vertex_bg_property };
  64. enum edge_bg_property_t { edge_bg_property };
  65. namespace boost
  66. {
  67. BOOST_INSTALL_PROPERTY(vertex, bg_property);
  68. BOOST_INSTALL_PROPERTY(edge, bg_property);
  69. }
  70. // To calculate distance, declare and construct a strategy with average earth radius
  71. boost::geometry::strategy::distance::haversine<double> const haversine(6372795.0);
  72. // Define properties for vertex
  73. template <typename Point>
  74. struct bg_vertex_property
  75. {
  76. bg_vertex_property()
  77. {
  78. boost::geometry::assign_zero(location);
  79. }
  80. bg_vertex_property(Point const& loc)
  81. : location(loc)
  82. {
  83. }
  84. Point location;
  85. };
  86. // Define properties for edge
  87. template <typename Linestring>
  88. struct bg_edge_property
  89. {
  90. bg_edge_property(Linestring const& line)
  91. : m_length(boost::geometry::length(line, haversine))
  92. , m_line(line)
  93. {
  94. }
  95. inline operator double() const
  96. {
  97. return m_length;
  98. }
  99. inline Linestring const& line() const
  100. {
  101. return m_line;
  102. }
  103. private :
  104. double m_length;
  105. Linestring m_line;
  106. };
  107. // Utility function to add a vertex to a graph. It might exist already. Then do not insert,
  108. // but return vertex descriptor back. It might not exist. Then add it (and return).
  109. // To efficiently handle this, a std::map is used.
  110. template <typename M, typename K, typename G>
  111. inline typename boost::graph_traits<G>::vertex_descriptor find_or_insert(M& map, K const& key, G& graph)
  112. {
  113. typename M::const_iterator it = map.find(key);
  114. if (it == map.end())
  115. {
  116. // Add a vertex to the graph
  117. typename boost::graph_traits<G>::vertex_descriptor new_vertex
  118. = boost::add_vertex(graph);
  119. // Set the property (= location)
  120. boost::put(boost::get(vertex_bg_property, graph), new_vertex,
  121. bg_vertex_property<typename M::key_type>(key));
  122. // Add to the map, using POINT as key
  123. map[key] = new_vertex;
  124. return new_vertex;
  125. }
  126. return it->second;
  127. }
  128. template
  129. <
  130. typename Graph,
  131. typename RoadTupleVector,
  132. typename CityTupleVector
  133. >
  134. void add_roads_and_connect_cities(Graph& graph,
  135. RoadTupleVector const& roads,
  136. CityTupleVector& cities)
  137. {
  138. typedef typename boost::range_value<RoadTupleVector>::type road_type;
  139. typedef typename boost::tuples::element<0, road_type>::type line_type;
  140. typedef typename boost::geometry::point_type<line_type>::type point_type;
  141. typedef typename boost::graph_traits<Graph>::vertex_descriptor vertex_type;
  142. // Define a map to be used during graph filling
  143. // Maps from point to vertex-id's
  144. typedef std::map<point_type, vertex_type, boost::geometry::less<point_type> > map_type;
  145. map_type map;
  146. // Fill the graph
  147. BOOST_FOREACH(road_type const& road, roads)
  148. {
  149. line_type const& line = road.template get<0>();
  150. // Find or add begin/end point of these line
  151. vertex_type from = find_or_insert(map, line.front(), graph);
  152. vertex_type to = find_or_insert(map, line.back(), graph);
  153. boost::add_edge(from, to, bg_edge_property<line_type>(line), graph);
  154. }
  155. // Find nearest graph vertex for each city, using the map
  156. typedef typename boost::range_value<CityTupleVector>::type city_type;
  157. BOOST_FOREACH(city_type& city, cities)
  158. {
  159. double min_distance = 1e300;
  160. for(typename map_type::const_iterator it = map.begin(); it != map.end(); ++it)
  161. {
  162. double dist = boost::geometry::distance(it->first, city.template get<0>());
  163. if (dist < min_distance)
  164. {
  165. min_distance = dist;
  166. // Set the vertex
  167. city.template get<2>() = it->second;
  168. }
  169. }
  170. }
  171. }
  172. template <typename Graph, typename Route>
  173. inline void add_edge_to_route(Graph const& graph,
  174. typename boost::graph_traits<Graph>::vertex_descriptor vertex1,
  175. typename boost::graph_traits<Graph>::vertex_descriptor vertex2,
  176. Route& route)
  177. {
  178. std::pair
  179. <
  180. typename boost::graph_traits<Graph>::edge_descriptor,
  181. bool
  182. > opt_edge = boost::edge(vertex1, vertex2, graph);
  183. if (opt_edge.second)
  184. {
  185. // Get properties of edge and of vertex
  186. bg_edge_property<Route> const& edge_prop =
  187. boost::get(boost::get(edge_bg_property, graph), opt_edge.first);
  188. bg_vertex_property<typename boost::geometry::point_type<Route>::type> const& vertex_prop =
  189. boost::get(boost::get(vertex_bg_property, graph), vertex2);
  190. // Depending on how edge connects to vertex, copy it forward or backward
  191. if (boost::geometry::equals(edge_prop.line().front(), vertex_prop.location))
  192. {
  193. std::copy(edge_prop.line().begin(), edge_prop.line().end(),
  194. std::back_inserter(route));
  195. }
  196. else
  197. {
  198. std::reverse_copy(edge_prop.line().begin(), edge_prop.line().end(),
  199. std::back_inserter(route));
  200. }
  201. }
  202. }
  203. template <typename Graph, typename Route>
  204. inline void build_route(Graph const& graph,
  205. std::vector<typename boost::graph_traits<Graph>::vertex_descriptor> const& predecessors,
  206. typename boost::graph_traits<Graph>::vertex_descriptor vertex1,
  207. typename boost::graph_traits<Graph>::vertex_descriptor vertex2,
  208. Route& route)
  209. {
  210. typedef typename boost::graph_traits<Graph>::vertex_descriptor vertex_type;
  211. vertex_type pred = predecessors[vertex2];
  212. add_edge_to_route(graph, vertex2, pred, route);
  213. while (pred != vertex1)
  214. {
  215. add_edge_to_route(graph, predecessors[pred], pred, route);
  216. pred = predecessors[pred];
  217. }
  218. }
  219. int main()
  220. {
  221. // Define a point in the Geographic coordinate system (currently Spherical)
  222. // (geographic calculations are in an extension; for sample it makes no difference)
  223. typedef boost::geometry::model::point
  224. <
  225. double, 2, boost::geometry::cs::spherical_equatorial<boost::geometry::degree>
  226. > point_type;
  227. typedef boost::geometry::model::linestring<point_type> line_type;
  228. // Define the graph, lateron containing the road network
  229. typedef boost::adjacency_list
  230. <
  231. boost::vecS, boost::vecS, boost::undirectedS
  232. , boost::property<vertex_bg_property_t, bg_vertex_property<point_type> >
  233. , boost::property<edge_bg_property_t, bg_edge_property<line_type> >
  234. > graph_type;
  235. typedef boost::graph_traits<graph_type>::vertex_descriptor vertex_type;
  236. // Init a bounding box, lateron used to define SVG map
  237. boost::geometry::model::box<point_type> box;
  238. boost::geometry::assign_inverse(box);
  239. // Read the cities
  240. typedef boost::tuple<point_type, std::string, vertex_type> city_type;
  241. std::vector<city_type> cities;
  242. read_wkt<point_type>("data/cities.wkt", cities, box);
  243. // Read the road network
  244. typedef boost::tuple<line_type, std::string> road_type;
  245. std::vector<road_type> roads;
  246. read_wkt<line_type>("data/roads.wkt", roads, box);
  247. graph_type graph;
  248. // Add roads and connect cities
  249. add_roads_and_connect_cities(graph, roads, cities);
  250. double const km = 1000.0;
  251. std::cout << "distances, all in KM" << std::endl
  252. << std::fixed << std::setprecision(0);
  253. // Main functionality: calculate shortest routes from/to all cities
  254. // For the first one, the complete route is stored as a linestring
  255. bool first = true;
  256. line_type route;
  257. int const n = boost::num_vertices(graph);
  258. BOOST_FOREACH(city_type const& city1, cities)
  259. {
  260. std::vector<vertex_type> predecessors(n);
  261. std::vector<double> costs(n);
  262. // Call Dijkstra (without named-parameter to be compatible with all VC)
  263. boost::dijkstra_shortest_paths(graph, city1.get<2>(),
  264. &predecessors[0], &costs[0],
  265. boost::get(edge_bg_property, graph),
  266. boost::get(boost::vertex_index, graph),
  267. std::less<double>(), std::plus<double>(),
  268. (std::numeric_limits<double>::max)(), double(),
  269. boost::dijkstra_visitor<boost::null_visitor>());
  270. BOOST_FOREACH(city_type const& city2, cities)
  271. {
  272. if (! boost::equals(city1.get<1>(), city2.get<1>()))
  273. {
  274. double distance = costs[city2.get<2>()] / km;
  275. double acof = boost::geometry::distance(city1.get<0>(), city2.get<0>(), haversine) / km;
  276. std::cout
  277. << std::setiosflags (std::ios_base::left) << std::setw(15)
  278. << city1.get<1>() << " - "
  279. << std::setiosflags (std::ios_base::left) << std::setw(15)
  280. << city2.get<1>()
  281. << " -> through the air: " << std::setw(4) << acof
  282. << " , over the road: " << std::setw(4) << distance
  283. << std::endl;
  284. if (first)
  285. {
  286. build_route(graph, predecessors,
  287. city1.get<2>(), city2.get<2>(),
  288. route);
  289. first = false;
  290. }
  291. }
  292. }
  293. }
  294. #if defined(HAVE_SVG)
  295. // Create the SVG
  296. std::ofstream stream("routes.svg");
  297. boost::geometry::svg_mapper<point_type> mapper(stream, 600, 600);
  298. // Map roads
  299. BOOST_FOREACH(road_type const& road, roads)
  300. {
  301. mapper.add(road.get<0>());
  302. }
  303. BOOST_FOREACH(road_type const& road, roads)
  304. {
  305. mapper.map(road.get<0>(),
  306. "stroke:rgb(128,128,128);stroke-width:1");
  307. }
  308. mapper.map(route,
  309. "stroke:rgb(0, 255, 0);stroke-width:6;opacity:0.5");
  310. // Map cities
  311. BOOST_FOREACH(city_type const& city, cities)
  312. {
  313. mapper.map(city.get<0>(),
  314. "fill:rgb(255,255,0);stroke:rgb(0,0,0);stroke-width:1");
  315. mapper.text(city.get<0>(), city.get<1>(),
  316. "fill:rgb(0,0,0);font-family:Arial;font-size:10px", 5, 5);
  317. }
  318. #endif
  319. return 0;
  320. }