PostGIS  2.1.10dev-r@@SVN_REVISION@@
Datum RASTER_addBandOutDB ( PG_FUNCTION_ARGS  )

Definition at line 5753 of file rt_pg.c.

References ovdump::band, ES_ERROR, ES_NONE, FALSE, window::gt, POSTGIS_RT_DEBUG, PT_END, rtrowdump::raster, rt_band_new_offline(), rt_raster_add_band(), rt_raster_deserialize(), rt_raster_destroy(), rt_raster_get_geotransform_matrix(), rt_raster_get_height(), rt_raster_get_num_bands(), rt_raster_get_srid(), rt_raster_get_width(), rt_raster_is_empty(), rt_raster_new(), rt_raster_same_alignment(), rt_raster_serialize(), rt_raster_set_geotransform_matrix(), rt_raster_set_srid(), rt_util_gdal_datatype_to_pixtype(), rt_util_gdal_open(), rt_util_gdal_register_all(), rt_util_gdal_sr_auth_info(), rt_raster_serialized_t::size, SRID_UNKNOWN, and TRUE.

5754 {
5755  rt_pgraster *pgraster = NULL;
5756  rt_pgraster *pgrtn = NULL;
5757 
5758  rt_raster raster = NULL;
5759  rt_band band = NULL;
5760  int numbands = 0;
5761  int dstnband = 1; /* 1-based */
5762  int appendband = FALSE;
5763  char *outdbfile = NULL;
5764  int *srcnband = NULL; /* 1-based */
5765  int numsrcnband = 0;
5766  int allbands = FALSE;
5767  int hasnodata = FALSE;
5768  double nodataval = 0.;
5769  uint16_t width = 0;
5770  uint16_t height = 0;
5771  char *authname = NULL;
5772  char *authcode = NULL;
5773 
5774  int i = 0;
5775  int j = 0;
5776 
5777  GDALDatasetH hdsOut;
5778  GDALRasterBandH hbandOut;
5779  GDALDataType gdpixtype;
5780 
5781  rt_pixtype pt = PT_END;
5782  double gt[6] = {0.};
5783  double ogt[6] = {0.};
5784  rt_raster _rast = NULL;
5785  int aligned = 0;
5786  int err = 0;
5787 
5788  /* destination raster */
5789  if (!PG_ARGISNULL(0)) {
5790  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
5791 
5792  /* raster */
5793  raster = rt_raster_deserialize(pgraster, FALSE);
5794  if (!raster) {
5795  PG_FREE_IF_COPY(pgraster, 0);
5796  elog(ERROR, "RASTER_addBandOutDB: Could not deserialize destination raster");
5797  PG_RETURN_NULL();
5798  }
5799 
5800  POSTGIS_RT_DEBUG(4, "destination raster isn't NULL");
5801  }
5802 
5803  /* destination band index (1) */
5804  if (!PG_ARGISNULL(1))
5805  dstnband = PG_GETARG_INT32(1);
5806  else
5807  appendband = TRUE;
5808 
5809  /* outdb file (2) */
5810  if (PG_ARGISNULL(2)) {
5811  elog(NOTICE, "Out-db raster file not provided. Returning original raster");
5812  if (pgraster != NULL) {
5813  rt_raster_destroy(raster);
5814  PG_RETURN_POINTER(pgraster);
5815  }
5816  else
5817  PG_RETURN_NULL();
5818  }
5819  else {
5820  outdbfile = text_to_cstring(PG_GETARG_TEXT_P(2));
5821  if (!strlen(outdbfile)) {
5822  elog(NOTICE, "Out-db raster file not provided. Returning original raster");
5823  if (pgraster != NULL) {
5824  rt_raster_destroy(raster);
5825  PG_RETURN_POINTER(pgraster);
5826  }
5827  else
5828  PG_RETURN_NULL();
5829  }
5830  }
5831 
5832  /* outdb band index (3) */
5833  if (!PG_ARGISNULL(3)) {
5834  ArrayType *array;
5835  Oid etype;
5836  Datum *e;
5837  bool *nulls;
5838 
5839  int16 typlen;
5840  bool typbyval;
5841  char typalign;
5842 
5843  allbands = FALSE;
5844 
5845  array = PG_GETARG_ARRAYTYPE_P(3);
5846  etype = ARR_ELEMTYPE(array);
5847  get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
5848 
5849  switch (etype) {
5850  case INT2OID:
5851  case INT4OID:
5852  break;
5853  default:
5854  if (pgraster != NULL) {
5855  rt_raster_destroy(raster);
5856  PG_FREE_IF_COPY(pgraster, 0);
5857  }
5858  elog(ERROR, "RASTER_addBandOutDB: Invalid data type for band indexes");
5859  PG_RETURN_NULL();
5860  break;
5861  }
5862 
5863  deconstruct_array(array, etype, typlen, typbyval, typalign, &e, &nulls, &numsrcnband);
5864 
5865  srcnband = palloc(sizeof(int) * numsrcnband);
5866  if (srcnband == NULL) {
5867  if (pgraster != NULL) {
5868  rt_raster_destroy(raster);
5869  PG_FREE_IF_COPY(pgraster, 0);
5870  }
5871  elog(ERROR, "RASTER_addBandOutDB: Could not allocate memory for band indexes");
5872  PG_RETURN_NULL();
5873  }
5874 
5875  for (i = 0, j = 0; i < numsrcnband; i++) {
5876  if (nulls[i]) continue;
5877 
5878  switch (etype) {
5879  case INT2OID:
5880  srcnband[j] = DatumGetInt16(e[i]);
5881  break;
5882  case INT4OID:
5883  srcnband[j] = DatumGetInt32(e[i]);
5884  break;
5885  }
5886  j++;
5887  }
5888 
5889  if (j < numsrcnband) {
5890  srcnband = repalloc(srcnband, sizeof(int) * j);
5891  if (srcnband == NULL) {
5892  if (pgraster != NULL) {
5893  rt_raster_destroy(raster);
5894  PG_FREE_IF_COPY(pgraster, 0);
5895  }
5896  elog(ERROR, "RASTER_addBandOutDB: Could not reallocate memory for band indexes");
5897  PG_RETURN_NULL();
5898  }
5899 
5900  numsrcnband = j;
5901  }
5902  }
5903  else
5904  allbands = TRUE;
5905 
5906  /* nodataval (4) */
5907  if (!PG_ARGISNULL(4)) {
5908  hasnodata = TRUE;
5909  nodataval = PG_GETARG_FLOAT8(4);
5910  }
5911  else
5912  hasnodata = FALSE;
5913 
5914  /* validate input */
5915 
5916  /* make sure dstnband is valid */
5917  if (raster != NULL) {
5918  numbands = rt_raster_get_num_bands(raster);
5919  if (!appendband) {
5920  if (dstnband < 1) {
5921  elog(NOTICE, "Invalid band index %d for adding bands. Using band index 1", dstnband);
5922  dstnband = 1;
5923  }
5924  else if (dstnband > numbands) {
5925  elog(NOTICE, "Invalid band index %d for adding bands. Using band index %d", dstnband, numbands);
5926  dstnband = numbands + 1;
5927  }
5928  }
5929  else
5930  dstnband = numbands + 1;
5931  }
5932 
5933  /* open outdb raster file */
5935  hdsOut = rt_util_gdal_open(outdbfile, GA_ReadOnly, 1);
5936  if (hdsOut == NULL) {
5937  if (pgraster != NULL) {
5938  rt_raster_destroy(raster);
5939  PG_FREE_IF_COPY(pgraster, 0);
5940  }
5941  elog(ERROR, "RASTER_addBandOutDB: Could not open out-db file with GDAL");
5942  PG_RETURN_NULL();
5943  }
5944 
5945  /* get offline raster's geotransform */
5946  if (GDALGetGeoTransform(hdsOut, ogt) != CE_None) {
5947  ogt[0] = 0;
5948  ogt[1] = 1;
5949  ogt[2] = 0;
5950  ogt[3] = 0;
5951  ogt[4] = 0;
5952  ogt[5] = -1;
5953  }
5954 
5955  /* raster doesn't exist, create it now */
5956  if (raster == NULL) {
5957  raster = rt_raster_new(GDALGetRasterXSize(hdsOut), GDALGetRasterYSize(hdsOut));
5958  if (rt_raster_is_empty(raster)) {
5959  elog(ERROR, "RASTER_addBandOutDB: Could not create new raster");
5960  PG_RETURN_NULL();
5961  }
5962  rt_raster_set_geotransform_matrix(raster, ogt);
5964 
5965  if (rt_util_gdal_sr_auth_info(hdsOut, &authname, &authcode) == ES_NONE) {
5966  if (
5967  authname != NULL &&
5968  strcmp(authname, "EPSG") == 0 &&
5969  authcode != NULL
5970  ) {
5971  rt_raster_set_srid(raster, atoi(authcode));
5972  }
5973  else
5974  elog(INFO, "Unknown SRS auth name and code from out-db file. Defaulting SRID of new raster to %d", SRID_UNKNOWN);
5975  }
5976  else
5977  elog(INFO, "Could not get SRS auth name and code from out-db file. Defaulting SRID of new raster to %d", SRID_UNKNOWN);
5978  }
5979 
5980  /* some raster info */
5981  width = rt_raster_get_width(raster);
5982  height = rt_raster_get_height(raster);
5983 
5984  /* are rasters aligned? */
5985  _rast = rt_raster_new(1, 1);
5987  rt_raster_set_srid(_rast, rt_raster_get_srid(raster));
5988  err = rt_raster_same_alignment(raster, _rast, &aligned, NULL);
5989  rt_raster_destroy(_rast);
5990 
5991  if (err != ES_NONE) {
5992  GDALClose(hdsOut);
5993  if (raster != NULL)
5994  rt_raster_destroy(raster);
5995  if (pgraster != NULL)
5996  PG_FREE_IF_COPY(pgraster, 0);
5997  elog(ERROR, "RASTER_addBandOutDB: Could not test alignment of out-db file");
5998  return ES_ERROR;
5999  }
6000  else if (!aligned)
6001  elog(WARNING, "The in-db representation of the out-db raster is not aligned. Band data may be incorrect");
6002 
6003  numbands = GDALGetRasterCount(hdsOut);
6004 
6005  /* build up srcnband */
6006  if (allbands) {
6007  numsrcnband = numbands;
6008  srcnband = palloc(sizeof(int) * numsrcnband);
6009  if (srcnband == NULL) {
6010  GDALClose(hdsOut);
6011  if (raster != NULL)
6012  rt_raster_destroy(raster);
6013  if (pgraster != NULL)
6014  PG_FREE_IF_COPY(pgraster, 0);
6015  elog(ERROR, "RASTER_addBandOutDB: Could not allocate memory for band indexes");
6016  PG_RETURN_NULL();
6017  }
6018 
6019  for (i = 0, j = 1; i < numsrcnband; i++, j++)
6020  srcnband[i] = j;
6021  }
6022 
6023  /* check band properties and add band */
6024  for (i = 0, j = dstnband - 1; i < numsrcnband; i++, j++) {
6025  /* valid index? */
6026  if (srcnband[i] < 1 || srcnband[i] > numbands) {
6027  elog(NOTICE, "Out-db file does not have a band at index %d. Returning original raster", srcnband[i]);
6028  GDALClose(hdsOut);
6029  if (raster != NULL)
6030  rt_raster_destroy(raster);
6031  if (pgraster != NULL)
6032  PG_RETURN_POINTER(pgraster);
6033  else
6034  PG_RETURN_NULL();
6035  }
6036 
6037  /* get outdb band */
6038  hbandOut = NULL;
6039  hbandOut = GDALGetRasterBand(hdsOut, srcnband[i]);
6040  if (NULL == hbandOut) {
6041  GDALClose(hdsOut);
6042  if (raster != NULL)
6043  rt_raster_destroy(raster);
6044  if (pgraster != NULL)
6045  PG_FREE_IF_COPY(pgraster, 0);
6046  elog(ERROR, "RASTER_addBandOutDB: Could not get band %d from GDAL dataset", srcnband[i]);
6047  PG_RETURN_NULL();
6048  }
6049 
6050  /* supported pixel type */
6051  gdpixtype = GDALGetRasterDataType(hbandOut);
6052  pt = rt_util_gdal_datatype_to_pixtype(gdpixtype);
6053  if (pt == PT_END) {
6054  elog(NOTICE, "Pixel type %s of band %d from GDAL dataset is not supported. Returning original raster", GDALGetDataTypeName(gdpixtype), srcnband[i]);
6055  GDALClose(hdsOut);
6056  if (raster != NULL)
6057  rt_raster_destroy(raster);
6058  if (pgraster != NULL)
6059  PG_RETURN_POINTER(pgraster);
6060  else
6061  PG_RETURN_NULL();
6062  }
6063 
6064  /* use out-db band's nodata value if nodataval not already set */
6065  if (!hasnodata)
6066  nodataval = GDALGetRasterNoDataValue(hbandOut, &hasnodata);
6067 
6068  /* add band */
6069  band = rt_band_new_offline(
6070  width, height,
6071  pt,
6072  hasnodata, nodataval,
6073  srcnband[i] - 1, outdbfile
6074  );
6075  if (band == NULL) {
6076  GDALClose(hdsOut);
6077  if (raster != NULL)
6078  rt_raster_destroy(raster);
6079  if (pgraster != NULL)
6080  PG_FREE_IF_COPY(pgraster, 0);
6081  elog(ERROR, "RASTER_addBandOutDB: Could not create new out-db band");
6082  PG_RETURN_NULL();
6083  }
6084 
6085  if (rt_raster_add_band(raster, band, j) < 0) {
6086  GDALClose(hdsOut);
6087  if (raster != NULL)
6088  rt_raster_destroy(raster);
6089  if (pgraster != NULL)
6090  PG_FREE_IF_COPY(pgraster, 0);
6091  elog(ERROR, "RASTER_addBandOutDB: Could not add new out-db band to raster");
6092  PG_RETURN_NULL();
6093  }
6094  }
6095 
6096  pgrtn = rt_raster_serialize(raster);
6097  rt_raster_destroy(raster);
6098  if (pgraster != NULL)
6099  PG_FREE_IF_COPY(pgraster, 0);
6100  if (!pgrtn)
6101  PG_RETURN_NULL();
6102 
6103  SET_VARSIZE(pgrtn, pgrtn->size);
6104  PG_RETURN_POINTER(pgrtn);
6105 }
int rt_raster_get_num_bands(rt_raster raster)
Definition: rt_api.c:5677
int rt_util_gdal_register_all(int force_register_all)
Definition: rt_api.c:445
void * rt_raster_serialize(rt_raster raster)
Return this raster in serialized form.
Definition: rt_api.c:8158
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_api.c:5387
tuple gt
Definition: window.py:79
Definition: rt_api.h:184
tuple band
Definition: ovdump.py:57
tuple raster
Be careful!! Zeros function's input parameter can be a (height x width) array, not (width x height): ...
Definition: rtrowdump.py:123
void rt_raster_set_srid(rt_raster raster, int32_t srid)
Set raster's SRID.
Definition: rt_api.c:5668
int32_t rt_raster_get_srid(rt_raster raster)
Get raster's SRID.
Definition: rt_api.c:5661
GDALDatasetH rt_util_gdal_open(const char *fn, GDALAccess fn_access, int shared)
Definition: rt_api.c:487
void rt_raster_get_geotransform_matrix(rt_raster raster, double *gt)
Get 6-element array of raster geotransform matrix.
Definition: rt_api.c:6005
int rt_raster_add_band(rt_raster raster, rt_band band, int index)
Add band data to a raster.
Definition: rt_api.c:5706
rt_pixtype rt_util_gdal_datatype_to_pixtype(GDALDataType gdt)
Convert GDALDataType to rt_pixtype.
Definition: rt_api.c:264
rt_pixtype
Definition: rt_api.h:172
#define POSTGIS_RT_DEBUG(level, msg)
Definition: rt_pg.h:58
rt_band rt_band_new_offline(uint16_t width, uint16_t height, rt_pixtype pixtype, uint32_t hasnodata, double nodataval, uint8_t bandNum, const char *path)
Create an out-db rt_band.
Definition: rt_api.c:1527
#define SRID_UNKNOWN
Unknown SRID value.
Definition: liblwgeom.h:154
void rt_raster_set_geotransform_matrix(rt_raster raster, double *gt)
Set raster's geotransform using 6-element array.
Definition: rt_api.c:6026
rt_errorstate rt_util_gdal_sr_auth_info(GDALDatasetH hds, char **authname, char **authcode)
Get auth name and code.
Definition: rt_api.c:381
uint16_t rt_raster_get_height(rt_raster raster)
Definition: rt_api.c:5434
#define FALSE
Definition: dbfopen.c:169
Struct definitions.
Definition: rt_api.h:2175
int rt_raster_is_empty(rt_raster raster)
Return TRUE if the raster is empty.
Definition: rt_api.c:8550
rt_raster rt_raster_deserialize(void *serialized, int header_only)
Return a raster from a serialized form.
Definition: rt_api.c:8350
rt_raster rt_raster_new(uint32_t width, uint32_t height)
Construct a raster with given dimensions.
Definition: rt_api.c:5353
rt_errorstate rt_raster_same_alignment(rt_raster rast1, rt_raster rast2, int *aligned, char **reason)
Definition: rt_api.c:12774
#define TRUE
Definition: dbfopen.c:170
uint16_t rt_raster_get_width(rt_raster raster)
Definition: rt_api.c:5426

Here is the call graph for this function: