Shallow Thoughts : tags : qgis

Making OsmAnd Overlays with QGIS (2024 Edition)

Several years ago I wrote about Making a Land Ownership Overlay with QGIS and Making Overlay Maps for OsmAnd. I've been using that land use overlay for years. But recently I needed to make several more overlays: land ownership for Utah for a hiking trip, one for the eclipse, and I wanted to refresh my New Mexico land ownership overlay since it was several years out of date. It turns out some things have changed, so here's an update, starting from the point where your intended overlay is loaded as a layer in QGIS.

Read more ...

A Smaller, Lighter Dataset: Clip Layers in QGIS

The dataset I used for mapping fire perimeters is huge: not surprising if it's all historic fires for the US. Classifying it in QGIS gave a warning, and operations were very slow. Here's how to clip a big dataset in QGIS to restrict it to a smaller geographic area.

Read more ...

Making a Land Ownership overlay: Categorized Styles in QGIS

Now that I know how to make a map overlay for OsmAnd, I wanted a land ownership overlay. When we're hiking, we often wonder whether we're on Forest Service, BLM, or NPS land, or private land, or Indian land. It's not easy to tell.

Finding Land Ownership Data

The first trick was finding the data. The New Mexico State Land Office has an interactive New Mexico Land Status map, but that's no help when walking around, and their downloadable GIS files only cover the lands administered by the state land office, which mostly doesn't include any areas where we hike. They do have some detailed PDF maps of New Mexico Lands if you have a printer capable of printing enormous pages, which most of us don't.

In theory I could download their 11" x 17" Land Status PDF, convert it to a raster file, and georeference it as I described in the earlier article; but since they obviously have the GIS data (used for the interactive map) I'd much rather download the data and save myself all that extra work.

Eventually I found New Mexico ownership data at UNM's RGIS page, which has an excellent collection of GIS data available for download. Click on Boundaries, then download Surface Land Ownership. It's available in a variety of formats; I chose the geojson format because I find it the most readable and the easiest to parse with Python, though ESRI shapefiles arguably might have been easier in QGIS.

Colorizing Polygons in QGIS

You can run qgis on a geojson file directly. When it loads it shows the boundaries, and you can use the Info tool to click on a polygon and see its metadata -- ownership might be BLM, DOE, FS, I, or whatever. But they're all the same color, so it's hard to get a sense of land ownership just clicking around.

To colorize the polygons differently, right-click on the layer name and choose Properties. For Style, choose Categorized. For Column, pick the attribute you want to use to choose colors: for this dataset, it's "own", for ownership.

Color ramp is initially set to random. Click Classify to generate an initial color ramp, then click Apply to see what it looks like on the map.

Then you can customize the colors by doubleclicking on specific color swatches. For instance, by unstated convention most maps show Forest Service land as green, BLM and Indian land as various shades of brown. Click Apply as you change colors, until you're happy with the result.

Exporting to GeoTIFF

Update: In 2024, some of these steps have changed; see Making OsmAnd Overlays with QGIS (2024 Edition).

You can export the colored layer to GeoTIFF using QGIS' confusing and poorly documented Print Composer. Create one with: Project > New Print Composer, which will open with a blank white canvas.

Zoom and pan in the QGIS window so the full extent of the image you want to export is visible. Then, in the Print Composer, Layout > Add Map. Click and drag in the blank canvas, going from one corner to the opposite corner, and some portion of the map should appear.

There doesn't seem to be any way to Print Composer to import your whole map automatically, or for you to control what portion of the map from the QGIS window will show up in the Print Composer when you drag. If you guess wrong and don't get all of your map, hit Delete, switch to the QGIS window and drag and/or zoom your map a little, then switch back to Print Composer and try adding it again.

You can also make adjustments by changing the Extents in the Item Properties tab, and clicking the Set to map canvas extent button in that tab will enlarge your extents to cover approximately what's currently showing in the QGIS window.

It's a fiddly process and there's not much control, but when you decide it's close enough, Composer > Export as Image... and choose TIFF format. (Print Composer offers both TIFF and TIF; I don't know if there's a difference. I only tried TIFF with two effs.) That should write a GeoTIFF format; to verify that, go to a terminal and run gdalinfo on the saved TIFF file and make sure it says it's GeoTIFF.

Load into OsmAnd

Finally, load the image into OsmAnd's tiles folder as discussed in the previous article, then bring up the Configure map menu and enable the overlay.

I found that the black lines dividing the various pieces of land are a bit thicker than I'd like. You can't get that super accurate "I'm standing with one foot in USFS land and the other foot in BLM land" feeling because of the thick black DMZ dividing them. But that's probably just as well: I suspect the data doesn't have pinpoint accuracy either. I'm sure there's a way to reduce the thickness of the black line or eliminate it entirely, but for now, I'm happy with what I have.

Update: Here's another, easier, way to show land use on OsmAnd using overlay tiles from the BLM (in the US): Adding BLM Land Use Maps to Osmand on Android. It isn't as general (you can only show something you can get from an online tiled source) and it updates in real-time, meaning it might use cellphone data rather than working entirely offline, but it's still a great option to know about.

Making Overlay Maps for OsmAnd on Linux

For many years I've wished I could take a raster map image, like a geology map, an old historical map, or a trail map, and overlay it onto the map shown in OsmAnd so I can use it on my phone while walking around. I've tried many times, but there are so many steps and I never found a method that worked.

Last week, the ever helpful Bart Eisenberg posted to the OsmAnd list a video he'd made: Displaying web-based maps with MAPC2MAPC: OsmAnd Maps & Navigation. Bart makes great videos ... but in this case, MAPC2MAPC turned out to be a Windows program so it's no help to a Linux user. Darn!

But seeing his steps laid out inspired me to try again, and gave me some useful terms for web searching. And this time I finally succeeded. I was also helped by a post to the OsmAnd list by A Thompson, How to get aerial image into offline use?, though I needed to change a few of the steps. (Note: click on any of the screenshots here to see a larger version.)

Georeference the Image Using QGIS

Update in Feb 2024: Several things have changed in QGIS georeferencing (the version I'm using now is 3.28.15-Firenze), so note the updated sections below.

The first step is to georeference the image: turn the plain raster image into a GeoTiff that has references showing where on Earth its corners are. It turns out there's an open source program that can do that, QGIS. Although it's poorly documented, it's fairly easy once you figure out the trick.

I started with the tutorial Georeferencing Basics, but it omits one important point, which I finally found in BBRHUFT's How to Georeference a map in QGIS. Step 11 is the key: the Coordinate Reference System (CRS) must be the same in the georeferencing window as it is in the main QGIS window. That sounds like a no-brainer, but in practice, the lists of possible CRSes shown in the two windows don't overlap, so unless you follow BBRHUFT's advice and type 3857 into the filter box in both windows, you'll likely end up with CRSes that don't match. It'll look like it's working, but the resulting GeoTiff will have coordinates nowhere near where they should be

Instead, follow BBRHUFT's advice and type 3857 into the filter box in both windows. The "WGS 84 / Pseudo Mercator" CRS will show up and you can use it in both places. Then the GeoTiff will come out in the right place.

If you're starting from a PDF, you may need to convert it to a raster format like PNG or JPG first. GIMP can do that.

So, the full QGIS steps are:

• Start qgis, and bring up Project > Project Properties.
• Click on CRS, type 3857 in the filter box, and make sure the CRS is set to WGS 84 / Pseudo Mercator.
• (optional) While you're in Project Properties, you may want to click on General and set Display Coordinates Using: Decimal degrees. By default qgis always uses "meters" (measured from where?) which seems like a bizarre and un-useful default.
• In Plugins > Manage and Install Plugins..., install two plugins: Quick Map Services and Georeferencer GDAL.
  Update: both OSM and the Georeferencer are now included by default, so you no longer have to install these plugins.
• Back in the main QGIS window, go to Web > QuickMapServices (that's one of the plugins you just installed) and pick a background map, probably something from OpenStreetMap. Once you find one that loads, navigate to the area you're interested in.
  Update: Sometimes, including an OpenStreetMap is as easy as clicking on the OpenStreetMap line in the upper left pane of the QGIS window. Other times there is no such line, and you have to use Web > QuickMapServices. I don't know how QGIS decides.
  
  
  (All the OpenStreetMap tile sets print "API Key Required" all over everything. You can get an API key from Thunderforest -- I use one for PyTopo -- but it won't help you here, because nobody seems to know how to get QGIS to use an API key, though I found dozens of threads with people asking. If you find a way, let me know.)
  
• Start the georeferencer: Raster > Georeferencer > Georeferencer.... It may prompt you for a CRS; if so, choose WGS 84 / Pseudo Mercator if it's available, and if it's not, type 3857 in the filter box to find it.
• In the Georeferencer window, File > Open Raster...
  Update: If it isn't showing the files you want to open, you may need to adjust the Files of type menu to make sure it includes your file type. When I ran it just now, I had two JPG files to georeference. The first time, it showed those files automatically. When I went to open the second file, it showed nothing about I had to change Files of type to JPEG JFIF before I could open the file.
  
  
• Select matching points. This part is time consuming but fun.
  
  Zoom in to your raster map in the Georeferencer window and pick a point that's easy to identify: ideally the intersection of two roads, or a place where a road crosses a river, or similar identifying features. Make sure you can identify the same point in both the Georeferencer window and the main QGIS window (see the yellow arrows in the diagram).
  
• Click on the point in the Georeferencer window, then click "from map canvas".
• In the main QGIS window, click on the corresponding point. You can zoom with the mousewheel and pan by middle mouse dragging. Repeat for a total of four to eight points.
  
  You'll see yellow boxes pop up in both windows showing the matching points. If you don't, there's something wrong; probably your CRSes don't match.
  
• In the Georeferencer window, call up Settings > Transformation Settings. Choose a Transformation type (most people seem to like Thin Plate Spline), a Resampling method (Lanczos is good), a target SRS (make sure it's 3857 WGS 84 / Pseudo Mercator to match the main window). Set Output raster to your desired filename, ending in .tiff. Be sure to check the "Load in QGIS when done" box at the bottom. Click OK.
  
• File > Start Georeferencing. It's quite fast, and will create your .tiff file when it's finished, and should load the TIFF as a new layer in the main QGIS window. Check it and make sure it seems to be lined up properly with geographic features. If you need to you can right-click on the layer and adjust its transparency.
• In the Georeferencer window, File > Save GCP Points as... in case you want to go back and do this again. Now you can close the Georeferencer window.
• Project > Save to save the QGIS project, so if you want to try again you don't have to repeat all the steps.
  
  Now you can close qgis. Whew!

Convert the GeoTiff to Map Tiles

The ultimate goal is to convert to OsmAnd's sqlite format, but there's no way to get there directly. First you have to convert it to map tiles in a format called mbtiles.

QGIS has a plug-in called QTiles but it didn't work for me: it briefly displayed a progress bar which then disappeared without creating any files. Fortunately, you can do the conversion much more easily with gdal_translate, which at least on Debian is part of the gdal-bin package.

gdal_translate filename.tiff filename.mbtiles

That will create tiles for a limited range of zoom levels (maybe only one zoom level). gdalinfo will tell you the zoom levels in the file. If you want to be able to zoom out and still see your overlay, you might want to add wider zoom levels, which you can do like this:

gdaladdo -r nearest filename.mbtiles 2 4 8 16

Incidentally, gdal can also create a directory of tiles suitable for a web slippy map, though you don't need that for OsmAnd. For that, use gdal2tiles, which on Debian is part of the python-gdal package:

mkdir tiles
gdal2tiles filename.tiff tiles

Not only does it create tiles, it also includes multiple HTML files you can use to display those tiles using the Leaflet, OpenLayers or Google Maps JavaScript libraries. Very cool!

Create the OsmAnd sqlite file

Tarwirdur has written a nice simple Python script to translate from mbtiles to OsmAnd sqlite: mbtiles2osmand.py. Download it then run

mbtiles2osmand.py filename.mbtiles filename.sqlitedb

So easy to use! Most of the other references I saw said to use Mobile Atlas Creator (MOBAC) and that looked a lot more complicated.

Incidentally, Bart's video says MAPC2MAPC calls the format "Locus/Rmaps/Galileo/OSMAND (sqlite)", which might be useful to know for web search purposes.

Install in OsmAnd

Once you have the .sqlitedb file, copy it to OsmAnd's tiles folder in whatever way you prefer. For me, that's adb push file.sqlitedb $androidSD/Android/data/net.osmand.plus/files/tiles where $androidSD is the /storage/whatever location of my device's SD card.

Then start OsmAnd and tap on the icon in the upper left for your current mode (car, bike, walking etc.) to bring up the Configure map menu. Scroll down to Overlay or Underlay map, enable one of those two and you should be able to choose your newly installed map.

You can adjust the overlay's transparency with a slider that's visible at the bottom of the map (the blue slider just above the distance scale), so you can see your overlay and the main map at the same time.

The overlay disappears if you zoom out too far, and I haven't yet figured out what controls that; I'm still working on those details.

Sure, this process is a lot of work. But the result is worth it. Check out the geologic layers we walked through on a portion of a recent hike in Rendija Canyon (our hike is the purple path).