Oups, I figured it would be of no use. I of course works on pure geometry, raster functions only sees postgis points and ints. Here few lines of the point table
qgis_id, st_astext(geom) AS geom,z,reflectance 1 POINT Z (2248.03170000003 20680.0229999999 49.02441796875) 49024 8110 2 POINT Z (2247.76560000001 20680.2519999999 49.0198828125) 49020 7890 3 POINT Z (2248.25289999996 20680.2479999999 49.02616015625) 49026 8404 Is there any function to work on pixel sets , (like read/write many at a time) Now the code : (I can make a test case out of this if you want) Thanks, Rémi-C /////////////////////////////// > > SELECT postgis_full_version(); > --POSTGIS="2.1.0 r11822" GEOS="3.5.0dev-CAPI-1.9.0 r3963" PROJ="Rel. > 4.8.0, 6 March 2012" GDAL="GDAL 1.10.0, released 2013/04/24" LIBXML="2.8.0" > RASTER > > > > > ----converting a patch to raster > --get patch extend > --create a raster with this extend > > --explode patch to points > --set pixel value to point intersecting it. > > > --preparing raster : > --creating the raster table > DROP TABLE IF EXISTS patch_to_raster; > CREATE TABLE patch_to_raster(rid serial primary key > , rast raster); > > -- > --DELETE from patch_to_raster > --inserting an empty raster > INSERT INTO patch_to_raster (rast) VALUES ( > ST_MakeEmptyRaster( > 50 > ,50 > ,upperleftx:=0 > ,upperlefty:=0 > ,scalex:=1 > , scaley:=1 > , skewx:=0 > ,skewy:=0 > , srid:=932011 > ) --srid of translated lambert 93 to match laser referential > ); > --setting the correct pixel size > UPDATE patch_to_raster SET rast = ST_SetScale(rast,0.02,0.02); > > --checking the content: > SELECT ST_Summary(rast) > FROM patch_to_raster; > --Raster of 50x50 pixels has 0 bands and extent of BOX(0 0,1 1) > > > --adding band > --first band is Z : float > --second band is reflectance : float > UPDATE patch_to_raster SET rast = ST_AddBand( --1'st band, Z > rast > ,pixeltype:='32BF' -- '32BUI' > , initialvalue:=NULL > , nodataval:=NULL > ); > UPDATE patch_to_raster SET rast = ST_AddBand( --2nd band, reflectance > rast > , pixeltype:='32BF' --'32BUI' > , initialvalue:=NULL > , nodataval:=NULL > ); > > > --checking the content: > SELECT ST_Summary(rast) FROM patch_to_raster; > --Raster of 50x50 pixels has 2 bands and extent of BOX(0 0,1 1) > --band 1 of pixtype 32BF is in-db with no NODATA value > --band 2 of pixtype 32BF is in-db with no NODATA value > > --getting a patch and creating a table with it > DROP TABLE IF EXISTS one_patch_into_points; > CREATE TABLE one_patch_into_points AS > WITH patch AS ( > SELECT gid,patch > FROM acquisition_tmob_012013.riegl_pcpatch_space > WHERE gid = 87263 > ) > ,point AS ( > SELECT PC_Explode(patch) AS point > FROM patch > ) > SELECT row_number() over() AS qgis_id, ST_SetSRID(point::geometry,932011) > AS geom, @(PC_Get(point, 'z')*1000)::int AS z, @((PC_Get(point, > 'reflectance')+50)*200)::int aS reflectance > FROM point; > --2746 points > > --checking > SELECT * > FROM one_patch_into_points > --1 01010000A0AB380E00DC7EFB3A1090A140CDFDD4780132D440A69BC42020834840 > 49024 8110 > --2 POINT Z (2247.76560000001 20680.2519999999 49.0198828125) 49020 7890 > > --updating the raster to set it on the patch place. > --we get the upper left coordinate of patch > WITH patch AS ( > SELECT gid,patch, ST_AsText(patch::geometry) > FROM acquisition_tmob_012013.riegl_pcpatch_space > WHERE gid = 87263 > ) > ,centroid AS ( > SELECT ST_Centroid(patch::geometry) AS centroid > FROM patch > ) > ,upperleft AS ( > SELECT round(ST_X(centroid)) -0.5 AS upper_left_x, > round(ST_y(centroid))-0.5 AS upper_left_y > FROM centroid > )--updating the raster with it > UPDATE patch_to_raster SET rast = > ST_SetUpperLeft(rast,upper_left_x,upper_left_y) > FROM upperleft; > > --checking the raster > SELECT ST_Summary(rast) FROM patch_to_raster; > --Raster of 50x50 pixels has 2 bands and extent of BOX(2247.5 > 20679.5,2248.5 20680.5) > --note : the patch bbox : POLYGON((2247.51 20679.50,2247.51 > 20680.49,2248.49 20680.49,2248.49 20679.50,2247.51 20679.50)), some pixels > will be empty > > > > --updating the raster with pixel value > --using postgis_addons function : > https://raw.github.com/pedrogit/postgisaddons/master/postgis_addons.sql > -- > UPDATE patch_to_raster SET rast = > ST_ExtractToRaster( > rast > ,band:=1 > ,schemaname:= 'test_raster' > ,tablename:= 'one_patch_into_points' > ,geomrastcolumnname:= 'geom' > ,valuecolumnname:='z' > ); > --5 sec : way too long! > UPDATE patch_to_raster SET rast = > ST_ExtractToRaster( > rast > ,band:=2 > ,schemaname:= 'test_raster' > ,tablename:= 'one_patch_into_points' > ,geomrastcolumnname:= 'geom' > ,valuecolumnname:='reflectance' > ,method:='MEAN_OF_VALUES_AT_PIXEL_CENTROID' > ); > --6 sec, way too long ! > > > --checking the value of pixel in both band (no shortcut to get several > bands at a time??) > SELECT s1 as x, s2 as y, ST_Value(rast, 1,s1, s2, > exclude_nodata_value:=true) AS value_band_1, ST_Value(rast, 2,s1, s2, > exclude_nodata_value:=true) AS value_band_2 > FROM generate_series(1,50) as s1, generate_series(1,50) as s2, > patch_to_raster; > --values are NULL! > > --checking the histogram > SELECT ST_Histogram( rast, 1 ) --, boolean exclude_nodata_value=true > FROM patch_to_raster; > --no line returned > > --checking result > SELECT ST_Summary(rast) FROM patch_to_raster; > --Raster of 50x50 pixels has 2 bands and extent of BOX(2247.5 > 20679.5,2248.5 20680.5) > -- band 1 of pixtype 32BF is in-db with NODATA value of > -3.40282346638529e+38 > -- band 2 of pixtype 32BF is in-db with NODATA value of > -3.40282346638529e+38 > //////////////////////////////// >
_______________________________________________ postgis-users mailing list postgis-users@lists.osgeo.org http://lists.osgeo.org/cgi-bin/mailman/listinfo/postgis-users