PostGIS  2.2.8dev-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 2165 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().

2165  {
2166  rt_raster rast = NULL;
2167  double gt[6] = {0};
2168  CPLErr cplerr;
2169  uint32_t width = 0;
2170  uint32_t height = 0;
2171  uint32_t numBands = 0;
2172  int i = 0;
2173  char *authname = NULL;
2174  char *authcode = NULL;
2175 
2176  GDALRasterBandH gdband = NULL;
2177  GDALDataType gdpixtype = GDT_Unknown;
2178  rt_band band;
2179  int32_t idx;
2180  rt_pixtype pt = PT_END;
2181  uint32_t ptlen = 0;
2182  int hasnodata = 0;
2183  double nodataval;
2184 
2185  int x;
2186  int y;
2187 
2188  int nXBlocks, nYBlocks;
2189  int nXBlockSize, nYBlockSize;
2190  int iXBlock, iYBlock;
2191  int nXValid, nYValid;
2192  int iY;
2193 
2194  uint8_t *values = NULL;
2195  uint32_t valueslen = 0;
2196  uint8_t *ptr = NULL;
2197 
2198  assert(NULL != ds);
2199 
2200  /* raster size */
2201  width = GDALGetRasterXSize(ds);
2202  height = GDALGetRasterYSize(ds);
2203  RASTER_DEBUGF(3, "Raster dimensions (width x height): %d x %d", width, height);
2204 
2205  /* create new raster */
2206  RASTER_DEBUG(3, "Creating new raster");
2207  rast = rt_raster_new(width, height);
2208  if (NULL == rast) {
2209  rterror("rt_raster_from_gdal_dataset: Out of memory allocating new raster");
2210  return NULL;
2211  }
2212  RASTER_DEBUGF(3, "Created raster dimensions (width x height): %d x %d", rast->width, rast->height);
2213 
2214  /* get raster attributes */
2215  cplerr = GDALGetGeoTransform(ds, gt);
2216  if (GDALGetGeoTransform(ds, gt) != CE_None) {
2217  RASTER_DEBUG(4, "Using default geotransform matrix (0, 1, 0, 0, 0, -1)");
2218  gt[0] = 0;
2219  gt[1] = 1;
2220  gt[2] = 0;
2221  gt[3] = 0;
2222  gt[4] = 0;
2223  gt[5] = -1;
2224  }
2225 
2226  /* apply raster attributes */
2228 
2229  RASTER_DEBUGF(3, "Raster geotransform (%f, %f, %f, %f, %f, %f)",
2230  gt[0], gt[1], gt[2], gt[3], gt[4], gt[5]);
2231 
2232  /* srid */
2233  if (rt_util_gdal_sr_auth_info(ds, &authname, &authcode) == ES_NONE) {
2234  if (
2235  authname != NULL &&
2236  strcmp(authname, "EPSG") == 0 &&
2237  authcode != NULL
2238  ) {
2239  rt_raster_set_srid(rast, atoi(authcode));
2240  RASTER_DEBUGF(3, "New raster's SRID = %d", rast->srid);
2241  }
2242 
2243  if (authname != NULL)
2244  rtdealloc(authname);
2245  if (authcode != NULL)
2246  rtdealloc(authcode);
2247  }
2248 
2249  numBands = GDALGetRasterCount(ds);
2250 
2251 #if POSTGIS_DEBUG_LEVEL > 3
2252  for (i = 1; i <= numBands; i++) {
2253  GDALRasterBandH _grb = NULL;
2254  double _min;
2255  double _max;
2256  double _mean;
2257  double _stddev;
2258 
2259  _grb = GDALGetRasterBand(ds, i);
2260  GDALComputeRasterStatistics(_grb, FALSE, &_min, &_max, &_mean, &_stddev, NULL, NULL);
2261  RASTER_DEBUGF(4, "GDAL Band %d stats: %f, %f, %f, %f", i, _min, _max, _mean, _stddev);
2262  }
2263 #endif
2264 
2265  /* copy bands */
2266  for (i = 1; i <= numBands; i++) {
2267  RASTER_DEBUGF(3, "Processing band %d of %d", i, numBands);
2268  gdband = NULL;
2269  gdband = GDALGetRasterBand(ds, i);
2270  if (NULL == gdband) {
2271  rterror("rt_raster_from_gdal_dataset: Could not get GDAL band");
2272  rt_raster_destroy(rast);
2273  return NULL;
2274  }
2275  RASTER_DEBUGF(4, "gdband @ %p", gdband);
2276 
2277  /* pixtype */
2278  gdpixtype = GDALGetRasterDataType(gdband);
2279  RASTER_DEBUGF(4, "gdpixtype, size = %s, %d", GDALGetDataTypeName(gdpixtype), GDALGetDataTypeSize(gdpixtype) / 8);
2280  pt = rt_util_gdal_datatype_to_pixtype(gdpixtype);
2281  if (pt == PT_END) {
2282  rterror("rt_raster_from_gdal_dataset: Unknown pixel type for GDAL band");
2283  rt_raster_destroy(rast);
2284  return NULL;
2285  }
2286  ptlen = rt_pixtype_size(pt);
2287 
2288  /* size: width and height */
2289  width = GDALGetRasterBandXSize(gdband);
2290  height = GDALGetRasterBandYSize(gdband);
2291  RASTER_DEBUGF(3, "GDAL band dimensions (width x height): %d x %d", width, height);
2292 
2293  /* nodata */
2294  nodataval = GDALGetRasterNoDataValue(gdband, &hasnodata);
2295  RASTER_DEBUGF(3, "(hasnodata, nodataval) = (%d, %f)", hasnodata, nodataval);
2296 
2297  /* create band object */
2299  rast, pt,
2300  (hasnodata ? nodataval : 0),
2301  hasnodata, nodataval, rt_raster_get_num_bands(rast)
2302  );
2303  if (idx < 0) {
2304  rterror("rt_raster_from_gdal_dataset: Could not allocate memory for raster band");
2305  rt_raster_destroy(rast);
2306  return NULL;
2307  }
2308  band = rt_raster_get_band(rast, idx);
2309  RASTER_DEBUGF(3, "Created band of dimension (width x height): %d x %d", band->width, band->height);
2310 
2311  /* this makes use of GDAL's "natural" blocks */
2312  GDALGetBlockSize(gdband, &nXBlockSize, &nYBlockSize);
2313  nXBlocks = (width + nXBlockSize - 1) / nXBlockSize;
2314  nYBlocks = (height + nYBlockSize - 1) / nYBlockSize;
2315  RASTER_DEBUGF(4, "(nXBlockSize, nYBlockSize) = (%d, %d)", nXBlockSize, nYBlockSize);
2316  RASTER_DEBUGF(4, "(nXBlocks, nYBlocks) = (%d, %d)", nXBlocks, nYBlocks);
2317 
2318  /* allocate memory for values */
2319  valueslen = ptlen * nXBlockSize * nYBlockSize;
2320  values = rtalloc(valueslen);
2321  if (values == NULL) {
2322  rterror("rt_raster_from_gdal_dataset: Could not allocate memory for GDAL band pixel values");
2323  rt_raster_destroy(rast);
2324  return NULL;
2325  }
2326  RASTER_DEBUGF(3, "values @ %p of length = %d", values, valueslen);
2327 
2328  for (iYBlock = 0; iYBlock < nYBlocks; iYBlock++) {
2329  for (iXBlock = 0; iXBlock < nXBlocks; iXBlock++) {
2330  x = iXBlock * nXBlockSize;
2331  y = iYBlock * nYBlockSize;
2332  RASTER_DEBUGF(4, "(iXBlock, iYBlock) = (%d, %d)", iXBlock, iYBlock);
2333  RASTER_DEBUGF(4, "(x, y) = (%d, %d)", x, y);
2334 
2335  memset(values, 0, valueslen);
2336 
2337  /* valid block width */
2338  if ((iXBlock + 1) * nXBlockSize > width)
2339  nXValid = width - (iXBlock * nXBlockSize);
2340  else
2341  nXValid = nXBlockSize;
2342 
2343  /* valid block height */
2344  if ((iYBlock + 1) * nYBlockSize > height)
2345  nYValid = height - (iYBlock * nYBlockSize);
2346  else
2347  nYValid = nYBlockSize;
2348 
2349  RASTER_DEBUGF(4, "(nXValid, nYValid) = (%d, %d)", nXValid, nYValid);
2350 
2351  cplerr = GDALRasterIO(
2352  gdband, GF_Read,
2353  x, y,
2354  nXValid, nYValid,
2355  values, nXValid, nYValid,
2356  gdpixtype,
2357  0, 0
2358  );
2359  if (cplerr != CE_None) {
2360  rterror("rt_raster_from_gdal_dataset: Could not get data from GDAL raster");
2361  rtdealloc(values);
2362  rt_raster_destroy(rast);
2363  return NULL;
2364  }
2365 
2366  /* if block width is same as raster width, shortcut */
2367  if (nXBlocks == 1 && nYBlockSize > 1 && nXValid == width) {
2368  x = 0;
2369  y = nYBlockSize * iYBlock;
2370 
2371  RASTER_DEBUGF(4, "Setting set of pixel lines at (%d, %d) for %d pixels", x, y, nXValid * nYValid);
2372  rt_band_set_pixel_line(band, x, y, values, nXValid * nYValid);
2373  }
2374  else {
2375  ptr = values;
2376  x = nXBlockSize * iXBlock;
2377  for (iY = 0; iY < nYValid; iY++) {
2378  y = iY + (nYBlockSize * iYBlock);
2379 
2380  RASTER_DEBUGF(4, "Setting pixel line at (%d, %d) for %d pixels", x, y, nXValid);
2381  rt_band_set_pixel_line(band, x, y, ptr, nXValid);
2382  ptr += (nXValid * ptlen);
2383  }
2384  }
2385  }
2386  }
2387 
2388  /* free memory */
2389  rtdealloc(values);
2390  }
2391 
2392  return rast;
2393 }
int32_t srid
Definition: librtcore.h:2263
band
Definition: ovdump.py:57
int rt_raster_get_num_bands(rt_raster raster)
Definition: rt_raster.c:372
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:2280
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:197
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:2265
ds
Definition: pixval.py:66
uint16_t width
Definition: librtcore.h:2264
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:311
uint16_t width
Definition: librtcore.h:2279
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:307
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:720
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
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: