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
51 /**
52 * Make a new raster with no bands
53 */
54 PG_FUNCTION_INFO_V1(RASTER_makeEmpty);
RASTER_makeEmpty(PG_FUNCTION_ARGS)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
117 /**
118 * Add band(s) to the given raster at the given position(s).
119 */
120 PG_FUNCTION_INFO_V1(RASTER_addBand);
RASTER_addBand(PG_FUNCTION_ARGS)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
340 /**
341 * Add bands from array of rasters to a destination raster
342 */
343 PG_FUNCTION_INFO_V1(RASTER_addBandRasterArray);
RASTER_addBandRasterArray(PG_FUNCTION_ARGS)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 /* source rasters' band index, 1-based */
387 if (!PG_ARGISNULL(2))
388 srcnband = PG_GETARG_INT32(2);
389 if (srcnband < 1) {
390 elog(NOTICE, "Invalid band index for source rasters (must be 1-based). Returning original raster");
391 if (raster != NULL) {
392 rt_raster_destroy(raster);
393 PG_RETURN_POINTER(pgraster);
394 }
395 else
396 PG_RETURN_NULL();
397 }
398 POSTGIS_RT_DEBUGF(4, "srcnband = %d", srcnband);
399
400 /* destination raster's band index, 1-based */
401 if (!PG_ARGISNULL(3)) {
402 dstnband = PG_GETARG_INT32(3);
403 appendband = FALSE;
404
405 if (dstnband < 1) {
406 elog(NOTICE, "Invalid band index for destination raster (must be 1-based). Returning original raster");
407 if (raster != NULL) {
408 rt_raster_destroy(raster);
409 PG_RETURN_POINTER(pgraster);
410 }
411 else
412 PG_RETURN_NULL();
413 }
414 }
415 else
416 appendband = TRUE;
417
418 /* additional processing of dstnband */
419 if (raster != NULL) {
420 dstnumbands = rt_raster_get_num_bands(raster);
421
422 if (dstnumbands < 1) {
423 appendband = TRUE;
424 dstnband = 1;
425 }
426 else if (appendband)
427 dstnband = dstnumbands + 1;
428 else if (dstnband > dstnumbands) {
429 elog(NOTICE, "Band index provided for destination raster is greater than the number of bands in the raster. Bands will be appended");
430 appendband = TRUE;
431 dstnband = dstnumbands + 1;
432 }
433 }
434 POSTGIS_RT_DEBUGF(4, "appendband = %d", appendband);
435 POSTGIS_RT_DEBUGF(4, "dstnband = %d", dstnband);
436
437 /* process set of source rasters */
438 POSTGIS_RT_DEBUG(3, "Processing array of source rasters");
439 array = PG_GETARG_ARRAYTYPE_P(1);
440 etype = ARR_ELEMTYPE(array);
441 get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
442
443 deconstruct_array(array, etype, typlen, typbyval, typalign, &e,
444 &nulls, &n);
445
446 /* decrement srcnband and dstnband by 1, now 0-based */
447 srcnband--;
448 dstnband--;
449 POSTGIS_RT_DEBUGF(4, "0-based nband (src, dst) = (%d, %d)", srcnband, dstnband);
450
451 /* time to copy bands */
452 for (i = 0; i < n; i++) {
453 if (nulls[i]) continue;
454 src = NULL;
455
456 pgsrc = (rt_pgraster *) PG_DETOAST_DATUM(e[i]);
457 src = rt_raster_deserialize(pgsrc, FALSE);
458 if (src == NULL) {
459 pfree(nulls);
460 pfree(e);
461 if (raster != NULL)
462 rt_raster_destroy(raster);
463 if (pgraster != NULL)
464 PG_FREE_IF_COPY(pgraster, 0);
465 elog(ERROR, "RASTER_addBandRasterArray: Could not deserialize source raster at index %d", i + 1);
466 PG_RETURN_NULL();
467 }
468
469 srcnumbands = rt_raster_get_num_bands(src);
470 POSTGIS_RT_DEBUGF(4, "source raster %d has %d bands", i + 1, srcnumbands);
471
472 /* band index isn't valid */
473 if (srcnband > srcnumbands - 1) {
474 elog(NOTICE, "Invalid band index for source raster at index %d. Returning original raster", i + 1);
475 pfree(nulls);
476 pfree(e);
477 rt_raster_destroy(src);
478 if (raster != NULL) {
479 rt_raster_destroy(raster);
480 PG_RETURN_POINTER(pgraster);
481 }
482 else
483 PG_RETURN_NULL();
484 }
485
486 /* destination raster is empty, new raster */
487 if (raster == NULL) {
488 uint32_t srcnbands[1] = {srcnband};
489
490 POSTGIS_RT_DEBUG(4, "empty destination raster, using rt_raster_from_band");
491
492 raster = rt_raster_from_band(src, srcnbands, 1);
493 rt_raster_destroy(src);
494 if (raster == NULL) {
495 pfree(nulls);
496 pfree(e);
497 if (pgraster != NULL)
498 PG_FREE_IF_COPY(pgraster, 0);
499 elog(ERROR, "RASTER_addBandRasterArray: Could not create raster from source raster at index %d", i + 1);
500 PG_RETURN_NULL();
501 }
502 }
503 /* copy band */
504 else {
505 rtn = rt_raster_copy_band(
506 raster, src,
507 srcnband, dstnband
508 );
509 rt_raster_destroy(src);
510
511 if (rtn == -1 || rt_raster_get_num_bands(raster) == dstnumbands) {
512 elog(NOTICE, "Could not add band from source raster at index %d to destination raster. Returning original raster", i + 1);
513 rt_raster_destroy(raster);
514 pfree(nulls);
515 pfree(e);
516 if (pgraster != NULL)
517 PG_RETURN_POINTER(pgraster);
518 else
519 PG_RETURN_NULL();
520 }
521 }
522
523 dstnband++;
524 dstnumbands++;
525 }
526
527 if (raster != NULL) {
528 pgrtn = rt_raster_serialize(raster);
529 rt_raster_destroy(raster);
530 if (pgraster != NULL)
531 PG_FREE_IF_COPY(pgraster, 0);
532 if (!pgrtn)
533 PG_RETURN_NULL();
534
535 SET_VARSIZE(pgrtn, pgrtn->size);
536 PG_RETURN_POINTER(pgrtn);
537 }
538
539 PG_RETURN_NULL();
540 }
541
542 /**
543 * Add out-db band to the given raster at the given position
544 */
545 PG_FUNCTION_INFO_V1(RASTER_addBandOutDB);
RASTER_addBandOutDB(PG_FUNCTION_ARGS)546 Datum RASTER_addBandOutDB(PG_FUNCTION_ARGS)
547 {
548 rt_pgraster *pgraster = NULL;
549 rt_pgraster *pgrtn = NULL;
550
551 rt_raster raster = NULL;
552 rt_band band = NULL;
553 int numbands = 0;
554 int dstnband = 1; /* 1-based */
555 int appendband = FALSE;
556 char *outdbfile = NULL;
557 int *srcnband = NULL; /* 1-based */
558 int numsrcnband = 0;
559 int allbands = FALSE;
560 int hasnodata = FALSE;
561 double nodataval = 0.;
562 uint16_t width = 0;
563 uint16_t height = 0;
564 char *authname = NULL;
565 char *authcode = NULL;
566
567 int i = 0;
568 int j = 0;
569
570 GDALDatasetH hdsOut;
571
572 double ogt[6] = {0.};
573 rt_raster _rast = NULL;
574 int aligned = 0;
575 int err = 0;
576
577 /* destination raster */
578 if (!PG_ARGISNULL(0)) {
579 pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
580
581 /* raster */
582 raster = rt_raster_deserialize(pgraster, FALSE);
583 if (!raster) {
584 PG_FREE_IF_COPY(pgraster, 0);
585 elog(ERROR, "RASTER_addBandOutDB: Cannot deserialize destination raster");
586 PG_RETURN_NULL();
587 }
588
589 POSTGIS_RT_DEBUG(4, "destination raster isn't NULL");
590 }
591
592 /* destination band index (1) */
593 if (!PG_ARGISNULL(1))
594 dstnband = PG_GETARG_INT32(1);
595 else
596 appendband = TRUE;
597
598 /* outdb file (2) */
599 if (PG_ARGISNULL(2)) {
600 elog(NOTICE, "Out-db raster file not provided. Returning original raster");
601 if (pgraster != NULL) {
602 rt_raster_destroy(raster);
603 PG_RETURN_POINTER(pgraster);
604 }
605 else
606 PG_RETURN_NULL();
607 }
608 else {
609 outdbfile = text_to_cstring(PG_GETARG_TEXT_P(2));
610 if (!strlen(outdbfile)) {
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 }
620
621 /* outdb band index (3) */
622 if (!PG_ARGISNULL(3)) {
623 ArrayType *array;
624 Oid etype;
625 Datum *e;
626 bool *nulls;
627
628 int16 typlen;
629 bool typbyval;
630 char typalign;
631
632 allbands = FALSE;
633
634 array = PG_GETARG_ARRAYTYPE_P(3);
635 etype = ARR_ELEMTYPE(array);
636 get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
637
638 switch (etype) {
639 case INT2OID:
640 case INT4OID:
641 break;
642 default:
643 if (pgraster != NULL) {
644 rt_raster_destroy(raster);
645 PG_FREE_IF_COPY(pgraster, 0);
646 }
647 elog(ERROR, "RASTER_addBandOutDB: Invalid data type for band indexes");
648 PG_RETURN_NULL();
649 break;
650 }
651
652 deconstruct_array(array, etype, typlen, typbyval, typalign, &e, &nulls, &numsrcnband);
653
654 srcnband = palloc(sizeof(int) * numsrcnband);
655 if (srcnband == NULL) {
656 if (pgraster != NULL) {
657 rt_raster_destroy(raster);
658 PG_FREE_IF_COPY(pgraster, 0);
659 }
660 elog(ERROR, "RASTER_addBandOutDB: Cannot allocate memory for band indexes");
661 PG_RETURN_NULL();
662 }
663
664 for (i = 0, j = 0; i < numsrcnband; i++) {
665 if (nulls[i]) continue;
666
667 switch (etype) {
668 case INT2OID:
669 srcnband[j] = DatumGetInt16(e[i]);
670 break;
671 case INT4OID:
672 srcnband[j] = DatumGetInt32(e[i]);
673 break;
674 }
675 j++;
676 }
677
678 if (j < numsrcnband) {
679 srcnband = repalloc(srcnband, sizeof(int) * j);
680 if (srcnband == NULL) {
681 if (pgraster != NULL) {
682 rt_raster_destroy(raster);
683 PG_FREE_IF_COPY(pgraster, 0);
684 }
685 elog(ERROR, "RASTER_addBandOutDB: Cannot reallocate memory for band indexes");
686 PG_RETURN_NULL();
687 }
688
689 numsrcnband = j;
690 }
691 }
692 else
693 allbands = TRUE;
694
695 /* nodataval (4) */
696 if (!PG_ARGISNULL(4)) {
697 hasnodata = TRUE;
698 nodataval = PG_GETARG_FLOAT8(4);
699 }
700 else
701 hasnodata = FALSE;
702
703 /* validate input */
704
705 /* make sure dstnband is valid */
706 if (raster != NULL) {
707 numbands = rt_raster_get_num_bands(raster);
708 if (!appendband) {
709 if (dstnband < 1) {
710 elog(NOTICE, "Invalid band index %d for adding bands. Using band index 1", dstnband);
711 dstnband = 1;
712 }
713 else if (numbands > 0 && dstnband > numbands) {
714 elog(NOTICE, "Invalid band index %d for adding bands. Using band index %d", dstnband, numbands);
715 dstnband = numbands + 1;
716 }
717 }
718 else
719 dstnband = numbands + 1;
720 }
721
722 /* open outdb raster file */
723 rt_util_gdal_register_all(0);
724 hdsOut = rt_util_gdal_open(outdbfile, GA_ReadOnly, 1);
725 if (hdsOut == NULL) {
726 if (pgraster != NULL) {
727 rt_raster_destroy(raster);
728 PG_FREE_IF_COPY(pgraster, 0);
729 }
730 elog(ERROR, "RASTER_addBandOutDB: Cannot open out-db file with GDAL");
731 PG_RETURN_NULL();
732 }
733
734 /* get offline raster's geotransform */
735 if (GDALGetGeoTransform(hdsOut, ogt) != CE_None) {
736 ogt[0] = 0;
737 ogt[1] = 1;
738 ogt[2] = 0;
739 ogt[3] = 0;
740 ogt[4] = 0;
741 ogt[5] = -1;
742 }
743
744 /* raster doesn't exist, create it now */
745 if (raster == NULL) {
746 raster = rt_raster_new(GDALGetRasterXSize(hdsOut), GDALGetRasterYSize(hdsOut));
747 if (rt_raster_is_empty(raster)) {
748 elog(ERROR, "RASTER_addBandOutDB: Cannot create new raster");
749 PG_RETURN_NULL();
750 }
751 rt_raster_set_geotransform_matrix(raster, ogt);
752
753 if (rt_util_gdal_sr_auth_info(hdsOut, &authname, &authcode) == ES_NONE) {
754 if (
755 authname != NULL &&
756 strcmp(authname, "EPSG") == 0 &&
757 authcode != NULL
758 ) {
759 rt_raster_set_srid(raster, atoi(authcode));
760 }
761 else
762 elog(INFO, "Unknown SRS auth name and code from out-db file. Defaulting SRID of new raster to %d", SRID_UNKNOWN);
763 }
764 else
765 elog(INFO, "Cannot get SRS auth name and code from out-db file. Defaulting SRID of new raster to %d", SRID_UNKNOWN);
766 }
767
768 /* some raster info */
769 width = rt_raster_get_width(raster);
770 height = rt_raster_get_height(raster);
771
772 /* are rasters aligned? */
773 _rast = rt_raster_new(1, 1);
774 rt_raster_set_geotransform_matrix(_rast, ogt);
775 rt_raster_set_srid(_rast, rt_raster_get_srid(raster));
776 err = rt_raster_same_alignment(raster, _rast, &aligned, NULL);
777 rt_raster_destroy(_rast);
778
779 if (err != ES_NONE) {
780 GDALClose(hdsOut);
781 if (raster != NULL)
782 rt_raster_destroy(raster);
783 if (pgraster != NULL)
784 PG_FREE_IF_COPY(pgraster, 0);
785 elog(ERROR, "RASTER_addBandOutDB: Cannot test alignment of out-db file");
786 return ES_ERROR;
787 }
788 else if (!aligned)
789 elog(WARNING, "The in-db representation of the out-db raster is not aligned. Band data may be incorrect");
790
791 /* build up srcnband */
792 if (allbands) {
793 numsrcnband = GDALGetRasterCount(hdsOut);
794 GDALClose(hdsOut);
795
796 srcnband = palloc(sizeof(int) * numsrcnband);
797 if (srcnband == NULL) {
798 if (raster != NULL)
799 rt_raster_destroy(raster);
800 if (pgraster != NULL)
801 PG_FREE_IF_COPY(pgraster, 0);
802 elog(ERROR, "RASTER_addBandOutDB: Cannot allocate memory for band indexes");
803 PG_RETURN_NULL();
804 }
805
806 for (i = 0, j = 1; i < numsrcnband; i++, j++)
807 srcnband[i] = j;
808 }
809 else
810 GDALClose(hdsOut);
811
812 /* add band */
813 for (i = 0, j = dstnband - 1; i < numsrcnband; i++, j++) {
814
815 /* create band with path */
816 band = rt_band_new_offline_from_path(
817 width, height,
818 hasnodata, nodataval,
819 srcnband[i], outdbfile,
820 FALSE
821 );
822 if (band == NULL) {
823 if (raster != NULL)
824 rt_raster_destroy(raster);
825 if (pgraster != NULL)
826 PG_FREE_IF_COPY(pgraster, 0);
827 elog(ERROR, "RASTER_addBandOutDB: Cannot create new out-db band");
828 PG_RETURN_NULL();
829 }
830
831 /* add band */
832 if (rt_raster_add_band(raster, band, j) < 0) {
833 if (raster != NULL)
834 rt_raster_destroy(raster);
835 if (pgraster != NULL)
836 PG_FREE_IF_COPY(pgraster, 0);
837 elog(ERROR, "RASTER_addBandOutDB: Cannot add new out-db band to raster");
838 PG_RETURN_NULL();
839 }
840 }
841
842 pgrtn = rt_raster_serialize(raster);
843 rt_raster_destroy(raster);
844 if (pgraster != NULL)
845 PG_FREE_IF_COPY(pgraster, 0);
846 if (!pgrtn)
847 PG_RETURN_NULL();
848
849 SET_VARSIZE(pgrtn, pgrtn->size);
850 PG_RETURN_POINTER(pgrtn);
851 }
852
853 /**
854 * Copy a band from one raster to another one at the given position.
855 */
856 PG_FUNCTION_INFO_V1(RASTER_copyBand);
RASTER_copyBand(PG_FUNCTION_ARGS)857 Datum RASTER_copyBand(PG_FUNCTION_ARGS)
858 {
859 rt_pgraster *pgto = NULL;
860 rt_pgraster *pgfrom = NULL;
861 rt_pgraster *pgrtn = NULL;
862 rt_raster torast = NULL;
863 rt_raster fromrast = NULL;
864 int toindex = 0;
865 int fromband = 0;
866 int oldtorastnumbands = 0;
867 int newtorastnumbands = 0;
868 int newbandindex = 0;
869
870 /* Deserialize torast */
871 if (PG_ARGISNULL(0)) PG_RETURN_NULL();
872 pgto = (rt_pgraster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
873
874 torast = rt_raster_deserialize(pgto, FALSE);
875 if (!torast) {
876 PG_FREE_IF_COPY(pgto, 0);
877 elog(ERROR, "RASTER_copyBand: Could not deserialize first raster");
878 PG_RETURN_NULL();
879 }
880
881 /* Deserialize fromrast */
882 if (!PG_ARGISNULL(1)) {
883 pgfrom = (rt_pgraster *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
884
885 fromrast = rt_raster_deserialize(pgfrom, FALSE);
886 if (!fromrast) {
887 rt_raster_destroy(torast);
888 PG_FREE_IF_COPY(pgfrom, 1);
889 PG_FREE_IF_COPY(pgto, 0);
890 elog(ERROR, "RASTER_copyBand: Could not deserialize second raster");
891 PG_RETURN_NULL();
892 }
893
894 oldtorastnumbands = rt_raster_get_num_bands(torast);
895
896 if (PG_ARGISNULL(2))
897 fromband = 1;
898 else
899 fromband = PG_GETARG_INT32(2);
900
901 if (PG_ARGISNULL(3))
902 toindex = oldtorastnumbands + 1;
903 else
904 toindex = PG_GETARG_INT32(3);
905
906 /* Copy band fromrast torast */
907 newbandindex = rt_raster_copy_band(
908 torast, fromrast,
909 fromband - 1, toindex - 1
910 );
911
912 newtorastnumbands = rt_raster_get_num_bands(torast);
913 if (newtorastnumbands == oldtorastnumbands || newbandindex == -1) {
914 elog(NOTICE, "RASTER_copyBand: Could not add band to raster. "
915 "Returning original raster."
916 );
917 }
918
919 rt_raster_destroy(fromrast);
920 PG_FREE_IF_COPY(pgfrom, 1);
921 }
922
923 /* Serialize and return torast */
924 pgrtn = rt_raster_serialize(torast);
925 rt_raster_destroy(torast);
926 PG_FREE_IF_COPY(pgto, 0);
927 if (!pgrtn) PG_RETURN_NULL();
928
929 SET_VARSIZE(pgrtn, pgrtn->size);
930 PG_RETURN_POINTER(pgrtn);
931 }
932
933 /**
934 * Break up a raster into smaller tiles. SRF function
935 */
936 PG_FUNCTION_INFO_V1(RASTER_tile);
RASTER_tile(PG_FUNCTION_ARGS)937 Datum RASTER_tile(PG_FUNCTION_ARGS)
938 {
939 FuncCallContext *funcctx;
940 int call_cntr;
941 int max_calls;
942 int i = 0;
943 int j = 0;
944
945 struct tile_arg_t {
946
947 struct {
948 rt_raster raster;
949 double gt[6];
950 int srid;
951 int width;
952 int height;
953 } raster;
954
955 struct {
956 int width;
957 int height;
958
959 int nx;
960 int ny;
961 } tile;
962
963 int numbands;
964 int *nbands;
965
966 struct {
967 int pad;
968 double hasnodata;
969 double nodataval;
970 } pad;
971 };
972 struct tile_arg_t *arg1 = NULL;
973 struct tile_arg_t *arg2 = NULL;
974
975 if (SRF_IS_FIRSTCALL()) {
976 MemoryContext oldcontext;
977 rt_pgraster *pgraster = NULL;
978 int numbands;
979
980 ArrayType *array;
981 Oid etype;
982 Datum *e;
983 bool *nulls;
984
985 int16 typlen;
986 bool typbyval;
987 char typalign;
988
989 POSTGIS_RT_DEBUG(2, "RASTER_tile: first call");
990
991 /* create a function context for cross-call persistence */
992 funcctx = SRF_FIRSTCALL_INIT();
993
994 /* switch to memory context appropriate for multiple function calls */
995 oldcontext = MemoryContextSwitchTo(funcctx->multi_call_memory_ctx);
996
997 /* Get input arguments */
998 if (PG_ARGISNULL(0)) {
999 MemoryContextSwitchTo(oldcontext);
1000 SRF_RETURN_DONE(funcctx);
1001 }
1002
1003 /* allocate arg1 */
1004 arg1 = palloc(sizeof(struct tile_arg_t));
1005 if (arg1 == NULL) {
1006 MemoryContextSwitchTo(oldcontext);
1007 elog(ERROR, "RASTER_tile: Could not allocate memory for arguments");
1008 SRF_RETURN_DONE(funcctx);
1009 }
1010
1011 pgraster = (rt_pgraster *) PG_DETOAST_DATUM_COPY(PG_GETARG_DATUM(0));
1012 arg1->raster.raster = rt_raster_deserialize(pgraster, FALSE);
1013 if (!arg1->raster.raster) {
1014 ereport(ERROR, (
1015 errcode(ERRCODE_OUT_OF_MEMORY),
1016 errmsg("Could not deserialize raster")
1017 ));
1018 pfree(arg1);
1019 PG_FREE_IF_COPY(pgraster, 0);
1020 MemoryContextSwitchTo(oldcontext);
1021 SRF_RETURN_DONE(funcctx);
1022 }
1023
1024 /* raster has bands */
1025 numbands = rt_raster_get_num_bands(arg1->raster.raster);
1026 /*
1027 if (!numbands) {
1028 elog(NOTICE, "Raster provided has no bands");
1029 rt_raster_destroy(arg1->raster.raster);
1030 pfree(arg1);
1031 PG_FREE_IF_COPY(pgraster, 0);
1032 MemoryContextSwitchTo(oldcontext);
1033 SRF_RETURN_DONE(funcctx);
1034 }
1035 */
1036
1037 /* width (1) */
1038 if (PG_ARGISNULL(1)) {
1039 elog(NOTICE, "Width cannot be NULL. Returning NULL");
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 arg1->tile.width = PG_GETARG_INT32(1);
1047 if (arg1->tile.width < 1) {
1048 elog(NOTICE, "Width must be greater than zero. Returning NULL");
1049 rt_raster_destroy(arg1->raster.raster);
1050 pfree(arg1);
1051 PG_FREE_IF_COPY(pgraster, 0);
1052 MemoryContextSwitchTo(oldcontext);
1053 SRF_RETURN_DONE(funcctx);
1054 }
1055
1056 /* height (2) */
1057 if (PG_ARGISNULL(2)) {
1058 elog(NOTICE, "Height cannot be NULL. Returning NULL");
1059 rt_raster_destroy(arg1->raster.raster);
1060 pfree(arg1);
1061 PG_FREE_IF_COPY(pgraster, 0);
1062 MemoryContextSwitchTo(oldcontext);
1063 SRF_RETURN_DONE(funcctx);
1064 }
1065 arg1->tile.height = PG_GETARG_INT32(2);
1066 if (arg1->tile.height < 1) {
1067 elog(NOTICE, "Height must be greater than zero. Returning NULL");
1068 rt_raster_destroy(arg1->raster.raster);
1069 pfree(arg1);
1070 PG_FREE_IF_COPY(pgraster, 0);
1071 MemoryContextSwitchTo(oldcontext);
1072 SRF_RETURN_DONE(funcctx);
1073 }
1074
1075 /* nband, array (3) */
1076 if (numbands && !PG_ARGISNULL(3)) {
1077 array = PG_GETARG_ARRAYTYPE_P(3);
1078 etype = ARR_ELEMTYPE(array);
1079 get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
1080
1081 switch (etype) {
1082 case INT2OID:
1083 case INT4OID:
1084 break;
1085 default:
1086 rt_raster_destroy(arg1->raster.raster);
1087 pfree(arg1);
1088 PG_FREE_IF_COPY(pgraster, 0);
1089 MemoryContextSwitchTo(oldcontext);
1090 elog(ERROR, "RASTER_tile: Invalid data type for band indexes");
1091 SRF_RETURN_DONE(funcctx);
1092 break;
1093 }
1094
1095 deconstruct_array(array, etype, typlen, typbyval, typalign, &e, &nulls, &(arg1->numbands));
1096
1097 arg1->nbands = palloc(sizeof(int) * arg1->numbands);
1098 if (arg1->nbands == NULL) {
1099 rt_raster_destroy(arg1->raster.raster);
1100 pfree(arg1);
1101 PG_FREE_IF_COPY(pgraster, 0);
1102 MemoryContextSwitchTo(oldcontext);
1103 elog(ERROR, "RASTER_tile: Could not allocate memory for band indexes");
1104 SRF_RETURN_DONE(funcctx);
1105 }
1106
1107 for (i = 0, j = 0; i < arg1->numbands; i++) {
1108 if (nulls[i]) continue;
1109
1110 switch (etype) {
1111 case INT2OID:
1112 arg1->nbands[j] = DatumGetInt16(e[i]) - 1;
1113 break;
1114 case INT4OID:
1115 arg1->nbands[j] = DatumGetInt32(e[i]) - 1;
1116 break;
1117 }
1118
1119 j++;
1120 }
1121
1122 if (j < arg1->numbands) {
1123 arg1->nbands = repalloc(arg1->nbands, sizeof(int) * j);
1124 if (arg1->nbands == NULL) {
1125 rt_raster_destroy(arg1->raster.raster);
1126 pfree(arg1);
1127 PG_FREE_IF_COPY(pgraster, 0);
1128 MemoryContextSwitchTo(oldcontext);
1129 elog(ERROR, "RASTER_tile: Could not reallocate memory for band indexes");
1130 SRF_RETURN_DONE(funcctx);
1131 }
1132
1133 arg1->numbands = j;
1134 }
1135
1136 /* validate nbands */
1137 for (i = 0; i < arg1->numbands; i++) {
1138 if (!rt_raster_has_band(arg1->raster.raster, arg1->nbands[i])) {
1139 elog(NOTICE, "Band at index %d not found in raster", arg1->nbands[i] + 1);
1140 rt_raster_destroy(arg1->raster.raster);
1141 pfree(arg1->nbands);
1142 pfree(arg1);
1143 PG_FREE_IF_COPY(pgraster, 0);
1144 MemoryContextSwitchTo(oldcontext);
1145 SRF_RETURN_DONE(funcctx);
1146 }
1147 }
1148 }
1149 else {
1150 arg1->numbands = numbands;
1151
1152 if (numbands) {
1153 arg1->nbands = palloc(sizeof(int) * arg1->numbands);
1154
1155 if (arg1->nbands == NULL) {
1156 rt_raster_destroy(arg1->raster.raster);
1157 pfree(arg1);
1158 PG_FREE_IF_COPY(pgraster, 0);
1159 MemoryContextSwitchTo(oldcontext);
1160 elog(ERROR, "RASTER_dumpValues: Could not allocate memory for pixel values");
1161 SRF_RETURN_DONE(funcctx);
1162 }
1163
1164 for (i = 0; i < arg1->numbands; i++) {
1165 arg1->nbands[i] = i;
1166 POSTGIS_RT_DEBUGF(4, "arg1->nbands[%d] = %d", arg1->nbands[i], i);
1167 }
1168 }
1169 }
1170
1171 /* pad (4) and padnodata (5) */
1172 if (!PG_ARGISNULL(4)) {
1173 arg1->pad.pad = PG_GETARG_BOOL(4) ? 1 : 0;
1174
1175 if (arg1->pad.pad && !PG_ARGISNULL(5)) {
1176 arg1->pad.hasnodata = 1;
1177 arg1->pad.nodataval = PG_GETARG_FLOAT8(5);
1178 }
1179 else {
1180 arg1->pad.hasnodata = 0;
1181 arg1->pad.nodataval = 0;
1182 }
1183 }
1184 else {
1185 arg1->pad.pad = 0;
1186 arg1->pad.hasnodata = 0;
1187 arg1->pad.nodataval = 0;
1188 }
1189
1190 /* store some additional metadata */
1191 arg1->raster.srid = rt_raster_get_srid(arg1->raster.raster);
1192 arg1->raster.width = rt_raster_get_width(arg1->raster.raster);
1193 arg1->raster.height = rt_raster_get_height(arg1->raster.raster);
1194 rt_raster_get_geotransform_matrix(arg1->raster.raster, arg1->raster.gt);
1195
1196 /* determine maximum number of tiles from raster */
1197 arg1->tile.nx = ceil(arg1->raster.width / (double) arg1->tile.width);
1198 arg1->tile.ny = ceil(arg1->raster.height / (double) arg1->tile.height);
1199 POSTGIS_RT_DEBUGF(4, "# of tiles (x, y) = (%d, %d)", arg1->tile.nx, arg1->tile.ny);
1200
1201 /* Store needed information */
1202 funcctx->user_fctx = arg1;
1203
1204 /* total number of tuples to be returned */
1205 funcctx->max_calls = (arg1->tile.nx * arg1->tile.ny);
1206
1207 MemoryContextSwitchTo(oldcontext);
1208 }
1209
1210 /* stuff done on every call of the function */
1211 funcctx = SRF_PERCALL_SETUP();
1212
1213 call_cntr = funcctx->call_cntr;
1214 max_calls = funcctx->max_calls;
1215 arg2 = funcctx->user_fctx;
1216
1217 /* do when there is more left to send */
1218 if (call_cntr < max_calls) {
1219 rt_pgraster *pgtile = NULL;
1220 rt_raster tile = NULL;
1221 rt_band _band = NULL;
1222 rt_band band = NULL;
1223 rt_pixtype pixtype = PT_END;
1224 int hasnodata = 0;
1225 double nodataval = 0;
1226 int width = 0;
1227 int height = 0;
1228
1229 int k = 0;
1230 int tx = 0;
1231 int ty = 0;
1232 int rx = 0;
1233 int ry = 0;
1234 int ex = 0; /* edge tile on right */
1235 int ey = 0; /* edge tile on bottom */
1236 double ulx = 0;
1237 double uly = 0;
1238 uint16_t len = 0;
1239 void *vals = NULL;
1240 uint16_t nvals;
1241
1242 POSTGIS_RT_DEBUGF(3, "call number %d", call_cntr);
1243
1244 /*
1245 find offset based upon tile #
1246
1247 0 1 2
1248 3 4 5
1249 6 7 8
1250 */
1251 ty = call_cntr / arg2->tile.nx;
1252 tx = call_cntr % arg2->tile.nx;
1253 POSTGIS_RT_DEBUGF(4, "tile (x, y) = (%d, %d)", tx, ty);
1254
1255 /* edge tile? only important if padding is false */
1256 if (!arg2->pad.pad) {
1257 if (ty + 1 == arg2->tile.ny)
1258 ey = 1;
1259 if (tx + 1 == arg2->tile.nx)
1260 ex = 1;
1261 }
1262
1263 /* upper-left of tile in raster coordinates */
1264 rx = tx * arg2->tile.width;
1265 ry = ty * arg2->tile.height;
1266 POSTGIS_RT_DEBUGF(4, "raster coordinates = %d, %d", rx, ry);
1267
1268 /* determine tile width and height */
1269 /* default to user-defined */
1270 width = arg2->tile.width;
1271 height = arg2->tile.height;
1272
1273 /* override user-defined if edge tile (only possible if padding is false */
1274 if (ex || ey) {
1275 /* right edge */
1276 if (ex)
1277 width = arg2->raster.width - rx;
1278 /* bottom edge */
1279 if (ey)
1280 height = arg2->raster.height - ry;
1281 }
1282
1283 /* create empty raster */
1284 tile = rt_raster_new(width, height);
1285 rt_raster_set_geotransform_matrix(tile, arg2->raster.gt);
1286 rt_raster_set_srid(tile, arg2->raster.srid);
1287
1288 /* upper-left of tile in spatial coordinates */
1289 if (rt_raster_cell_to_geopoint(arg2->raster.raster, rx, ry, &ulx, &uly, arg2->raster.gt) != ES_NONE) {
1290 rt_raster_destroy(tile);
1291 rt_raster_destroy(arg2->raster.raster);
1292 if (arg2->numbands) pfree(arg2->nbands);
1293 pfree(arg2);
1294 elog(ERROR, "RASTER_tile: Could not compute the coordinates of the upper-left corner of the output tile");
1295 SRF_RETURN_DONE(funcctx);
1296 }
1297 rt_raster_set_offsets(tile, ulx, uly);
1298 POSTGIS_RT_DEBUGF(4, "spatial coordinates = %f, %f", ulx, uly);
1299
1300 /* compute length of pixel line to read */
1301 len = arg2->tile.width;
1302 if (rx + arg2->tile.width >= arg2->raster.width)
1303 len = arg2->raster.width - rx;
1304 POSTGIS_RT_DEBUGF(3, "read line len = %d", len);
1305
1306 /* copy bands to tile */
1307 for (i = 0; i < arg2->numbands; i++) {
1308 POSTGIS_RT_DEBUGF(4, "copying band %d to tile %d", arg2->nbands[i], call_cntr);
1309
1310 _band = rt_raster_get_band(arg2->raster.raster, arg2->nbands[i]);
1311 if (_band == NULL) {
1312 int nband = arg2->nbands[i] + 1;
1313 rt_raster_destroy(tile);
1314 rt_raster_destroy(arg2->raster.raster);
1315 pfree(arg2->nbands);
1316 pfree(arg2);
1317 elog(ERROR, "RASTER_tile: Could not get band %d from source raster", nband);
1318 SRF_RETURN_DONE(funcctx);
1319 }
1320
1321 pixtype = rt_band_get_pixtype(_band);
1322 hasnodata = rt_band_get_hasnodata_flag(_band);
1323 if (hasnodata)
1324 rt_band_get_nodata(_band, &nodataval);
1325 else if (arg2->pad.pad && arg2->pad.hasnodata) {
1326 hasnodata = 1;
1327 nodataval = arg2->pad.nodataval;
1328 }
1329 else
1330 nodataval = rt_band_get_min_value(_band);
1331
1332 /* inline band */
1333 if (!rt_band_is_offline(_band)) {
1334 if (rt_raster_generate_new_band(tile, pixtype, nodataval, hasnodata, nodataval, i) < 0) {
1335 rt_raster_destroy(tile);
1336 rt_raster_destroy(arg2->raster.raster);
1337 pfree(arg2->nbands);
1338 pfree(arg2);
1339 elog(ERROR, "RASTER_tile: Could not add new band to output tile");
1340 SRF_RETURN_DONE(funcctx);
1341 }
1342 band = rt_raster_get_band(tile, i);
1343 if (band == NULL) {
1344 rt_raster_destroy(tile);
1345 rt_raster_destroy(arg2->raster.raster);
1346 pfree(arg2->nbands);
1347 pfree(arg2);
1348 elog(ERROR, "RASTER_tile: Could not get newly added band from output tile");
1349 SRF_RETURN_DONE(funcctx);
1350 }
1351
1352 /* if isnodata, set flag and continue */
1353 if (rt_band_get_isnodata_flag(_band)) {
1354 rt_band_set_isnodata_flag(band, 1);
1355 continue;
1356 }
1357
1358 /* copy data */
1359 for (j = 0; j < arg2->tile.height; j++) {
1360 k = ry + j;
1361
1362 if (k >= arg2->raster.height) {
1363 POSTGIS_RT_DEBUGF(4, "row %d is beyond extent of source raster. skipping", k);
1364 continue;
1365 }
1366
1367 POSTGIS_RT_DEBUGF(4, "getting pixel line %d, %d for %d pixels", rx, k, len);
1368 if (rt_band_get_pixel_line(_band, rx, k, len, &vals, &nvals) != ES_NONE) {
1369 rt_raster_destroy(tile);
1370 rt_raster_destroy(arg2->raster.raster);
1371 pfree(arg2->nbands);
1372 pfree(arg2);
1373 elog(ERROR, "RASTER_tile: Could not get pixel line from source raster");
1374 SRF_RETURN_DONE(funcctx);
1375 }
1376
1377 if (nvals && rt_band_set_pixel_line(band, 0, j, vals, nvals) != ES_NONE) {
1378 rt_raster_destroy(tile);
1379 rt_raster_destroy(arg2->raster.raster);
1380 pfree(arg2->nbands);
1381 pfree(arg2);
1382 elog(ERROR, "RASTER_tile: Could not set pixel line of output tile");
1383 SRF_RETURN_DONE(funcctx);
1384 }
1385 }
1386 }
1387 /* offline */
1388 else {
1389 uint8_t bandnum = 0;
1390 rt_band_get_ext_band_num(_band, &bandnum);
1391
1392 band = rt_band_new_offline(
1393 width, height,
1394 pixtype,
1395 hasnodata, nodataval,
1396 bandnum, rt_band_get_ext_path(_band)
1397 );
1398
1399 if (band == NULL) {
1400 rt_raster_destroy(tile);
1401 rt_raster_destroy(arg2->raster.raster);
1402 pfree(arg2->nbands);
1403 pfree(arg2);
1404 elog(ERROR, "RASTER_tile: Could not create new offline band for output tile");
1405 SRF_RETURN_DONE(funcctx);
1406 }
1407
1408 if (rt_raster_add_band(tile, band, i) < 0) {
1409 rt_band_destroy(band);
1410 rt_raster_destroy(tile);
1411 rt_raster_destroy(arg2->raster.raster);
1412 pfree(arg2->nbands);
1413 pfree(arg2);
1414 elog(ERROR, "RASTER_tile: Could not add new offline band to output tile");
1415 SRF_RETURN_DONE(funcctx);
1416 }
1417 }
1418 }
1419
1420 pgtile = rt_raster_serialize(tile);
1421 rt_raster_destroy(tile);
1422 if (!pgtile) {
1423 rt_raster_destroy(arg2->raster.raster);
1424 if (arg2->numbands) pfree(arg2->nbands);
1425 pfree(arg2);
1426 SRF_RETURN_DONE(funcctx);
1427 }
1428
1429 SET_VARSIZE(pgtile, pgtile->size);
1430 SRF_RETURN_NEXT(funcctx, PointerGetDatum(pgtile));
1431 }
1432 /* do when there is no more left */
1433 else {
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
1441 /**
1442 * Return new raster from selected bands of existing raster through ST_Band.
1443 * second argument is an array of band numbers (1 based)
1444 */
1445 PG_FUNCTION_INFO_V1(RASTER_band);
RASTER_band(PG_FUNCTION_ARGS)1446 Datum RASTER_band(PG_FUNCTION_ARGS)
1447 {
1448 rt_pgraster *pgraster;
1449 rt_pgraster *pgrast;
1450 rt_raster raster;
1451 rt_raster rast;
1452
1453 bool skip = FALSE;
1454 ArrayType *array;
1455 Oid etype;
1456 Datum *e;
1457 bool *nulls;
1458 int16 typlen;
1459 bool typbyval;
1460 char typalign;
1461
1462 uint32_t numBands;
1463 uint32_t *bandNums;
1464 uint32 idx = 0;
1465 int n;
1466 int i = 0;
1467 int j = 0;
1468
1469 if (PG_ARGISNULL(0))
1470 PG_RETURN_NULL();
1471 pgraster = (rt_pgraster *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
1472
1473 raster = rt_raster_deserialize(pgraster, FALSE);
1474 if (!raster) {
1475 PG_FREE_IF_COPY(pgraster, 0);
1476 elog(ERROR, "RASTER_band: Could not deserialize raster");
1477 PG_RETURN_NULL();
1478 }
1479
1480 /* process bandNums */
1481 if (PG_ARGISNULL(1)) {
1482 elog(NOTICE, "Band number(s) not provided. Returning original raster");
1483 skip = TRUE;
1484 }
1485 if (!skip) {
1486 numBands = rt_raster_get_num_bands(raster);
1487
1488 array = PG_GETARG_ARRAYTYPE_P(1);
1489 etype = ARR_ELEMTYPE(array);
1490 get_typlenbyvalalign(etype, &typlen, &typbyval, &typalign);
1491
1492 switch (etype) {
1493 case INT2OID:
1494 case INT4OID:
1495 break;
1496 default:
1497 rt_raster_destroy(raster);
1498 PG_FREE_IF_COPY(pgraster, 0);
1499 elog(ERROR, "RASTER_band: Invalid data type for band number(s)");
1500 PG_RETURN_NULL();
1501 break;
1502 }
1503
1504 deconstruct_array(array, etype, typlen, typbyval, typalign, &e,
1505 &nulls, &n);
1506
1507 bandNums = palloc(sizeof(uint32_t) * n);
1508 for (i = 0, j = 0; i < n; i++) {
1509 if (nulls[i]) continue;
1510
1511 switch (etype) {
1512 case INT2OID:
1513 idx = (uint32_t) DatumGetInt16(e[i]);
1514 break;
1515 case INT4OID:
1516 idx = (uint32_t) DatumGetInt32(e[i]);
1517 break;
1518 }
1519
1520 POSTGIS_RT_DEBUGF(3, "band idx (before): %d", idx);
1521 if (idx > numBands || idx < 1) {
1522 elog(NOTICE, "Invalid band index (must use 1-based). Returning original raster");
1523 skip = TRUE;
1524 break;
1525 }
1526
1527 bandNums[j] = idx - 1;
1528 POSTGIS_RT_DEBUGF(3, "bandNums[%d] = %d", j, bandNums[j]);
1529 j++;
1530 }
1531
1532 if (skip || j < 1) {
1533 pfree(bandNums);
1534 skip = TRUE;
1535 }
1536 }
1537
1538 if (!skip) {
1539 rast = rt_raster_from_band(raster, bandNums, j);
1540 pfree(bandNums);
1541 rt_raster_destroy(raster);
1542 PG_FREE_IF_COPY(pgraster, 0);
1543 if (!rast) {
1544 elog(ERROR, "RASTER_band: Could not create new raster");
1545 PG_RETURN_NULL();
1546 }
1547
1548 pgrast = rt_raster_serialize(rast);
1549 rt_raster_destroy(rast);
1550
1551 if (!pgrast)
1552 PG_RETURN_NULL();
1553
1554 SET_VARSIZE(pgrast, pgrast->size);
1555 PG_RETURN_POINTER(pgrast);
1556 }
1557
1558 PG_RETURN_POINTER(pgraster);
1559 }
1560