[postgis-users] Help with ST_CLIP

Regina Obe lr at pcorp.us
Fri Mar 12 15:21:41 PST 2021


Simon,

How big is the area you are trying to do?  

Couple of thoughts of what to try
1) ST_Clip is a bit finicky in that if a pixel is only partially covered it is not considered to be on. So what I usually do is just buffer my clipping geometry  a tincy bit like 1 meter like below to grab any pixels only partially overlapped.

2) Also you should only be returning the part that intersects.  You are missing a WHERE clause unless you have just one raser tile, but still good to have a WHERE anyway as a precaution.  You'll know immediately if you get no rows they don't even intersect.
You can speed up the intersect check by adding a functional index geometry(geog) to your geography table

3) I tend to have last arg be crop = true (you had false).  If you set to false, you end up with the same raster size as your original but filled up with nulls there were was data.  Generally not something you want.


Try this:
 create table gbr100x100clip
 as
SELECT ST_Clip(g.rast, 1, ST_Buffer(m.geog,1)::geometry, true) as rast   
from 
 gbr_marine_park_boundary as m
	INNER JOIN        gbr100x100 as g ON ST_Intersects(g.rast, m.geog::geometry)
;



>> This works but the value (depth) of the whole (slipped) rasters is NULL.
>>
>> I expected to see depth values but this returns null values.
>>
>> SELECT (stats).*
>>   FROM (SELECT ST_SummaryStats(rast, 1) As stats FROM gbr100x100clip) 
>> As foo  order by max desc;
>>

> -----Original Message-----
> From: postgis-users [mailto:postgis-users-bounces at lists.osgeo.org] On
> Behalf Of Simon G Greener
> Sent: Friday, March 12, 2021 4:56 PM
> To: postgis-users at lists.osgeo.org
> Subject: Re: [postgis-users] Help with ST_CLIP
> 
> Stephen,
> 
> A type 4326 not 4236.
> 
>  > select st_astext(geog::geometry) from gbr_marine_park_boundary;
> 
> The numbers look correct after the cast.
> 
> Even though the data is extensive should I consider transforming to a
> projected SRID?
> 
> Geographic rasters are not supported right? This is the reason for the cast to
> call the ST_Clip function.
> 
> regards
> 
> Simon
> 
> On 13/03/2021 5:42 am, Stephen Woodbridge wrote:
> > Hi Simon,
> >
> > A geography column in 4236 sounds strange. Check the coordinate values
> > of the cast polygon.
> >
> > select st_astext(geog::geometry) from gbr_marine_park_boundary;
> >
> > and see if they look like what you expect.
> >
> > -Steve W
> >
> > On 3/12/2021 1:37 AM, Simon SPDBA Greener wrote:
> >> "POSTGIS=""2.4.3 r16312"" PGSQL=""100"" GEOS=""3.6.2-CAPI-1.10.2
> >> 4d2925d6"" PROJ=""Rel. 4.9.3, 15 August 2016"" GDAL=""GDAL 2.2.3,
> >> released 2017/11/20"" LIBXML=""2.9.4"" LIBJSON=""0.12.1""
> >> LIBPROTOBUF=""1.2.1"" RASTER"
> >>
> >> I have a polygon that defines the boundary of the Great Barrier Reef.
> >>
> >> A) That boundary is stored as a vector geography column (4326) in a
> >> table.
> >> B) I have loaded a raster dataset of depths (4326) and checked its
> >> statistics.
> >>     The data loaded correctly and I can see valid depth values.
> >>
> >> Now I want to clip the data in B with the polygon in A.
> >>
> >> (I can't use a raster function with a geography object so I cast it
> >> to geometry.)
> >>
> >> create table gbr100x100clip
> >> as
> >> select ST_Clip(g.rast, 1, m.geog::geometry, false) as rast   from
> >> gbr_marine_park_boundary as m,        gbr100x100 as g;
> >>
> >> This works but the value (depth) of the whole (slipped) rasters is NULL.
> >>
> >> I expected to see depth values but this returns null values.
> >>
> >> SELECT (stats).*
> >>   FROM (SELECT ST_SummaryStats(rast, 1) As stats FROM gbr100x100clip)
> >> As foo  order by max desc;
> >>
> >> Can anyone give me an idea of what I am doing wrong?
> >>
> >> regards
> >> Simon
> >>
> >> _______________________________________________
> >> postgis-users mailing list
> >> postgis-users at lists.osgeo.org
> >> https://lists.osgeo.org/mailman/listinfo/postgis-users
> >
> >
> _______________________________________________
> postgis-users mailing list
> postgis-users at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/postgis-users



More information about the postgis-users mailing list