PostGIS  2.5.7dev-r@@SVN_REVISION@@
rtpg_create.c
Go to the documentation of this file.
1 /*
2  *
3  * WKTRaster - Raster Types for PostGIS
4  * http://trac.osgeo.org/postgis/wiki/WKTRaster
5  *
6  * Copyright (C) 2011-2013 Regents of the University of California
7  * <bkpark@ucdavis.edu>
8  * Copyright (C) 2010-2011 Jorge Arevalo <jorge.arevalo@deimos-space.com>
9  * Copyright (C) 2010-2011 David Zwarg <dzwarg@azavea.com>
10  * Copyright (C) 2009-2011 Pierre Racine <pierre.racine@sbf.ulaval.ca>
11  * Copyright (C) 2009-2011 Mateusz Loskot <mateusz@loskot.net>
12  * Copyright (C) 2008-2009 Sandro Santilli <strk@kbt.io>
13  *
14  * This program is free software; you can redistribute it and/or
15  * modify it under the terms of the GNU General Public License
16  * as published by the Free Software Foundation; either version 2
17  * of the License, or (at your option) any later version.
18  *
19  * This program is distributed in the hope that it will be useful,
20  * but WITHOUT ANY WARRANTY; without even the implied warranty of
21  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22  * GNU General Public License for more details.
23  *
24  * You should have received a copy of the GNU General Public License
25  * along with this program; if not, write to the Free Software Foundation,
26  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
27  *
28  */
29 
30 #include <postgres.h>
31 #include <fmgr.h>
32 #include <funcapi.h>
33 #include <utils/builtins.h> /* for text_to_cstring() */
34 #include "utils/lsyscache.h" /* for get_typlenbyvalalign */
35 #include "utils/array.h" /* for ArrayType */
36 #include "catalog/pg_type.h" /* for INT2OID, INT4OID, FLOAT4OID, FLOAT8OID and TEXTOID */
37 
38 #include "rtpostgis.h"
39 
40 /* Raster and band creation */
41 Datum RASTER_makeEmpty(PG_FUNCTION_ARGS);
42 Datum RASTER_addBand(PG_FUNCTION_ARGS);
43 Datum RASTER_addBandRasterArray(PG_FUNCTION_ARGS);
44 Datum RASTER_addBandOutDB(PG_FUNCTION_ARGS);
45 Datum RASTER_copyBand(PG_FUNCTION_ARGS);
46 Datum RASTER_tile(PG_FUNCTION_ARGS);
47 
48 /* create new raster from existing raster's bands */
49 Datum RASTER_band(PG_FUNCTION_ARGS);
50 
54 PG_FUNCTION_INFO_V1(RASTER_makeEmpty);
55 Datum RASTER_makeEmpty(PG_FUNCTION_ARGS)
56 {
57  uint16 width = 0, height = 0;
58  double ipx = 0, ipy = 0, scalex = 0, scaley = 0, skewx = 0, skewy = 0;
59  int32_t srid = SRID_UNKNOWN;
60  rt_pgraster *pgraster = NULL;
61  rt_raster raster;
62 
63  if (PG_NARGS() < 9) {
64  elog(ERROR, "RASTER_makeEmpty: ST_MakeEmptyRaster requires 9 args");
65  PG_RETURN_NULL();
66  }
67 
68  if (!PG_ARGISNULL(0))
69  width = PG_GETARG_UINT16(0);
70 
71  if (!PG_ARGISNULL(1))
72  height = PG_GETARG_UINT16(1);
73 
74  if (!PG_ARGISNULL(2))
75  ipx = PG_GETARG_FLOAT8(2);
76 
77  if (!PG_ARGISNULL(3))
78  ipy = PG_GETARG_FLOAT8(3);
79 
80  if (!PG_ARGISNULL(4))
81  scalex = PG_GETARG_FLOAT8(4);
82 
83  if (!PG_ARGISNULL(5))
84  scaley = PG_GETARG_FLOAT8(5);
85 
86  if (!PG_ARGISNULL(6))
87  skewx = PG_GETARG_FLOAT8(6);
88 
89  if (!PG_ARGISNULL(7))
90  skewy = PG_GETARG_FLOAT8(7);
91 
92  if (!PG_ARGISNULL(8))
93  srid = PG_GETARG_INT32(8);
94 
95  POSTGIS_RT_DEBUGF(4, "%dx%d, ip:%g,%g, scale:%g,%g, skew:%g,%g srid:%d",
96  width, height, ipx, ipy, scalex, scaley,
97  skewx, skewy, srid);
98 
99  raster = rt_raster_new(width, height);
100  if (raster == NULL)
101  PG_RETURN_NULL(); /* error was supposedly printed already */
102 
103  rt_raster_set_scale(raster, scalex, scaley);
104  rt_raster_set_offsets(raster, ipx, ipy);
105  rt_raster_set_skews(raster, skewx, skewy);
106  rt_raster_set_srid(raster, srid);
107 
108  pgraster = rt_raster_serialize(raster);
109  rt_raster_destroy(raster);
110  if (!pgraster)
111  PG_RETURN_NULL();
112 
113  SET_VARSIZE(pgraster, pgraster->size);
114  PG_RETURN_POINTER(pgraster);
115 }
116 
120 PG_FUNCTION_INFO_V1(RASTER_addBand);
121 Datum RASTER_addBand(PG_FUNCTION_ARGS)
122 {
123  rt_pgraster *pgraster = NULL;
124  rt_pgraster *pgrtn = NULL;
125  rt_raster raster = NULL;
126  int bandindex = 0;
127  int maxbandindex = 0;
128  int numbands = 0;
129  int lastnumbands = 0;
130 
131  text *text_pixtype = NULL;
132  char *char_pixtype = NULL;
133 
134  struct addbandarg {
135  int index;
136  bool append;
137  rt_pixtype pixtype;
138  double initialvalue;
139  bool hasnodata;
140  double nodatavalue;
141  };
142  struct addbandarg *arg = NULL;
143 
144  ArrayType *array;
145  Oid etype;
146  Datum *e;
147  bool *nulls;
148  int16 typlen;
149  bool typbyval;
150  char typalign;
151  int n = 0;
152 
153  HeapTupleHeader tup;
154  bool isnull;
155  Datum tupv;
156 
157  int i = 0;
158 
159  /* pgraster is null, return null */
160  if (PG_ARGISNULL(0))
161  PG_RETURN_NULL();
162  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
163 
164  /* raster */
165  raster = rt_raster_deserialize(pgraster, FALSE);
166  if (!raster) {
167  PG_FREE_IF_COPY(pgraster, 0);
168  elog(ERROR, "RASTER_addBand: Could not deserialize raster");
169  PG_RETURN_NULL();
170  }
171 
172  /* process set of addbandarg */
173  POSTGIS_RT_DEBUG(3, "Processing Arg 1 (addbandargset)");
174  array = PG_GETARG_ARRAYTYPE_P(1);
175  etype = ARR_ELEMTYPE(array);
176  get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
177 
178  deconstruct_array(array, etype, typlen, typbyval, typalign, &e,
179  &nulls, &n);
180 
181  if (!n) {
182  PG_FREE_IF_COPY(pgraster, 0);
183  elog(ERROR, "RASTER_addBand: Invalid argument for addbandargset");
184  PG_RETURN_NULL();
185  }
186 
187  /* allocate addbandarg */
188  arg = (struct addbandarg *) palloc(sizeof(struct addbandarg) * n);
189  if (arg == NULL) {
190  rt_raster_destroy(raster);
191  PG_FREE_IF_COPY(pgraster, 0);
192  elog(ERROR, "RASTER_addBand: Could not allocate memory for addbandarg");
193  PG_RETURN_NULL();
194  }
195 
196  /*
197  process each element of addbandargset
198  each element is the index of where to add the new band,
199  new band's pixeltype, the new band's initial value and
200  the new band's NODATA value if NOT NULL
201  */
202  for (i = 0; i < n; i++) {
203  if (nulls[i]) continue;
204 
205  POSTGIS_RT_DEBUGF(4, "Processing addbandarg at index %d", i);
206 
207  /* each element is a tuple */
208  tup = (HeapTupleHeader) DatumGetPointer(e[i]);
209  if (NULL == tup) {
210  pfree(arg);
211  rt_raster_destroy(raster);
212  PG_FREE_IF_COPY(pgraster, 0);
213  elog(ERROR, "RASTER_addBand: Invalid argument for addbandargset");
214  PG_RETURN_NULL();
215  }
216 
217  /* new band index, 1-based */
218  arg[i].index = 0;
219  arg[i].append = TRUE;
220  tupv = GetAttributeByName(tup, "index", &isnull);
221  if (!isnull) {
222  arg[i].index = DatumGetInt32(tupv);
223  arg[i].append = FALSE;
224  }
225 
226  /* for now, only check that band index is 1-based */
227  if (!arg[i].append && arg[i].index < 1) {
228  pfree(arg);
229  rt_raster_destroy(raster);
230  PG_FREE_IF_COPY(pgraster, 0);
231  elog(ERROR, "RASTER_addBand: Invalid argument for addbandargset. Invalid band index (must be 1-based) for addbandarg of index %d", i);
232  PG_RETURN_NULL();
233  }
234 
235  /* new band pixeltype */
236  arg[i].pixtype = PT_END;
237  tupv = GetAttributeByName(tup, "pixeltype", &isnull);
238  if (isnull) {
239  pfree(arg);
240  rt_raster_destroy(raster);
241  PG_FREE_IF_COPY(pgraster, 0);
242  elog(ERROR, "RASTER_addBand: Invalid argument for addbandargset. Pixel type cannot be NULL for addbandarg of index %d", i);
243  PG_RETURN_NULL();
244  }
245  text_pixtype = (text *) DatumGetPointer(tupv);
246  if (text_pixtype == NULL) {
247  pfree(arg);
248  rt_raster_destroy(raster);
249  PG_FREE_IF_COPY(pgraster, 0);
250  elog(ERROR, "RASTER_addBand: Invalid argument for addbandargset. Pixel type cannot be NULL for addbandarg of index %d", i);
251  PG_RETURN_NULL();
252  }
253  char_pixtype = text_to_cstring(text_pixtype);
254 
255  arg[i].pixtype = rt_pixtype_index_from_name(char_pixtype);
256  pfree(char_pixtype);
257  if (arg[i].pixtype == PT_END) {
258  pfree(arg);
259  rt_raster_destroy(raster);
260  PG_FREE_IF_COPY(pgraster, 0);
261  elog(ERROR, "RASTER_addBand: Invalid argument for addbandargset. Invalid pixel type for addbandarg of index %d", i);
262  PG_RETURN_NULL();
263  }
264 
265  /* new band initialvalue */
266  arg[i].initialvalue = 0;
267  tupv = GetAttributeByName(tup, "initialvalue", &isnull);
268  if (!isnull)
269  arg[i].initialvalue = DatumGetFloat8(tupv);
270 
271  /* new band NODATA value */
272  arg[i].hasnodata = FALSE;
273  arg[i].nodatavalue = 0;
274  tupv = GetAttributeByName(tup, "nodataval", &isnull);
275  if (!isnull) {
276  arg[i].hasnodata = TRUE;
277  arg[i].nodatavalue = DatumGetFloat8(tupv);
278  }
279  }
280 
281  /* add new bands to raster */
282  lastnumbands = rt_raster_get_num_bands(raster);
283  for (i = 0; i < n; i++) {
284  if (nulls[i]) continue;
285 
286  POSTGIS_RT_DEBUGF(3, "%d bands in old raster", lastnumbands);
287  maxbandindex = lastnumbands + 1;
288 
289  /* check that new band's index doesn't exceed maxbandindex */
290  if (!arg[i].append) {
291  if (arg[i].index > maxbandindex) {
292  elog(NOTICE, "Band index for addbandarg of index %d exceeds possible value. Adding band at index %d", i, maxbandindex);
293  arg[i].index = maxbandindex;
294  }
295  }
296  /* append, so use maxbandindex */
297  else
298  arg[i].index = maxbandindex;
299 
300  POSTGIS_RT_DEBUGF(4, "new band (index, pixtype, initialvalue, hasnodata, nodatavalue) = (%d, %s, %f, %s, %f)",
301  arg[i].index,
302  rt_pixtype_name(arg[i].pixtype),
303  arg[i].initialvalue,
304  arg[i].hasnodata ? "TRUE" : "FALSE",
305  arg[i].nodatavalue
306  );
307 
308  bandindex = rt_raster_generate_new_band(
309  raster,
310  arg[i].pixtype, arg[i].initialvalue,
311  arg[i].hasnodata, arg[i].nodatavalue,
312  arg[i].index - 1
313  );
314 
315  numbands = rt_raster_get_num_bands(raster);
316  if (numbands == lastnumbands || bandindex == -1) {
317  pfree(arg);
318  rt_raster_destroy(raster);
319  PG_FREE_IF_COPY(pgraster, 0);
320  elog(ERROR, "RASTER_addBand: Could not add band defined by addbandarg of index %d to raster", i);
321  PG_RETURN_NULL();
322  }
323 
324  lastnumbands = numbands;
325  POSTGIS_RT_DEBUGF(3, "%d bands in new raster", lastnumbands);
326  }
327 
328  pfree(arg);
329 
330  pgrtn = rt_raster_serialize(raster);
331  rt_raster_destroy(raster);
332  PG_FREE_IF_COPY(pgraster, 0);
333  if (!pgrtn)
334  PG_RETURN_NULL();
335 
336  SET_VARSIZE(pgrtn, pgrtn->size);
337  PG_RETURN_POINTER(pgrtn);
338 }
339 
343 PG_FUNCTION_INFO_V1(RASTER_addBandRasterArray);
344 Datum RASTER_addBandRasterArray(PG_FUNCTION_ARGS)
345 {
346  rt_pgraster *pgraster = NULL;
347  rt_pgraster *pgsrc = NULL;
348  rt_pgraster *pgrtn = NULL;
349 
350  rt_raster raster = NULL;
351  rt_raster src = NULL;
352 
353  int srcnband = 1;
354  bool appendband = FALSE;
355  int dstnband = 1;
356  int srcnumbands = 0;
357  int dstnumbands = 0;
358 
359  ArrayType *array;
360  Oid etype;
361  Datum *e;
362  bool *nulls;
363  int16 typlen;
364  bool typbyval;
365  char typalign;
366  int n = 0;
367 
368  int rtn = 0;
369  int i = 0;
370 
371  /* destination raster */
372  if (!PG_ARGISNULL(0)) {
373  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
374 
375  /* raster */
376  raster = rt_raster_deserialize(pgraster, FALSE);
377  if (!raster) {
378  PG_FREE_IF_COPY(pgraster, 0);
379  elog(ERROR, "RASTER_addBandRasterArray: Could not deserialize destination raster");
380  PG_RETURN_NULL();
381  }
382 
383  POSTGIS_RT_DEBUG(4, "destination raster isn't NULL");
384  }
385 
386  if (PG_ARGISNULL(1))
387  {
388  if (raster != NULL)
389  {
390  rt_raster_destroy(raster);
391  PG_RETURN_POINTER(pgraster);
392  }
393  else
394  PG_RETURN_NULL();
395  }
396 
397  /* source rasters' band index, 1-based */
398  if (!PG_ARGISNULL(2))
399  srcnband = PG_GETARG_INT32(2);
400  if (srcnband < 1) {
401  elog(NOTICE, "Invalid band index for source rasters (must be 1-based). Returning original raster");
402  if (raster != NULL) {
403  rt_raster_destroy(raster);
404  PG_RETURN_POINTER(pgraster);
405  }
406  else
407  PG_RETURN_NULL();
408  }
409  POSTGIS_RT_DEBUGF(4, "srcnband = %d", srcnband);
410 
411  /* destination raster's band index, 1-based */
412  if (!PG_ARGISNULL(3)) {
413  dstnband = PG_GETARG_INT32(3);
414  appendband = FALSE;
415 
416  if (dstnband < 1) {
417  elog(NOTICE, "Invalid band index for destination raster (must be 1-based). Returning original raster");
418  if (raster != NULL) {
419  rt_raster_destroy(raster);
420  PG_RETURN_POINTER(pgraster);
421  }
422  else
423  PG_RETURN_NULL();
424  }
425  }
426  else
427  appendband = TRUE;
428 
429  /* additional processing of dstnband */
430  if (raster != NULL) {
431  dstnumbands = rt_raster_get_num_bands(raster);
432 
433  if (dstnumbands < 1) {
434  appendband = TRUE;
435  dstnband = 1;
436  }
437  else if (appendband)
438  dstnband = dstnumbands + 1;
439  else if (dstnband > dstnumbands) {
440  elog(NOTICE, "Band index provided for destination raster is greater than the number of bands in the raster. Bands will be appended");
441  appendband = TRUE;
442  dstnband = dstnumbands + 1;
443  }
444  }
445  POSTGIS_RT_DEBUGF(4, "appendband = %d", appendband);
446  POSTGIS_RT_DEBUGF(4, "dstnband = %d", dstnband);
447 
448  /* process set of source rasters */
449  POSTGIS_RT_DEBUG(3, "Processing array of source rasters");
450  array = PG_GETARG_ARRAYTYPE_P(1);
451  etype = ARR_ELEMTYPE(array);
452  get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
453 
454  deconstruct_array(array, etype, typlen, typbyval, typalign, &e,
455  &nulls, &n);
456 
457  /* decrement srcnband and dstnband by 1, now 0-based */
458  srcnband--;
459  dstnband--;
460  POSTGIS_RT_DEBUGF(4, "0-based nband (src, dst) = (%d, %d)", srcnband, dstnband);
461 
462  /* time to copy bands */
463  for (i = 0; i < n; i++) {
464  if (nulls[i]) continue;
465  src = NULL;
466 
467  pgsrc = (rt_pgraster *) PG_DETOAST_DATUM(e[i]);
468  src = rt_raster_deserialize(pgsrc, FALSE);
469  if (src == NULL) {
470  pfree(nulls);
471  pfree(e);
472  if (raster != NULL)
473  rt_raster_destroy(raster);
474  if (pgraster != NULL)
475  PG_FREE_IF_COPY(pgraster, 0);
476  elog(ERROR, "RASTER_addBandRasterArray: Could not deserialize source raster at index %d", i + 1);
477  PG_RETURN_NULL();
478  }
479 
480  srcnumbands = rt_raster_get_num_bands(src);
481  POSTGIS_RT_DEBUGF(4, "source raster %d has %d bands", i + 1, srcnumbands);
482 
483  /* band index isn't valid */
484  if (srcnband > srcnumbands - 1) {
485  elog(NOTICE, "Invalid band index for source raster at index %d. Returning original raster", i + 1);
486  pfree(nulls);
487  pfree(e);
488  rt_raster_destroy(src);
489  if (raster != NULL) {
490  rt_raster_destroy(raster);
491  PG_RETURN_POINTER(pgraster);
492  }
493  else
494  PG_RETURN_NULL();
495  }
496 
497  /* destination raster is empty, new raster */
498  if (raster == NULL) {
499  uint32_t srcnbands[1] = {srcnband};
500 
501  POSTGIS_RT_DEBUG(4, "empty destination raster, using rt_raster_from_band");
502 
503  raster = rt_raster_from_band(src, srcnbands, 1);
504  rt_raster_destroy(src);
505  if (raster == NULL) {
506  pfree(nulls);
507  pfree(e);
508  if (pgraster != NULL)
509  PG_FREE_IF_COPY(pgraster, 0);
510  elog(ERROR, "RASTER_addBandRasterArray: Could not create raster from source raster at index %d", i + 1);
511  PG_RETURN_NULL();
512  }
513  }
514  /* copy band */
515  else {
516  rtn = rt_raster_copy_band(
517  raster, src,
518  srcnband, dstnband
519  );
520  rt_raster_destroy(src);
521 
522  if (rtn == -1 || rt_raster_get_num_bands(raster) == dstnumbands) {
523  elog(NOTICE, "Could not add band from source raster at index %d to destination raster. Returning original raster", i + 1);
524  rt_raster_destroy(raster);
525  pfree(nulls);
526  pfree(e);
527  if (pgraster != NULL)
528  PG_RETURN_POINTER(pgraster);
529  else
530  PG_RETURN_NULL();
531  }
532  }
533 
534  dstnband++;
535  dstnumbands++;
536  }
537 
538  if (raster != NULL) {
539  pgrtn = rt_raster_serialize(raster);
540  rt_raster_destroy(raster);
541  if (pgraster != NULL)
542  PG_FREE_IF_COPY(pgraster, 0);
543  if (!pgrtn)
544  PG_RETURN_NULL();
545 
546  SET_VARSIZE(pgrtn, pgrtn->size);
547  PG_RETURN_POINTER(pgrtn);
548  }
549 
550  PG_RETURN_NULL();
551 }
552 
556 PG_FUNCTION_INFO_V1(RASTER_addBandOutDB);
557 Datum RASTER_addBandOutDB(PG_FUNCTION_ARGS)
558 {
559  rt_pgraster *pgraster = NULL;
560  rt_pgraster *pgrtn = NULL;
561 
562  rt_raster raster = NULL;
563  rt_band band = NULL;
564  int numbands = 0;
565  int dstnband = 1; /* 1-based */
566  int appendband = FALSE;
567  char *outdbfile = NULL;
568  int *srcnband = NULL; /* 1-based */
569  int numsrcnband = 0;
570  int allbands = FALSE;
571  int hasnodata = FALSE;
572  double nodataval = 0.;
573  uint16_t width = 0;
574  uint16_t height = 0;
575  char *authname = NULL;
576  char *authcode = NULL;
577 
578  int i = 0;
579  int j = 0;
580 
581  GDALDatasetH hdsOut;
582 
583  double ogt[6] = {0.};
584  rt_raster _rast = NULL;
585  int aligned = 0;
586  int err = 0;
587 
588  /* destination raster */
589  if (!PG_ARGISNULL(0)) {
590  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
591 
592  /* raster */
593  raster = rt_raster_deserialize(pgraster, FALSE);
594  if (!raster) {
595  PG_FREE_IF_COPY(pgraster, 0);
596  elog(ERROR, "RASTER_addBandOutDB: Cannot deserialize destination raster");
597  PG_RETURN_NULL();
598  }
599 
600  POSTGIS_RT_DEBUG(4, "destination raster isn't NULL");
601  }
602 
603  /* destination band index (1) */
604  if (!PG_ARGISNULL(1))
605  dstnband = PG_GETARG_INT32(1);
606  else
607  appendband = TRUE;
608 
609  /* outdb file (2) */
610  if (PG_ARGISNULL(2)) {
611  elog(NOTICE, "Out-db raster file not provided. Returning original raster");
612  if (pgraster != NULL) {
613  rt_raster_destroy(raster);
614  PG_RETURN_POINTER(pgraster);
615  }
616  else
617  PG_RETURN_NULL();
618  }
619  else {
620  outdbfile = text_to_cstring(PG_GETARG_TEXT_P(2));
621  if (!strlen(outdbfile)) {
622  elog(NOTICE, "Out-db raster file not provided. Returning original raster");
623  if (pgraster != NULL) {
624  rt_raster_destroy(raster);
625  PG_RETURN_POINTER(pgraster);
626  }
627  else
628  PG_RETURN_NULL();
629  }
630  }
631 
632  /* outdb band index (3) */
633  if (!PG_ARGISNULL(3)) {
634  ArrayType *array;
635  Oid etype;
636  Datum *e;
637  bool *nulls;
638 
639  int16 typlen;
640  bool typbyval;
641  char typalign;
642 
643  allbands = FALSE;
644 
645  array = PG_GETARG_ARRAYTYPE_P(3);
646  etype = ARR_ELEMTYPE(array);
647  get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
648 
649  switch (etype) {
650  case INT2OID:
651  case INT4OID:
652  break;
653  default:
654  if (pgraster != NULL) {
655  rt_raster_destroy(raster);
656  PG_FREE_IF_COPY(pgraster, 0);
657  }
658  elog(ERROR, "RASTER_addBandOutDB: Invalid data type for band indexes");
659  PG_RETURN_NULL();
660  break;
661  }
662 
663  deconstruct_array(array, etype, typlen, typbyval, typalign, &e, &nulls, &numsrcnband);
664 
665  srcnband = palloc(sizeof(int) * numsrcnband);
666  if (srcnband == NULL) {
667  if (pgraster != NULL) {
668  rt_raster_destroy(raster);
669  PG_FREE_IF_COPY(pgraster, 0);
670  }
671  elog(ERROR, "RASTER_addBandOutDB: Cannot allocate memory for band indexes");
672  PG_RETURN_NULL();
673  }
674 
675  for (i = 0, j = 0; i < numsrcnband; i++) {
676  if (nulls[i]) continue;
677 
678  switch (etype) {
679  case INT2OID:
680  srcnband[j] = DatumGetInt16(e[i]);
681  break;
682  case INT4OID:
683  srcnband[j] = DatumGetInt32(e[i]);
684  break;
685  }
686  j++;
687  }
688 
689  if (j < numsrcnband) {
690  srcnband = repalloc(srcnband, sizeof(int) * j);
691  if (srcnband == NULL) {
692  if (pgraster != NULL) {
693  rt_raster_destroy(raster);
694  PG_FREE_IF_COPY(pgraster, 0);
695  }
696  elog(ERROR, "RASTER_addBandOutDB: Cannot reallocate memory for band indexes");
697  PG_RETURN_NULL();
698  }
699 
700  numsrcnband = j;
701  }
702  }
703  else
704  allbands = TRUE;
705 
706  /* nodataval (4) */
707  if (!PG_ARGISNULL(4)) {
708  hasnodata = TRUE;
709  nodataval = PG_GETARG_FLOAT8(4);
710  }
711  else
712  hasnodata = FALSE;
713 
714  /* validate input */
715 
716  /* make sure dstnband is valid */
717  if (raster != NULL) {
718  numbands = rt_raster_get_num_bands(raster);
719  if (!appendband) {
720  if (dstnband < 1) {
721  elog(NOTICE, "Invalid band index %d for adding bands. Using band index 1", dstnband);
722  dstnband = 1;
723  }
724  else if (numbands > 0 && dstnband > numbands) {
725  elog(NOTICE, "Invalid band index %d for adding bands. Using band index %d", dstnband, numbands);
726  dstnband = numbands + 1;
727  }
728  }
729  else
730  dstnband = numbands + 1;
731  }
732 
733  /* open outdb raster file */
734  rt_util_gdal_register_all(0);
735  hdsOut = rt_util_gdal_open(outdbfile, GA_ReadOnly, 1);
736  if (hdsOut == NULL) {
737  if (pgraster != NULL) {
738  rt_raster_destroy(raster);
739  PG_FREE_IF_COPY(pgraster, 0);
740  }
741  elog(ERROR, "RASTER_addBandOutDB: Cannot open out-db file with GDAL");
742  PG_RETURN_NULL();
743  }
744 
745  /* get offline raster's geotransform */
746  if (GDALGetGeoTransform(hdsOut, ogt) != CE_None) {
747  ogt[0] = 0;
748  ogt[1] = 1;
749  ogt[2] = 0;
750  ogt[3] = 0;
751  ogt[4] = 0;
752  ogt[5] = -1;
753  }
754 
755  /* raster doesn't exist, create it now */
756  if (raster == NULL) {
757  raster = rt_raster_new(GDALGetRasterXSize(hdsOut), GDALGetRasterYSize(hdsOut));
758  if (rt_raster_is_empty(raster)) {
759  elog(ERROR, "RASTER_addBandOutDB: Cannot create new raster");
760  PG_RETURN_NULL();
761  }
762  rt_raster_set_geotransform_matrix(raster, ogt);
763 
764  if (rt_util_gdal_sr_auth_info(hdsOut, &authname, &authcode) == ES_NONE) {
765  if (
766  authname != NULL &&
767  strcmp(authname, "EPSG") == 0 &&
768  authcode != NULL
769  ) {
770  rt_raster_set_srid(raster, atoi(authcode));
771  }
772  else
773  elog(INFO, "Unknown SRS auth name and code from out-db file. Defaulting SRID of new raster to %d", SRID_UNKNOWN);
774  }
775  else
776  elog(INFO, "Cannot get SRS auth name and code from out-db file. Defaulting SRID of new raster to %d", SRID_UNKNOWN);
777  }
778 
779  /* some raster info */
780  width = rt_raster_get_width(raster);
781  height = rt_raster_get_height(raster);
782 
783  /* are rasters aligned? */
784  _rast = rt_raster_new(1, 1);
785  rt_raster_set_geotransform_matrix(_rast, ogt);
786  rt_raster_set_srid(_rast, rt_raster_get_srid(raster));
787  err = rt_raster_same_alignment(raster, _rast, &aligned, NULL);
788  rt_raster_destroy(_rast);
789 
790  if (err != ES_NONE) {
791  GDALClose(hdsOut);
792  if (raster != NULL)
793  rt_raster_destroy(raster);
794  if (pgraster != NULL)
795  PG_FREE_IF_COPY(pgraster, 0);
796  elog(ERROR, "RASTER_addBandOutDB: Cannot test alignment of out-db file");
797  return ES_ERROR;
798  }
799  else if (!aligned)
800  elog(WARNING, "The in-db representation of the out-db raster is not aligned. Band data may be incorrect");
801 
802  /* build up srcnband */
803  if (allbands) {
804  numsrcnband = GDALGetRasterCount(hdsOut);
805  GDALClose(hdsOut);
806 
807  srcnband = palloc(sizeof(int) * numsrcnband);
808  if (srcnband == NULL) {
809  if (raster != NULL)
810  rt_raster_destroy(raster);
811  if (pgraster != NULL)
812  PG_FREE_IF_COPY(pgraster, 0);
813  elog(ERROR, "RASTER_addBandOutDB: Cannot allocate memory for band indexes");
814  PG_RETURN_NULL();
815  }
816 
817  for (i = 0, j = 1; i < numsrcnband; i++, j++)
818  srcnband[i] = j;
819  }
820  else
821  GDALClose(hdsOut);
822 
823  /* add band */
824  for (i = 0, j = dstnband - 1; i < numsrcnband; i++, j++) {
825 
826  /* create band with path */
827  band = rt_band_new_offline_from_path(
828  width, height,
829  hasnodata, nodataval,
830  srcnband[i], outdbfile,
831  FALSE
832  );
833  if (band == NULL) {
834  if (raster != NULL)
835  rt_raster_destroy(raster);
836  if (pgraster != NULL)
837  PG_FREE_IF_COPY(pgraster, 0);
838  elog(ERROR, "RASTER_addBandOutDB: Cannot create new out-db band");
839  PG_RETURN_NULL();
840  }
841 
842  /* add band */
843  if (rt_raster_add_band(raster, band, j) < 0) {
844  if (raster != NULL)
845  rt_raster_destroy(raster);
846  if (pgraster != NULL)
847  PG_FREE_IF_COPY(pgraster, 0);
848  elog(ERROR, "RASTER_addBandOutDB: Cannot add new out-db band to raster");
849  PG_RETURN_NULL();
850  }
851  }
852 
853  pgrtn = rt_raster_serialize(raster);
854  rt_raster_destroy(raster);
855  if (pgraster != NULL)
856  PG_FREE_IF_COPY(pgraster, 0);
857  if (!pgrtn)
858  PG_RETURN_NULL();
859 
860  SET_VARSIZE(pgrtn, pgrtn->size);
861  PG_RETURN_POINTER(pgrtn);
862 }
863 
867 PG_FUNCTION_INFO_V1(RASTER_copyBand);
868 Datum RASTER_copyBand(PG_FUNCTION_ARGS)
869 {
870  rt_pgraster *pgto = NULL;
871  rt_pgraster *pgfrom = NULL;
872  rt_pgraster *pgrtn = NULL;
873  rt_raster torast = NULL;
874  rt_raster fromrast = NULL;
875  int toindex = 0;
876  int fromband = 0;
877  int oldtorastnumbands = 0;
878  int newtorastnumbands = 0;
879  int newbandindex = 0;
880 
881  /* Deserialize torast */
882  if (PG_ARGISNULL(0)) PG_RETURN_NULL();
883  pgto = (rt_pgraster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
884 
885  torast = rt_raster_deserialize(pgto, FALSE);
886  if (!torast) {
887  PG_FREE_IF_COPY(pgto, 0);
888  elog(ERROR, "RASTER_copyBand: Could not deserialize first raster");
889  PG_RETURN_NULL();
890  }
891 
892  /* Deserialize fromrast */
893  if (!PG_ARGISNULL(1)) {
894  pgfrom = (rt_pgraster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
895 
896  fromrast = rt_raster_deserialize(pgfrom, FALSE);
897  if (!fromrast) {
898  rt_raster_destroy(torast);
899  PG_FREE_IF_COPY(pgfrom, 1);
900  PG_FREE_IF_COPY(pgto, 0);
901  elog(ERROR, "RASTER_copyBand: Could not deserialize second raster");
902  PG_RETURN_NULL();
903  }
904 
905  oldtorastnumbands = rt_raster_get_num_bands(torast);
906 
907  if (PG_ARGISNULL(2))
908  fromband = 1;
909  else
910  fromband = PG_GETARG_INT32(2);
911 
912  if (PG_ARGISNULL(3))
913  toindex = oldtorastnumbands + 1;
914  else
915  toindex = PG_GETARG_INT32(3);
916 
917  /* Copy band fromrast torast */
918  newbandindex = rt_raster_copy_band(
919  torast, fromrast,
920  fromband - 1, toindex - 1
921  );
922 
923  newtorastnumbands = rt_raster_get_num_bands(torast);
924  if (newtorastnumbands == oldtorastnumbands || newbandindex == -1) {
925  elog(NOTICE, "RASTER_copyBand: Could not add band to raster. "
926  "Returning original raster."
927  );
928  }
929 
930  rt_raster_destroy(fromrast);
931  PG_FREE_IF_COPY(pgfrom, 1);
932  }
933 
934  /* Serialize and return torast */
935  pgrtn = rt_raster_serialize(torast);
936  rt_raster_destroy(torast);
937  PG_FREE_IF_COPY(pgto, 0);
938  if (!pgrtn) PG_RETURN_NULL();
939 
940  SET_VARSIZE(pgrtn, pgrtn->size);
941  PG_RETURN_POINTER(pgrtn);
942 }
943 
947 PG_FUNCTION_INFO_V1(RASTER_tile);
948 Datum RASTER_tile(PG_FUNCTION_ARGS)
949 {
950  FuncCallContext *funcctx;
951  int call_cntr;
952  int max_calls;
953  int i = 0;
954  int j = 0;
955 
956  struct tile_arg_t {
957 
958  struct {
959  rt_raster raster;
960  double gt[6];
961  int srid;
962  int width;
963  int height;
964  } raster;
965 
966  struct {
967  int width;
968  int height;
969 
970  int nx;
971  int ny;
972  } tile;
973 
974  int numbands;
975  int *nbands;
976 
977  struct {
978  int pad;
979  double hasnodata;
980  double nodataval;
981  } pad;
982  };
983  struct tile_arg_t *arg1 = NULL;
984  struct tile_arg_t *arg2 = NULL;
985 
986  if (SRF_IS_FIRSTCALL()) {
987  MemoryContext oldcontext;
988  rt_pgraster *pgraster = NULL;
989  int numbands;
990 
991  ArrayType *array;
992  Oid etype;
993  Datum *e;
994  bool *nulls;
995 
996  int16 typlen;
997  bool typbyval;
998  char typalign;
999 
1000  POSTGIS_RT_DEBUG(2, "RASTER_tile: first call");
1001 
1002  /* create a function context for cross-call persistence */
1003  funcctx = SRF_FIRSTCALL_INIT();
1004 
1005  /* switch to memory context appropriate for multiple function calls */
1006  oldcontext = MemoryContextSwitchTo(funcctx->multi_call_memory_ctx);
1007 
1008  /* Get input arguments */
1009  if (PG_ARGISNULL(0)) {
1010  MemoryContextSwitchTo(oldcontext);
1011  SRF_RETURN_DONE(funcctx);
1012  }
1013 
1014  /* allocate arg1 */
1015  arg1 = palloc(sizeof(struct tile_arg_t));
1016  if (arg1 == NULL) {
1017  MemoryContextSwitchTo(oldcontext);
1018  elog(ERROR, "RASTER_tile: Could not allocate memory for arguments");
1019  SRF_RETURN_DONE(funcctx);
1020  }
1021 
1022  pgraster = (rt_pgraster *) PG_DETOAST_DATUM_COPY(PG_GETARG_DATUM(0));
1023  arg1->raster.raster = rt_raster_deserialize(pgraster, FALSE);
1024  if (!arg1->raster.raster) {
1025  ereport(ERROR, (
1026  errcode(ERRCODE_OUT_OF_MEMORY),
1027  errmsg("Could not deserialize raster")
1028  ));
1029  pfree(arg1);
1030  PG_FREE_IF_COPY(pgraster, 0);
1031  MemoryContextSwitchTo(oldcontext);
1032  SRF_RETURN_DONE(funcctx);
1033  }
1034 
1035  /* raster has bands */
1036  numbands = rt_raster_get_num_bands(arg1->raster.raster);
1037  /*
1038  if (!numbands) {
1039  elog(NOTICE, "Raster provided has no bands");
1040  rt_raster_destroy(arg1->raster.raster);
1041  pfree(arg1);
1042  PG_FREE_IF_COPY(pgraster, 0);
1043  MemoryContextSwitchTo(oldcontext);
1044  SRF_RETURN_DONE(funcctx);
1045  }
1046  */
1047 
1048  /* width (1) */
1049  if (PG_ARGISNULL(1)) {
1050  elog(NOTICE, "Width cannot be NULL. Returning NULL");
1051  rt_raster_destroy(arg1->raster.raster);
1052  pfree(arg1);
1053  PG_FREE_IF_COPY(pgraster, 0);
1054  MemoryContextSwitchTo(oldcontext);
1055  SRF_RETURN_DONE(funcctx);
1056  }
1057  arg1->tile.width = PG_GETARG_INT32(1);
1058  if (arg1->tile.width < 1) {
1059  elog(NOTICE, "Width must be greater than zero. Returning NULL");
1060  rt_raster_destroy(arg1->raster.raster);
1061  pfree(arg1);
1062  PG_FREE_IF_COPY(pgraster, 0);
1063  MemoryContextSwitchTo(oldcontext);
1064  SRF_RETURN_DONE(funcctx);
1065  }
1066 
1067  /* height (2) */
1068  if (PG_ARGISNULL(2)) {
1069  elog(NOTICE, "Height cannot be NULL. Returning NULL");
1070  rt_raster_destroy(arg1->raster.raster);
1071  pfree(arg1);
1072  PG_FREE_IF_COPY(pgraster, 0);
1073  MemoryContextSwitchTo(oldcontext);
1074  SRF_RETURN_DONE(funcctx);
1075  }
1076  arg1->tile.height = PG_GETARG_INT32(2);
1077  if (arg1->tile.height < 1) {
1078  elog(NOTICE, "Height must be greater than zero. Returning NULL");
1079  rt_raster_destroy(arg1->raster.raster);
1080  pfree(arg1);
1081  PG_FREE_IF_COPY(pgraster, 0);
1082  MemoryContextSwitchTo(oldcontext);
1083  SRF_RETURN_DONE(funcctx);
1084  }
1085 
1086  /* nband, array (3) */
1087  if (numbands && !PG_ARGISNULL(3)) {
1088  array = PG_GETARG_ARRAYTYPE_P(3);
1089  etype = ARR_ELEMTYPE(array);
1090  get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
1091 
1092  switch (etype) {
1093  case INT2OID:
1094  case INT4OID:
1095  break;
1096  default:
1097  rt_raster_destroy(arg1->raster.raster);
1098  pfree(arg1);
1099  PG_FREE_IF_COPY(pgraster, 0);
1100  MemoryContextSwitchTo(oldcontext);
1101  elog(ERROR, "RASTER_tile: Invalid data type for band indexes");
1102  SRF_RETURN_DONE(funcctx);
1103  break;
1104  }
1105 
1106  deconstruct_array(array, etype, typlen, typbyval, typalign, &e, &nulls, &(arg1->numbands));
1107 
1108  arg1->nbands = palloc(sizeof(int) * arg1->numbands);
1109  if (arg1->nbands == NULL) {
1110  rt_raster_destroy(arg1->raster.raster);
1111  pfree(arg1);
1112  PG_FREE_IF_COPY(pgraster, 0);
1113  MemoryContextSwitchTo(oldcontext);
1114  elog(ERROR, "RASTER_tile: Could not allocate memory for band indexes");
1115  SRF_RETURN_DONE(funcctx);
1116  }
1117 
1118  for (i = 0, j = 0; i < arg1->numbands; i++) {
1119  if (nulls[i]) continue;
1120 
1121  switch (etype) {
1122  case INT2OID:
1123  arg1->nbands[j] = DatumGetInt16(e[i]) - 1;
1124  break;
1125  case INT4OID:
1126  arg1->nbands[j] = DatumGetInt32(e[i]) - 1;
1127  break;
1128  }
1129 
1130  j++;
1131  }
1132 
1133  if (j < arg1->numbands) {
1134  arg1->nbands = repalloc(arg1->nbands, sizeof(int) * j);
1135  if (arg1->nbands == NULL) {
1136  rt_raster_destroy(arg1->raster.raster);
1137  pfree(arg1);
1138  PG_FREE_IF_COPY(pgraster, 0);
1139  MemoryContextSwitchTo(oldcontext);
1140  elog(ERROR, "RASTER_tile: Could not reallocate memory for band indexes");
1141  SRF_RETURN_DONE(funcctx);
1142  }
1143 
1144  arg1->numbands = j;
1145  }
1146 
1147  /* validate nbands */
1148  for (i = 0; i < arg1->numbands; i++) {
1149  if (!rt_raster_has_band(arg1->raster.raster, arg1->nbands[i])) {
1150  elog(NOTICE, "Band at index %d not found in raster", arg1->nbands[i] + 1);
1151  rt_raster_destroy(arg1->raster.raster);
1152  pfree(arg1->nbands);
1153  pfree(arg1);
1154  PG_FREE_IF_COPY(pgraster, 0);
1155  MemoryContextSwitchTo(oldcontext);
1156  SRF_RETURN_DONE(funcctx);
1157  }
1158  }
1159  }
1160  else {
1161  arg1->numbands = numbands;
1162 
1163  if (numbands) {
1164  arg1->nbands = palloc(sizeof(int) * arg1->numbands);
1165 
1166  if (arg1->nbands == NULL) {
1167  rt_raster_destroy(arg1->raster.raster);
1168  pfree(arg1);
1169  PG_FREE_IF_COPY(pgraster, 0);
1170  MemoryContextSwitchTo(oldcontext);
1171  elog(ERROR, "RASTER_dumpValues: Could not allocate memory for pixel values");
1172  SRF_RETURN_DONE(funcctx);
1173  }
1174 
1175  for (i = 0; i < arg1->numbands; i++) {
1176  arg1->nbands[i] = i;
1177  POSTGIS_RT_DEBUGF(4, "arg1->nbands[%d] = %d", arg1->nbands[i], i);
1178  }
1179  }
1180  }
1181 
1182  /* pad (4) and padnodata (5) */
1183  if (!PG_ARGISNULL(4)) {
1184  arg1->pad.pad = PG_GETARG_BOOL(4) ? 1 : 0;
1185 
1186  if (arg1->pad.pad && !PG_ARGISNULL(5)) {
1187  arg1->pad.hasnodata = 1;
1188  arg1->pad.nodataval = PG_GETARG_FLOAT8(5);
1189  }
1190  else {
1191  arg1->pad.hasnodata = 0;
1192  arg1->pad.nodataval = 0;
1193  }
1194  }
1195  else {
1196  arg1->pad.pad = 0;
1197  arg1->pad.hasnodata = 0;
1198  arg1->pad.nodataval = 0;
1199  }
1200 
1201  /* store some additional metadata */
1202  arg1->raster.srid = rt_raster_get_srid(arg1->raster.raster);
1203  arg1->raster.width = rt_raster_get_width(arg1->raster.raster);
1204  arg1->raster.height = rt_raster_get_height(arg1->raster.raster);
1205  rt_raster_get_geotransform_matrix(arg1->raster.raster, arg1->raster.gt);
1206 
1207  /* determine maximum number of tiles from raster */
1208  arg1->tile.nx = ceil(arg1->raster.width / (double) arg1->tile.width);
1209  arg1->tile.ny = ceil(arg1->raster.height / (double) arg1->tile.height);
1210  POSTGIS_RT_DEBUGF(4, "# of tiles (x, y) = (%d, %d)", arg1->tile.nx, arg1->tile.ny);
1211 
1212  /* Store needed information */
1213  funcctx->user_fctx = arg1;
1214 
1215  /* total number of tuples to be returned */
1216  funcctx->max_calls = (arg1->tile.nx * arg1->tile.ny);
1217 
1218  MemoryContextSwitchTo(oldcontext);
1219  }
1220 
1221  /* stuff done on every call of the function */
1222  funcctx = SRF_PERCALL_SETUP();
1223 
1224  call_cntr = funcctx->call_cntr;
1225  max_calls = funcctx->max_calls;
1226  arg2 = funcctx->user_fctx;
1227 
1228  /* do when there is more left to send */
1229  if (call_cntr < max_calls) {
1230  rt_pgraster *pgtile = NULL;
1231  rt_raster tile = NULL;
1232  rt_band _band = NULL;
1233  rt_band band = NULL;
1234  rt_pixtype pixtype = PT_END;
1235  int hasnodata = 0;
1236  double nodataval = 0;
1237  int width = 0;
1238  int height = 0;
1239 
1240  int k = 0;
1241  int tx = 0;
1242  int ty = 0;
1243  int rx = 0;
1244  int ry = 0;
1245  int ex = 0; /* edge tile on right */
1246  int ey = 0; /* edge tile on bottom */
1247  double ulx = 0;
1248  double uly = 0;
1249  uint16_t len = 0;
1250  void *vals = NULL;
1251  uint16_t nvals;
1252 
1253  POSTGIS_RT_DEBUGF(3, "call number %d", call_cntr);
1254 
1255  /*
1256  find offset based upon tile #
1257 
1258  0 1 2
1259  3 4 5
1260  6 7 8
1261  */
1262  ty = call_cntr / arg2->tile.nx;
1263  tx = call_cntr % arg2->tile.nx;
1264  POSTGIS_RT_DEBUGF(4, "tile (x, y) = (%d, %d)", tx, ty);
1265 
1266  /* edge tile? only important if padding is false */
1267  if (!arg2->pad.pad) {
1268  if (ty + 1 == arg2->tile.ny)
1269  ey = 1;
1270  if (tx + 1 == arg2->tile.nx)
1271  ex = 1;
1272  }
1273 
1274  /* upper-left of tile in raster coordinates */
1275  rx = tx * arg2->tile.width;
1276  ry = ty * arg2->tile.height;
1277  POSTGIS_RT_DEBUGF(4, "raster coordinates = %d, %d", rx, ry);
1278 
1279  /* determine tile width and height */
1280  /* default to user-defined */
1281  width = arg2->tile.width;
1282  height = arg2->tile.height;
1283 
1284  /* override user-defined if edge tile (only possible if padding is false */
1285  if (ex || ey) {
1286  /* right edge */
1287  if (ex)
1288  width = arg2->raster.width - rx;
1289  /* bottom edge */
1290  if (ey)
1291  height = arg2->raster.height - ry;
1292  }
1293 
1294  /* create empty raster */
1295  tile = rt_raster_new(width, height);
1296  rt_raster_set_geotransform_matrix(tile, arg2->raster.gt);
1297  rt_raster_set_srid(tile, arg2->raster.srid);
1298 
1299  /* upper-left of tile in spatial coordinates */
1300  if (rt_raster_cell_to_geopoint(arg2->raster.raster, rx, ry, &ulx, &uly, arg2->raster.gt) != ES_NONE) {
1301  rt_raster_destroy(tile);
1302  rt_raster_destroy(arg2->raster.raster);
1303  if (arg2->numbands) pfree(arg2->nbands);
1304  pfree(arg2);
1305  elog(ERROR, "RASTER_tile: Could not compute the coordinates of the upper-left corner of the output tile");
1306  SRF_RETURN_DONE(funcctx);
1307  }
1308  rt_raster_set_offsets(tile, ulx, uly);
1309  POSTGIS_RT_DEBUGF(4, "spatial coordinates = %f, %f", ulx, uly);
1310 
1311  /* compute length of pixel line to read */
1312  len = arg2->tile.width;
1313  if (rx + arg2->tile.width >= arg2->raster.width)
1314  len = arg2->raster.width - rx;
1315  POSTGIS_RT_DEBUGF(3, "read line len = %d", len);
1316 
1317  /* copy bands to tile */
1318  for (i = 0; i < arg2->numbands; i++) {
1319  POSTGIS_RT_DEBUGF(4, "copying band %d to tile %d", arg2->nbands[i], call_cntr);
1320 
1321  _band = rt_raster_get_band(arg2->raster.raster, arg2->nbands[i]);
1322  if (_band == NULL) {
1323  int nband = arg2->nbands[i] + 1;
1324  rt_raster_destroy(tile);
1325  rt_raster_destroy(arg2->raster.raster);
1326  pfree(arg2->nbands);
1327  pfree(arg2);
1328  elog(ERROR, "RASTER_tile: Could not get band %d from source raster", nband);
1329  SRF_RETURN_DONE(funcctx);
1330  }
1331 
1332  pixtype = rt_band_get_pixtype(_band);
1333  hasnodata = rt_band_get_hasnodata_flag(_band);
1334  if (hasnodata)
1335  rt_band_get_nodata(_band, &nodataval);
1336  else if (arg2->pad.pad && arg2->pad.hasnodata) {
1337  hasnodata = 1;
1338  nodataval = arg2->pad.nodataval;
1339  }
1340  else
1341  nodataval = rt_band_get_min_value(_band);
1342 
1343  /* inline band */
1344  if (!rt_band_is_offline(_band)) {
1345  if (rt_raster_generate_new_band(tile, pixtype, nodataval, hasnodata, nodataval, i) < 0) {
1346  rt_raster_destroy(tile);
1347  rt_raster_destroy(arg2->raster.raster);
1348  pfree(arg2->nbands);
1349  pfree(arg2);
1350  elog(ERROR, "RASTER_tile: Could not add new band to output tile");
1351  SRF_RETURN_DONE(funcctx);
1352  }
1353  band = rt_raster_get_band(tile, i);
1354  if (band == NULL) {
1355  rt_raster_destroy(tile);
1356  rt_raster_destroy(arg2->raster.raster);
1357  pfree(arg2->nbands);
1358  pfree(arg2);
1359  elog(ERROR, "RASTER_tile: Could not get newly added band from output tile");
1360  SRF_RETURN_DONE(funcctx);
1361  }
1362 
1363  /* if isnodata, set flag and continue */
1364  if (rt_band_get_isnodata_flag(_band)) {
1365  rt_band_set_isnodata_flag(band, 1);
1366  continue;
1367  }
1368 
1369  /* copy data */
1370  for (j = 0; j < arg2->tile.height; j++) {
1371  k = ry + j;
1372 
1373  if (k >= arg2->raster.height) {
1374  POSTGIS_RT_DEBUGF(4, "row %d is beyond extent of source raster. skipping", k);
1375  continue;
1376  }
1377 
1378  POSTGIS_RT_DEBUGF(4, "getting pixel line %d, %d for %d pixels", rx, k, len);
1379  if (rt_band_get_pixel_line(_band, rx, k, len, &vals, &nvals) != ES_NONE) {
1380  rt_raster_destroy(tile);
1381  rt_raster_destroy(arg2->raster.raster);
1382  pfree(arg2->nbands);
1383  pfree(arg2);
1384  elog(ERROR, "RASTER_tile: Could not get pixel line from source raster");
1385  SRF_RETURN_DONE(funcctx);
1386  }
1387 
1388  if (nvals && rt_band_set_pixel_line(band, 0, j, vals, nvals) != ES_NONE) {
1389  rt_raster_destroy(tile);
1390  rt_raster_destroy(arg2->raster.raster);
1391  pfree(arg2->nbands);
1392  pfree(arg2);
1393  elog(ERROR, "RASTER_tile: Could not set pixel line of output tile");
1394  SRF_RETURN_DONE(funcctx);
1395  }
1396  }
1397  }
1398  /* offline */
1399  else {
1400  uint8_t bandnum = 0;
1401  rt_band_get_ext_band_num(_band, &bandnum);
1402 
1403  band = rt_band_new_offline(
1404  width, height,
1405  pixtype,
1406  hasnodata, nodataval,
1407  bandnum, rt_band_get_ext_path(_band)
1408  );
1409 
1410  if (band == NULL) {
1411  rt_raster_destroy(tile);
1412  rt_raster_destroy(arg2->raster.raster);
1413  pfree(arg2->nbands);
1414  pfree(arg2);
1415  elog(ERROR, "RASTER_tile: Could not create new offline band for output tile");
1416  SRF_RETURN_DONE(funcctx);
1417  }
1418 
1419  if (rt_raster_add_band(tile, band, i) < 0) {
1420  rt_band_destroy(band);
1421  rt_raster_destroy(tile);
1422  rt_raster_destroy(arg2->raster.raster);
1423  pfree(arg2->nbands);
1424  pfree(arg2);
1425  elog(ERROR, "RASTER_tile: Could not add new offline band to output tile");
1426  SRF_RETURN_DONE(funcctx);
1427  }
1428  }
1429  }
1430 
1431  pgtile = rt_raster_serialize(tile);
1432  rt_raster_destroy(tile);
1433  if (!pgtile) {
1434  rt_raster_destroy(arg2->raster.raster);
1435  if (arg2->numbands) pfree(arg2->nbands);
1436  pfree(arg2);
1437  SRF_RETURN_DONE(funcctx);
1438  }
1439 
1440  SET_VARSIZE(pgtile, pgtile->size);
1441  SRF_RETURN_NEXT(funcctx, PointerGetDatum(pgtile));
1442  }
1443  /* do when there is no more left */
1444  else {
1445  rt_raster_destroy(arg2->raster.raster);
1446  if (arg2->numbands) pfree(arg2->nbands);
1447  pfree(arg2);
1448  SRF_RETURN_DONE(funcctx);
1449  }
1450 }
1451 
1456 PG_FUNCTION_INFO_V1(RASTER_band);
1457 Datum RASTER_band(PG_FUNCTION_ARGS)
1458 {
1459  rt_pgraster *pgraster;
1460  rt_pgraster *pgrast;
1461  rt_raster raster;
1462  rt_raster rast;
1463 
1464  bool skip = FALSE;
1465  ArrayType *array;
1466  Oid etype;
1467  Datum *e;
1468  bool *nulls;
1469  int16 typlen;
1470  bool typbyval;
1471  char typalign;
1472 
1473  uint32_t numBands;
1474  uint32_t *bandNums;
1475  uint32 idx = 0;
1476  int n;
1477  int i = 0;
1478  int j = 0;
1479 
1480  if (PG_ARGISNULL(0))
1481  PG_RETURN_NULL();
1482  pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
1483 
1484  raster = rt_raster_deserialize(pgraster, FALSE);
1485  if (!raster) {
1486  PG_FREE_IF_COPY(pgraster, 0);
1487  elog(ERROR, "RASTER_band: Could not deserialize raster");
1488  PG_RETURN_NULL();
1489  }
1490 
1491  /* process bandNums */
1492  if (PG_ARGISNULL(1)) {
1493  elog(NOTICE, "Band number(s) not provided. Returning original raster");
1494  skip = TRUE;
1495  }
1496  if (!skip) {
1497  numBands = rt_raster_get_num_bands(raster);
1498 
1499  array = PG_GETARG_ARRAYTYPE_P(1);
1500  etype = ARR_ELEMTYPE(array);
1501  get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
1502 
1503  switch (etype) {
1504  case INT2OID:
1505  case INT4OID:
1506  break;
1507  default:
1508  rt_raster_destroy(raster);
1509  PG_FREE_IF_COPY(pgraster, 0);
1510  elog(ERROR, "RASTER_band: Invalid data type for band number(s)");
1511  PG_RETURN_NULL();
1512  break;
1513  }
1514 
1515  deconstruct_array(array, etype, typlen, typbyval, typalign, &e,
1516  &nulls, &n);
1517 
1518  bandNums = palloc(sizeof(uint32_t) * n);
1519  for (i = 0, j = 0; i < n; i++) {
1520  if (nulls[i]) continue;
1521 
1522  switch (etype) {
1523  case INT2OID:
1524  idx = (uint32_t) DatumGetInt16(e[i]);
1525  break;
1526  case INT4OID:
1527  idx = (uint32_t) DatumGetInt32(e[i]);
1528  break;
1529  }
1530 
1531  POSTGIS_RT_DEBUGF(3, "band idx (before): %d", idx);
1532  if (idx > numBands || idx < 1) {
1533  elog(NOTICE, "Invalid band index (must use 1-based). Returning original raster");
1534  skip = TRUE;
1535  break;
1536  }
1537 
1538  bandNums[j] = idx - 1;
1539  POSTGIS_RT_DEBUGF(3, "bandNums[%d] = %d", j, bandNums[j]);
1540  j++;
1541  }
1542 
1543  if (skip || j < 1) {
1544  pfree(bandNums);
1545  skip = TRUE;
1546  }
1547  }
1548 
1549  if (!skip) {
1550  rast = rt_raster_from_band(raster, bandNums, j);
1551  pfree(bandNums);
1552  rt_raster_destroy(raster);
1553  PG_FREE_IF_COPY(pgraster, 0);
1554  if (!rast) {
1555  elog(ERROR, "RASTER_band: Could not create new raster");
1556  PG_RETURN_NULL();
1557  }
1558 
1559  pgrast = rt_raster_serialize(rast);
1560  rt_raster_destroy(rast);
1561 
1562  if (!pgrast)
1563  PG_RETURN_NULL();
1564 
1565  SET_VARSIZE(pgrast, pgrast->size);
1566  PG_RETURN_POINTER(pgrast);
1567  }
1568 
1569  PG_RETURN_POINTER(pgraster);
1570 }
TRUE
#define TRUE
Definition: dbfopen.c:169
FALSE
#define FALSE
Definition: dbfopen.c:168
SRID_UNKNOWN
#define SRID_UNKNOWN
Unknown SRID value.
Definition: liblwgeom.h:188
rt_band_set_isnodata_flag
rt_errorstate rt_band_set_isnodata_flag(rt_band band, int flag)
Set isnodata flag value.
Definition: rt_band.c:695
rt_raster_cell_to_geopoint
rt_errorstate rt_raster_cell_to_geopoint(rt_raster raster, double xr, double yr, double *xw, double *yw, double *gt)
Convert an xr, yr raster point to an xw, yw point on map.
Definition: rt_raster.c:755
rt_raster_get_srid
int32_t rt_raster_get_srid(rt_raster raster)
Get raster's SRID.
Definition: rt_raster.c:356
rt_util_gdal_register_all
int rt_util_gdal_register_all(int force_register_all)
Definition: rt_util.c:336
rt_raster_generate_new_band
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_band_get_ext_path
const char * rt_band_get_ext_path(rt_band band)
Return band's external path (only valid when rt_band_is_offline returns non-zero).
Definition: rt_band.c:363
rt_raster_set_geotransform_matrix
void rt_raster_set_geotransform_matrix(rt_raster raster, double *gt)
Set raster's geotransform using 6-element array.
Definition: rt_raster.c:727
rt_raster_set_scale
void rt_raster_set_scale(rt_raster raster, double scaleX, double scaleY)
Set scale in projection units.
Definition: rt_raster.c:137
rt_raster_add_band
int rt_raster_add_band(rt_raster raster, rt_band band, int index)
Add band data to a raster.
Definition: rt_raster.c:405
rt_band_get_hasnodata_flag
int rt_band_get_hasnodata_flag(rt_band band)
Get hasnodata flag value.
Definition: rt_band.c:674
rt_pixtype_index_from_name
rt_pixtype rt_pixtype_index_from_name(const char *pixname)
Definition: rt_pixel.c:80
rt_raster_destroy
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_raster.c:82
rt_pixtype
rt_pixtype
Definition: librtcore.h:185
PT_END
@ PT_END
Definition: librtcore.h:197
rt_raster_set_skews
void rt_raster_set_skews(rt_raster raster, double skewX, double skewY)
Set skews about the X and Y axis.
Definition: rt_raster.c:168
rt_band_get_isnodata_flag
int rt_band_get_isnodata_flag(rt_band band)
Get isnodata flag value.
Definition: rt_band.c:714
rt_raster_new
rt_raster rt_raster_new(uint32_t width, uint32_t height)
Construct a raster with given dimensions.
Definition: rt_raster.c:48
rt_raster_has_band
int rt_raster_has_band(rt_raster raster, int nband)
Return TRUE if the raster has a band of this number.
Definition: rt_raster.c:1347
rt_raster_serialize
void * rt_raster_serialize(rt_raster raster)
Return this raster in serialized form.
Definition: rt_serialize.c:521
rt_band_get_min_value
double rt_band_get_min_value(rt_band band)
Returns the minimal possible value for the band according to the pixel type.
Definition: rt_band.c:1745
rt_util_gdal_open
GDALDatasetH rt_util_gdal_open(const char *fn, GDALAccess fn_access, int shared)
Definition: rt_util.c:381
ES_NONE
@ ES_NONE
Definition: librtcore.h:180
ES_ERROR
@ ES_ERROR
Definition: librtcore.h:181
rt_band_set_pixel_line
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:853
rt_band_destroy
void rt_band_destroy(rt_band band)
Destroy a raster band.
Definition: rt_band.c:340
rt_raster_get_num_bands
uint16_t rt_raster_get_num_bands(rt_raster raster)
Definition: rt_raster.c:372
rt_raster_get_height
uint16_t rt_raster_get_height(rt_raster raster)
Definition: rt_raster.c:129
rt_band_new_offline
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_band.c:124
rt_band_get_ext_band_num
rt_errorstate rt_band_get_ext_band_num(rt_band band, uint8_t *bandnum)
Return bands' external band number (only valid when rt_band_is_offline returns non-zero).
Definition: rt_band.c:376
rt_raster_set_srid
void rt_raster_set_srid(rt_raster raster, int32_t srid)
Set raster's SRID.
Definition: rt_raster.c:363
rt_raster_from_band
rt_raster rt_raster_from_band(rt_raster raster, uint32_t *bandNums, int count)
Construct a new rt_raster from an existing rt_raster and an array of band numbers.
Definition: rt_raster.c:1435
rt_band_get_nodata
rt_errorstate rt_band_get_nodata(rt_band band, double *nodata)
Get NODATA value.
Definition: rt_band.c:1730
rt_band_get_pixtype
rt_pixtype rt_band_get_pixtype(rt_band band)
Return pixeltype of this band.
Definition: rt_band.c:631
rt_pixtype_name
const char * rt_pixtype_name(rt_pixtype pixtype)
Definition: rt_pixel.c:110
rt_band_new_offline_from_path
rt_band rt_band_new_offline_from_path(uint16_t width, uint16_t height, int hasnodata, double nodataval, uint8_t bandNum, const char *path, int force)
Create an out-db rt_band from path.
Definition: rt_band.c:199
rt_raster_get_width
uint16_t rt_raster_get_width(rt_raster raster)
Definition: rt_raster.c:121
rt_util_gdal_sr_auth_info
rt_errorstate rt_util_gdal_sr_auth_info(GDALDatasetH hds, char **authname, char **authcode)
Get auth name and code.
Definition: rt_util.c:270
rt_band_is_offline
int rt_band_is_offline(rt_band band)
Return non-zero if the given band data is on the filesystem.
Definition: rt_band.c:329
rt_raster_copy_band
int rt_raster_copy_band(rt_raster torast, rt_raster fromrast, int fromindex, int toindex)
Copy one band from one raster to another.
Definition: rt_raster.c:1370
rt_raster_same_alignment
rt_errorstate rt_raster_same_alignment(rt_raster rast1, rt_raster rast2, int *aligned, char **reason)
Definition: rt_spatial_relationship.c:45
rt_raster_get_geotransform_matrix
void rt_raster_get_geotransform_matrix(rt_raster raster, double *gt)
Get 6-element array of raster geotransform matrix.
Definition: rt_raster.c:706
rt_raster_set_offsets
void rt_raster_set_offsets(rt_raster raster, double x, double y)
Set insertion points in projection units.
Definition: rt_raster.c:199
rt_raster_deserialize
rt_raster rt_raster_deserialize(void *serialized, int header_only)
Return a raster from a serialized form.
Definition: rt_serialize.c:725
rt_band_get_pixel_line
rt_errorstate rt_band_get_pixel_line(rt_band band, int x, int y, uint16_t len, void **vals, uint16_t *nvals)
Get values of multiple pixels.
Definition: rt_band.c:1137
rt_raster_is_empty
int rt_raster_is_empty(rt_raster raster)
Return TRUE if the raster is empty.
Definition: rt_raster.c:1334
rt_raster_get_band
rt_band rt_raster_get_band(rt_raster raster, int bandNum)
Return Nth band, or NULL if unavailable.
Definition: rt_raster.c:381
ovdump.band
band
Definition: ovdump.py:57
pixval.nband
nband
Definition: pixval.py:52
rtrowdump.raster
raster
Be careful!! Zeros function's input parameter can be a (height x width) array, not (width x height): ...
Definition: rtrowdump.py:121
window.gt
gt
Definition: window.py:77
text_to_cstring
char * text_to_cstring(const text *textptr)
PG_FUNCTION_INFO_V1
PG_FUNCTION_INFO_V1(RASTER_makeEmpty)
Make a new raster with no bands.
RASTER_tile
Datum RASTER_tile(PG_FUNCTION_ARGS)
Definition: rtpg_create.c:948
RASTER_addBand
Datum RASTER_addBand(PG_FUNCTION_ARGS)
Definition: rtpg_create.c:121
RASTER_band
Datum RASTER_band(PG_FUNCTION_ARGS)
Definition: rtpg_create.c:1457
RASTER_makeEmpty
Datum RASTER_makeEmpty(PG_FUNCTION_ARGS)
Definition: rtpg_create.c:55
RASTER_copyBand
Datum RASTER_copyBand(PG_FUNCTION_ARGS)
Definition: rtpg_create.c:868
RASTER_addBandOutDB
Datum RASTER_addBandOutDB(PG_FUNCTION_ARGS)
Definition: rtpg_create.c:557
RASTER_addBandRasterArray
Datum RASTER_addBandRasterArray(PG_FUNCTION_ARGS)
Definition: rtpg_create.c:344
POSTGIS_RT_DEBUG
#define POSTGIS_RT_DEBUG(level, msg)
Definition: rtpostgis.h:61
POSTGIS_RT_DEBUGF
#define POSTGIS_RT_DEBUGF(level, msg,...)
Definition: rtpostgis.h:65
rt_band_t
Definition: librtcore.h:2313
rt_raster_serialized_t::size
uint32_t size
Definition: librtcore.h:2252
rt_raster_serialized_t
Struct definitions.
Definition: librtcore.h:2250
rt_raster_t::srid
int32_t srid
Definition: librtcore.h:2300
rt_raster_t
Definition: librtcore.h:2284
uint32_t
unsigned int uint32_t
Definition: uthash.h:78
uint8_t
unsigned char uint8_t
Definition: uthash.h:79