OSM and speed bump map in navigators

    imageSince starting from the beginning is not interesting, I will start from the end.
    We did it. We got speed bumps from OpenStreetMap data, crossed them with the terrible commercial program Navitel, made a web viewer of these speed bumps, and the interface for adding them for beginners at http://latlon.org/tc/ . And they even wrote a small press release , a link to which can be sent to friends and acquaintances, motorists.
    But something special can be said for the Habr: how it all works inside, and how it was done.

    The main rule of OpenStreetMap is Have fun. Everything described was done not for the sake of money or some other goodies, and not even for its own convenience.

    On one of the Minsk mini OSMovoks (mini, because nowhere exceptThey don’t announce the IRC channel ) I was awarded the task: “Do you want to do a useful job? get the lodgers from the osmosis into the navitel. " Then everything was somehow forgotten, and after a couple of months I still remembered, thanks to the heroic work of comrade stas_symba from the forum w3bsit3-dns.com. For several months now, he alone had been collecting and uploading updates on speed bumps in Belarus. “But after all, such data is easier to collect immediately in OSM,” I thought, and I sat down for export.

    Format


    Since I don’t have a navigator, Navitel or a car, I had to do everything blindly, reading forums and thinking with my brains. To start and understand the scale of the disaster - in brief about what we managed to get about the format of the data from the forums.

    The cops.txt file (and / or speedcam.txt) is a regular csv file, the first line contains the header, and all subsequent ones contain data. Title:

    idx,x,y,type,speed,dirtype,direction
    where:
    • idx - unique object number;
    • x, y - coordinates, in the projection longitude-latitude EPSG: 4326;
    • speed - allowed speed for the object;
    • dirtype - type of direction of action - in 1 or 2 directions;
    • direction - direction of the action;


    Tool selection

    Tools for solving a problem are usually chosen from what is at hand. At hand were a lot of scripts for threading OSM XML , and the PostGIS database . Despite the fact that it’s more usual to work with stream scripts, the last field of the speedcam file hinted that you would have to work with the geometry of objects, which PostGIS is famous for.

    To start, on the home server using osm2pgsql, the usual base for Mapnik was imported . The first request was written quickly and simply:
    select * from planet_osm_point where traffic_calming is not NULL;
    In response to this, many fields were returned to me, including geometry in the WKB format , encoded in hex. Not good at all.

    I had to get into PostGIS mahual, in which the functions ST_X and ST_Y were found. It seemed that it was necessary. Rewriting: Suddenly? Expectedly, but we need to make an explanation. In web mapping, two projections are widely used: for transmission and display to the user - latitude-longitude, it is also EPSG: 4326 , and “Google”, it is also “Mercator on a sphere”, it is EPSG: 3857, it is also EPSG: 900913, she is EPSG: 3785 (why are there so many codes? a long history of disputes between large corporations, hobbyists and registrars, worthy of a separate post). The projection is good in that the transition from it to 4326 in any language takes two lines of mathematics from the strength . It is it, metric, and not in angular coordinates, that is used when displaying maps on the screen. And to facilitate calculations, it is in it that the Mapnik database is stored. Well okay. Re-projecting.
    gis=> select st_x(way), st_y(way) from planet_osm_point where traffic_calming is not NULL limit 1;
    st_x | st_y
    ------------------+------------------
    3085590.21426068 | 7159526.18035388
    (1 запись)







    gis=> select ST_X(transform(p.way,94326)) as X, ST_Y(transform(p.way,94326)) as Y from planet_osm_point p where traffic_calming is not NULL limit 1;
    x | y
    ------------+------------
    27.7183285 | 53.9438334
    (1 запись)


    What's next? Speed.

    Nobody in OSM puts speed on speed bumps. We must figure out where to get it. We ask the first person that came across (it turned out to be wildMan ), we get the answer - according to the SDA, the allowed speed is 60 in the city, 90 outside the city. And there may also be a speed limit on the road itself, where it will be entered in the corresponding tag - maxspeed. Well, at the same time, you can get a sign of one-sidedness from the oneway tag - in vain, perhaps, in the format, the field is reserved for it?

    Combine with the table of polygons - it stores the administrative boundaries, and the table of lines - it stores the lines of roads.

    select ST_X(transform(p.way,94326)) as X,
    ST_Y(transform(p.way,94326)) as Y, '102' as TYPE,
    case when l.maxspeed is not NULL then l.maxspeed else case when t.admin_level = '8' then '60' else '90' end end as SPEED,
    case when l.oneway = 'yes' then '1' else '2' end as DirType

    from (planet_osm_point p
    join planet_osm_line l on (l.highway is not NULL and ST_DWithin(p.way,l.way,.1)))
    LEFT OUTER JOIN planet_osm_polygon t on (t.admin_level = '8' and ST_Within(p.way, t.way))
    where p.traffic_calming is not NULL;


    Using the standard SQL join, we compare all three tables - for each speed bump we find in a couple a road line no further than ten centimeters from it, and the polygon of the city surrounding it, if possible. Next, we try to fill in the speed field using the maxspeed tag of the line, and if we did not fill it out, we look to see if the dot fell into the administrative border of the city: if yes, put 60, if not - 90.

    idx - serial number - it was a temptation to take from osm_id . But then it became clear that the policeman, who is at the junction of the two lines, will meet in both lines on both sides, and his id will turn out to be the same in different directions, in connection with which the best option seemed to be to number them all in order.

    Fear and lack of documentation


    Since Navitel is a commercial product, and it does not have either the source code or the normal demo version, the last figure - direction was left without comment. No hints about how it is considered were found, except that this is clearly a figure in degrees.
    “Well, since all the cards are in the Mercator projection, and it preserves angles and directions , then they are counted in it!”

    To calculate the azimuth, PostGIS requires two points. Obviously, they should lie on the same road on which the speed bump is installed, at some distance on both sides of it.
    The ST_Line_Locate_Point method allowed to find the length along the line from its beginning to the specified point, and ST_Line_Interpolate_Point - on the contrary, to find the point along the length along the line.

    But psql draws a pseudographics tablet in response to my requests when I need a CSV! Is it possible to write a wrapper?
    It turned out that everything is much simpler. The COPY function allows you to deliver the result of the query immediately to the CSV - just what you need. Here is such a little cat joy.

    I drive out the first result, I find the owner of Navitel, I force him to check ... and oops. Directions are not parallel to the roads. The thought “really ...” arises, and the paws themselves write the reprojection in 4326 when calculating the azimuth.

    Yes, the valiant Navitel developers consider the direction in a projection that does not preserve direction. Honor to them, praise and honor.

    Quasi-tot


    As a result, the request was just such a monster, extracting from OpenStreetMap data what they initially did not seem to contain - speed bumps in Navitel format. You can write a request that gets speed cameras from them as homework. :) But Kotyara would not be Kotyara if his cat ears did not stick out from everywhere . Therefore, as a warm-up, in addition to a five-minute unloading of updated speed bumps and cameras for the CIS , a rendering of all this disgrace into tiles was organized so that those who do not have a navigator with Navitel could admire it. Closer to the weekend comrade andrewsh
    CREATE TEMP SEQUENCE idx ;
    COPY (

    select nextval('idx') as idx,
    ST_X(transform(p.way,94326)) as X,
    ST_Y(transform(p.way,94326)) as Y, '102' as TYPE,
    case when l.maxspeed is not NULL then l.maxspeed else case when t.admin_level = '8' then '60' else '90' end end as SPEED,
    case when l.oneway = 'yes' then '1' else '2' end as DirType,
    floor(ST_Azimuth(
    transform(ST_Line_Interpolate_Point(l.way,
    case when (ST_Line_Locate_Point(l.way,p.way)-0.01 < 0) then 0 else ST_Line_Locate_Point(l.way,p.way)-0.0000001 end),94326),
    transform(ST_Line_Interpolate_Point(l.way,
    case when (ST_Line_Locate_Point(l.way,p.way)+0.01 > 1) then 1 else ST_Line_Locate_Point(l.way,p.way)+0.0000001 end)
    ,94326))/(2*pi())*360) as Direction
    from (planet_osm_point p
    join planet_osm_line l on (l.highway is not NULL and ST_DWithin(p.way,l.way,.1)))
    LEFT OUTER JOIN planet_osm_polygon t on (t.admin_level = '8' and ST_Within(p.way, t.way))
    where p.traffic_calming is not NULL
    ) TO STDOUT WITH CSV HEADER;




    (I think he also has something to tell if you ask him :) I wrote a convenient interface for reporting new speed bumps without registering - latlon.org/tc , and it was decided to put all the services on the public . I hope you will like it.

    Have fun! ;)

    Also popular now: