Union rasters previously cut by bounding box

Question

How do I union rasters, where the original data is grouped by an attribute without running into an issue where the rasters that should be unioned, aren't the same sizes?

I'm trying to generate a single PNG for the users' current bounding box. Cutting the required geometries works great and I end up with a few clipped geoms that match the bounding box exactly. The geometries are grouped by an attribute int and coloured accordingly.

However after this, I can't seem to get the rasters to join up into a single image.

---

Some detail

Bounds is drawn in blue

--bounds
select 
   st_geomfromtext(
        'POLYGON((18.48287 -33.95910, 18.49523 -33.95910,18.49523 -33.96829,18.48287 -33.96829,18.48287 -33.95910))
       ')
  ) as bounding_box

Data

--bounds
select data.wkb_geometry, data.attibute
from (
         values ('POLYGON((18.4893249 -33.9693479,18.487338 -33.975303,18.4836563 -33.9722486,18.4810635 -33.9694083,18.4804401 -33.9701518,18.4788584 -33.969387,18.4769812 -33.9711783,18.4758629 -33.9699897,18.47497 -33.96881,18.474341 -33.9686,18.4763888 -33.9636919,18.4771354 -33.9638799,18.4780189 -33.964247,18.4793621 -33.9646786,18.4798761 -33.964869,18.480315 -33.9650772,18.4808655 -33.9654135,18.482002 -33.9662287,18.4893249 -33.9693479))',1),
                ('POLYGON((18.4874775 -33.9753787,18.4895958 -33.9693783,18.4911785 -33.9672392,18.4912114 -33.9667627,18.4911048 -33.9664196,18.4895065 -33.9645728,18.4884931 -33.96334,18.4879777 -33.9624253,18.4876715 -33.9620845,18.4863472 -33.9608744,18.4858276 -33.959958,18.4856109 -33.959617,18.4855713 -33.959227,18.4856073 -33.9587142,18.4859062 -33.9558208,18.4861392 -33.953389,18.4862619 -33.9529841,18.4865074 -33.9523433,18.4870734 -33.9526725,18.4884629 -33.9534733,18.4896968 -33.9542027,18.4909845 -33.955003,18.4919826 -33.95558,18.4933138 -33.9563781,18.4936203 -33.9565488,18.4937229 -33.9567017,18.4938761 -33.9568949,18.4939786 -33.9570856,18.4936162 -33.9585412,18.4932268 -33.9600814,18.4933172 -33.9602339,18.4935444 -33.9602731,18.4937287 -33.9603389,18.4936413 -33.9608175,18.4937095 -33.9615053,18.4948511 -33.9614817,18.4956911 -33.9636467,18.4954551 -33.9660477,18.4944031 -33.9682377,18.4932861 -33.9725157,18.4924889 -33.9766531,18.4911568 -33.9783253,18.4905335 -33.9781892,18.4874775 -33.9753787))',1),
                ('POLYGON((18.478434 -33.959864,18.479759 -33.959663,18.48191 -33.959769,18.485405 -33.959819,18.4860366 -33.9609406,18.487591 -33.962258,18.487962 -33.962984,18.489653 -33.964911,18.4906361 -33.9661385,18.4909099 -33.9666043,18.49098 -33.96699,18.4893785 -33.9692589,18.4857733 -33.9677527,18.4820824 -33.9661842,18.481005 -33.965369,18.4803673 -33.9650016,18.4796866 -33.9647409,18.4764746 -33.9635707,18.476956 -33.962229,18.47781 -33.961372,18.478434 -33.959864))',2)
     ) as data (wkb_geometry, attibute);

Now the original polygons are grouped by attribute and unioned. Generating the PNGs works great however the result is not unioned into a single image. Infact the result is three individual PNGs.

The result above is missing the green polygon in the previous image, which I need to be included to create a final single PNG with all clipped polygons, regardless of attribute.

ignoring the fact that the above blue polygon is missing, unintentional when the screenshot was taken

This is where I run into the issue where if I attempt to union the rasters, I get the error

ST_Union(raster): ERROR: rt_raster_from_two_rasters: The two rasters provided do not have the same alignment

Ideally I guess would be to create a reference raster (?) as to help keep context of where the geometries should be unioned.

I'm not sure if this is the solution as I can imagine you'd end up with a randomly generating a result where attribute=1 would overlap attribute=2 and vice versa, depending on the union order of the geoms.

This is rough. Anyway here is the query.

Query

with bbox as (select st_setsrid( st_geomfromtext(
                           'POLYGON((18.48287 -33.95910,18.49523 -33.95910,18.49523 -33.96829,18.48287 -33.96829,18.48287 -33.95910))'), 4326) as bounding_box
    ), clipped as (
        SELECT ST_union(st_intersection( wkb_geometry,bbox.bounding_box))
                   as wkb_geometry,attribute
        from data_table,
             bbox where st_intersects(wkb_geometry, bbox.bounding_box) group by attribute
    ), rasters as (
    select ST_AsRaster(
                clipped.wkb_geometry,
                256, 256, ARRAY ['8BUI', '8BUI', '8BUI'],
                (CASE
                     WHEN attribute = 1 THEN ARRAY [66,135,245]
                     WHEN attribute = 2 THEN ARRAY [48, 59, 99]
                    --...
                     ELSE ARRAY [147,150,150] END),
                ARRAY [0,0,0]
               ) as rast
    from clipped
)
SELECT ST_AsPNG(rasters.rast) as tile
from rasters;

--------------------------------------------------------------------------------

UPDATE

Image generation works with the answer from [Pierre Racine] however I'm struggling to generate low quality pngs where the sides aren't cut off (see second image with darker background

This results in tiles that dont fit their div exactly.

The only solution for now is generating high resolution images, takes up to 30 sec, which reduces these slivers to mostly unnoticeable. However not ideal.

Answer 1

You have to align all the rasters created with ST_AsRaster() with additional parameters (scalex, scaley, gridx, gridy) so that you can then ST_Union() them together into a single raster before ST_AsPNG() it.

In the example below I have set both pixel sizes (scalex and scaley) to 0.05 and both pixel alignment (gridx and gridy) to 0 when calling ST_AsRaster(). This correspond to variant 4 of ST_AsRaster() in the doc.

WITH bbox AS (
  SELECT ST_SetSRID(ST_GeomFromText('POLYGON((18.48287 -33.95910,18.49523 
  -33.95910,18.49523 -33.96829,18.48287 -33.96829,18.48287 -33.95910))'), 4326) 
  AS bounding_box
), clipped AS (
  SELECT ST_Union(ST_Intersection(wkb_geometry, bbox.bounding_box)) AS wkb_geometry, 
         attribute
  FROM data_table, bbox 
  WHERE ST_Intersects(wkb_geometry, bbox.bounding_box) 
  GROUP BY attribute
), rasters AS (
  SELECT ST_AsRaster(clipped.wkb_geometry, 0.05, -0.05, 0, 0, 
                     ARRAY['8BUI', '8BUI', '8BUI'], 
                     (CASE WHEN attribute = 1 THEN ARRAY [66, 135, 245] 
                           WHEN attribute = 2 THEN ARRAY[48, 59, 99]
                           --...
                           ELSE ARRAY[147,150,150] END),
                     ARRAY[0,0,0]) AS rast
  FROM clipped
)
SELECT ST_AsPNG(ST_Union(rasters.rast)) as tile
FROM rasters;