PostGIS  2.5.0dev-r@@SVN_REVISION@@

◆ rt_raster_from_gdal_dataset()

rt_raster rt_raster_from_gdal_dataset ( GDALDatasetH  ds)

Return a raster from a GDAL dataset.

Parameters
ds: the GDAL dataset to convert to a raster
Returns
raster or NULL

Definition at line 2186 of file rt_raster.c.

References ovdump::band, ES_NONE, FALSE, window::gt, rt_raster_t::height, rt_band_t::height, PT_END, rtpixdump::rast, RASTER_DEBUG, RASTER_DEBUGF, rt_band_set_pixel_line(), rt_pixtype_size(), rt_raster_destroy(), rt_raster_generate_new_band(), rt_raster_get_band(), rt_raster_get_num_bands(), rt_raster_new(), rt_raster_set_geotransform_matrix(), rt_raster_set_srid(), rt_util_gdal_datatype_to_pixtype(), rt_util_gdal_sr_auth_info(), rtalloc(), rtdealloc(), rterror(), rt_raster_t::srid, rt_raster_t::width, rt_band_t::width, pixval::x, and pixval::y.

Referenced by build_overview(), convert_raster(), RASTER_fromGDALRaster(), rt_band_load_offline_data(), rt_raster_gdal_rasterize(), rt_raster_gdal_warp(), and test_gdal_to_raster().

2186  {
2187  rt_raster rast = NULL;
2188  double gt[6] = {0};
2189  CPLErr cplerr;
2190  uint32_t width = 0;
2191  uint32_t height = 0;
2192  uint32_t numBands = 0;
2193  uint32_t i = 0;
2194  char *authname = NULL;
2195  char *authcode = NULL;
2196 
2197  GDALRasterBandH gdband = NULL;
2198  GDALDataType gdpixtype = GDT_Unknown;
2199  rt_band band;
2200  int32_t idx;
2201  rt_pixtype pt = PT_END;
2202  uint32_t ptlen = 0;
2203  int hasnodata = 0;
2204  double nodataval;
2205 
2206  int x;
2207  int y;
2208 
2209  uint32_t nXBlocks, nYBlocks;
2210  int nXBlockSize, nYBlockSize;
2211  uint32_t iXBlock, iYBlock;
2212  uint32_t nXValid, nYValid;
2213  uint32_t iY;
2214 
2215  uint8_t *values = NULL;
2216  uint32_t valueslen = 0;
2217  uint8_t *ptr = NULL;
2218 
2219  assert(NULL != ds);
2220 
2221  /* raster size */
2222  width = GDALGetRasterXSize(ds);
2223  height = GDALGetRasterYSize(ds);
2224  RASTER_DEBUGF(3, "Raster dimensions (width x height): %d x %d", width, height);
2225 
2226  /* create new raster */
2227  RASTER_DEBUG(3, "Creating new raster");
2228  rast = rt_raster_new(width, height);
2229  if (NULL == rast) {
2230  rterror("rt_raster_from_gdal_dataset: Out of memory allocating new raster");
2231  return NULL;
2232  }
2233  RASTER_DEBUGF(3, "Created raster dimensions (width x height): %d x %d", rast->width, rast->height);
2234 
2235  /* get raster attributes */
2236  cplerr = GDALGetGeoTransform(ds, gt);
2237  if (GDALGetGeoTransform(ds, gt) != CE_None) {
2238  RASTER_DEBUG(4, "Using default geotransform matrix (0, 1, 0, 0, 0, -1)");
2239  gt[0] = 0;
2240  gt[1] = 1;
2241  gt[2] = 0;
2242  gt[3] = 0;
2243  gt[4] = 0;
2244  gt[5] = -1;
2245  }
2246 
2247  /* apply raster attributes */
2249 
2250  RASTER_DEBUGF(3, "Raster geotransform (%f, %f, %f, %f, %f, %f)",
2251  gt[0], gt[1], gt[2], gt[3], gt[4], gt[5]);
2252 
2253  /* srid */
2254  if (rt_util_gdal_sr_auth_info(ds, &authname, &authcode) == ES_NONE) {
2255  if (
2256  authname != NULL &&
2257  strcmp(authname, "EPSG") == 0 &&
2258  authcode != NULL
2259  ) {
2260  rt_raster_set_srid(rast, atoi(authcode));
2261  RASTER_DEBUGF(3, "New raster's SRID = %d", rast->srid);
2262  }
2263 
2264  if (authname != NULL)
2265  rtdealloc(authname);
2266  if (authcode != NULL)
2267  rtdealloc(authcode);
2268  }
2269 
2270  numBands = GDALGetRasterCount(ds);
2271 
2272 #if POSTGIS_DEBUG_LEVEL > 3
2273  for (i = 1; i <= numBands; i++) {
2274  GDALRasterBandH _grb = NULL;
2275  double _min;
2276  double _max;
2277  double _mean;
2278  double _stddev;
2279 
2280  _grb = GDALGetRasterBand(ds, i);
2281  GDALComputeRasterStatistics(_grb, FALSE, &_min, &_max, &_mean, &_stddev, NULL, NULL);
2282  RASTER_DEBUGF(4, "GDAL Band %d stats: %f, %f, %f, %f", i, _min, _max, _mean, _stddev);
2283  }
2284 #endif
2285 
2286  /* copy bands */
2287  for (i = 1; i <= numBands; i++) {
2288  RASTER_DEBUGF(3, "Processing band %d of %d", i, numBands);
2289  gdband = NULL;
2290  gdband = GDALGetRasterBand(ds, i);
2291  if (NULL == gdband) {
2292  rterror("rt_raster_from_gdal_dataset: Could not get GDAL band");
2293  rt_raster_destroy(rast);
2294  return NULL;
2295  }
2296  RASTER_DEBUGF(4, "gdband @ %p", gdband);
2297 
2298  /* pixtype */
2299  gdpixtype = GDALGetRasterDataType(gdband);
2300  RASTER_DEBUGF(4, "gdpixtype, size = %s, %d", GDALGetDataTypeName(gdpixtype), GDALGetDataTypeSize(gdpixtype) / 8);
2301  pt = rt_util_gdal_datatype_to_pixtype(gdpixtype);
2302  if (pt == PT_END) {
2303  rterror("rt_raster_from_gdal_dataset: Unknown pixel type for GDAL band");
2304  rt_raster_destroy(rast);
2305  return NULL;
2306  }
2307  ptlen = rt_pixtype_size(pt);
2308 
2309  /* size: width and height */
2310  width = GDALGetRasterBandXSize(gdband);
2311  height = GDALGetRasterBandYSize(gdband);
2312  RASTER_DEBUGF(3, "GDAL band dimensions (width x height): %d x %d", width, height);
2313 
2314  /* nodata */
2315  nodataval = GDALGetRasterNoDataValue(gdband, &hasnodata);
2316  RASTER_DEBUGF(3, "(hasnodata, nodataval) = (%d, %f)", hasnodata, nodataval);
2317 
2318  /* create band object */
2320  rast, pt,
2321  (hasnodata ? nodataval : 0),
2322  hasnodata, nodataval, rt_raster_get_num_bands(rast)
2323  );
2324  if (idx < 0) {
2325  rterror("rt_raster_from_gdal_dataset: Could not allocate memory for raster band");
2326  rt_raster_destroy(rast);
2327  return NULL;
2328  }
2329  band = rt_raster_get_band(rast, idx);
2330  RASTER_DEBUGF(3, "Created band of dimension (width x height): %d x %d", band->width, band->height);
2331 
2332  /* this makes use of GDAL's "natural" blocks */
2333  GDALGetBlockSize(gdband, &nXBlockSize, &nYBlockSize);
2334  nXBlocks = (width + nXBlockSize - 1) / nXBlockSize;
2335  nYBlocks = (height + nYBlockSize - 1) / nYBlockSize;
2336  RASTER_DEBUGF(4, "(nXBlockSize, nYBlockSize) = (%d, %d)", nXBlockSize, nYBlockSize);
2337  RASTER_DEBUGF(4, "(nXBlocks, nYBlocks) = (%d, %d)", nXBlocks, nYBlocks);
2338 
2339  /* allocate memory for values */
2340  valueslen = ptlen * nXBlockSize * nYBlockSize;
2341  values = rtalloc(valueslen);
2342  if (values == NULL) {
2343  rterror("rt_raster_from_gdal_dataset: Could not allocate memory for GDAL band pixel values");
2344  rt_raster_destroy(rast);
2345  return NULL;
2346  }
2347  RASTER_DEBUGF(3, "values @ %p of length = %d", values, valueslen);
2348 
2349  for (iYBlock = 0; iYBlock < nYBlocks; iYBlock++) {
2350  for (iXBlock = 0; iXBlock < nXBlocks; iXBlock++) {
2351  x = iXBlock * nXBlockSize;
2352  y = iYBlock * nYBlockSize;
2353  RASTER_DEBUGF(4, "(iXBlock, iYBlock) = (%d, %d)", iXBlock, iYBlock);
2354  RASTER_DEBUGF(4, "(x, y) = (%d, %d)", x, y);
2355 
2356  memset(values, 0, valueslen);
2357 
2358  /* valid block width */
2359  if ((iXBlock + 1) * nXBlockSize > width)
2360  nXValid = width - (iXBlock * nXBlockSize);
2361  else
2362  nXValid = nXBlockSize;
2363 
2364  /* valid block height */
2365  if ((iYBlock + 1) * nYBlockSize > height)
2366  nYValid = height - (iYBlock * nYBlockSize);
2367  else
2368  nYValid = nYBlockSize;
2369 
2370  RASTER_DEBUGF(4, "(nXValid, nYValid) = (%d, %d)", nXValid, nYValid);
2371 
2372  cplerr = GDALRasterIO(
2373  gdband, GF_Read,
2374  x, y,
2375  nXValid, nYValid,
2376  values, nXValid, nYValid,
2377  gdpixtype,
2378  0, 0
2379  );
2380  if (cplerr != CE_None) {
2381  rterror("rt_raster_from_gdal_dataset: Could not get data from GDAL raster");
2382  rtdealloc(values);
2383  rt_raster_destroy(rast);
2384  return NULL;
2385  }
2386 
2387  /* if block width is same as raster width, shortcut */
2388  if (nXBlocks == 1 && nYBlockSize > 1 && nXValid == width) {
2389  x = 0;
2390  y = nYBlockSize * iYBlock;
2391 
2392  RASTER_DEBUGF(4, "Setting set of pixel lines at (%d, %d) for %d pixels", x, y, nXValid * nYValid);
2393  rt_band_set_pixel_line(band, x, y, values, nXValid * nYValid);
2394  }
2395  else {
2396  ptr = values;
2397  x = nXBlockSize * iXBlock;
2398  for (iY = 0; iY < nYValid; iY++) {
2399  y = iY + (nYBlockSize * iYBlock);
2400 
2401  RASTER_DEBUGF(4, "Setting pixel line at (%d, %d) for %d pixels", x, y, nXValid);
2402  rt_band_set_pixel_line(band, x, y, ptr, nXValid);
2403  ptr += (nXValid * ptlen);
2404  }
2405  }
2406  }
2407  }
2408 
2409  /* free memory */
2410  rtdealloc(values);
2411  }
2412 
2413  return rast;
2414 }
int32_t srid
Definition: librtcore.h:2280
band
Definition: ovdump.py:57
void rterror(const char *fmt,...)
Wrappers used for reporting errors and info.
Definition: rt_context.c:199
void * rtalloc(size_t size)
Wrappers used for managing memory.
Definition: rt_context.c:171
uint16_t height
Definition: librtcore.h:2297
void rt_raster_set_srid(rt_raster raster, int32_t srid)
Set raster&#39;s SRID.
Definition: rt_raster.c:363
gt
Definition: window.py:77
rt_pixtype
Definition: librtcore.h:185
uint16_t rt_raster_get_num_bands(rt_raster raster)
Definition: rt_raster.c:372
void rt_raster_set_geotransform_matrix(rt_raster raster, double *gt)
Set raster&#39;s geotransform using 6-element array.
Definition: rt_raster.c:727
uint16_t height
Definition: librtcore.h:2282
ds
Definition: pixval.py:66
unsigned int uint32_t
Definition: uthash.h:78
uint16_t width
Definition: librtcore.h:2281
rt_pixtype rt_util_gdal_datatype_to_pixtype(GDALDataType gdt)
Convert GDALDataType to rt_pixtype.
Definition: rt_util.c:153
#define RASTER_DEBUGF(level, msg,...)
Definition: librtcore.h:299
uint16_t width
Definition: librtcore.h:2296
int rt_pixtype_size(rt_pixtype pixtype)
Return size in bytes of a value in the given pixtype.
Definition: rt_pixel.c:39
rt_raster rt_raster_new(uint32_t width, uint32_t height)
Construct a raster with given dimensions.
Definition: rt_raster.c:48
void rtdealloc(void *mem)
Definition: rt_context.c:186
#define FALSE
Definition: dbfopen.c:168
#define RASTER_DEBUG(level, msg)
Definition: librtcore.h:295
rt_errorstate rt_band_set_pixel_line(rt_band band, int x, int y, void *vals, uint32_t len)
Set values of multiple pixels.
Definition: rt_band.c:805
int rt_raster_generate_new_band(rt_raster raster, rt_pixtype pixtype, double initialvalue, uint32_t hasnodata, double nodatavalue, int index)
Generate a new inline band and add it to a raster.
Definition: rt_raster.c:485
rt_errorstate rt_util_gdal_sr_auth_info(GDALDatasetH hds, char **authname, char **authcode)
Get auth name and code.
Definition: rt_util.c:270
unsigned char uint8_t
Definition: uthash.h:79
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_raster.c:82
rt_band rt_raster_get_band(rt_raster raster, int n)
Return Nth band, or NULL if unavailable.
Definition: rt_raster.c:381
Here is the call graph for this function:
Here is the caller graph for this function: