<html><body><div style="color:#000; background-color:#fff; font-family:Courier New, courier, monaco, monospace, sans-serif;font-size:12pt"><div><span>Sorry there was an error in that code it should be:</span></div><div><span><br></span></div><div><span>import rpy2.robjects as robjects<div></div>
<div>import sys</div>
<div>#Create R stuff and load libraries</div>
<div>R = robjects.rs</div>
<div>#R.library('sp')</div>
<div>R.library('gstat')</div>
<div>R.library('raster')</div>
<div>R.library('RPostgreSQL')</div>
<div>R("drv <- dbDriver('PostgreSQL')")</div>
<div>R("con <- dbConnect(drv, dbname='progsci', user='postgres')")</div>
<div>R("rs <- dbSendQuery(con, 'SELECT * FROM crimes WHERE crimetype='Burglary')")</div>
<div>R("data <- fetch(rs, n=-1)")</div>
<div>R("coordinates(data) <- ~x+y")</div>
<div>R("x.range <- as.integer(range(e@coords[,1]))")</div>
<div>R("y.range <- as.integer(range(e@coords[,2]))")</div>
<div>R("grd <- expand.grid(x=seq(from=x.range[1], to=x.range[2], bx=25), y=seq(from=y.range[1], to=y.range[2], by=25))")</div>
<div>R("coordinates(grd) <- ~x+y")</div>
<div>R("gridded(grd) <- TRUE")</div>
<div>R("idw.out <- idw(value~1, e, grd, idp=2.5)")</div>
<div>spplottext = "spplot(idw.out, \"var1.pred\", col.regions=topo.colors(20), pretty=TRUE)"</div>
<div></div>
<div>print(spplottext)</div>
<div>res = R(spplottext)</div>
<div>print(res)</div>
<div>sys.exit(0)</div>
<div></div></span></div><div><br><blockquote style="border-left: 2px solid rgb(16, 16, 255); margin-left: 5px; margin-top: 5px; padding-left: 5px;"> <div style="font-size: 12pt; font-family: 'Courier New', courier, monaco, monospace, sans-serif; "> <div style="font-size: 12pt; font-family: 'times new roman', 'new york', times, serif; "> <font size="2" face="Arial"> <hr size="1"> <b><span style="font-weight:bold;">From:</span></b> Phil James <borntopedal@yahoo.co.uk><br> <b><span style="font-weight: bold;">To:</span></b> PostGIS Users Discussion <postgis-users@postgis.refractions.net> <br> <b><span style="font-weight: bold;">Sent:</span></b> Tuesday, 6 December 2011, 7:58<br> <b><span style="font-weight: bold;">Subject:</span></b> Re: [postgis-users] spatial distribution maps (heat maps?)<br> </font> <br><div id="yiv2134649023"><div><div style="color: rgb(0, 0, 0); background-color: rgb(255, 255, 255); font-size: 12pt; font-family: 'Courier
New', courier, monaco, monospace, sans-serif; "><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><span>I would use something that is designed to do this analysis - R would be an obvious choice but there are others - using Python you can connect to Postgis to grab the data and then rpy to run R commands from python - we have used this configuration successfully using Django and it works fast enough for the web generating an image - this could of course be optimised to cache outputs if performance is an issue.</span></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><span><br></span></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><span>A synopsis of spatial functions in R is shown below:</span></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><span><br></span></div><div style="font-size: 12pt;
font-family: courier, monaco, monospace, sans-serif; "><span><a rel="nofollow" target="_blank" href="http://cran.r-project.org/web/views/Spatial.html">http://cran.r-project.org/web/views/Spatial.html</a><br></span></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><br></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; ">Some sample code to run IDW below (just a hack but you get the idea).</div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><br></div><div><div>import rpy2.robjects as robjects</div><div>import sys</div><div>#Create R stuff and load libraries</div><div>R =
robjects.r</div><div>#R.library('sp')</div><div>R.library('gstat')</div><div>R.library('raster')</div><div>R.library('RPostgreSQL')</div><div>R("drv <- dbDriver('PostgreSQL')")</div><div>R("con <- dbConnect(drv, dbname='geoprocessing', user='postgres')")</div><div>R("rs <- dbSendQuery(con, 'SELECT x, y, value FROM randompts')")</div><div>R("na <- fetch(rs, n=-1)")</div><div>R("e <- na.omit(na)")</div><div>R("coordinates(e) <- ~x+y")</div><div>R("x.range <- as.integer(range(e@coords[,1]))")</div><div>R("y.range <- as.integer(range(e@coords[,2]))")</div><div>R("grd <- expand.grid(x=seq(from=x.range[1], to=x.range[2], bx=25), y=seq(from=y.range[1], to=y.range[2], by=25))")</div><div>R("coordinates(grd) <- ~x+y")</div><div>R("gridded(grd) <- TRUE")</div><div>R("idw.out <- idw(value~1, e, grd, idp=2.5)")</div><div>spplottext = "spplot(idw.out, \"var1.pred\", col.regions=topo.colors(20),
pretty=TRUE)"</div><div><br></div><div>print(spplottext)</div><div>res = R(spplottext)</div><div>print(res)</div><div>sys.exit(0)</div></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><span><br></span></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><span>I am sure there are other solutions along the same lines.</span></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><span><br></span></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><span>Phil</span></div><div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "><blockquote style="border-left:2px solid rgb(16, 16, 255);margin-left:5px;margin-top:5px;padding-left:5px;"> <div style="font-size: 12pt; font-family: courier, monaco, monospace, sans-serif; "> <div style="font-size: 12pt; font-family: times, serif; "> <font
size="2" face="Arial"> <hr size="1"> <b><span style="font-weight:bold;">From:</span></b> Puneet Kishor <punk.kish@gmail.com><br> <b><span style="font-weight:bold;">To:</span></b> PostGIS Users Discussion <postgis-users@postgis.refractions.net> <br> <b><span style="font-weight:bold;">Sent:</span></b> Tuesday, 6 December 2011, 1:54<br> <b><span style="font-weight:bold;">Subject:</span></b> Re: [postgis-users] spatial distribution maps (heat maps?)<br> </font> <br>Steve,<br><br>Thanks for your reply.<br><br><br>On Dec 5, 2011, at 6:52 PM, Stephen Woodbridge wrote:<br><br>> On 12/5/2011 6:45 PM, Puneet Kishor wrote:<br>>> I am looking to create some heatmap-ish or DEM-ish spatial<br>>> distribution visualizations. Any suggestions on how to go about doing<br>>> so given a lon-lat containing
geom column in PostGIS?<br>> <br>> ..<br>> <br>> Create an empty raster. Then select your points and for each points and for each one "add" the sample dot centered on your location to to the raster. By "add", I mean you need to add the value of the pixels in the dot to the pixels in the raster. So raster=1 + dot=3 => raster=4 and limit it at maximum saturation. This should work for multiple channels also. This should give you a heat map.<br>> <br><br><br>Reading the above made me realize that I should have rephrased my question -- I don't want to create images on the server side. I realize now that what I really want to do is to do spatial clustering on the server side and then send the data to the browser. I wrote my own very naive clustering routine in Perl, and also tried Algorithm::KMeans [1]. This kind of analysis allows me to create a summary of my data that I can then plot on the client (see image at [2]).<br><br>Of course, my
algorithm is way too naive, and waaaaay too slow, although I *can* compute the summary and cache it using Storable.<br><br>So, here is my rephrased question -- I am looking to do spatial clustering on my Pg data. The added complication is that I do not have access to WKTRaster.<br><br><br>[1] <a rel="nofollow" target="_blank" href="http://search.cpan.org/~avikak/Algorithm-KMeans-1.30/lib/Algorithm/KMeans.pm">http://search.cpan.org/~avikak/Algorithm-KMeans-1.30/lib/Algorithm/KMeans.pm</a><br>[2] <a rel="nofollow" target="_blank" href="http://dl.dropbox.com/u/3526821/occurrences.png">http://dl.dropbox.com/u/3526821/occurrences.png</a><br>_______________________________________________<br>postgis-users mailing list<br><a rel="nofollow" ymailto="mailto:postgis-users@postgis.refractions.net" target="_blank" href="mailto:postgis-users@postgis.refractions.net">postgis-users@postgis.refractions.net</a><br><a rel="nofollow" target="_blank"
href="http://postgis.refractions.net/mailman/listinfo/postgis-users">http://postgis.refractions.net/mailman/listinfo/postgis-users</a><br><br><br> </div> </div> </blockquote></div> </div></div></div><br>_______________________________________________<br>postgis-users mailing list<br><a ymailto="mailto:postgis-users@postgis.refractions.net" href="mailto:postgis-users@postgis.refractions.net">postgis-users@postgis.refractions.net</a><br><a href="http://postgis.refractions.net/mailman/listinfo/postgis-users" target="_blank">http://postgis.refractions.net/mailman/listinfo/postgis-users</a><br><br><br> </div> </div> </blockquote></div> </div></body></html>