[postgis-users] PostGIS Raster efficient enough to get height given coords from a large DEM?

[Stefan Keller]

Hi Volkmar

2011/12/4 Volkmar Herbst <volkmar.herbst at gmx.de>:
> Yes, we also do use plr for kriging of soil data, which works very fast, but
> there are many different variants of it. It helps to know about geo
> statistics and get the right type and formula for the kriging. It is also
> possible to do the following steps:
> 1. create your grid (in plr)
> 2. calculate kriging (in plr)
> 3. create contours (in plr)
> 4. create contour polygons (I have done it in a postgis function but it
> should be also possible in plr)
> I abandoned the raster data type as kriging output to postgis as I could not
> create the postgis raster from the plr raster in an acceptable time. The
> ratser postgis function st_setvalue which works on each pixel takes too long
> and there seem to be no other way to set raster values programmatically,
> that’s why I switched to contours.

That's obvious to me that ST_SetValue is slow compared to a setter
method on a raw 2D array.
If you need fast setters you have to do this 1. in memory and 2. using
a raw 2D array and a programming (or scripting) language.
In my case I only need fast access and analysis functions.

Yours, Stefan

> We do need array to pixel value or so-
> see
> http://trac.osgeo.org/postgis/ticket/1228
> Hope that helps
> Volkmar
>
>
>
> -----Ursprüngliche Nachricht-----
> Von: postgis-users-bounces at postgis.refractions.net
> [mailto:postgis-users-bounces at postgis.refractions.net] Im Auftrag von Bborie
> Park
> Gesendet: Sonntag, 4. Dezember 2011 13:46
> An: PostGIS Users Discussion
> Betreff: Re: [postgis-users] PostGIS Raster efficient enough to get height
> given coords from a large DEM?
>
> On Sat, Dec 3, 2011 at 7:16 PM, Stefan Keller <sfkeller at gmail.com> wrote:
>> Hi Bborie
>>
>> Many thanks for your response.
>> To summarize: My goal is to get the height (in meters) given an
>> arbitrary lat/lon coordinate from the latest free digital elevation
>> model (like mentioned here
>>
> http://slashgeo.org/2011/12/01/New-Free-Elevation-Dataset-named-Global-Multi
> -resolution-Terrain-Elevation-Data-2010-GMTE
>> ).
>>
>> The approach using PostGIS to access raster data is
>> 1st. to convert an ASCII GRID text file (an ESRI format) to GeoTIFF
>> (not tiled nor "pyramide"-ized).
>> Then (2nd) to load it into PostGIS, like described here
>> http://trac.osgeo.org/postgis/wiki/WKTRasterTutorial01 .
>>
>
> I don't remember but I think GDAL supports ASCII GRID so you don't
> need an intermediary file.  The loader uses GDAL to convert them into
> the PostGIS Raster format.
>
>> I still have some basic questions about this interesting approach:
>>
>> Regarding my question of "rules of thumb" in order to set the tile
>> size you answered:
>>> Generally speaking, the smaller the tile size, the faster the operation.
>  In
>>> my testing of my rasters, 200x200 is slower than 100x100 which is slower
>>> than 50x50.  The smaller tile size results in more rows in the table BUT
> the
>>> performance dramatically increases.  I can't remember of the top of my
> head
>>> but the performance gain was linear or better.
>>
>> Thats a good hint.
>> Q1: But are there no calculcation rules on how to estimate an ideal
>> absolute size (in pixels)?
>>
>
> Nope.  Testing is the only way as there are other possible factors
> involved.  And even my estimates above are only guidelines.  I know
> I'm going to have to test for my production database servers before I
> know what works best for them.  In a pinch, I'd say you'd be happy
> between 30 - 70.
>
>> Regarding a query given lat/lon you proposed:
>>> Assuming you test to see what the ideal (or "good-enough") tile size
>>> is for the elevation rasterset, you should see good performance using
>>>
>>> ST_Value(rast, ST_SetSRID(ST_MakePoint(long, lat), srid))
>>>
>>> You could eliminate the ST_SetSRID and ST_MakePoint by passing a
>>> correct EWKT instead
>>>
>>> ST_Value(rast, 'YOUREWKT'::geometry)
>>
>> Q2: Given the raster has a resolution (grid) of - say at best - 30m,
>> does this function really interpolate, given the input is an arbitrary
>> coordinate (I can't find a hint to that in the docs
>>
> (http://postgis.refractions.net/documentation/manual-svn/RT_ST_Value.html))?
>>
>
> There is no interpolation.  It literally is a pin stuck through the raster.
>
>> Q3: If yes, I assume its linear. How could I proceed to get a result
>> using non-linear interpolation (like Kriging)?
>>
>
> I don't know if you can using only PostGIS raster.  If you're able to
> use PL/R as well, you should be able to tap into their kriging tools.
>
> -bborie
>
> --
> Bborie Park
> Programmer
> Center for Vectorborne Diseases
> UC Davis
> 530-752-8380
> bkpark at ucdavis.edu
> _______________________________________________
> postgis-users mailing list
> postgis-users at postgis.refractions.net
> http://postgis.refractions.net/mailman/listinfo/postgis-users
>
> _______________________________________________
> postgis-users mailing list
> postgis-users at postgis.refractions.net
> http://postgis.refractions.net/mailman/listinfo/postgis-users