PostGIS  2.5.0rc1-r@@SVN_REVISION@@

◆ build_overview()

static int build_overview ( int  idx,
RTLOADERCFG config,
RASTERINFO info,
uint32_t  ovx,
STRINGBUFFER tileset,
STRINGBUFFER buffer 
)
static

Definition at line 1353 of file raster2pgsql.c.

References _, append_stringbuffer(), raster_loader_config::copy_statements, rasterinfo_t::dim, FALSE, raster_loader_config::file_column, raster_loader_config::file_column_name, rasterinfo_t::gdalbandtype, window::gt, rasterinfo_t::gt, rasterinfo_t::hasnodata, if(), insert_records(), stringbuffer_t::length, MAXOVFACTOR, MINOVFACTOR, rasterinfo_t::nband, rasterinfo_t::nband_count, rasterinfo_t::nodataval, raster_loader_config::out_srid, raster_loader_config::overview, raster_loader_config::overview_count, raster_loader_config::overview_table, raster_loader_config::pad_tile, rtpixdump::rast, raster_loader_config::raster_column, raster_destroy(), raster_loader_config::rt_file, raster_loader_config::rt_filename, rt_raster_from_gdal_dataset(), rt_raster_set_srid(), rt_raster_to_hexwkb(), rtdealloc_stringbuffer(), rterror(), raster_loader_config::schema, rasterinfo_t::srid, rasterinfo_t::srs, and raster_loader_config::tile_size.

Referenced by process_rasters().

1353  {
1354  GDALDatasetH hdsSrc;
1355  VRTDatasetH hdsOv;
1356  VRTSourcedRasterBandH hbandOv;
1357  double gtOv[6] = {0.};
1358  int dimOv[2] = {0};
1359 
1360  uint32_t j = 0;
1361  int factor;
1362  const char *ovtable = NULL;
1363 
1364  VRTDatasetH hdsDst;
1365  VRTSourcedRasterBandH hbandDst;
1366  int tile_size[2] = {0};
1367  int _tile_size[2] = {0};
1368  int ntiles[2] = {1, 1};
1369  int xtile = 0;
1370  int ytile = 0;
1371  double gt[6] = {0.};
1372 
1373  rt_raster rast = NULL;
1374  char *hex;
1375  uint32_t hexlen = 0;
1376 
1377  hdsSrc = GDALOpenShared(config->rt_file[idx], GA_ReadOnly);
1378  if (hdsSrc == NULL) {
1379  rterror(_("build_overview: Could not open raster: %s"), config->rt_file[idx]);
1380  return 0;
1381  }
1382 
1383  /* working copy of geotransform matrix */
1384  memcpy(gtOv, info->gt, sizeof(double) * 6);
1385 
1386  if (ovx >= config->overview_count) {
1387  rterror(_("build_overview: Invalid overview index: %d"), ovx);
1388  return 0;
1389  }
1390  factor = config->overview[ovx];
1391  ovtable = (const char *) config->overview_table[ovx];
1392 
1393  /* factor must be within valid range */
1394  if (factor < MINOVFACTOR || factor > MAXOVFACTOR) {
1395  rterror(_("build_overview: Overview factor %d is not between %d and %d"), factor, MINOVFACTOR, MAXOVFACTOR);
1396  return 0;
1397  }
1398 
1399  dimOv[0] = (int) (info->dim[0] + (factor / 2)) / factor;
1400  dimOv[1] = (int) (info->dim[1] + (factor / 2)) / factor;
1401 
1402  /* create VRT dataset */
1403  hdsOv = VRTCreate(dimOv[0], dimOv[1]);
1404  /*
1405  GDALSetDescription(hdsOv, "/tmp/ov.vrt");
1406  */
1407  GDALSetProjection(hdsOv, info->srs);
1408 
1409  /* adjust scale */
1410  gtOv[1] *= factor;
1411  gtOv[5] *= factor;
1412 
1413  GDALSetGeoTransform(hdsOv, gtOv);
1414 
1415  /* add bands as simple sources */
1416  for (j = 0; j < info->nband_count; j++) {
1417  GDALAddBand(hdsOv, info->gdalbandtype[j], NULL);
1418  hbandOv = (VRTSourcedRasterBandH) GDALGetRasterBand(hdsOv, j + 1);
1419 
1420  if (info->hasnodata[j])
1421  GDALSetRasterNoDataValue(hbandOv, info->nodataval[j]);
1422 
1423  VRTAddSimpleSource(
1424  hbandOv, GDALGetRasterBand(hdsSrc, info->nband[j]),
1425  0, 0,
1426  info->dim[0], info->dim[1],
1427  0, 0,
1428  dimOv[0], dimOv[1],
1429  "near", VRT_NODATA_UNSET
1430  );
1431  }
1432 
1433  /* make sure VRT reflects all changes */
1434  VRTFlushCache(hdsOv);
1435 
1436  /* decide on tile size */
1437  if (!config->tile_size[0])
1438  tile_size[0] = dimOv[0];
1439  else
1440  tile_size[0] = config->tile_size[0];
1441  if (!config->tile_size[1])
1442  tile_size[1] = dimOv[1];
1443  else
1444  tile_size[1] = config->tile_size[1];
1445 
1446  /* number of tiles */
1447  if (
1448  tile_size[0] != dimOv[0] &&
1449  tile_size[1] != dimOv[1]
1450  ) {
1451  ntiles[0] = (dimOv[0] + tile_size[0] - 1) / tile_size[0];
1452  ntiles[1] = (dimOv[1] + tile_size[1] - 1) / tile_size[1];
1453  }
1454 
1455  /* working copy of geotransform matrix */
1456  memcpy(gt, gtOv, sizeof(double) * 6);
1457 
1458  /* tile overview */
1459  /* each tile is a VRT with constraints set for just the data required for the tile */
1460  for (ytile = 0; ytile < ntiles[1]; ytile++) {
1461 
1462  /* edge y tile */
1463  if (!config->pad_tile && ntiles[1] > 1 && (ytile + 1) == ntiles[1])
1464  _tile_size[1] = dimOv[1] - (ytile * tile_size[1]);
1465  else
1466  _tile_size[1] = tile_size[1];
1467 
1468  for (xtile = 0; xtile < ntiles[0]; xtile++) {
1469  /*
1470  char fn[100];
1471  sprintf(fn, "/tmp/ovtile%d.vrt", (ytile * ntiles[0]) + xtile);
1472  */
1473 
1474  /* edge x tile */
1475  if (!config->pad_tile && ntiles[0] > 1 && (xtile + 1) == ntiles[0])
1476  _tile_size[0] = dimOv[0] - (xtile * tile_size[0]);
1477  else
1478  _tile_size[0] = tile_size[0];
1479 
1480  /* compute tile's upper-left corner */
1481  GDALApplyGeoTransform(
1482  gtOv,
1483  xtile * tile_size[0], ytile * tile_size[1],
1484  &(gt[0]), &(gt[3])
1485  );
1486 
1487  /* create VRT dataset */
1488  hdsDst = VRTCreate(_tile_size[0], _tile_size[1]);
1489  /*
1490  GDALSetDescription(hdsDst, fn);
1491  */
1492  GDALSetProjection(hdsDst, info->srs);
1493  GDALSetGeoTransform(hdsDst, gt);
1494 
1495  /* add bands as simple sources */
1496  for (j = 0; j < info->nband_count; j++) {
1497  GDALAddBand(hdsDst, info->gdalbandtype[j], NULL);
1498  hbandDst = (VRTSourcedRasterBandH) GDALGetRasterBand(hdsDst, j + 1);
1499 
1500  if (info->hasnodata[j])
1501  GDALSetRasterNoDataValue(hbandDst, info->nodataval[j]);
1502 
1503  VRTAddSimpleSource(
1504  hbandDst, GDALGetRasterBand(hdsOv, j + 1),
1505  xtile * tile_size[0], ytile * tile_size[1],
1506  _tile_size[0], _tile_size[1],
1507  0, 0,
1508  _tile_size[0], _tile_size[1],
1509  "near", VRT_NODATA_UNSET
1510  );
1511  }
1512 
1513  /* make sure VRT reflects all changes */
1514  VRTFlushCache(hdsDst);
1515 
1516  /* convert VRT dataset to rt_raster */
1517  rast = rt_raster_from_gdal_dataset(hdsDst);
1518  if (rast == NULL) {
1519  rterror(_("build_overview: Could not convert VRT dataset to PostGIS raster"));
1520  GDALClose(hdsDst);
1521  return 0;
1522  }
1523 
1524  /* set srid if provided */
1525  rt_raster_set_srid(rast, info->srid);
1526 
1527  /* convert rt_raster to hexwkb */
1528  hex = rt_raster_to_hexwkb(rast, FALSE, &hexlen);
1529  raster_destroy(rast);
1530 
1531  if (hex == NULL) {
1532  rterror(_("build_overview: Could not convert PostGIS raster to hex WKB"));
1533  GDALClose(hdsDst);
1534  return 0;
1535  }
1536 
1537  /* add hexwkb to tileset */
1538  append_stringbuffer(tileset, hex);
1539 
1540  GDALClose(hdsDst);
1541 
1542  /* flush if tileset gets too big */
1543  if (tileset->length > 10) {
1544  if (!insert_records(
1545  config->schema, ovtable, config->raster_column,
1546  (config->file_column ? config->rt_filename[idx] : NULL), config->file_column_name,
1547  config->copy_statements, config->out_srid,
1548  tileset, buffer
1549  )) {
1550  rterror(_("build_overview: Could not convert raster tiles into INSERT or COPY statements"));
1551  GDALClose(hdsSrc);
1552  return 0;
1553  }
1554 
1555  rtdealloc_stringbuffer(tileset, 0);
1556  }
1557  }
1558  }
1559 
1560  GDALClose(hdsOv);
1561  GDALClose(hdsSrc);
1562  return 1;
1563 }
GDALDataType * gdalbandtype
Definition: raster2pgsql.h:179
static void rtdealloc_stringbuffer(STRINGBUFFER *buffer, int freebuffer)
Definition: raster2pgsql.c:777
#define _(String)
Definition: shpcommon.h:24
void rterror(const char *fmt,...)
Wrappers used for reporting errors and info.
Definition: rt_context.c:199
uint32_t length
Definition: raster2pgsql.h:196
static int append_stringbuffer(STRINGBUFFER *buffer, const char *str)
Definition: raster2pgsql.c:810
gt
Definition: window.py:77
unsigned int uint32_t
Definition: uthash.h:78
static int insert_records(const char *schema, const char *table, const char *column, const char *filename, const char *file_column_name, int copy_statements, int out_srid, STRINGBUFFER *tileset, STRINGBUFFER *buffer)
Definition: raster2pgsql.c:880
#define MINOVFACTOR
Definition: raster2pgsql.h:58
char * rt_raster_to_hexwkb(rt_raster raster, int outasin, uint32_t *hexwkbsize)
Return this raster in HEXWKB form (null-terminated hex)
Definition: rt_wkb.c:679
uint32_t dim[2]
Definition: raster2pgsql.h:172
void rt_raster_set_srid(rt_raster raster, int32_t srid)
Set raster&#39;s SRID.
Definition: rt_raster.c:363
#define FALSE
Definition: dbfopen.c:168
static void raster_destroy(rt_raster raster)
Definition: raster2pgsql.c:76
rt_raster rt_raster_from_gdal_dataset(GDALDatasetH ds)
Return a raster from a GDAL dataset.
Definition: rt_raster.c:2182
double * nodataval
Definition: raster2pgsql.h:185
#define MAXOVFACTOR
Definition: raster2pgsql.h:59
uint32_t nband_count
Definition: raster2pgsql.h:176
double gt[6]
Definition: raster2pgsql.h:188
if(!(yy_init))
Here is the call graph for this function:
Here is the caller graph for this function: