Category Archives: commandline

One SRTMGL1 GeoTIFF to rule them all

NOTE: In up to date GDAL things might be easier AND I messed up with --config and -co below. Poke me to update it please.

So about half a year ago Lukas Martinelli asked about a global SRTM GeoTIFF. The SRTM elevation data is usually shared in many small tiles which can be ideal for some cases but annoying for others. I like downloading and processing big files so I took that as a challenge. It’s probably some mental thing. I never finished this blog post back then. It’s probably another mental thing. 8) Read on if you are interested in some GNU coreutils hackery as well as GDAL magick.

Here is how I did it. Endless thanks as usual to Even Rouault who fixed GDAL bugs and gave great hints. I learned more about GeoTIFFs and GDAL than I should have needed to know.

Downloading the source files

The source files are available at (warning, visiting this URL will make your browser cry and potentially render your system unresponsive). To download them you must create an NASA Earthdata Login at :\


                         U.S. GOVERNMENT COMPUTER

This US Government computer is for authorized users only.  By accessing this
system you are consenting to complete monitoring with no expectation of privacy.
Unauthorized access or use may subject you to disciplinary action and criminal

OMGOMGOMG. I am sure they would have used <blink> if accessibility guidelines allowed it. As I am not sure what their Terms of Service are, I will not give you a copy’n’paste ready line to download them all. wget or aria2 work well. You should get 14297 files with a total size of 98G.

Inspecting the ZIP files and preparing for GDAL’s vsizip

Each of those ZIP files has just a single file “.hgt” inside. GDAL specifically has support for the “SRTM HGT Format“, so we know it can read those files. We just need to extract them. (If you clicked that link you see that with GDAL 2.2 it will support directly loading the data from the ZIPs, that’s just Even being awesome.)

We don’t want to uncompress all those files just because we want to build a GeoTIFF from them later, do we? Luckily GDAL has a vsizip thingie which allows it to read files inside zip archives. Wicked! For this we need the “deep” paths to the files though, for example

The files inside the archives are always simply the same filename minus the “.SRTMGL1” and the .zip extension. Perfect, now we know all the files from inside the ZIPs! Right? Nope. Unfortunately some (17) of the archives do include the “.SRTMGL1” bit in their files, for example N38E051.SRTMGL1.hgt… >:(

That’s why we need to look into each zip file and determine the filename inside. (Alternatively *you* could extract them all after all, you with your fancy, huge SSD.)

We can simply use `unzip` with the `-t` switch to make it show what’s inside (and as a “free” benefit it will also check the file’s integrity for us). `tee` is used here to save the output to a file. This will take a while as it needs to read through all the files…

$ unzip -t "*" | tee unzip-t.log


testing: N29E000.hgt OK
No errors detected in compressed data of

testing: N27E033.hgt OK
No errors detected in compressed data of

testing: N45E066.hgt OK
No errors detected in compressed data of

14297 archives were successfully processed.
14297 archives were successfully processed.

Yay, no errors!

In the output we see the path to the archive and we see the name of the file inside. With some sed and grep we can easily construct “deep” paths.

First we remove the status lines, the blank lines and the very last full status line with `grep`:

$ cat unzip-t.log | grep -v -e "No errors" -e '^$' -e 'successfully processed'

testing: N12E044.hgt OK
testing: N29W002.hgt OK
testing: N45E116.hgt OK

Then we concatenate every consecutive two lines into one line with `paste` (I LOVE THIS 1 TRICK!):

$ cat unzip-t.log | grep -v -e "No errors" -e '^$' -e 'successfully processed' | paste - -

Archive: testing: N12E044.hgt OK
Archive: testing: N29W002.hgt OK
Archive: testing: N45E116.hgt OK

Finally we strip things away with sed (I don’t care that you would do it differently, this is how I quickly did it with my flair of hammering) and direct the output into a new file called `zips`:

$ cat unzip-t.log | grep -v -e "No errors" -e '^$' -e 'successfully processed' | paste - - | sed -e 's#Archive:[ ]*##' -e 's#\t[ ]*testing: #/#' -e 's#[ ]*OK##' > zips

Yay, we have all the actual inside-zip paths now!

If you are curious, you can try gdalinfo with those now. You need to prepend /vsizip/ to the path to make it read inside zip files.

$ gdalinfo /vsizip/

Driver: SRTMHGT/SRTMHGT File Format
Files: /vsizip/
Size is 3601, 3601
Coordinate System is:
        SPHEROID["WGS 84",6378137,298.257223563,
Origin = (43.999861111111109,13.000138888888889)
Pixel Size = (0.000277777777778,-0.000277777777778)
Corner Coordinates:
Upper Left  (  43.9998611,  13.0001389) ( 43d59'59.50"E, 13d 0' 0.50"N)
Lower Left  (  43.9998611,  11.9998611) ( 43d59'59.50"E, 11d59'59.50"N)
Upper Right (  45.0001389,  13.0001389) ( 45d 0' 0.50"E, 13d 0' 0.50"N)
Lower Right (  45.0001389,  11.9998611) ( 45d 0' 0.50"E, 11d59'59.50"N)
Center      (  44.5000000,  12.5000000) ( 44d30' 0.00"E, 12d30' 0.00"N)
Band 1 Block=3601x1 Type=Int16, ColorInterp=Undefined
  NoData Value=-32768
  Unit Type: m

Hooray! GDAL reads the hgt file from inside its zip!

We need to prepend all the paths with `/vsizip/` so let’s do that:

$ sed 's#^#/vsizip/##' zips > vsizips


Unfortunately those 17 misnamed files from earlier need special treatment… This is what my notes say, not sure what it was supposed to do. If you are recreating this all, please just ask me and I will look at it again. For now, let’s just pretend that this leads to a file called `hgtfiles` in which all the paths are perfect and all the files are perfect.

grep -v 'SRTMGL1.hgt$' vsizips > vsizipswithoutmisnamedfiles
mkdir misnamedfiles
cd misnamedfiles/
nano ../listofmisnamedfiles # insert the paths to those 17 zips here #TODO
while read filename; do unzip "../${filename}"; done < ../listofmisnamedfiles
rename 's/SRTMGL1.//' *.hgt
cd ..
find misnamedfiles/ -type f > listofmisnamedfileshgt
cat listofmisnamedfileshgt vsizipswithoutmisnamedfiles > hgtfiles

As I said above though, go with a recent GDAL and this is all much easier. Even even included a check for the different filenames inside, how can you not like that guy!

Building a Virtual Raster Table

Virtual Raster Tables (VRT) are some kind of files that pretend to be rasters. They are awesome. Here we use a VRT that simply turns all the small rasters we have into one big ass raster.

Ok, ready to create a Virtual Raster Table!

$ time gdalbuildvrt -input_file_list hgtfiles srtmgl1.003.vrt

0...10...20...30...40...50...60...70...80...90...100 - done.

real 0m22.679s
user 0m19.250s
sys 0m3.105s


You could go ahead and use this for your work/leisure now. But remember, it is tens of gigabytes of data so if you do not use it at a 1:1 scale things will not be fun and might fry your cat. You want overviews/pyramids.

Turning the VRT into a GeoTIFF (Optional)

Let’s make a HUGEGEOTIFF because that’s cool! You don’t have to, instead you could build overviews for the .vrt file.

We want it quick to read and small so I used DEFLATE, TILED and the horizontal predictor. I ran this on a weak i7 with 2G of RAM and can’t remember what the worst bottleneck was. Probably CPU.

$ time gdal_translate -co NUM_THREADS=ALL_CPUS --config PREDICTOR 2 -co COMPRESS=DEFLATE -co TILED=YES -co BIGTIFF=YES srtmgl1.003.vrt srtmgl1.003.vrt.tif

Input file size is 1296001, 417601
0...10...20...30...40...50...60...70...80...90...100 - done.

real 231m57.961s
user 503m23.044s
sys 3m14.600s

If you are courageous you can load that file in your GIS now. But again, unless you watch it at a 1:1 scale or something close to that it WILL not be much fun and potentially expose weaknesses in your GIS and fill up your system’s memory and crash and you had no unsaved work open, didn’t you?

Building Overviews

Overviews aka pyramids are usually about 1/3 the size of the full image. If they are not, you probably used different compression settings. This was the step that I expected to be just boring to wait for, but it turned out the most problematic. I tried all available free tools but none worked properly with an input this big. With GDAL we found a workaround after a while.

gdaladdo has problems building multiple overview levels with a file this big… The trick is to built the overviews sequentially, not in one invocation. The best solution at the moment was building overviews for overviews for overviews and so on until you reach a comfortable size. Something like:

gdaladdo -ro srtmgl1.003.tif 2
gdaladdo -ro srtmgl1.003.tif.ovr 2
gdaladdo -ro srtmgl1.003.tif.ovr.ovr 2
gdaladdo -ro srtmgl1.003.tif.ovr.ovr.ovr 2
gdaladdo -ro srtmgl1.003.tif.ovr.ovr.ovr.ovr 2 4 8 16 32 64 128 256 512

I used the following parameters:


Apparently gdaladdo automatically makes them tiled, which is good. I used to find out.

Creating overviews took just about 3 hours with this weird trick. So, in total the processing just takes half a day.


Enjoy! I impose no license/copyright/whatever on this derivative work.


Hillshading? Slope? You do it!

Found a better, faster way? You blog it!

We could add the overviews to the main image with `gdal_translate srtmgl1.003.tif srtmgl1.003.withovr.tif -co COPY_SRC_OVERVIEWS=YES [other options]` says Even.

ffmpeg on raspbian / Raspberry Pi

Since is a bit messy, here is how you can compile ffmpeg with x264 on raspbian. Changes are building in your home directory, getting just a shallow git clone and building with all CPU cores. Also no unnecessary sudo…

Read the comments below!

# In a directory of your choosing (I used ~/ffmpeg):

# build and install x264
git clone --depth 1 git://
cd x264
./configure --host=arm-unknown-linux-gnueabi --enable-static --disable-opencl
make -j 4
sudo make install
# build and make ffmpeg
git clone --depth=1 git://
cd ffmpeg
./configure --arch=armel --target-os=linux --enable-gpl --enable-libx264 --enable-nonfree
make -j4
sudo make install

Read the comments below!

Hopefully someone, somewhere will provide a repository for this kind of stuff some day.

It takes just 25 minutes on a Raspberry Pi 3. Not hours or days like some old internet sources on old Raspis say.

In case you are wondering v4l2 should work with this.

Properly splitting a file at specific intervals with ffmpeg

Ever wanted to split a media file (video, audio, both) into segments of 10 minutes or something like that? The internet is full of terrible hacks and shitty Stack Overflow answers for this. So here is how you easily, properly split a file into same-length segments with ffmpeg.

-f segment -segment_time SECONDS fileprefix%04d.ext

Done. segment_time takes seconds as argument. For example:

ffmpeg -i recording.opus -c:a libvorbis -f segment -segment_time 3600 recording_%04d.ogg


ffmpeg -i huge.wav -c:a copy -f segment -segment_time 600 huge-%04d.flac

Not so hard, is it?

Let’s Encrypt with Lighttpd

Save the page locally and use a clean browser to open the HTML file.
Create your keys and everything.
Cat domain.key and chained.pem into a new pemfile.
Grab the Let’s Encrypt Authority X1 (IdenTrust cross-signed) pemfile from
Put those two files in a safe place on your server.
Don’t forget to set proper permissions on the files. = "/path/to/lets-encrypt-x1-cross-signed.pem" # the Let's Encrypt certificate
ssl.pemfile = "/path/to/pemfile.pem" # your pemfile

Restart lighttpd.

Converting coordinates with cs2cs

In the spirit of Blog Little Things and after reading WEBKID’s You Might Not Need QGIS earlier, here you go:

So you have thousands or millions of coordinates in “the wrong” or “that stupid” coordinate reference system and want to convert them to “the one you need”? Easily done with proj‘s cs2cs tool. You feed it a list of coordinates and tell it how you want them transformed and put out.

Let’s say you have these crazy accurate super coordinates in EPSG:4326 (WGS84) in a file called MyStupidCoordinates.txt (for example from copying columns out of Excel).

9.92510775592794303 53.63833616755876932
9.89715616602172332 53.61639299685715798
9.91541274879985579 53.63589289755918799
9.91922922611725966 53.63415202989637010
9.92072211107566382 53.63179675033637750
9.89998015840083490 53.62284810375095390
9.90427723676020832 53.60740624866674153
9.95012889485460583 53.64563499608360075
9.89590860878436196 53.62979225409544171
9.92944379841924452 53.60061248161031955

and you want them in the much superior and wonderfully metric EPSG:25832 (ETRS89 / UTM zone 32N). You simply use +init=sourceCRS +to +init=targetCRS. if you have no idea what CRS your coordinates are in, enjoy 5 minutes with and you will either know or you might want to take a walk (don’t forget to try them flipped though).

cs2cs +init=epsg:4326 +to +init=epsg:25832 MyStupidCoordinates.txt

and it will print

561164.00 5943681.64 0.00
559346.79 5941216.81 0.00
560526.52 5943401.54 0.00
560781.36 5943211.12 0.00
560883.46 5942950.37 0.00
559524.51 5941937.29 0.00
559830.54 5940223.00 0.00
562807.40 5944515.43 0.00
559245.50 5942706.43 0.00
561505.49 5939488.65 0.00

You probably want more decimal places though, since your input coordinates were accurate down to the attodegree. For this, you can specify the output format with -f

cs2cs +init=epsg:4326 +to +init=epsg:25832 -f "%.17f" MyStupidCoordinates.txt

561164.00100000295788050 5943681.64279999211430550 0.00000000000000000
559346.79200000269338489 5941216.80899999570101500 0.00000000000000000
560526.52400000276975334 5943401.53899999428540468 0.00000000000000000
560781.36300000245682895 5943211.12199999392032623 0.00000000000000000
560883.46400000224821270 5942950.37499999348074198 0.00000000000000000
559524.51200000231619924 5941937.29399999510496855 0.00000000000000000
559830.54200000269338489 5940223.00299999490380287 0.00000000000000000
562807.39800000295508653 5944515.43399999290704727 0.00000000000000000
559245.50000000221189111 5942706.42899999395012856 0.00000000000000000
561505.48800000257324427 5939488.65499999187886715 0.00000000000000000


You could direct these transformed coordinates into a new file called MyCoolCoordinates.txt by adding a redirection of its output:

cs2cs +init=epsg:4326 +to +init=epsg:25832 MyStupidCoordinates.txt > MyCoolCoordinates.txt

You can find out more about cs2cs by reading its manpage:

man cs2cs

PS: You can handle that Z coordinate, right?

Merge GML files into one SQLite or Spatialite file

For example the buildings in Hamburg, Germany:

for file in *.xml
 if [ -f "${layer}".spatialite ]
  ogr2ogr -f "SQLite" -update -append "${layer}".spatialite "${file}" "${layer}" -dsco SPATIALITE=YES
  ogr2ogr -a_srs EPSG:25832 -f "SQLite" "${layer}".spatialite "${file}" "${layer}" -dsco SPATIALITE=YES

Remove the -dsco SPATIALITE=YES and change the output filename for SQLite. QGIS can work with both.

$ bash AX_Gebaeude

Be aware that Spatialite is much more sensitive to geometric problems. You might get things like

ERROR 1: sqlite3_step() failed:
ax_gebaeude.GEOMETRY violates Geometry constraint [geom-type or SRID not allowed] (19)
ERROR 1: ROLLBACK transaction failed: cannot rollback - no transaction is active
ERROR 1: Unable to write feature 1712 from layer AX_Gebaeude.

ERROR 1: Terminating translation prematurely after failed translation of layer AX_Gebaeude (use -skipfailures to skip errors)

but on the other hand, you get spatial indexing which makes queries or high zoom interaction much quicker.

Be aware that if you try to merge files into a Shapefile and fields are getting truncated, those fields will only be filled with data for the first file you merge. On the later files OGR will try to match the input field names to the merged file’s fieldnames, notice the difference and discard them. If you still want to convert to Shapefiles, check out the -fieldTypeToString IntegerList,StringList options.

A GeoTIFF of the German 2011 census population density

tl;dr: GMT is documented for people who use it since the 80s.

update 2: You can use gmtinfo to calculate the extents, using -I- will make it output the -R parameter! xyz2grd $(gmtinfo -I- *.xyz) ...

update: xyz2grd supports GeoTIFF via its GDAL driver (now?)! For example -Gfile.tif=gd:GTiff. See eg So the way below is overly complicated. I do not know how to apply the advanced settings of GDAL though like compression and prediction etc.


The Statistische Ämter des Bundes und der Länder offer a 100 meter grid of Germany’s population density: (110MB). Datensatzbeschreibung_Bevoelkerung_100m_Gitter.xlsx provides additional information.

Let’s turn that dataset into a GeoTIFF so we can use it in our GIS. We will use free and open-source tools from GMT and GDAL. GDAL loves to interpolate values but our data is discrete/regular. We do not want any kind of interpolation. So xyz2grd from GMT is the best choice for turning the xyz data into a “continuous” GIS format (tell me if not).

Getting to know the data

Inside the zip is a 1.3GB file Zensus_Bevoelkerung_100m-Gitter.csv with about 36 million lines.



Datensatzbeschreibung_Bevoelkerung_100m_Gitter.xlsx says the coordinates are in ETRS89-LAEA Europe – EPSG:3035.

First we need to find out the geographic extends of the data, you could use your favourite cli tools for that, I wrote a quick .vrt file and used ogrinfo on that:

$ cat Zensus_Bevoelkerung_100m-Gitter.csv.vrt

 <OGRVRTLayer name="Zensus_Bevoelkerung_100m-Gitter">
  <GeometryField encoding="PointFromColumns" x="x_mp_100m" y="y_mp_100m" />

$ ogrinfo -al Zensus_Bevoelkerung_100m-Gitter.csv.vrt

INFO: Open of `Zensus_Bevoelkerung_100m-Gitter.csv.vrt’ using driver `VRT’ successful.

Layer name: Zensus_Bevoelkerung_100m-Gitter
Geometry: Point
Feature Count: 35785840
Extent: (4031350.000000, 2684050.000000) – (4672550.000000, 3551450.000000)

Creating a netcdf/grd file from the xyz data with xyz2grd

xyz2grd wants xyz, nothing else. The Gitter_ID_100m column is redundant in any case, you can calculate it yourself from the x and y fields if needed. So first let’s convert it to a “x y z” format with awk. The separator is ;

awk 'FS=";" {print $2" "$3" "$4}' Zensus_Bevoelkerung_100m-Gitter.csv >

Now we can write our xyz2grd commandline.

We have our extends:

We know the spacing is 100 units (meters):

There is one header line:

And of course we know that we want a “classic” netcdf4 chunk size, whatever that means. We knew that right away, not after googling helplessly for an hour and eventually finding the hint on some mailing list. Not knowing this might have lead to QGIS only seeing NaN values for z, R’s ncdf/raster saying “Error in substr(w, 1, 3) : invalid multibyte string at ‘<89>HDF” and GDAL “0ERROR 1: nBlockYSize = 130, only 1 supported when reading bottom-up dataset”.

The resulting commandline:
xyz2grd -Vl -R4031350/4672550/2684050/3551450 -I100 -h1 --IO_NC4_CHUNK_SIZE=c -GZensus_Bevoelkerung_100m-Gitter.cdf

You can inspect the file with grdinfo and gdalinfo now if you want.

Creating a GeoTIFF from the netcdf/grd file

Let’s turn it into a GeoTIFF with gdal_translate. We will need a bunch of commandline parameters.

The spatial reference system is EPSG:3035, so:
-a_srs EPSG:3035

Values of -1 mean “no data”:
-a_nodata -1

TIFF is uncompressed by default, we want good lossless compression:

The resulting commandline:
gdal_translate -co COMPRESS=DEFLATE -a_srs EPSG:3035 -a_nodata -1 Zensus_Bevoelkerung_100m-Gitter.cdf Zensus_Bevoelkerung_100m-Gitter.tif

The resulting file is about 8 Megabytes and should work in any reasonable GIS. Have fun!

TODO: What license is this now?


“Commandline File Hosts”

These are the “terminal to URL” file hosts I know. Please tell me about others. The speed tests are not scientific at all, my servers might be uncapable of better speeds because of saturation or bad routing.

curl --upload-file "01_Name_Game_(Intro).mp3"

Maximum Filesize: 5 Gigabytes
Expiration after: 2 weeks
Serves error 500 for non-existing files. Does not serve direct downloads anymore. Bleh.
Download URLs are modeled after your original filename (“special” characters are replaced): upload: ~30-110 Megabit/s download: ~240 Megabit/s

curl -F "file=@01_Name_Game_(Intro).mp3"

Maximum Filesize: 20 Gigabytes
Expiration after: 4 hours
Requires you to visit their website to get a key (URL) assigned.
Download URLs are random and but it serves original filename with content-disposition: upload: ~60-75 Megabit/s download: 100 Megabit/s

curl --upload-file 01_Name_Game_\(Intro\).mp3

Maximum Filesize: 50 Megabytes
Expiration after: 6 months
They might increase filesize and expiration if there is enough interest, as of now you can not “sign-up” for those.
Download URLs are unreadable, filenames are lost: upload: ~8 Megabit/s (until it errored :) )
Requires sign-up.