Friday, May 20, 2016

In-Browser filtering of raster pixels in Leaflet, Part 3

In my last posting, I mentioned described the process of creating a three-band TIFF and then slicing it into three-band PNG tiles for display in Leaflet. The tiles are not visually impressive, being almost wholly black, but they have the proper R, G, and B codes corresponding to the location's PHZ, PCP, and GEZ. (if you don't know what those are, go back and read that post from a few days ago)

So now the last step: transforming these tiles based on a set of arbitrary RGB codes determined by the user. Let's go through the working source code for the L.TileLayer.PixelFilter and I'll describe how it works.

First off, we need two RGBA codes so we can color the pixels as hit or miss, and we need the list of [r,g,b] trios that would be considered a match. The initialize method accepts these as constructor params, and the various set functions allow this to be set later as well.

(The decision of what list of RGB codes should be passed into setPixelCodes() is a matter of business logic: selectors, clicks, lists... outside the scope of this discussion. See the demo code for a trivial example.)

Next, we need to intercept the map tile when it's in the browser. Thus the
tileload event handler set up in the constructor, which runs the tile through
applyFiltersToTile() to do the real work.

applyFiltersToTile() is the interesting part, which only took me a couple of hours sitting down with some HTML5 Canvas tutorials. Let's dissect it one piece at a time:
  • "copy the image data onto a canvas" The first paragraph creates a Canvas of the same width and height as the image, and copies the data into it. The IMG object itself won't let us access raw pixel data, but once it's in a Canvas we can.
  • "target imagedata" We then use createImageData() to construct a new and empty byte sequence, which can later be assigned back into a Canvas using putImageData() This is effectively a straight list of bytes, and as we write to it, we need to always write 4 bytes at a time for every pixel: R, G, B, and A.
  • "generate the list of integers" Before we start looping over pixels, a performance optimization. Arrays have the indexOf() method so we can determine whether a given value is in an array, but it only works on primitive values and not three-element arrays. The lodash library has findIndex() which would work, but it means a new dependency... and also the performance was not so great (256 x 256 = 65536 pixels per tile). So I cheat, and translate the list of desired
    pixelCodes into a simple list of integers so that indexOf() can work after all.
  • "iterate over the pixels" Now we start looping over the pixels in our Canvas and doing the actual substitution. Again, we step over the bytes in fours (R, G, B, A) and we would assign them into the output imagedata in fours as well. Each pixel is translated into an integer, and if that integer appears in the pixelcodes list of integers, it's a hit. Since indexOf() is a native function it's pretty swift, and since our pixelcodes list tends to be very short (10-50 values) the usually-avoided loop-in-loop is actually quite fast.
  • "push a R, a G, and a B" Whatever the results, we push a new R, G, B, A set of 4 bytes onto the output imagedata.
  •  "write the image back" And here we go: write the imagedata sequence back into the Canvas to replace the old data, then reassign the img.src to this Canvas's base64 representation of the newly-created imagedata. Now the visible tile isn't based on a HTTP file at all, but from a base64-encoded image in memory.
The only gotcha I found, was that upon calling setPixelCodes() the tiles would flash to their proper color and then instantly all turn into the no-match color. It worked... then would un-work itself?

The tileload event wraps the IMG element's load event. This means that when I assigned the img.src to replace the tile's visible representation... this was itself a tileload event! It would call applyFiltersToTile() again, and this time the RGB codes didn't match anything on my filters so all pixels were no-match. Worse, the assignment of img.src was again a tileload event, so it was in an infinite loop of processing the tile.

Thus the already_pixel_swapped flag. After processing the IMG element, this flag is set and
applyFiltersToTile() will skip out on subsequent runs on that IMG element. If we need to change pixel codes via setPixelCodes() that calls layer.redraw() which empties out all these old IMG elements anyway, replacing them with fresh new ones that do not have already_pixel_swapped set.

So yeah, it was an educational day. Not only did we achieve what the client asked for (the more complex version, not the simplified case I presented last week) and as a neat reusable package, but the performance is near instant.


Tuesday, May 17, 2016

In-Browser filtering of raster pixels in Leaflet, Part 2

A few days back I posted a link to L.TileLayer.PixelFilter This is a Leaflet TileLayer extension which rewrites the tiles after they have loaded, comparing each pixel against a set of pixel-codes and replacing the pixel with either a "matched" or "not matched" color. It's pretty useful for generating dynamic masks and highlights, if your back-end data is a raster and not vector data.

I had said that the client's request was to display plant hardiness zones, and that I was using the Unique Values colors as the RGB codes to match against zone codes. I lied. The reality was much more complicated than that, but would have been distracting from the end and the result. So here's the longer story.

The Rasters and the Requirement


The client has 3 rasters: Plant Hardiness Zones (PHZ), Precipitation classifications (PCP), Global Ecoregions (GEZ). Each of these is a worldwide dataset, and each has about 20 discrete integer values numbered 10 to 30ish.

The user selects multiple combination of these three factors, e.g. these two:
  • Zone 7 (PHZ=7)
  • Rainfall 20-30 inches per year (PCP=3)
  • Temperate oceanic forest (GEZ=31)

  • Zone 7 (PHZ=7)
  • Rainfall 20-30 inches per year (PCP=3)
  • Temperate continental forest (GEZ=32)
The user would then see any areas which match any of these "criteria trios" highlighted on the map. The idea is that a plant that's invasive in these areas, would also be invasive in the other areas highlighted on the map.

Fortunately, the rasters are not intended to be visible and do not need to be visually pleasing. They need to be either colored (the pixel matches any of the trios) or else transparent (not a match).

Raster To Tiles To Leaflet


Three variables and a need for one raster... sounds like we could have a three-band raster! I used ArcMap's Composite Bands tool to merge the three rasters into a three-band raster. Voila, one TIFF with three-part pixels such as (7, 3, 31) and (7, 3, 32) I just have to keep straight that R is PHZ, G is PCP, and B is GEZ.

Second, I needed to slice this up into map tiles. But it's very important that:
  • The tiles be in a format that preserves R, G, and B exactly. Something GIF or JPEG would be entirely unsuitable since GIF picks a 256-color palette, and JPEG is lossy and fudges together colors. TIFF on the other hand is not viewable in browsers as a map tile. But PNG is just perfect: RGB by default, and viewable in browsers.
  • The tile-slicing process (I picked gdal2tiles.py) must also preserve exact RGB values. The default average resampler uses interpolation, so a pixel in between a PCP=3 and PCP=2 would get PCP=2.5 and that's no good! Use the nearest neighbor resampling algorithm, which guarantees to use an existing pixel value.
Slicing up the TIFF into web-read PNGs was pretty straightforward:
gdal2tiles.py -z 0-7 -r near -t_srs epsg:3857 threeband.tif tiles

Point Leaflet at it, and I get some nearly-black tiles loading and displaying just as I expected. Well, after I remembered that gdal2tiles.py generates tiles in TMS numbering scheme, so I had to set tms:true in my L.TileLayer constructor.

I downloaded a few of the tiles and opened them up in GIMP and sure enough, that shade of black is actually (7,3,31) The tiles are not meant to be visually attractive... but those preserved RGB codes are beautiful to me.


Client-Side Pixel Detection and Swapping


That's the topic of my next posting: now that I have RGB tiles, how can I intercept them, compare them against a set of RGB codes, and appropriately color/clear pixels...? Turns out it was easier and more elegant than I had imagined.

Saturday, May 14, 2016

In-Browser filtering of raster pixels in Leaflet, Part 1

Last week I was asked to do something I haven't done before:

A client has a raster of values, this being World Plant Hardiness Zones. Each pixel is a code number, which ranges from 1 to 13. If you're into gardening, you've seen the US variation on this on the back of seed packets.

The client wants the user to select one or more hardiness zones, and have the appropriate pixels highlighted on the map. For reasons I'll cover in a later blog post, this needs to be done using the raster's RGB codes and not through some of the simpler techniques such as GeoJSON polygons with varying opacity.

I ended up writing this Leaflet TileLayer extension. It takes a list of RGB codes, and then rewrites the tile images in-canvas so as to have only two colors: those which matched and those which did not.
http://greeninfo-network.github.io/L.TileLayer.PixelFilter/
https://github.com/GreenInfo-Network/L.TileLayer.PixelFilter/

And there's our basic need filled: near-instant pixel-filtering into "selected" and "non-selected" categories, with appropriate color fill.

If you're in a hurry, there's the link to the demo and to the library. But the technical process involved, and the application of it for the client turned out to be very interesting indeed, so they will form my next few postings.