2

Due to tetraphobia in Taiwan we are not allowed to have streets whose names end in "4".

My mission is to reconstruct where 14th Street would have gone.

missing contour

I shall use GDAL/OGR, pretending the streets are contours!

Here's my data,

$ cat e.csv WKT "LINESTRING Z (120.71991 24.17197 12,120.72250 24.17178 12,120.72420 24.17151 12)" "LINESTRING Z (120.72420 24.17225 13,120.72251 24.17254 13,120.71991 24.17265 13)" "LINESTRING Z (120.71991 24.17409 15,120.72420 24.17399 15)" "LINESTRING Z (120.72420 24.17504 16,120.71991 24.17497 16)"

On it, I used gdal_rasterize and then gdal_contour -fl 14, but the results were so embarrassing that I dare not mention the exact command line I used, and I'm starting to get tetraphobia too.

Improve this question
6
  • 2
    Perhaps it is as simple as inserting a line in your input data??? It is anybody's guess where the missing street should be, but you could always plot it at the mid-point between "13" and "15". Who is to say you will be wrong? Frankly, this doesn't look like a question with merit. Commented Feb 1 at 9:57
  • 2
    Whaaat? This is one super fun question! But I see no easy way with just GDAL commandline tools. Are using QGIS or Python code options? Commented Feb 1 at 12:18
  • 2
    Create a polygon between 13 and 15 and generate a skeleton with some tool like postgis.net/docs/ST_StraightSkeleton.html. Doing all that with GDAL command line tools feels complicated. Commented Feb 1 at 15:21
  • @LeighBettenay yes I want the computer to plot it at the midpoint, and maybe 14.25, 14.50... 16.25 supplementary contour lines too someday (with -i .25). Commented Feb 2 at 1:38
  • @bugmenot123 thanks but I don't want to use any GUIs (Qgis) and well, learned too much Perl so too late to use Python. Commented Feb 2 at 1:40

2 Answers 2

Reset to default
1

Here's your missing contour,

$ gdal_rasterize e.csv f.tif -q -3d -ts 256 256 $ gdal_fillnodata f.tif g.tif -q $ gdal_contour g.tif h.csv -q -a name -fl 14 -lco GEOMETRY=AS_WKT $ ogr2ogr k.csv h.csv -simplify .001 -lco GEOMETRY=AS_WKT \ -dim XYZ -zfield name -sql 'SELECT name FROM h' $ cat k.csv $ WKT,name $ "LINESTRING Z (120.7242 24.1731230229664 14,120.71991 24.1733715227135 14)","14"

contour 14

I hope you are happy.

I see. I forgot gdal_fillnodata. How simple-minded of me.

What about the other contours? I hope they are doing OK.

As a matter of fact, no. They kind of got blown around a little. You know, like in "high winds".

But just paste the one new contour I gave you along with your other old contours, and nobody will know the difference!

https://gdal.org/en/latest/api/gdal_alg.html suggests for sparse point data you might try https://gdal.org/en/latest/programs/gdal_grid.html in your spare time instead of gdal_fillnodata.

Got it! Thanks.

Hey, by the way, I don't need that final "14" at back.

But the SELECT needs name else say bye bye to all the 14's.

Improve this answer
0

Here is a an additional totally different approach, that doesn't disrupt the other contours very much at all.

The key for me to be able to use gdal_grid successfully was to get out of longitude and latitude, and into a rectangular system first.

$ ogr2ogr f.csv e.csv -lco GEOMETRY=AS_WKT \ -sql 'SELECT GEOMETRY AS WKT FROM e' \ -s_srs EPSG:4326 -t_srs EPSG:3826 -dialect SQLite $ gdal_grid -a linear f.csv g.tif $ gdal_contour -a Name g.tif h.csv \ -lco GEOMETRY=AS_WKT -i 1 $ ogr2ogr h.kml h.csv -f LIBKML -sql \ 'SELECT Name FROM h' \ -s_srs EPSG:3826 -t_srs EPSG:4326

Street contours

What about Delaunay Triangulation?

Delaunay Triangulation is included in the gdal_grid -a linear step.

Anyway, we now have a general way to draw contours, starting with just points or lines.

Improve this answer

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.