Hello, I've created a WKTRasterBand-to-GDALRasterBand conversor, in the frame of WKT Raster project. I'm not sure if the results are correct. You can get the code from here https://svn.osgeo.org/postgis/spike/wktraster/rt_core/rt_api.c, but basically, the steps are:
1.- Register GDAL/OGR drivers 2.- Create empty OGR in-memory Datasource: OGR_Dr_CreateDataSource(OGRGetDriverByName("Memory"),"", NULL); 3.- Create empty GDAL Dataset, taking width and height from WKT Raster Band ("band" variable): GDALCreate(GDALGetDriverByName("MEM"), "",rt_band_get_width(ctx, band), rt_band_get_height(ctx, band), 0, GDT_Byte, NULL); 4.- Add default geotransform to GDAL empty dataset: {0.0, 1.0, 0.0, 0.0, 0.0,1.0 }; (I can't get geotransform from WKT Raster Band) 5.- Add GDAL band, getting data from WKT Raster Band: sprintf(szDataPointer, "DATAPOINTER=%p", rt_band_get_data(ctx, band)); apszOptions[0] = szDataPointer; GDALAddBand(memdataset, nPixelType, apszOptions) (pixel type got from WKT Raster Band too) 6.- Create empty OGR layer, to store polygonized raster: OGR_DS_CreateLayer(memdatasource, "Polygonized layer", NULL, wkbPolygon, NULL); 7.- Get GDAL RasterBand from GDAL Dataset, and polygonize it gdal_band = GDALGetRasterBand(memdataset, 1); GDALPolygonize(gdal_band, NULL, hLayer, -1, NULL, NULL, NULL); 8.- Get all the features from layer and export them to WKT format: for(j = 0; j < nFeatureCount; j++) { hFeature = OGR_L_GetFeature(hLayer, j); hGeom = OGR_F_GetGeometryRef(hFeature); OGR_G_ExportToWkt(hGeom, &pszSrcText); } So, at the end, pszSrcText has each of the WKT polygons from original raster. My problems are: 1.- I'm testing this code with really simple raster, and I'm confused. For example: * 2x2 32 bits signed int raster, with "1" in all cells. The result: POLYGON((0 0, 0 2, 2 2, 2 0, 0 0)). I was expecting POLYGON((0 0, 0 1, 1 1, 1 0, 0 0)) * the raster from http://trac.osgeo.org/postgis/attachment/wiki/WKTRaster/SpecificationWorking01/WKTRasterEnvelopeConvexHullAndShape.gif. The result: POLYGON ((3 1,3 2,2 2,2 3,1 3,1 6,2 6,2 7,3 7,3 8,5 8,5 6,3 6,3 3,4 3,5 3,5 1,3 1)) POLYGON ((3 3,3 6,6 6,6 3,3 3)) POLYGON ((5 1,5 3,6 3,6 6,5 6,5 8,6 8,6 7,7 7,7 6,8 6,8 3,7 3,7 2,6 2,6 1,5 1)) POLYGON ((0 0,0 9,9 9,9 0,0 0),(6 7,6 8,3 8,3 7,2 7,2 6,1 6,1 3,2 3,2 2,3 2,3 1,6 1,6 2,7 2,7 3,8 3,8 6,7 6,7 7,6 7)) Something strange here? Am I doing things right? 2.- I'd like to store the pixel value in each polygon. Is a good idea to create a field in the layer to store pixel value? 3.- From a GDAL RasterBand, I can fetch the owning dataset handle (GDALGetBandDataset).I think it could be useful to get the WKT Raster handle from a WKT RasterBand too... Am I right? Sorry, I know this post is really long. Thanks for reading it, anyway :-) Best regards, Jorge http://www.gis4free.org/blog
_______________________________________________ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev