Wednesday, June 12, 2019

Drawing Lines in .NET Core: Comparing ImageSharp and System.Drawing.Common

This is a follow-up post from an early experiment I did almost 3 years ago: https://build-failed.blogspot.com/2016/08/creating-simple-tileserver-with-net.html

At the time, I was trying to setup a simple Tile Server in C# that ran in .NET Core in Win, Mac and Linux, generating images on-the-fly. As I didn't find any C# lib that supported drawing lines and such, I built it myself (see above blog-post for details) over a lib called "ImageProcessor", which allowed me to set individual pixels.

Now, as I'm doing another project that also requires drawing stuff I decided to see how much the landscape has evolved over the last couple of years.

Spoiler alert, it changed significantly. As a summary:
  • ImageProcessor is now marked as "in soft archive mode", with focus shifting to another library called ImageSharp, from the same author
  • ImageSharp now supports Drawing functionality and actually seems very feature rich on that regard. Also, it's fully managed, not relying on native OS calls. If any of you worked with GDI+ in the past, you'll probably remember the mem leaks and thread-"unsafety" of it, particularly if trying to use it on the web.
  • Meanwhile, Microsoft released a Nuget package called "Microsoft.Drawing.Common", which is part of the Windows Compatibility Pack, aiming to help developers migrate their existing .NET code to .NET Core
  • As opposed to ImageSharp, Microsoft.Drawing.Common acts a a bridge to OS specific logic. In Windows it relies on GDI+, on Linux requires installing libgdiplus (from Mono).
I'll quote Scott Hanselman on this (https://www.hanselman.com/blog/HowDoYouUseSystemDrawingInNETCore.aspx):
There's lots of great options for image processing on .NET Core now! It's important to understand that this System.Drawing layer is great for existing System.Drawing code, but you probably shouldn't write NEW image management code with it. Instead, consider one of the great other open source options.
With that said, I'm very interested to understand how both libs compare, thus I've setup a very, very simple test. I'm not going to test the functional differences between both libs and I'm going to focus on the use-case that's most relevant to me: drawing lines.

All the relevant code is available on the Gitrepo: https://github.com/pmcxs/core-linedrawing-benchmark/

The following options can parametrized: number of lines being generated, image size and line width.

Ok, let's start:

Basic Functionality

Starting with a simple test, just to compare (visually) how both results fare:

Number of Lines: 10
Image Size: 800 px
Line Width: 10 px
(default zoom)
First the good news: both work :)

I did find a small difference on the output though. Although the images are very similar, both libraries handle corners differently (at least on the default behavior).

System.Drawing creates a strange protrusion on sharp edges. Zooming in on the previous image:
(zoomed in)
That seems really strange. Increasing the line thickness the effect looks even worse.

Number of Lines: 10
Image Size: 800 px
Line Width: 50 px

(default zoom)
Yeah, even on the default zoom it looks awful with System.Drawing. Basically all corners converge to a single pixel. Thus, the larger the line, the worse the effect gets.

Performance

The code I've placed above already outputs the duration of both approaches. Unfortunately (although probably expected) System.Drawing is still faster than ImageSharp. Versions being tested:
  • SixLabors.ImageSharp.Drawing: 1.0.0-dev002737
  • System.Drawing.Common: 4.6.0-preview5.19224.8
Before starting, I actually found this page: https://docs.sixlabors.com/articles/ImageSharp/GettingStarted.html and it mentions a couple of potential issues that might cause performance problems:
A few troubleshooting steps to try:
  • Make sure your code runs on 64bit! Older .NET Framework versions are using the legacy runtime on 32 bits, having no built-in SIMD support.

I can confirm my test isn't affect by those problems: That flag returns true and I'm running in 64 bits.

I've varied the number of lines (keeping the line width to 10px and image size to 800x800) and the results are as follows:
line width: 10px
I was slightly surprised by the jump from 10 to 100 lines in ImageSharp, particularly as from 100 lines to 1000 its performance was almost the same.

I then did the same test, but increasing the line width from 10px to 100px:
line width: 100px
I don't understand how ImageSharp is faster when drawing 1000 lines vs the previous test with just 100 lines:
  • 100 lines (10 pixel width): 1315 ms
  • 1000 lines (100 pixel width): 776 ms
The ImageSharp performance seems to get worse the thinner the lines are, which I can try to confirm that by going to an extreme of 1px lines.

But, lets make things even more interesting. I'll include my own custom implementation from 3 years ago to the mix. It didn't compile at first, but was easy enough to update, including using the new recommendation for a Span class to set pixels.

The results were interesting. My custom implementation is at its peak the thinnest the line is, hence it gets some really, really strong results. I'm tempted to say there's a bug or an unintended side-effect from the current implementation in ImageSharp. I'll do some additional experiments before submitting a issue.

I had to change the scale to be logarithmic, as otherwise the results would be hidden because of the strange behavior in ImageSharp, as it takes over 2 minutes to render 1000 lines.
line width: 1px
For actual numbers, with 1000 lines of 1px (from slowest to fastest):
  • ImageSharp (Linux): 140000 ms
  • ImageSharp (Windows): 127000 ms
  • System.Drawing (Linux): 132 ms
  • Custom (Linux): 74 ms
  • System.Drawing (Windows): 64 ms
  • Custom (Windows): 62 ms
Increasing the width to 10px the results aren't as strong, but still quite good:
  • ImageSharp (Linux): 1530 ms
  • ImageSharp (Windows): 1311 ms
  • Custom (Linux): 146 ms
  • Custom (Windows): 127 ms
  • System.Drawing (Linux): 123 ms
  • System.Drawing (Windows): 122 ms
There's a catch though:
  • This is literally the only use case I've built. No polygons, beziers, text, etc
  • Also, my "corner handling" logic is really crappy
Regardless, there could be some potential there, thus I might revisit this topic later on. I've included my custom logic on the repo I've mentioned above: https://github.com/pmcxs/core-linedrawing-benchmark/tree/master/src/ImageSharpCustomDrawing

Sunday, June 2, 2019

Building a lightweight api (ASP.NET Core 2.0) to find the geo region that contains a given Coordinate

Context: Given a set of geographic areas I needed an API that would, for a supplied coordinate, return the area that contained it.

So, it needed to:

- load a set of features from a Shapefile
- provide a Web API that receives as input a latitude and longitude
- return the feature that contains that point



Conceptually this is pretty standard, and I expected to have multiple implementations in c# readily available on the web. Interestingly enough, that's not the case (or I probably didn't search properly).

I've decided to build a simple service to address the requirements stated above. It's very simple and available at: https://github.com/pmcxs/region-locator. Its README provides some info on how to setup and use.

How does it work? Breakdown below:

Saturday, May 25, 2019

Playing with Mapbox Vector Tiles (Part 3 - Using Mabox GL)


On this post I'm going to pick-up what I've described on my previous two posts and creating a demo of Mapbox Vector Tiles integrated with Mapbox GL.

I'll be more or less recreating what I did on this experiment: http://psousa.net/demos/maps-and-boardgames-part-2/demo3.html


That old experiment used:
  • TileMill to create the raster tiles
  • UTFGrid tiles to provide meta-information on the hexagons
  • Leaflet as the mapping technology
This new one uses:
  • Tipecanoe to generate the vector tiles (from geojson sources)
  • Mapbox GL JS as the mapping technology
The tricky bit is making sure that the helicopters snap to the displayed hexagons, ideally leveraging the vector data that's present on the tiles.

With that said, I'm actually going to start with the end-result, followed by the break-down on how it's built.


I didn't try to replicate the experience at 100%, but it's close enough.

Live demo at: http://psousa.net/demos/vector/part3/

How this was done:

Thursday, April 18, 2019

Creating a large Geotif with Forest Coverage for the whole World

For a pet-project that I'm making I was trying to find accurate forest coverage for the whole World.

A raster file seemed more adequate and I wanted something like this, but with a much higher resolution (and ideally already georeferenced)


I found the perfect data-source from Global Forest Watch: https://www.globalforestwatch.org/

Global Forest Watch provide an incredible dataset called "Tree Cover (2000)" that has a 30x30m resolution which includes the density of tree canopy coverage overall.

It's too good to be true, right?

Well, in a sense yes. The main problem is that it's just too much data and you can't download the image as a whole.

Alternatively, they provide you an interactive map where you can download each section separately, at: http://earthenginepartners.appspot.com/science-2013-global-forest/download_v1.6.html

This consists of 504 (36x14) images, already georeferenced. For example, if you download the highlighted square above you'll get the the following picture:
https://storage.googleapis.com/earthenginepartners-hansen/GFC-2018-v1.6/Hansen_GFC-2018-v1.6_treecover2000_50N_010W.tif
It's "just" 218MB, hence you can somehow imagine the size of the whole lot. Should be massive.

So, three challenges:
  1. How to download all images
  2. How to merge them together to a single file
  3. (Optional, but recommended) Reducing the resolution a bit to make it more manageable 

1. How to Download all images

Well, doing it manually is definitively an option, although it's probably easier to do it programmatically.
import ssl
import urllib.request

ctx = ssl.create_default_context()
ctx.check_hostname = False
ctx.verify_mode = ssl.CERT_NONE

sections = []

for i in range(8,-6,-1):
    for j in range(-18,18):
        sections.append(f'{abs(i)*10}{"N" if i >= 0 else "S"}_{str(abs(j)*10).zfill(3)}{"E" if j >= 0 else "W"}')

for section in sections:
    url = 'https://storage.googleapis.com/earthenginepartners-hansen/GFC-2018-v1.6/' + \
            f'Hansen_GFC-2018-v1.6_treecover2000_{section}.tif'

    with urllib.request.urlopen(url, context=ctx) as u, open(f"{section}.tif", 'wb') as f:
        f.write(u.read())
The code above, in Python 3.x, iterates all the grid squares, prepares the proper download url and downloads the image.

As the https certificate isn't valid you need to turn off the ssl checks, hence the code at the beginning.

2. How to merge them together to a single file

It's actually quite simple, but you'll need GDAL for that, hence you'll need to install it first.

gdal_merge is incredibly simple to use:

gdal_merge.py -o output-file-name.tif file1.tif file2.tif fileN.tif

Adding to those parameters I would suggest compressing the output, as otherwise an already large file could become horrendously huge.

gdal_merge.py -o output-file-name.tif <files> -co COMPRESS=DEFLATE

And that's it. I'll show how this all ties together on the Python script in the end, but you can "easily" do it manually if you concatenate the 504 file names on this command.

3. Reducing the resolution a bit to make it more manageable 

As I've mentioned, the source images combined result in lots and lots of GBs, which I currently don't have available on my machine. Hence, I've reduced the resolution of each image.

Please note that this isn't simply a resolution change on a Graphics Software, as it needs to preserve the geospatial information. Again, GDAL to the rescue, now using the gdalwarp command:
gdalwarp -tr 0.0025 -0.0025 file.tif

The first two parameters represent the pixel size. From running the command gdalinfo on any of the original tifs I can see that the original pixel size is:

Pixel Size = (0.0002500000000000,-0.0002500000000000)

Empirically I've decided to keep 1/10th of the original precision, hence I've supplied the aforementioned values (0.0025 -0.0025)

As before, I would suggest compressing the content
gdalwarp -tr 0.0025 -0.0025 file.tif -co COMPRESS=DEFLATE

You do lose some quality, but it's a trade-off. If you have plenty of RAM + Disk Space you can keep an higher resolution.

Original
1/10th of resolution
Final script

The following Python 3 does everything in one go. The interesting bit is that I change the resolution of each individual tile before merging the complete map. The script also cleans up after itself, only leaving the final tif file, named "treecover2000.tif"
import ssl
import urllib.request
import os

ctx = ssl.create_default_context()
ctx.check_hostname = False
ctx.verify_mode = ssl.CERT_NONE
extension = ".small.tif"
sections = []

for i in range(8,-6,-1):
    for j in range(-18,18):
        sections.append(f'{abs(i)*10}{"N" if i >= 0 else "S"}_{str(abs(j)*10).zfill(3)}{"E" if j >= 0 else "W"}')

for section in sections:
    print(f'Downloading section {section}')
    url = 'https://storage.googleapis.com/earthenginepartners-hansen/GFC-2018-v1.6/' + \
            f'Hansen_GFC-2018-v1.6_treecover2000_{section}.tif'

    with urllib.request.urlopen(url, context=ctx) as u, open(f"{section}.tif", 'wb') as f:
        f.write(u.read())

    os.system(f'gdalwarp -tr 0.0025 -0.0025 -overwrite {section}.tif {section}{extension} -co COMPRESS=DEFLATE')
    os.system(f'rm {section}.tif')

os.system(f'gdal_merge.py -o treecover2000.tif { (extension + " ").join(sections)}{extension} -co COMPRESS=DEFLATE')
os.system(f'rm *{extension}')
The "treecover2000.tif" ends-up with 751MB and looks AWESOME. Zooming in on Portugal, Spain and a bit of France

Tuesday, March 7, 2017

Playing with Mapbox Vector Tiles (Part 2 - Generating custom tiles)

Part 1. End-to-End experiment
Part 2. Generating custom tiles
Part 3. Using Mabox GL

On my previous post I've played around with various tools on the Mapbox vector tiles ecosystem. On this post I'm taking it further: generating my own vector tiles.

There are various options to generate vector tiles:
  • using a server that generates the vector tiles dynamically based on other data (such as PostGIS)
  • programatically generating the protobuf (rolling your own code or using existing libs)
  • using a tool that receives various inputs and outputs the corresponding vector tiles, which can then be used as shown on my previous post
Most options for these approaches are properly identified at this page: https://github.com/mapbox/awesome-vector-tiles

On this post I'm going to focus on this third option, particularly using a tool from Mapbox called Tippecanoe. I'm not an expert on any of these alternatives (I'm learning as I'm writing this post) but Tippecanoe seems incredibly robust, including a great deal of customisation options when generating the tiles.

Tuesday, February 28, 2017

Playing with Mapbox Vector Tiles (Part 1 - End-to-end experiment)


So, what are vector tiles?

For those that are not familiar with this approach, it provides an alternative to raster tiles where the main difference is that vector tiles provide data instead of a rendered image, although still using a similar xyz tiling structure.

This provides some benefits, such as the ability to have meta-data associated with the tile and being able to generate different presentations for the same information.

One of the problems is that there isn't an official standard for vector tiles. Regardless, with time, a particular vector tile format seems to have gained more traction than all others: the Mapbox Vector Tiles.

Mapbox has provided a consistent ecosystem around it and various providers have started to support it and, on this post, I'm going to play around with some of the existing technology around Mapbox Vector Tiles.

Tuesday, October 25, 2016

Splitting vector and raster files in QGIS (Part 2 - Creating a proper QGIS plugin)

Part 1. Splitting vector and raster files in QGIS
Part 2. Creating a proper QGIS plugin

On my previous post I've created a python script which could be executed inside QGIS in order to split the various layers onto smaller, more manageable chunks.

I've been wanting to create a proper QGIS plugin to package that logic for some time but I've been postponing it as I was expecting it to be slightly tricky, mostly due to the fact that I've never done a QGIS plugin before and anticipated that it could quite challenging.

Truth be told, it was actually quite simple. Due to a mix of awesome online documentation and some starter tools this was mostly a straightforward process.

Let me start by showing the end-result and them I'll explain how it was built.

Tuesday, August 2, 2016

Creating a simple TileServer with .NET Core 1.0 (Part 2 - Drawing Lines)

Part 1. Setting up the the project
Part 2. Improving drawing logic

On my last post I've setup an asp.net core project which outputs simple maps tiles dynamically generated using the ImageProcessor library (which now supports .NET Core).

As I've mentioned that lib doesn't yet include drawing functionality, so on this post I'll try to address that as a pre-requisite to being able to draw proper map tiles. For now I'll focus on drawing lines, including support for variable width and anti-aliasing, doing some benchmarking along the way to make sure that the performance is adequate.

By-the-way, since my first post the proper 1.0 release has been launched, so I'm also updating the code to reflect the latest bits :)

Lots of stuff to do so let's get started:

Friday, June 10, 2016

Creating a simple TileServer with .NET Core RC2 (Part 1 - Setting up the project)

Part 1. Setting up the the project
Part 2. Improving drawing logic

As most of you might have heard Microsoft has recently released .Net core RC2. Although it's still subject to lots of changes I think now is a good time to get on the bandwagon, particularly with the various improvements that have been done.

On this series I'm going to do a cross-platform tile-server that generates tiles with hexagons displayable on a map. As I'm just learning .NET core this will be a learning exercise for me and I'll post all of the steps that I've done in order to achieve the end-result.

Sunday, June 14, 2015

Game-development Log (14. CDN Improvements)


I haven't been really active on the development of this project as I've been occupied with other stuff. Regardless, I'm now focused on returning to active development on it.

Interestingly enough what actually triggered my return was a couple of Azure announcements last week on the CDN front. First, some context:

Although I like Azure quite a lot its CDN offering has been quite lacking to say the least. The community was quite vocal in requesting some essential features but Microsoft neglected to provide any updates or expected delivery dates, as seen here: http://feedback.azure.com/forums/169397-cdn

For example, the top voted feature request was the ability to force content to be refreshed, which is, IMHO, completely essential for a CDN offering:

http://feedback.azure.com/forums/169397-cdn/suggestions/556307-ability-to-force-the-cdn-to-refresh-any-cached-con

The feature was requested 5 years ago, eventually marked as "planned", and no further update was provided by Microsoft, similarly to the other requested features.

Well, all of this until last week, when Microsoft apparently woke up.

First they've provided feedback on most of the feature requests and provided an expectation around release dates. Not ideal (as most of these features are late by a few years) but positive nevertheless.

Then, the icing on the top of the cake was this post:

https://azure.microsoft.com/blog/2015/06/04/announcing-custom-origin-support-for-azure-cdn/

Basically Microsoft shipped three awesome features for the CDN:

I'll just copy&paste from that post:
Custom Origins Supported
Azure CDN can now be used with any origin. Previously, Azure CDN only supported a limited set of Azure Services (i.e. Web Apps, Storage, Cloud Services and Media Services) and you only had the ability to create a CDN endpoint for an Azure Service that was in your Azure Subscription. With this recent update, you can now create a CDN endpoint for any origin you like. This includes the ability to create an origin in your own data center, an origin provided by third party cloud providers, etc. and gives you the flexibility to use any origin you like with Azure CDN!
Multiple CDN Endpoints with the Same Origin
Several of you may have tried to create multiple CDN endpoints for the same origin and found this wasn’t possible due to restrictions. We have now removed the restrictions and you now have the ability to create multiple endpoints for the same origin URL. This provides you a more refined control over content management and can be used to improve performance as multiple host names can be used to access assets from the same origin.
Save Content in any Origin Folder
Previously, when you created a CDN endpoint for cloud services you were required to use “/cdn/” as the default origin path. For example, if the path for your cloud service washttp://strasbourg.cloudapp.net you were required to use http://strasbourg.cloudapp.net/cdn/ as the root path to get content from your origin when you created a CDN endpoint. This restriction has been removed and you can store content in any folder. Using the previous example, you can now use http://strasbourg.cloudapp.net/as the root path to get content from your origin.
These might seem minor changes but let me explain how they positively affect this project:

Sunday, May 3, 2015

Using Genetic Algorithms to solve the Traveling Salesman Problem on Bing Maps

A couple of years ago I was really into Genetic Algorithms and Ant Colony Systems, mostly focusing on solving known NP-Complex problems such as the TRP (Traveling Salesman Problem) and the VRP (Vehicle Routing Problem).

As I have a couple of interesting use-cases that could benefit from these types of algorithms what better way to refresh my knowledge than making a simple mapping experiment solving the TRP problem?

Here's a short summary of these concepts:
  • TRP - Optimization problem that tries to find the shortest route that passes on all supplied points
  • Genetic Algorithm
    • There's a population were each element represents a solution to the problem
    • The algorithm progresses through various iterations, called generations
    • On each generation the various elements of the population mate and create new elements
    • The fittest elements survive and the weakest die
It obviously has much more going on than that, including stuff like roulette selection, mutation, elitism, etc.

I'll create a map where the user can input waypoints and the algorithm will find the shortest path.

Saturday, March 21, 2015

Improve tile-loading at the browser

Slippy Map, such as Bing Maps or Google Maps, is composed of multiple tiles. Each tile, typically a 256x256 image, is individually fetched from the server.

As displays are now supporting incredibly high resolutions, this means that tons of server requests will be required to fill a single map view. For example, on my laptop with retina display, opening a single full screen map will result in about 48 individual tiles being requested simultaneously.

This is a problem as by default web browsers will limit the number of active connections for each domain. This value varies per browser, but we're talking about an average of 6 concurrent downloads per domain, which is quite low. So, assuming all tiles are served from the same domain, lots of throttling will occur.

So, how to cope with this?

1. Tile Size

If you control the tile generation process a "simple" option will be to generate bigger tiles, hence reducing the number of requests. Bing Maps, for instance, supports setting different sizes for the tiles (reference).

For example, setting the tile size to be 512 instead of 256:

var MM = Microsoft.Maps;
var map = new MM.Map(document.getElementById("mapDiv"), {
    center: new MM.Location(45.0, 10),
    zoom: 5,
    credentials:"your key here"});

var tileSource = new MM.TileSource({
    width: 512,
    height: 512,
    uriConstructor:  function(tile) {
        return "images/square512.png";              
    }});

    var tileLayer = new MM.TileLayer({ mercator: tileSource});
    map.entities.push(tileLayer);
In this particular case we're talking about 15 tiles being requested (albeit each one being bigger), which is a big difference from the previous 48.

2. Serve tiles from different urls

A technique called domain sharding can be also be used, on which different domains are used to fetch the same information, thus bypassing the "same-domain" browser limitation.

A good example of this can be seed on Bing Maps, as it's using this technique to speed up serving the tiles.

Taking a look at the web-traffic for the base tiles we can see 4 different hostnames being used:

  • https://t0.ssl.ak.dynamic.tiles.virtualearth.net
  • https://t1.ssl.ak.dynamic.tiles.virtualearth.net
  • https://t2.ssl.ak.dynamic.tiles.virtualearth.net
  • https://t3.ssl.ak.dynamic.tiles.virtualearth.net

The corresponding hostname is determined by the last digit of the quadkey that identifies the tile. For example, tile 0331 will use t1, tile 0330 will use t0, and so on.

Monday, March 16, 2015

Game-development Log (13. Polishing the experience)


Although I haven't been really active on my blog I've been playing a lot with this project adding tons of new stuff. I'm going to detail the various elements that I've updated/implemented during the last month:
  • Pre-Generating additional Zoom Level images
  • Representing Altitude
  • Unit Energy
  • Unit Direction
  • Unit LOD Icons
  • Infantry Unit Type
  • Movement Restriction
  • Coordinate display on higher zoom levels

Tuesday, February 3, 2015

Displaying 3d objects on Bing Maps (using Three.js)

On my previous post I've made a couple of experiments on displaying 2d content on top of Bing Maps using Pixi. Performance was top-notch, particularly if the browser supported WebGL (fallbacking to Canvas otherwise).

This time I'm trying to take the WebGL experiment even further by adding 3d content on top of a map. For that I'm going to use the most popular WebGL 3D lib out there: Three.js

Let me start by showing the end-result.


I've placed some boxes on top of a map. Although the map by itself is 2D the boxes are rendered in 3D using WebGL. Thus, as further away from the screen center (the vanishing point) the more pronounced the depth effect will be.

Tuesday, January 6, 2015

Displaying WebGL on Bing Maps (using Pixi.js)

Something that has been on my backlog for some time is trying to mix Bing Maps and WebGL, similarly to what I've done for an "old" Google Maps experiment.

That previous demo was done on top of a Google maps sample, hence just requiring some small tweaks and improvements. Also, was very low-level and not really practical to adapt to more "real-world" usage, as it required programming the shaders, computing the transformation matrixes, etc.
Thus, I was trying to find a alternative WebGL JS lib that was:
  • Fast
  • Easy to use, albeit still providing some low-level control, namely on primitives drawing
After some research I ended up with two candidates:
IvanK Lib is pretty good (and fast) but Pixi.js takes the cake with tons of functionality and a large community using it.

I'm going to enumerate the various experiments I did, showing a sample page for each.

Thursday, January 1, 2015

Game-development Log (12. Scaling up the loader process)


During the last month I've been mostly refactoring the loader process. As it was I had lots of trouble loading larger areas. I've detailed the problem (and the solution) on a previous post. I've also optimized the loading process and I'm able to load the whole world in about 1 day, already including pre-generating the lower zoom level images.

Regardless, I've haven't yet imported the whole world as everything is still work in progress. I'm planning to add some additional stuff to the map, like mountains, regions, cities and overall improved aesthetics.

I've currently imported to Azure a rectangle that goes from coordinate [15E 60N] to [15W 30N]. Basically this area:


If you zoom in inside this area the Bing Maps tiles are replaced with my own (after zoom level 7). For example, zooming in into London:


Or a random area on Norway:

Or Béchar in Algeria:

The most important thing to note is that those hexagons are not simply cosmetic. All information is being pushed to client-side so that the units will behave differently depending on the terrain. For example, a tank can only cross a river through a bridge:


I've also added:
Deserts
Two levels of forests: dense and sparse

Anyway, you can play around with it at: https://site-win.azurewebsites.net/



Thursday, December 25, 2014

Processing GeoTiff files in .NET (without using GDAL)

I've recently created a process that is able to import geoferenced raster data, namely GeoTiff files... with a catch:
As I didn't find a proper lib to load GeoTiff files I was converting the source files to an ASCII Gridded XYZ format prior to importing. Despite the fancy name it's simply an ASCII file where each line contains:
<LONGITUDE>     <LATITUDE>      <VALUE1>   [ <VALUE 2>]   [<VALUE 3>]
Each line on this file corresponds to a pixel on a raster image, thus it's incredibly easy to parse in C#. The code would be something like this:
foreach (var line in File.ReadLines(filePath))
{
    if (line == null)
    {
        continue;
    }

    var components = line.Split(
        new[] { ' ' }, 
        StringSplitOptions.RemoveEmptyEntries);

    var longitude = double.Parse(components[0]);
    var latitude = double.Parse(components[1]);

    IEnumerable<string> values = components.Skip(2);

    //... process each data item

    yield return dataItem;
}

Nice and easy. But this has two disadvantages:

First, it requires an extra step to convert to XYZ from Geotiff, although easily done with the GDAL:
gdal_translate -of XYZ <input.tif> <output.xyz>
Another disadvantage is that the XYZ becomes much larger than its TIFF counterpart. Depending on the source file the difference could be something like 10 GB vs 500 MB (yes, 50x bigger)

So, I would really like to be able to process a GeoTiff file directly in c#. One obvious candidate is using GDAL, particularly its C# bindings. Although these bindings work relatively well they rely on an unmanaged GDAL installation. For me this is a problem as I want to host this on Azure WebSites, where installing stuff is not an option.
So, I embraced the challenge of being able to process a GeoTIFF file in a fully managed way.


Wednesday, December 10, 2014

Splitting vector and raster files in QGIS programmatically

On the context of the game I'm developing I have to load gigantic .shp/.tif files and parse them to hexagons.

The problem is that this operation is executed in memory (because I need to consolidate/merge data) and the process fails with large files. I had two options:
  • Change the process so that it would store the preliminary results on a separate datastore
  • Instead of loading full shapefiles or raster images, split them into smaller chunks
I did try the first option and, although working, became so much slower (like 100x) that I had to discard it (or at least shelve it while I searched for a better alternative).

I also tried the smaller chunks approach and started by creating "packages" per Country, manually picking data from sources such as Geofabrik.

But this posed two problems:
    • Very tedious work, particularly as there are hundreds of countries.
    • Wouldn't work for larger countries, hitting the memory roadblock as well.
So I opted to split the files in a grid like manner. I decided to use QGIS as it provides all the required tooling.  

Doing this splitting is quite easy to do manually:

1. Open the shapefile (For this example I'm using Countries data from NaturalEarth).


2. Generate a grid.
  • Open the toolbox

  • Choose Geoalgorithms > Vector > Creation > Create graticule
  • Set the desired size (in degrees per rectangle) and set Grid type to Rectangle (polygon)

  • A grid is created over the map (I've changed the opacity of the grid layer to show the map beneath).

3. Select one of the rectangles



4. Perform an intersection between both layers
  • Choose Vector > Geoprocessing Tools > Intersect

  • Choose the base layer as input and the rectangle as the intersect layer, using the "Use only selected features" option.
  • After executing a new layer is created with the intersecting polygons


Now I would just need to repeat this for all the 648 grid items x number of layers to process. I'm assuming this would take about 1 minute per each one and about 10 layers. So, approximately 108 hours non-stop... Not going to happen :). So, what about automating this inside QGIS?

QGIS does provide a command-line/editor/plugins to programmatically leverage its functionalities. Unfortunately for me, it's in Python, which I had never used before. Regardless, the optimistic in me jumped at the opportunity to learn something new.

So, here it is, a sample Python script for QGIS that basically mimics the manual steps I did above:
  • Generates a grid (size hardcoded)
  • Iterates the various tiles in the grid
    • Iterates all layers currently visible (both raster and vector)
      • Outputs the intersection between the layer and the tile
      • Creates a subfolder (harcoded) with all the intersected layers for each tile.
Complete code (update: there's a newer implementation on the bottom)

import processing
import os

#Create a GRID
result = processing.runalg("qgis:creategrid", 1, 360, 180, 10, 10, 0, 0, "epsg:4326", None)

#Add it to the canvas
gridLayer = iface.addVectorLayer(result.get("OUTPUT"),"grid","ogr")

#Iterate ever square on the grid
i = 0
for square in gridLayer.getFeatures():
    i = i + 1
    
    #Create a new layer
    newSquareLayer = iface.addVectorLayer("Polygon?crs=epsg:4326", "temporary_polygon_" + str(i), "memory")        
    provider = newSquareLayer.dataProvider()

    #feature that simply holds one square
    newSquare = QgsFeature()
    newSquare.setGeometry( QgsGeometry.fromPolygon(square.geometry().asPolygon()))
    provider.addFeatures( [ newSquare ] )

    #Make sure the target folder exists
    folder = "c:\\temp\\grid\\grid_" + str(i) 
    if not os.path.exists(folder):
        os.makedirs(folder)

    #iterate the various layers except the grid
    for mapLayer in iface.mapCanvas().layers():
        
        layerType = mapLayer.type()
        layerName = mapLayer.name()
        intersectionName = "intersection_" + layerName + "_" + str(i)

        #vector layers and raster layers are processed differently
        if layerType == QgsMapLayer.VectorLayer and layerName != "grid":
        
            #Calculate the intersection between the specific grid rectangle and the layer
            intersection = processing.runalg("qgis:intersection", mapLayer, newSquareLayer, None)

            iface.addVectorLayer(intersection.get("OUTPUT"),intersectionName,"ogr")
            
            #create a shapefile for this new intersection layer on the filesystem. 
            #A separate folder will be added for each square
            intersectionLayer = QgsMapLayerRegistry.instance().mapLayersByName(intersectionName)[0]
            QgsVectorFileWriter.writeAsVectorFormat(
                        intersectionLayer, 
                        folder + "\\" + layerName + ".shp", 
                        "utf-8", 
                        QgsCoordinateReferenceSystem(4326), 
                        "ESRI Shapefile")

            #remove the intersection layer from the canvas
            QgsMapLayerRegistry.instance().removeMapLayers( [intersectionLayer.id()] )
                
        elif layerType == QgsMapLayer.RasterLayer:
            
            #Calculate the intersection between the specific grid rectangle and the raster layer
            intersection = processing.runalg('saga:clipgridwithpolygon', mapLayer, newSquareLayer, None)
            
            #add the intersection to the map
            iface.addRasterLayer(intersection.get("OUTPUT"), intersectionName)

            #export to file
            intersectionLayer = QgsMapLayerRegistry.instance().mapLayersByName(intersectionName)[0]
            
            pipe = QgsRasterPipe()
            provider = intersectionLayer.dataProvider()
            pipe.set(provider.clone())
            
            rasterWriter = QgsRasterFileWriter(folder + "\\" + layerName + ".tif")
            xSize = provider.xSize()
            ySize = provider.ySize()
            
            rasterWriter.writeRaster(pipe, xSize, ySize, provider.extent(), provider.crs())
            
            #remove the intersection layer from the canvas
            QgsMapLayerRegistry.instance().removeMapLayers( [intersectionLayer.id()] )
            
        else:
            print "layer type not supported"
    
    #Now that all the intersections have been calculated remove the new square
    print "Removing temporary grid item " + str(i) + "(" + newSquareLayer.id() + ")"
    QgsMapLayerRegistry.instance().removeMapLayers( [newSquareLayer.id()] )
    
#Remove the grid
QgsMapLayerRegistry.instance().removeMapLayers( [gridLayer.id()] )
To use this script:
  • Open the Python console:
  1. Open the editor and copy&paste the script there
  2. Save the script
  3. Execute it

Could take a couple of minutes for large files, particularly the raster ones, but at least it's automatic.


Update (11/12/2014)

Although the above script worked as planned, the import process didn't, particularly on the data between the generated tiles. As there's no overlap on the clipped data, after processing small gaps would appear on the end-result. Thus, I've created a brand new implementation that is a little bit more polished and supports a new "buffer" parameter. This allows the tiles to overlap slightly like:
Also, the grid is now created programmatically without using "creategrid" function, which also allowed me to use a more logical X,Y naming for the tiles.

The new code is:
import processing
import os

#######  PARAMS  #######

originX = -180
originY = 90

stepX = 10
stepY = 10

width =  360
height = 180

iterationsX = width / stepX
iterationsY = height / stepY

buffer = 1

j = 0
i = 0

targetBaseFolder = "C:\\temp\\grid"

#######  MAIN   #######

for i in xrange(0,iterationsX):
    
    for j in xrange(0,iterationsY):
        
        tileId = str(i) + "_" + str(j)
        
        folder = targetBaseFolder + "\\" + tileId
        
        if not os.path.exists(folder):
            os.makedirs(folder)
        
        print "Processing tile " + tileId

        minX = (originX + i * stepX) - buffer
        maxY = (originY - j * stepY) + buffer
        maxX = (minX + stepX) + buffer
        minY = (maxY - stepY) -  buffer

        wkt = "POLYGON ((" + str(minX) + " " + str(maxY)+ ", " + str(maxX) + " " + str(maxY) + ", " + str(maxX) + " " + str(minY) + ", " + str(minX) + " " + str(minY) + ", " + str(minX) + " " + str(maxY) + "))"

        tileLayer = iface.addVectorLayer("Polygon?crs=epsg:4326", "tile", "memory")
        provider = tileLayer.dataProvider()
        tileFeature = QgsFeature()

        tileFeature.setGeometry(QgsGeometry.fromWkt(wkt))
        provider.addFeatures( [ tileFeature ] )
        
        for mapLayer in iface.mapCanvas().layers():
            
            layerType = mapLayer.type()
            layerName = mapLayer.name()
            intersectionName = "intersection_" + layerName + "_" + tileId

            #vector layers and raster layers are processed differently
            if layerType == QgsMapLayer.VectorLayer and layerName != "tile":
            
                #Calculate the intersection between the specific grid rectangle and the layer
                intersection = processing.runalg("qgis:intersection", mapLayer, tileLayer, None)

                iface.addVectorLayer(intersection.get("OUTPUT"),intersectionName,"ogr")
                
                #create a shapefile for this new intersection layer on the filesystem. 
                #A separate folder will be added for each square
                intersectionLayer = QgsMapLayerRegistry.instance().mapLayersByName(intersectionName)[0]
                QgsVectorFileWriter.writeAsVectorFormat(
                        intersectionLayer, 
                        folder + "\\" + layerName + ".shp", 
                        "utf-8", QgsCoordinateReferenceSystem(4326), 
                        "ESRI Shapefile")

                #remove the intersection layer from the canvas
                QgsMapLayerRegistry.instance().removeMapLayers( [intersectionLayer.id()] )
                    
            elif layerType == QgsMapLayer.RasterLayer:
                
                #Calculate the intersection between the specific grid rectangle and the raster layer
                intersection = processing.runalg('saga:clipgridwithpolygon', mapLayer, tileLayer, None)
                
                #add the intersection to the map
                iface.addRasterLayer(intersection.get("OUTPUT"), intersectionName)

                #export to file
                intersectionLayer = QgsMapLayerRegistry.instance().mapLayersByName(intersectionName)[0]
                
                pipe = QgsRasterPipe()
                provider = intersectionLayer.dataProvider()
                pipe.set(provider.clone())
                
                rasterWriter = QgsRasterFileWriter(folder + "\\" + layerName + ".tif")
                xSize = provider.xSize()
                ySize = provider.ySize()
                
                rasterWriter.writeRaster(pipe, xSize, ySize, provider.extent(), provider.crs())
                
                #remove the intersection layer from the canvas
                QgsMapLayerRegistry.instance().removeMapLayers( [intersectionLayer.id()] )
                
            else:
                print "layer type not supported"
        
        #remove the temporary tile
        QgsMapLayerRegistry.instance().removeMapLayers( [tileLayer.id()] )

Eventually I'll create a proper QGIS plugin for this that allows setting input parameters without all that hardcoded logic.

Tuesday, November 11, 2014

Game-development Log (11. Persistence and Authentication)


Work-in-progress: https://site-win.azurewebsites.net

On this new iteration the units now have owners. A user can only move his own units (represented in blue) and won't be able to change the other player's units (represented in red).

Authentication

I've used the default template that comes with ASP.NET MVC 5. It includes most of the logic to create and manage users, integrate with Facebook/Google/Twitter, including persistence support with Entity Framework.

I did have to tweak it but most was actually pretty simple. Here's the end-result:
  • I've added register/log in buttons at the top.

  • After Log-in/Register the username is displayed, included a dropdown menu with additional options. The "Profile" button opens a new page to link between local and external accounts and to change the password.


The aesthetic still needs tons of work, but at least most of the functionality is there. Currently I only allow local users, Facebook and Google.

Persistence

The unit movement is now persisted on a SQL Server database at Azure, using Entity Framework. To validate it just move an unit and refresh the browser.


Add Units

I've added a "Developer" menu. Currently it allows units to be created on the map.




HTTPS (and CDN)

As I've added authentication to the site I've changed everything from HTTP to HTTPS. Unfortunately, although Azure's CDN supports HTTPS endpoints, it gets much slower than its HTTP counterpart. Also, I was getting random 503 responses.

So, for now, I've removed my usage of Azure's CDN altogether, either pointing to the blob storage for the static image tiles or to the Tile-Service webAPI for the dynamic ones.

I really love Azure, but its CDN really sucks, particularly comparing against the likes of Amazon, Cloudfront, Akamai, Level3, etc.

Sunday, October 26, 2014

Game-development Log (10. Real-time with SignalR)


Work-in-progress: https://site-win.azurewebsites.net/ 

On this iteration I'm adding real-time notifications to the unit movement using the (awesome) SignalR library.

The idea for this will be simple: if a player moves a unit on his browser window all the others players will see that unit move in realtime on their browsers.

I've created an ultra-simple video demoing it. I've opened three separate browser windows: one with Chrome, Firefox and Safari. The movement of the units should be synchronised in real-time across the various browsers, regardless of the one that triggered the movement.


So, how was this setup? Easy as pie.

1. Add SignalR NuGet package
install-package Microsoft.AspNet.SignalR

2. Initialise SignalR on the Startup routine of ASP.NET MVC 5 (on previous versions it's slightly different)
public partial class Startup
{
    public void Configuration(IAppBuilder app)
    {
        app.MapSignalR();
    }
}
3. Create a Hub
My hub will be ultra simple. It will simply receive a movement notification and share it with all the other SignalR clients as is.
public class UnitHub : Hub
{
    public void SetPosition(dynamic unitdata)
    {
        Clients.Others.unitUpdated(unitdata);
    }
}

4. Create Javascript code both to send a SignalR message when the player moves a unit or when a notification of another player moving a unit has been received.
//code to receive notification from server
unitUpdated = function (unit) {
 // (...)
};

// Code to submit to server
server.setPosition(unit);
5. (Optional) If you're using Azure, you should enable WebSockets on your website.
And that's it.

Please note that I haven't yet added persistence nor support for different users. Thus, everyone will be moving the same units. Anyway, I don't have that many page views on my blog for it to be a problem :)

So, what's next?
- Adding a database to persist the units position
- Adding users, each with his own units.
- Authentication with Facebook, Twitter, Google, Windows, Email