During the years, people keep asking me the same question: can you store (and query) raster data in your Carto (former CartoDB) account? How?

Well, the answer is yes. And this is more a PostGIS Raster issue than a Carto issue. Let’s go.

Disclaimer: This code was written couple of years ago. It may be outdated, but still works.

Store raster data

First thing you need to do is easy: just drop your georeferenced raster file into the Dataset page of your account. The raster will be imported, and you’ll see something like this

Do you see that foto_pnoa table with gray shading name? That means you cannot click on it to see the table’s data. But don’t worry: your data is safe in the table. Don’t trust me? Just enter another table and run this query (of course, using your own table’s name)

select st_summarystats(the_raster_webmercator, 1) as stats from foto_pnoa

You should see the stats of your raster data now.

View your raster data over a map

Now, to the fun part. We’ll code a bit of JavaScript to put that data over a Carto map. We’ll even allow some calculation with the data (average raster value within a region). So, go for it.

The HTML

Pretty simple and straightforward. This is the relevant part of the HTML. But don’t worry, I’ll provide a link to the complete example at the end of the post.

<div class="header">
<h1>PostGIS Raster test</h1>
<h2>Draw a figure and click on it to see the avg raster value</h2>
</div>
<div id="map"></div>

The JS

We’ll create a Leaflet map, using:

Relevant parts of the code here (again, you’ll get the complete code later)

Create the map and the controls

We create a Leaflet map and the draw controls

// Create map
var map = new L.Map('map', {
zoomControl: true,
drawnControl: true,
center: [37.383333, -5.983333],
zoom: 11
});

// Add CartoDB basemaps
L.tileLayer('http://{s}.basemaps.cartocdn.com/light_all/{z}/{x}/{y}.png', {
attribution: '<a href="http://cartodb.com">CartoDB</a> © 2014',
maxZoom: 18
}).addTo(map);

// Add drawn controls
var drawnItems = new L.FeatureGroup();
map.addLayer(drawnItems);
var drawControl = new L.Control.Draw({
position: 'bottomleft',
draw: {
polyline: false,// Turns off this drawing tool
marker: false,
polygon: false,

rectangle: {
shapeOptions: {
color: '#a63b55'
},
showArea: true
},

circle: {
shapeOptions: {
color: '#662d91'
},

showArea: true
}

},
edit: {
featureGroup: drawnItems
}
});
map.addControl(drawControl);

What can I do with those controls? Let’s see

Handle draw actions

Whenever we draw a figure:

  1. Using the figure’s coords, we build a PostGIS geometry by calling ST_MakeBox2D, if we drawn a rectangle, or ST_Buffer, if we drawn a circle.
  2. Run a query to check the average raster value within that geometry and show the value

Check it out in the next snippet

map.on('draw:created', function (e) {
var type = e.layerType,
layer = e.layer;

var pol_pgis = null;

switch(type) {

// Create a Rectangle geometry in PostGIS
case 'rectangle':
var coords = layer.getLatLngs();

var southWest = L.latLng(coords[1].lat, coords[1].lng);
var northEast = L.latLng(coords[3].lat, coords[3].lng);

var pol_pgis = "st_transform(ST_SetSRID(ST_MakeBox2D(ST_Point(" +
coords[1].lng + ", " + coords[1].lat + "),ST_Point(" +
coords[3].lng + "," + coords[3].lat + ")),4326), 3857)";

break;

// Create a circle geometry in PostGIS
case 'circle':
var center = layer.getLatLng();

var pol_pgis = "st_transform(geometry(st_buffer(geography(st_setsrid(st_point(" +
center.lng + ", " + center.lat + "), 4326)), " + layer.getRadius() + ")),3857)";

break;

case 'polygon':
console.log(layer.toGeoJSON());

}

if (pol_pgis) {
q = "SELECT avg((stats).mean) as m from (select st_summarystats(the_raster_webmercator, 1) as stats from foto_pnoa where st_intersects(the_raster_webmercator, " + pol_pgis +")) as foo";

console.log("QUERY: " + q);

var sql = new cartodb.SQL({user: 'libregis'});
sql.execute(q)

.done(function(data) {
if (data.rows && data.rows.length > 0)
layer.bindPopup("Average raster value inside the " + type + ": " + data.rows[0].m);

else
layer.bindPopup("Could not get avg value!");
})

.error(function(errors) {
layer.bindPopup("Could not get avg value!");
})
}

else {
layer.bindPopup("Could not get avg value!");
}

drawnItems.addLayer(layer);
});

Show the raster tiles

The final touch. We put the raster tiles over the map using the Maps API, using AJAX (old fashioned way, I know…)

var config = {
"version": "1.3.1",
"layers": [
{
"type": "cartodb",
"options": {
"sql": "select * from foto_pnoa",
"cartocss": "#foto_pnoa {raster-opacity: 0.5;}",
"cartocss_version": "2.3.0",
"geom_column": "the_raster_webmercator",
"geom_type": "raster"
}
}
]
};

var request = new XMLHttpRequest();
request.open('POST', currentEndpoint(), true);
request.setRequestHeader('Content-Type', 'application/json; charset=UTF-8');
request.onload = function() {
if (this.status >= 200 && this.status < 400){             var layergroup = JSON.parse(this.response);             var tilesEndpoint = currentEndpoint() + '/' + layergroup.layergroupid + '/{z}/{x}/{y}.png';             var protocol = 'https:' == document.location.protocol ? 'https' : 'http';             if (layergroup.cdn_url && layergroup.cdn_url[protocol]) {                 var domain = layergroup.cdn_url[protocol];                 if ('http' === protocol) {                     domain = '{s}.' + domain;                 }                 tilesEndpoint = protocol + '://' + domain + '/' + currentUser() + '/api/v1/map/' + layergroup.layergroupid + '/{z}/{x}/{y}.png';             }             rasterLayer = L.tileLayer(tilesEndpoint, {                 maxZoom: 18             }).addTo(map);         } else {             throw 'Error calling server: Error ' + this.status + ' -> ' + this.response;
}
};
request.send(JSON.stringify(config));

That’s it. You can check the final result in the next codepen

You can grab the code here

Credits

The very first version of the code was made by Raul Ochoa

Advertisements