Followup: Density Mapping

[Bill Binko]

Hi everyone.

A few weeks ago, I asked a question about density mapping, and got
several responses, including several very helpful pointers.

I thought I would write back with where the progress I've made, along
with another request for feedback.

Here are three images that show some of the density maps I've been able
to make on the fly:

http://m.mapshine.com/images/dazzlehX6EfG.png
http://m.mapshine.com/images/dazzlesSqxTl.png
http://m.mapshine.com/images/dazzlePjFMQF.png

I'm sure there are statistics professors out there who will roll their
eyes at my approach, (and please, if I'm making any horrific mistakes,
please let me know) but the proof-of-concept I developed is performing
quite well in terms of speed.

The proof-of-concept works like this (in PHP/Mapscript):

1) Read a Standard Deviation in geographic units (1/3 mile for the
images above)
2) Generate a "stamp" image that is a gradient that starts at 100%
opaque at the center and "fades out" to 0% opaque at the center using
the Parzen Window Gaussian Kernel
equation:http://en.wikipedia.org/wiki/Parzen_window where "x" is
actually the distance from the center of the stamp
3) Create a new image to lay over the Mapserver generate map with the
same size, and "ImageAlphaBlending()" set to true.
4) Query the map layer using "whichShapes()" to get the points on the
layer we're interested in
5) Iterate through all the points on the layer and draw the "stamp"
image into the new image.  Because of GD's Alpha Blending, these
"stamps" build up effectively performing the Summation in the Parzen
Window document above.
6) Walk through the entire image (one point at a time) and map the 0-127
alpha value of each point onto table of images (just to make it look
good :)

I have thought about how this type of functionality could make it back
into the core Mapserver, and I actually think it could be fairly simple.
However, this is obviously a proof of concept, and I'm not pushing for
this... I just wanted to update everyone.

To take this same approach with the C Mapserver code, we'd need to do
the following:

1) Implement semi-transparent symbols per danmo's suggestion (right now,
we only have alpha transparency at the layer level, not the symbol
level)
2) Setup a hook so that user code can create the png "stamp" after the
layer's requested extent is known and before any drawing is done
3) Setup a hook where the user can get access to the raw layer image for
post-processing.

I am not 100% sure that the approach I've taken is completely valid.
However, I am quite sure that is far faster than the alternatives I
faced which were a) spawning out to GRASS/R to generate a large raster
or b) doing a distance calculation + a gaussian equation for each pair
of pixel/Point in the extent!  

I look forward to your responses.

Bill

-- 
Bill Binko <bill at binko.net>