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 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
553 /**
554 * Add out-db band to the given raster at the given position
555 */
556 PG_FUNCTION_INFO_V1(RASTER_addBandOutDB);
RASTER_addBandOutDB(PG_FUNCTION_ARGS)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
864 /**
865 * Copy a band from one raster to another one at the given position.
866 */
867 PG_FUNCTION_INFO_V1(RASTER_copyBand);
RASTER_copyBand(PG_FUNCTION_ARGS)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
944 /**
945 * Break up a raster into smaller tiles. SRF function
946 */
947 PG_FUNCTION_INFO_V1(RASTER_tile);
RASTER_tile(PG_FUNCTION_ARGS)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 int32_t 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
1452 /**
1453 * Return new raster from selected bands of existing raster through ST_Band.
1454 * second argument is an array of band numbers (1 based)
1455 */
1456 PG_FUNCTION_INFO_V1(RASTER_band);
RASTER_band(PG_FUNCTION_ARGS)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 }
1571