PostGIS  2.4.9dev-r@@SVN_REVISION@@
rt_mapalgebra.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 "librtcore.h"
31 #include "librtcore_internal.h"
32 
33 /******************************************************************************
34 * rt_band_reclass()
35 ******************************************************************************/
36 
49 rt_band
50 rt_band_reclass(
51  rt_band srcband, rt_pixtype pixtype,
52  uint32_t hasnodata, double nodataval,
53  rt_reclassexpr *exprset, int exprcount
54 ) {
55  rt_band band = NULL;
56  uint32_t width = 0;
57  uint32_t height = 0;
58  int numval = 0;
59  int memsize = 0;
60  void *mem = NULL;
61  uint32_t src_hasnodata = 0;
62  double src_nodataval = 0.0;
63  int isnodata = 0;
64 
65  int rtn;
66  uint32_t x;
67  uint32_t y;
68  int i;
69  double or = 0;
70  double ov = 0;
71  double nr = 0;
72  double nv = 0;
73  int do_nv = 0;
74  rt_reclassexpr expr = NULL;
75 
76  assert(NULL != srcband);
77  assert(NULL != exprset && exprcount > 0);
78  RASTER_DEBUGF(4, "exprcount = %d", exprcount);
79  RASTER_DEBUGF(4, "exprset @ %p", exprset);
80 
81  /* source nodata */
82  src_hasnodata = rt_band_get_hasnodata_flag(srcband);
83  if (src_hasnodata)
84  rt_band_get_nodata(srcband, &src_nodataval);
85 
86  /* size of memory block to allocate */
87  width = rt_band_get_width(srcband);
88  height = rt_band_get_height(srcband);
89  numval = width * height;
90  memsize = rt_pixtype_size(pixtype) * numval;
91  mem = (int *) rtalloc(memsize);
92  if (!mem) {
93  rterror("rt_band_reclass: Could not allocate memory for band");
94  return 0;
95  }
96 
97  /* initialize to zero */
98  if (!hasnodata) {
99  memset(mem, 0, memsize);
100  }
101  /* initialize to nodataval */
102  else {
103  int32_t checkvalint = 0;
104  uint32_t checkvaluint = 0;
105  double checkvaldouble = 0;
106  float checkvalfloat = 0;
107 
108  switch (pixtype) {
109  case PT_1BB:
110  {
111  uint8_t *ptr = mem;
112  uint8_t clamped_initval = rt_util_clamp_to_1BB(nodataval);
113  for (i = 0; i < numval; i++)
114  ptr[i] = clamped_initval;
115  checkvalint = ptr[0];
116  break;
117  }
118  case PT_2BUI:
119  {
120  uint8_t *ptr = mem;
121  uint8_t clamped_initval = rt_util_clamp_to_2BUI(nodataval);
122  for (i = 0; i < numval; i++)
123  ptr[i] = clamped_initval;
124  checkvalint = ptr[0];
125  break;
126  }
127  case PT_4BUI:
128  {
129  uint8_t *ptr = mem;
130  uint8_t clamped_initval = rt_util_clamp_to_4BUI(nodataval);
131  for (i = 0; i < numval; i++)
132  ptr[i] = clamped_initval;
133  checkvalint = ptr[0];
134  break;
135  }
136  case PT_8BSI:
137  {
138  int8_t *ptr = mem;
139  int8_t clamped_initval = rt_util_clamp_to_8BSI(nodataval);
140  for (i = 0; i < numval; i++)
141  ptr[i] = clamped_initval;
142  checkvalint = ptr[0];
143  break;
144  }
145  case PT_8BUI:
146  {
147  uint8_t *ptr = mem;
148  uint8_t clamped_initval = rt_util_clamp_to_8BUI(nodataval);
149  for (i = 0; i < numval; i++)
150  ptr[i] = clamped_initval;
151  checkvalint = ptr[0];
152  break;
153  }
154  case PT_16BSI:
155  {
156  int16_t *ptr = mem;
157  int16_t clamped_initval = rt_util_clamp_to_16BSI(nodataval);
158  for (i = 0; i < numval; i++)
159  ptr[i] = clamped_initval;
160  checkvalint = ptr[0];
161  break;
162  }
163  case PT_16BUI:
164  {
165  uint16_t *ptr = mem;
166  uint16_t clamped_initval = rt_util_clamp_to_16BUI(nodataval);
167  for (i = 0; i < numval; i++)
168  ptr[i] = clamped_initval;
169  checkvalint = ptr[0];
170  break;
171  }
172  case PT_32BSI:
173  {
174  int32_t *ptr = mem;
175  int32_t clamped_initval = rt_util_clamp_to_32BSI(nodataval);
176  for (i = 0; i < numval; i++)
177  ptr[i] = clamped_initval;
178  checkvalint = ptr[0];
179  break;
180  }
181  case PT_32BUI:
182  {
183  uint32_t *ptr = mem;
184  uint32_t clamped_initval = rt_util_clamp_to_32BUI(nodataval);
185  for (i = 0; i < numval; i++)
186  ptr[i] = clamped_initval;
187  checkvaluint = ptr[0];
188  break;
189  }
190  case PT_32BF:
191  {
192  float *ptr = mem;
193  float clamped_initval = rt_util_clamp_to_32F(nodataval);
194  for (i = 0; i < numval; i++)
195  ptr[i] = clamped_initval;
196  checkvalfloat = ptr[0];
197  break;
198  }
199  case PT_64BF:
200  {
201  double *ptr = mem;
202  for (i = 0; i < numval; i++)
203  ptr[i] = nodataval;
204  checkvaldouble = ptr[0];
205  break;
206  }
207  default:
208  {
209  rterror("rt_band_reclass: Unknown pixeltype %d", pixtype);
210  rtdealloc(mem);
211  return 0;
212  }
213  }
214 
215  /* Overflow checking */
216  rt_util_dbl_trunc_warning(
217  nodataval,
218  checkvalint, checkvaluint,
219  checkvalfloat, checkvaldouble,
220  pixtype
221  );
222  }
223  RASTER_DEBUGF(3, "rt_band_reclass: width = %d height = %d", width, height);
224 
225  band = rt_band_new_inline(width, height, pixtype, hasnodata, nodataval, mem);
226  if (!band) {
227  rterror("rt_band_reclass: Could not create new band");
228  rtdealloc(mem);
229  return 0;
230  }
231  rt_band_set_ownsdata_flag(band, 1); /* we DO own this data!!! */
232  RASTER_DEBUGF(3, "rt_band_reclass: new band @ %p", band);
233 
234  for (x = 0; x < width; x++) {
235  for (y = 0; y < height; y++) {
236  rtn = rt_band_get_pixel(srcband, x, y, &ov, &isnodata);
237 
238  /* error getting value, skip */
239  if (rtn != ES_NONE) {
240  RASTER_DEBUGF(3, "Cannot get value at %d, %d", x, y);
241  continue;
242  }
243  RASTER_DEBUGF(4, "(x, y, ov, isnodata) = (%d, %d, %f, %d)", x, y, ov, isnodata);
244 
245  do {
246  do_nv = 0;
247 
248  /* no data*/
249  if (hasnodata && isnodata) {
250  do_nv = 1;
251  break;
252  }
253 
254  for (i = 0; i < exprcount; i++) {
255  expr = exprset[i];
256 
257  /* ov matches min and max*/
258  if (
259  FLT_EQ(expr->src.min, ov) &&
260  FLT_EQ(expr->src.max, ov)
261  ) {
262  do_nv = 1;
263  break;
264  }
265 
266  /* process min */
267  if ((
268  expr->src.exc_min && (
269  expr->src.min > ov ||
270  FLT_EQ(expr->src.min, ov)
271  )) || (
272  expr->src.inc_min && (
273  expr->src.min < ov ||
274  FLT_EQ(expr->src.min, ov)
275  )) || (
276  expr->src.min < ov
277  )) {
278  /* process max */
279  if ((
280  expr->src.exc_max && (
281  ov > expr->src.max ||
282  FLT_EQ(expr->src.max, ov)
283  )) || (
284  expr->src.inc_max && (
285  ov < expr->src.max ||
286  FLT_EQ(expr->src.max, ov)
287  )) || (
288  ov < expr->src.max
289  )) {
290  do_nv = 1;
291  break;
292  }
293  }
294  }
295  }
296  while (0);
297 
298  /* no expression matched, do not continue */
299  if (!do_nv) continue;
300  RASTER_DEBUGF(3, "Using exprset[%d] unless NODATA", i);
301 
302  /* converting a value from one range to another range
303  OldRange = (OldMax - OldMin)
304  NewRange = (NewMax - NewMin)
305  NewValue = (((OldValue - OldMin) * NewRange) / OldRange) + NewMin
306  */
307 
308  /* NODATA */
309  if (hasnodata && isnodata) {
310  nv = nodataval;
311  }
312  /*
313  "src" min and max is the same, prevent division by zero
314  set nv to "dst" min, which should be the same as "dst" max
315  */
316  else if (FLT_EQ(expr->src.max, expr->src.min)) {
317  nv = expr->dst.min;
318  }
319  else {
320  or = expr->src.max - expr->src.min;
321  nr = expr->dst.max - expr->dst.min;
322  nv = (((ov - expr->src.min) * nr) / or) + expr->dst.min;
323 
324  /* if dst range is from high to low */
325  if (expr->dst.min > expr->dst.max) {
326  if (nv > expr->dst.min)
327  nv = expr->dst.min;
328  else if (nv < expr->dst.max)
329  nv = expr->dst.max;
330  }
331  /* if dst range is from low to high */
332  else {
333  if (nv < expr->dst.min)
334  nv = expr->dst.min;
335  else if (nv > expr->dst.max)
336  nv = expr->dst.max;
337  }
338  }
339 
340  /* round the value for integers */
341  switch (pixtype) {
342  case PT_1BB:
343  case PT_2BUI:
344  case PT_4BUI:
345  case PT_8BSI:
346  case PT_8BUI:
347  case PT_16BSI:
348  case PT_16BUI:
349  case PT_32BSI:
350  case PT_32BUI:
351  nv = round(nv);
352  break;
353  default:
354  break;
355  }
356 
357  RASTER_DEBUGF(4, "(%d, %d) ov: %f or: %f - %f nr: %f - %f nv: %f"
358  , x
359  , y
360  , ov
361  , (NULL != expr) ? expr->src.min : 0
362  , (NULL != expr) ? expr->src.max : 0
363  , (NULL != expr) ? expr->dst.min : 0
364  , (NULL != expr) ? expr->dst.max : 0
365  , nv
366  );
367  if (rt_band_set_pixel(band, x, y, nv, NULL) != ES_NONE) {
368  rterror("rt_band_reclass: Could not assign value to new band");
369  rt_band_destroy(band);
370  rtdealloc(mem);
371  return 0;
372  }
373 
374  expr = NULL;
375  }
376  }
377 
378  return band;
379 }
380 
381 /******************************************************************************
382 * rt_raster_iterator()
383 ******************************************************************************/
384 
385 typedef struct _rti_iterator_arg_t* _rti_iterator_arg;
386 struct _rti_iterator_arg_t {
387  int count;
388 
389  rt_raster *raster;
390  int *isempty;
391  double **offset;
392  int *width;
393  int *height;
394 
395  struct {
396  rt_band *rtband;
397  int *hasnodata;
398  int *isnodata;
399  double *nodataval;
400  double *minval;
401  } band;
402 
403  struct {
404  uint16_t x;
405  uint16_t y;
406  } distance;
407 
408  struct {
409  uint32_t rows;
410  uint32_t columns;
411  } dimension;
412 
413  struct {
414  double **values;
415  int **nodata;
416  } empty;
417 
418  rt_iterator_arg arg;
419 };
420 
421 static _rti_iterator_arg
422 _rti_iterator_arg_init() {
423  _rti_iterator_arg _param;
424 
425  _param = rtalloc(sizeof(struct _rti_iterator_arg_t));
426  if (_param == NULL) {
427  rterror("_rti_iterator_arg_init: Could not allocate memory for _rti_iterator_arg");
428  return NULL;
429  }
430 
431  _param->count = 0;
432 
433  _param->raster = NULL;
434  _param->isempty = NULL;
435  _param->offset = NULL;
436  _param->width = NULL;
437  _param->height = NULL;
438 
439  _param->band.rtband = NULL;
440  _param->band.hasnodata = NULL;
441  _param->band.isnodata = NULL;
442  _param->band.nodataval = NULL;
443  _param->band.minval = NULL;
444 
445  _param->distance.x = 0;
446  _param->distance.y = 0;
447 
448  _param->dimension.rows = 0;
449  _param->dimension.columns = 0;
450 
451  _param->empty.values = NULL;
452  _param->empty.nodata = NULL;
453 
454  _param->arg = NULL;
455 
456  return _param;
457 }
458 
459 static void
460 _rti_iterator_arg_destroy(_rti_iterator_arg _param) {
461  int i = 0;
462 
463  if (_param->raster != NULL)
464  rtdealloc(_param->raster);
465  if (_param->isempty != NULL)
466  rtdealloc(_param->isempty);
467  if (_param->width != NULL)
468  rtdealloc(_param->width);
469  if (_param->height != NULL)
470  rtdealloc(_param->height);
471 
472  if (_param->band.rtband != NULL)
473  rtdealloc(_param->band.rtband);
474  if (_param->band.hasnodata != NULL)
475  rtdealloc(_param->band.hasnodata);
476  if (_param->band.isnodata != NULL)
477  rtdealloc(_param->band.isnodata);
478  if (_param->band.nodataval != NULL)
479  rtdealloc(_param->band.nodataval);
480  if (_param->band.minval != NULL)
481  rtdealloc(_param->band.minval);
482 
483  if (_param->offset != NULL) {
484  for (i = 0; i < _param->count; i++) {
485  if (_param->offset[i] == NULL)
486  continue;
487  rtdealloc(_param->offset[i]);
488  }
489  rtdealloc(_param->offset);
490  }
491 
492  if (_param->empty.values != NULL) {
493  for (i = 0; i < _param->dimension.rows; i++) {
494  if (_param->empty.values[i] == NULL)
495  continue;
496  rtdealloc(_param->empty.values[i]);
497  }
498  rtdealloc(_param->empty.values);
499  }
500  if (_param->empty.nodata != NULL) {
501  for (i = 0; i < _param->dimension.rows; i++) {
502  if (_param->empty.nodata[i] == NULL)
503  continue;
504  rtdealloc(_param->empty.nodata[i]);
505  }
506  rtdealloc(_param->empty.nodata);
507  }
508 
509  if (_param->arg != NULL) {
510  if (_param->arg->values != NULL)
511  rtdealloc(_param->arg->values);
512  if (_param->arg->nodata != NULL)
513  rtdealloc(_param->arg->nodata);
514  if (_param->arg->src_pixel != NULL) {
515  for (i = 0; i < _param->count; i++) {
516  if (_param->arg->src_pixel[i] == NULL)
517  continue;
518  rtdealloc(_param->arg->src_pixel[i]);
519  }
520 
521  rtdealloc(_param->arg->src_pixel);
522  }
523 
524  rtdealloc(_param->arg);
525  }
526 
527  rtdealloc(_param);
528 }
529 
530 static int
531 _rti_iterator_arg_populate(
532  _rti_iterator_arg _param,
533  rt_iterator itrset, uint16_t itrcount,
534  uint16_t distancex, uint16_t distancey,
535  int *allnull, int *allempty
536 ) {
537  int i = 0;
538  int hasband = 0;
539 
540  _param->count = itrcount;
541  _param->distance.x = distancex;
542  _param->distance.y = distancey;
543  _param->dimension.columns = distancex * 2 + 1;
544  _param->dimension.rows = distancey * 2 + 1;
545 
546  /* allocate memory for children */
547  _param->raster = rtalloc(sizeof(rt_raster) * itrcount);
548  _param->isempty = rtalloc(sizeof(int) * itrcount);
549  _param->width = rtalloc(sizeof(int) * itrcount);
550  _param->height = rtalloc(sizeof(int) * itrcount);
551 
552  _param->offset = rtalloc(sizeof(double *) * itrcount);
553 
554  _param->band.rtband = rtalloc(sizeof(rt_band) * itrcount);
555  _param->band.hasnodata = rtalloc(sizeof(int) * itrcount);
556  _param->band.isnodata = rtalloc(sizeof(int) * itrcount);
557  _param->band.nodataval = rtalloc(sizeof(double) * itrcount);
558  _param->band.minval = rtalloc(sizeof(double) * itrcount);
559 
560  if (
561  _param->raster == NULL ||
562  _param->isempty == NULL ||
563  _param->width == NULL ||
564  _param->height == NULL ||
565  _param->offset == NULL ||
566  _param->band.rtband == NULL ||
567  _param->band.hasnodata == NULL ||
568  _param->band.isnodata == NULL ||
569  _param->band.nodataval == NULL ||
570  _param->band.minval == NULL
571  ) {
572  rterror("_rti_iterator_arg_populate: Could not allocate memory for children of _rti_iterator_arg");
573  return 0;
574  }
575 
576  *allnull = 0;
577  *allempty = 0;
578 
579  /*
580  check input rasters
581  not empty, band # is valid
582  copy raster pointers and set flags
583  */
584  for (i = 0; i < itrcount; i++) {
585  /* initialize elements */
586  _param->raster[i] = NULL;
587  _param->isempty[i] = 0;
588  _param->width[i] = 0;
589  _param->height[i] = 0;
590 
591  _param->offset[i] = NULL;
592 
593  _param->band.rtband[i] = NULL;
594  _param->band.hasnodata[i] = 0;
595  _param->band.isnodata[i] = 0;
596  _param->band.nodataval[i] = 0;
597  _param->band.minval[i] = 0;
598 
599  /* set isempty */
600  if (itrset[i].raster == NULL) {
601  _param->isempty[i] = 1;
602 
603  (*allnull)++;
604  (*allempty)++;
605 
606  continue;
607  }
608  else if (rt_raster_is_empty(itrset[i].raster)) {
609  _param->isempty[i] = 1;
610 
611  (*allempty)++;
612 
613  continue;
614  }
615 
616  /* check band number */
617  hasband = rt_raster_has_band(itrset[i].raster, itrset[i].nband);
618  if (!hasband) {
619  if (!itrset[i].nbnodata) {
620  rterror("_rti_iterator_arg_populate: Band %d not found for raster %d", itrset[i].nband, i);
621  return 0;
622  }
623  else {
624  RASTER_DEBUGF(4, "Band %d not found for raster %d. Using NODATA", itrset[i].nband, i);
625  }
626  }
627 
628  _param->raster[i] = itrset[i].raster;
629  if (hasband) {
630  _param->band.rtband[i] = rt_raster_get_band(itrset[i].raster, itrset[i].nband);
631  if (_param->band.rtband[i] == NULL) {
632  rterror("_rti_iterator_arg_populate: Could not get band %d for raster %d", itrset[i].nband, i);
633  return 0;
634  }
635 
636  /* hasnodata */
637  _param->band.hasnodata[i] = rt_band_get_hasnodata_flag(_param->band.rtband[i]);
638 
639  /* hasnodata = TRUE */
640  if (_param->band.hasnodata[i]) {
641  /* nodataval */
642  rt_band_get_nodata(_param->band.rtband[i], &(_param->band.nodataval[i]));
643 
644  /* isnodata */
645  _param->band.isnodata[i] = rt_band_get_isnodata_flag(_param->band.rtband[i]);
646  }
647  /* hasnodata = FALSE */
648  else {
649  /* minval */
650  _param->band.minval[i] = rt_band_get_min_value(_param->band.rtband[i]);
651  }
652  }
653 
654  /* width, height */
655  _param->width[i] = rt_raster_get_width(_param->raster[i]);
656  _param->height[i] = rt_raster_get_height(_param->raster[i]);
657 
658  /* init offset */
659  _param->offset[i] = rtalloc(sizeof(double) * 2);
660  if (_param->offset[i] == NULL) {
661  rterror("_rti_iterator_arg_populate: Could not allocate memory for offsets");
662  return 0;
663  }
664  }
665 
666  return 1;
667 }
668 
669 static int
670 _rti_iterator_arg_empty_init(_rti_iterator_arg _param) {
671  int x = 0;
672  int y = 0;
673 
674  _param->empty.values = rtalloc(sizeof(double *) * _param->dimension.rows);
675  _param->empty.nodata = rtalloc(sizeof(int *) * _param->dimension.rows);
676  if (_param->empty.values == NULL || _param->empty.nodata == NULL) {
677  rterror("_rti_iterator_arg_empty_init: Could not allocate memory for empty values and NODATA");
678  return 0;
679  }
680 
681  for (y = 0; y < _param->dimension.rows; y++) {
682  _param->empty.values[y] = rtalloc(sizeof(double) * _param->dimension.columns);
683  _param->empty.nodata[y] = rtalloc(sizeof(int) * _param->dimension.columns);
684 
685  if (_param->empty.values[y] == NULL || _param->empty.nodata[y] == NULL) {
686  rterror("_rti_iterator_arg_empty_init: Could not allocate memory for elements of empty values and NODATA");
687  return 0;
688  }
689 
690  for (x = 0; x < _param->dimension.columns; x++) {
691  _param->empty.values[y][x] = 0;
692  _param->empty.nodata[y][x] = 1;
693  }
694  }
695 
696  return 1;
697 }
698 
699 static int
700 _rti_iterator_arg_callback_init(_rti_iterator_arg _param) {
701  int i = 0;
702 
703  _param->arg = rtalloc(sizeof(struct rt_iterator_arg_t));
704  if (_param->arg == NULL) {
705  rterror("_rti_iterator_arg_callback_init: Could not allocate memory for rt_iterator_arg");
706  return 0;
707  }
708 
709  _param->arg->values = NULL;
710  _param->arg->nodata = NULL;
711  _param->arg->src_pixel = NULL;
712 
713  /* initialize argument components */
714  _param->arg->values = rtalloc(sizeof(double **) * _param->count);
715  _param->arg->nodata = rtalloc(sizeof(int **) * _param->count);
716  _param->arg->src_pixel = rtalloc(sizeof(int *) * _param->count);
717  if (_param->arg->values == NULL || _param->arg->nodata == NULL || _param->arg->src_pixel == NULL) {
718  rterror("_rti_iterator_arg_callback_init: Could not allocate memory for element of rt_iterator_arg");
719  return 0;
720  }
721  memset(_param->arg->values, 0, sizeof(double **) * _param->count);
722  memset(_param->arg->nodata, 0, sizeof(int **) * _param->count);
723 
724  /* initialize pos */
725  for (i = 0; i < _param->count; i++) {
726 
727  _param->arg->src_pixel[i] = rtalloc(sizeof(int) * 2);
728  if (_param->arg->src_pixel[i] == NULL) {
729  rterror("_rti_iterator_arg_callback_init: Could not allocate memory for position elements of rt_iterator_arg");
730  return 0;
731  }
732  memset(_param->arg->src_pixel[i], 0, sizeof(int) * 2);
733  }
734 
735  _param->arg->rasters = _param->count;
736  _param->arg->rows = _param->dimension.rows;
737  _param->arg->columns = _param->dimension.columns;
738 
739  _param->arg->dst_pixel[0] = 0;
740  _param->arg->dst_pixel[1] = 0;
741 
742  return 1;
743 }
744 
745 static void
746 _rti_iterator_arg_callback_clean(_rti_iterator_arg _param) {
747  int i = 0;
748  int y = 0;
749 
750  for (i = 0; i < _param->count; i++) {
751  RASTER_DEBUGF(5, "empty at @ %p", _param->empty.values);
752  RASTER_DEBUGF(5, "values at @ %p", _param->arg->values[i]);
753 
754  if (_param->arg->values[i] == _param->empty.values) {
755  _param->arg->values[i] = NULL;
756  _param->arg->nodata[i] = NULL;
757 
758  continue;
759  }
760 
761  for (y = 0; y < _param->dimension.rows; y++) {
762  rtdealloc(_param->arg->values[i][y]);
763  rtdealloc(_param->arg->nodata[i][y]);
764  }
765 
766  rtdealloc(_param->arg->values[i]);
767  rtdealloc(_param->arg->nodata[i]);
768 
769  _param->arg->values[i] = NULL;
770  _param->arg->nodata[i] = NULL;
771  }
772 }
773 
811 rt_errorstate
812 rt_raster_iterator(
813  rt_iterator itrset, uint16_t itrcount,
814  rt_extenttype extenttype, rt_raster customextent,
815  rt_pixtype pixtype,
816  uint8_t hasnodata, double nodataval,
817  uint16_t distancex, uint16_t distancey,
818  rt_mask mask,
819  void *userarg,
820  int (*callback)(
821  rt_iterator_arg arg,
822  void *userarg,
823  double *value,
824  int *nodata
825  ),
826  rt_raster *rtnraster
827 ) {
828  /* output raster */
829  rt_raster rtnrast = NULL;
830  /* output raster's band */
831  rt_band rtnband = NULL;
832 
833  /* working raster */
834  rt_raster rast = NULL;
835 
836  _rti_iterator_arg _param = NULL;
837  int allnull = 0;
838  int allempty = 0;
839  int aligned = 0;
840  double offset[4] = {0.};
841  rt_pixel npixels;
842 
843  int i = 0;
844  int status = 0;
845  int inextent = 0;
846  int x = 0;
847  int y = 0;
848  int _x = 0;
849  int _y = 0;
850 
851  int _width = 0;
852  int _height = 0;
853 
854  double minval;
855  double value;
856  int isnodata;
857  int nodata;
858 
859  RASTER_DEBUG(3, "Starting...");
860 
861  assert(itrset != NULL && itrcount > 0);
862  assert(rtnraster != NULL);
863 
864  /* init rtnraster to NULL */
865  *rtnraster = NULL;
866 
867  /* check that callback function is not NULL */
868  if (callback == NULL) {
869  rterror("rt_raster_iterator: Callback function not provided");
870  return ES_ERROR;
871  }
872 
873  /* check that custom extent is provided if extenttype = ET_CUSTOM */
874  if (extenttype == ET_CUSTOM && rt_raster_is_empty(customextent)) {
875  rterror("rt_raster_iterator: Custom extent cannot be empty if extent type is ET_CUSTOM");
876  return ES_ERROR;
877  }
878 
879  /* check that pixtype != PT_END */
880  if (pixtype == PT_END) {
881  rterror("rt_raster_iterator: Pixel type cannot be PT_END");
882  return ES_ERROR;
883  }
884 
885  /* initialize _param */
886  if ((_param = _rti_iterator_arg_init()) == NULL) {
887  rterror("rt_raster_iterator: Could not initialize internal variables");
888  return ES_ERROR;
889  }
890 
891  /* fill _param */
892  if (!_rti_iterator_arg_populate(_param, itrset, itrcount, distancex, distancey, &allnull, &allempty)) {
893  rterror("rt_raster_iterator: Could not populate for internal variables");
894  _rti_iterator_arg_destroy(_param);
895  return ES_ERROR;
896  }
897 
898  /* shortcut if all null, return NULL */
899  if (allnull == itrcount) {
900  RASTER_DEBUG(3, "all rasters are NULL, returning NULL");
901 
902  _rti_iterator_arg_destroy(_param);
903 
904  return ES_NONE;
905  }
906  /* shortcut if all empty, return empty raster */
907  else if (allempty == itrcount) {
908  RASTER_DEBUG(3, "all rasters are empty, returning empty raster");
909 
910  _rti_iterator_arg_destroy(_param);
911 
912  rtnrast = rt_raster_new(0, 0);
913  if (rtnrast == NULL) {
914  rterror("rt_raster_iterator: Could not create empty raster");
915  return ES_ERROR;
916  }
917  rt_raster_set_scale(rtnrast, 0, 0);
918 
919  *rtnraster = rtnrast;
920  return ES_NONE;
921  }
922 
923  /* check that all rasters are aligned */
924  RASTER_DEBUG(3, "checking alignment of all rasters");
925  rast = NULL;
926 
927  /* find raster to use as reference */
928  /* use custom if provided */
929  if (extenttype == ET_CUSTOM) {
930  RASTER_DEBUG(4, "using custom extent as reference raster");
931  rast = customextent;
932  }
933  /* use first valid one in _param->raster */
934  else {
935  for (i = 0; i < itrcount; i++) {
936  if (!_param->isempty[i]) {
937  RASTER_DEBUGF(4, "using raster at index %d as reference raster", i);
938  rast = _param->raster[i];
939  break;
940  }
941  }
942  }
943 
944  /* no rasters found, SHOULD NEVER BE HERE! */
945  if (rast == NULL) {
946  rterror("rt_raster_iterator: Could not find reference raster to use for alignment tests");
947 
948  _rti_iterator_arg_destroy(_param);
949 
950  return ES_ERROR;
951  }
952 
953  do {
954  aligned = 1;
955 
956  /* check custom first if set. also skip if rasters are the same */
957  if (extenttype == ET_CUSTOM && rast != customextent) {
958  if (rt_raster_same_alignment(rast, customextent, &aligned, NULL) != ES_NONE) {
959  rterror("rt_raster_iterator: Could not test for alignment between reference raster and custom extent");
960 
961  _rti_iterator_arg_destroy(_param);
962 
963  return ES_ERROR;
964  }
965 
966  RASTER_DEBUGF(5, "custom extent alignment: %d", aligned);
967  if (!aligned)
968  break;
969  }
970 
971  for (i = 0; i < itrcount; i++) {
972  /* skip NULL rasters and if rasters are the same */
973  if (_param->isempty[i] || rast == _param->raster[i])
974  continue;
975 
976  if (rt_raster_same_alignment(rast, _param->raster[i], &aligned, NULL) != ES_NONE) {
977  rterror("rt_raster_iterator: Could not test for alignment between reference raster and raster %d", i);
978 
979  _rti_iterator_arg_destroy(_param);
980 
981  return ES_ERROR;
982  }
983  RASTER_DEBUGF(5, "raster at index %d alignment: %d", i, aligned);
984 
985  /* abort checking since a raster isn't aligned */
986  if (!aligned)
987  break;
988  }
989  }
990  while (0);
991 
992  /* not aligned, error */
993  if (!aligned) {
994  rterror("rt_raster_iterator: The set of rasters provided (custom extent included, if appropriate) do not have the same alignment");
995 
996  _rti_iterator_arg_destroy(_param);
997 
998  return ES_ERROR;
999  }
1000 
1001  /* use extenttype to build output raster (no bands though) */
1002  i = -1;
1003  switch (extenttype) {
1004  case ET_INTERSECTION:
1005  case ET_UNION:
1006  /* make copy of first "real" raster */
1007  rtnrast = rtalloc(sizeof(struct rt_raster_t));
1008  if (rtnrast == NULL) {
1009  rterror("rt_raster_iterator: Could not allocate memory for output raster");
1010 
1011  _rti_iterator_arg_destroy(_param);
1012 
1013  return ES_ERROR;
1014  }
1015 
1016  for (i = 0; i < itrcount; i++) {
1017  if (!_param->isempty[i]) {
1018  memcpy(rtnrast, _param->raster[i], sizeof(struct rt_raster_serialized_t));
1019  break;
1020  }
1021  }
1022  rtnrast->numBands = 0;
1023  rtnrast->bands = NULL;
1024 
1025  /* get extent of output raster */
1026  rast = NULL;
1027  for (i = i + 1; i < itrcount; i++) {
1028  if (_param->isempty[i])
1029  continue;
1030 
1031  status = rt_raster_from_two_rasters(rtnrast, _param->raster[i], extenttype, &rast, NULL);
1032  rtdealloc(rtnrast);
1033 
1034  if (rast == NULL || status != ES_NONE) {
1035  rterror("rt_raster_iterator: Could not compute %s extent of rasters",
1036  extenttype == ET_UNION ? "union" : "intersection"
1037  );
1038 
1039  _rti_iterator_arg_destroy(_param);
1040 
1041  return ES_ERROR;
1042  }
1043  else if (rt_raster_is_empty(rast)) {
1044  rtinfo("rt_raster_iterator: Computed raster for %s extent is empty",
1045  extenttype == ET_UNION ? "union" : "intersection"
1046  );
1047 
1048  _rti_iterator_arg_destroy(_param);
1049 
1050  *rtnraster = rast;
1051  return ES_NONE;
1052  }
1053 
1054  rtnrast = rast;
1055  rast = NULL;
1056  }
1057 
1058  break;
1059  /*
1060  first, second and last have similar checks
1061  and continue into custom
1062  */
1063  case ET_FIRST:
1064  i = 0;
1065  case ET_SECOND:
1066  if (i < 0) {
1067  if (itrcount < 2)
1068  i = 0;
1069  else
1070  i = 1;
1071  }
1072  case ET_LAST:
1073  if (i < 0) i = itrcount - 1;
1074 
1075  /* input raster is null, return NULL */
1076  if (_param->raster[i] == NULL) {
1077  RASTER_DEBUGF(3, "returning NULL as %s raster is NULL and extent type is ET_%s",
1078  (i == 0 ? "first" : (i == 1 ? "second" : "last")),
1079  (i == 0 ? "FIRST" : (i == 1 ? "SECOND" : "LAST"))
1080  );
1081 
1082  _rti_iterator_arg_destroy(_param);
1083 
1084  return ES_NONE;
1085  }
1086  /* input raster is empty, return empty raster */
1087  else if (_param->isempty[i]) {
1088  RASTER_DEBUGF(3, "returning empty raster as %s raster is empty and extent type is ET_%s",
1089  (i == 0 ? "first" : (i == 1 ? "second" : "last")),
1090  (i == 0 ? "FIRST" : (i == 1 ? "SECOND" : "LAST"))
1091  );
1092 
1093  _rti_iterator_arg_destroy(_param);
1094 
1095  rtnrast = rt_raster_new(0, 0);
1096  if (rtnrast == NULL) {
1097  rterror("rt_raster_iterator: Could not create empty raster");
1098  return ES_ERROR;
1099  }
1100  rt_raster_set_scale(rtnrast, 0, 0);
1101 
1102  *rtnraster = rtnrast;
1103  return ES_NONE;
1104  }
1105  /* copy the custom extent raster */
1106  case ET_CUSTOM:
1107  rtnrast = rtalloc(sizeof(struct rt_raster_t));
1108  if (rtnrast == NULL) {
1109  rterror("rt_raster_iterator: Could not allocate memory for output raster");
1110 
1111  _rti_iterator_arg_destroy(_param);
1112 
1113  return ES_ERROR;
1114  }
1115 
1116  switch (extenttype) {
1117  case ET_CUSTOM:
1118  memcpy(rtnrast, customextent, sizeof(struct rt_raster_serialized_t));
1119  break;
1120  /* first, second, last */
1121  default:
1122  memcpy(rtnrast, _param->raster[i], sizeof(struct rt_raster_serialized_t));
1123  break;
1124  }
1125  rtnrast->numBands = 0;
1126  rtnrast->bands = NULL;
1127  break;
1128  }
1129 
1130  _width = rt_raster_get_width(rtnrast);
1131  _height = rt_raster_get_height(rtnrast);
1132 
1133  RASTER_DEBUGF(4, "rtnrast (width, height, ulx, uly, scalex, scaley, skewx, skewy, srid) = (%d, %d, %f, %f, %f, %f, %f, %f, %d)",
1134  _width,
1135  _height,
1136  rt_raster_get_x_offset(rtnrast),
1137  rt_raster_get_y_offset(rtnrast),
1138  rt_raster_get_x_scale(rtnrast),
1139  rt_raster_get_y_scale(rtnrast),
1140  rt_raster_get_x_skew(rtnrast),
1141  rt_raster_get_y_skew(rtnrast),
1142  rt_raster_get_srid(rtnrast)
1143  );
1144 
1145  /* init values and NODATA for use with empty rasters */
1146  if (!_rti_iterator_arg_empty_init(_param)) {
1147  rterror("rt_raster_iterator: Could not initialize empty values and NODATA");
1148 
1149  _rti_iterator_arg_destroy(_param);
1150  rt_raster_destroy(rtnrast);
1151 
1152  return ES_ERROR;
1153  }
1154 
1155  /* create output band */
1156  if (rt_raster_generate_new_band(
1157  rtnrast,
1158  pixtype,
1159  nodataval,
1160  hasnodata, nodataval,
1161  0
1162  ) < 0) {
1163  rterror("rt_raster_iterator: Could not add new band to output raster");
1164 
1165  _rti_iterator_arg_destroy(_param);
1166  rt_raster_destroy(rtnrast);
1167 
1168  return ES_ERROR;
1169  }
1170 
1171  /* get output band */
1172  rtnband = rt_raster_get_band(rtnrast, 0);
1173  if (rtnband == NULL) {
1174  rterror("rt_raster_iterator: Could not get new band from output raster");
1175 
1176  _rti_iterator_arg_destroy(_param);
1177  rt_raster_destroy(rtnrast);
1178 
1179  return ES_ERROR;
1180  }
1181 
1182  /* output band's minimum value */
1183  minval = rt_band_get_min_value(rtnband);
1184 
1185  /* initialize argument for callback function */
1186  if (!_rti_iterator_arg_callback_init(_param)) {
1187  rterror("rt_raster_iterator: Could not initialize callback function argument");
1188 
1189  _rti_iterator_arg_destroy(_param);
1190  rt_band_destroy(rtnband);
1191  rt_raster_destroy(rtnrast);
1192 
1193  return ES_ERROR;
1194  }
1195 
1196  /* fill _param->offset */
1197  for (i = 0; i < itrcount; i++) {
1198  if (_param->isempty[i])
1199  continue;
1200 
1201  status = rt_raster_from_two_rasters(rtnrast, _param->raster[i], ET_FIRST, &rast, offset);
1202  rtdealloc(rast);
1203  if (status != ES_NONE) {
1204  rterror("rt_raster_iterator: Could not compute raster offsets");
1205 
1206  _rti_iterator_arg_destroy(_param);
1207  rt_band_destroy(rtnband);
1208  rt_raster_destroy(rtnrast);
1209 
1210  return ES_ERROR;
1211  }
1212 
1213  _param->offset[i][0] = offset[2];
1214  _param->offset[i][1] = offset[3];
1215  RASTER_DEBUGF(4, "rast %d offset: %f %f", i, offset[2], offset[3]);
1216  }
1217 
1218  /* loop over each pixel (POI) of output raster */
1219  /* _x,_y are for output raster */
1220  /* x,y are for input raster */
1221  for (_y = 0; _y < _height; _y++) {
1222  for (_x = 0; _x < _width; _x++) {
1223  RASTER_DEBUGF(4, "iterating output pixel (x, y) = (%d, %d)", _x, _y);
1224  _param->arg->dst_pixel[0] = _x;
1225  _param->arg->dst_pixel[1] = _y;
1226 
1227  /* loop through each input raster */
1228  for (i = 0; i < itrcount; i++) {
1229  RASTER_DEBUGF(4, "raster %d", i);
1230 
1231  /*
1232  empty raster
1233  OR band does not exist and flag set to use NODATA
1234  OR band is NODATA
1235  */
1236  if (
1237  _param->isempty[i] ||
1238  (_param->band.rtband[i] == NULL && itrset[i].nbnodata) ||
1239  _param->band.isnodata[i]
1240  ) {
1241  RASTER_DEBUG(4, "empty raster, band does not exist or band is NODATA. using empty values and NODATA");
1242 
1243  x = _x;
1244  y = _y;
1245 
1246  _param->arg->values[i] = _param->empty.values;
1247  _param->arg->nodata[i] = _param->empty.nodata;
1248 
1249  continue;
1250  }
1251 
1252  /* input raster's X,Y */
1253  x = _x - (int) _param->offset[i][0];
1254  y = _y - (int) _param->offset[i][1];
1255  RASTER_DEBUGF(4, "source pixel (x, y) = (%d, %d)", x, y);
1256 
1257  _param->arg->src_pixel[i][0] = x;
1258  _param->arg->src_pixel[i][1] = y;
1259 
1260  /* neighborhood */
1261  npixels = NULL;
1262  status = 0;
1263  if (distancex > 0 && distancey > 0) {
1264  RASTER_DEBUG(4, "getting neighborhood");
1265 
1266  status = rt_band_get_nearest_pixel(
1267  _param->band.rtband[i],
1268  x, y,
1269  distancex, distancey,
1270  1,
1271  &npixels
1272  );
1273  if (status < 0) {
1274  rterror("rt_raster_iterator: Could not get pixel neighborhood");
1275 
1276  _rti_iterator_arg_destroy(_param);
1277  rt_band_destroy(rtnband);
1278  rt_raster_destroy(rtnrast);
1279 
1280  return ES_ERROR;
1281  }
1282  }
1283 
1284  /* get value of POI */
1285  /* get pixel's value */
1286  if (
1287  (x >= 0 && x < _param->width[i]) &&
1288  (y >= 0 && y < _param->height[i])
1289  ) {
1290  RASTER_DEBUG(4, "getting value of POI");
1291  if (rt_band_get_pixel(
1292  _param->band.rtband[i],
1293  x, y,
1294  &value,
1295  &isnodata
1296  ) != ES_NONE) {
1297  rterror("rt_raster_iterator: Could not get the pixel value of band");
1298 
1299  _rti_iterator_arg_destroy(_param);
1300  rt_band_destroy(rtnband);
1301  rt_raster_destroy(rtnrast);
1302 
1303  return ES_ERROR;
1304  }
1305  inextent = 1;
1306  }
1307  /* outside band extent, set to NODATA */
1308  else {
1309  RASTER_DEBUG(4, "Outside band extent, setting value to NODATA");
1310  /* has NODATA, use NODATA */
1311  if (_param->band.hasnodata[i])
1312  value = _param->band.nodataval[i];
1313  /* no NODATA, use min possible value */
1314  else
1315  value = _param->band.minval[i];
1316 
1317  inextent = 0;
1318  isnodata = 1;
1319  }
1320 
1321  /* add pixel to neighborhood */
1322  status++;
1323  if (status > 1)
1324  npixels = (rt_pixel) rtrealloc(npixels, sizeof(struct rt_pixel_t) * status);
1325  else
1326  npixels = (rt_pixel) rtalloc(sizeof(struct rt_pixel_t));
1327 
1328  if (npixels == NULL) {
1329  rterror("rt_raster_iterator: Could not reallocate memory for neighborhood");
1330 
1331  _rti_iterator_arg_destroy(_param);
1332  rt_band_destroy(rtnband);
1333  rt_raster_destroy(rtnrast);
1334 
1335  return ES_ERROR;
1336  }
1337 
1338  npixels[status - 1].x = x;
1339  npixels[status - 1].y = y;
1340  npixels[status - 1].nodata = 1;
1341  npixels[status - 1].value = value;
1342 
1343  /* set nodata flag */
1344  if ((!_param->band.hasnodata[i] && inextent) || !isnodata) {
1345  npixels[status - 1].nodata = 0;
1346  }
1347  RASTER_DEBUGF(4, "value, nodata: %f, %d", value, npixels[status - 1].nodata);
1348 
1349  /* convert set of rt_pixel to 2D array */
1350  status = rt_pixel_set_to_array(
1351  npixels, status, mask,
1352  x, y,
1353  distancex, distancey,
1354  &(_param->arg->values[i]),
1355  &(_param->arg->nodata[i]),
1356  NULL, NULL
1357  );
1358  rtdealloc(npixels);
1359  if (status != ES_NONE) {
1360  rterror("rt_raster_iterator: Could not create 2D array of neighborhood");
1361 
1362  _rti_iterator_arg_destroy(_param);
1363  rt_band_destroy(rtnband);
1364  rt_raster_destroy(rtnrast);
1365 
1366  return ES_ERROR;
1367  }
1368  }
1369 
1370  /* callback */
1371  RASTER_DEBUG(4, "calling callback function");
1372  value = 0;
1373  nodata = 0;
1374  status = callback(_param->arg, userarg, &value, &nodata);
1375 
1376  /* free memory from callback */
1377  _rti_iterator_arg_callback_clean(_param);
1378 
1379  /* handle callback status */
1380  if (status == 0) {
1381  rterror("rt_raster_iterator: Callback function returned an error");
1382 
1383  _rti_iterator_arg_destroy(_param);
1384  rt_band_destroy(rtnband);
1385  rt_raster_destroy(rtnrast);
1386 
1387  return ES_ERROR;
1388  }
1389 
1390  /* burn value to pixel */
1391  status = 0;
1392  if (!nodata) {
1393  status = rt_band_set_pixel(rtnband, _x, _y, value, NULL);
1394  RASTER_DEBUGF(4, "burning pixel (%d, %d) with value: %f", _x, _y, value);
1395  }
1396  else if (!hasnodata) {
1397  status = rt_band_set_pixel(rtnband, _x, _y, minval, NULL);
1398  RASTER_DEBUGF(4, "burning pixel (%d, %d) with minval: %f", _x, _y, minval);
1399  }
1400  else {
1401  RASTER_DEBUGF(4, "NOT burning pixel (%d, %d)", _x, _y);
1402  }
1403  if (status != ES_NONE) {
1404  rterror("rt_raster_iterator: Could not set pixel value");
1405 
1406  _rti_iterator_arg_destroy(_param);
1407  rt_band_destroy(rtnband);
1408  rt_raster_destroy(rtnrast);
1409 
1410  return ES_ERROR;
1411  }
1412  }
1413  }
1414 
1415  /* lots of cleanup */
1416  _rti_iterator_arg_destroy(_param);
1417 
1418  *rtnraster = rtnrast;
1419  return ES_NONE;
1420 }
1421 
1422 /******************************************************************************
1423 * rt_raster_colormap()
1424 ******************************************************************************/
1425 
1426 typedef struct _rti_colormap_arg_t* _rti_colormap_arg;
1427 struct _rti_colormap_arg_t {
1428  rt_raster raster;
1429  rt_band band;
1430 
1431  rt_colormap_entry nodataentry;
1432  int hasnodata;
1433  double nodataval;
1434 
1435  int nexpr;
1436  rt_reclassexpr *expr;
1437 
1438  int npos;
1439  uint16_t *pos;
1440 
1441 };
1442 
1443 static _rti_colormap_arg
1444 _rti_colormap_arg_init(rt_raster raster) {
1445  _rti_colormap_arg arg = NULL;
1446 
1447  arg = rtalloc(sizeof(struct _rti_colormap_arg_t));
1448  if (arg == NULL) {
1449  rterror("_rti_colormap_arg_init: Could not allocate memory for _rti_color_arg");
1450  return NULL;
1451  }
1452 
1453  arg->band = NULL;
1454  arg->nodataentry = NULL;
1455  arg->hasnodata = 0;
1456  arg->nodataval = 0;
1457 
1458  if (raster == NULL)
1459  arg->raster = NULL;
1460  /* raster provided */
1461  else {
1462  arg->raster = rt_raster_clone(raster, 0);
1463  if (arg->raster == NULL) {
1464  rterror("_rti_colormap_arg_init: Could not create output raster");
1465  return NULL;
1466  }
1467  }
1468 
1469  arg->nexpr = 0;
1470  arg->expr = NULL;
1471 
1472  arg->npos = 0;
1473  arg->pos = NULL;
1474 
1475  return arg;
1476 }
1477 
1478 static void
1479 _rti_colormap_arg_destroy(_rti_colormap_arg arg) {
1480  int i = 0;
1481 
1482  if (arg->raster != NULL) {
1483  rt_band band = NULL;
1484 
1485  for (i = rt_raster_get_num_bands(arg->raster) - 1; i >= 0; i--) {
1486  band = rt_raster_get_band(arg->raster, i);
1487  if (band != NULL)
1488  rt_band_destroy(band);
1489  }
1490 
1491  rt_raster_destroy(arg->raster);
1492  }
1493 
1494  if (arg->nexpr) {
1495  for (i = 0; i < arg->nexpr; i++) {
1496  if (arg->expr[i] != NULL)
1497  rtdealloc(arg->expr[i]);
1498  }
1499  rtdealloc(arg->expr);
1500  }
1501 
1502  if (arg->npos)
1503  rtdealloc(arg->pos);
1504 
1505  rtdealloc(arg);
1506  arg = NULL;
1507 }
1508 
1520 rt_raster rt_raster_colormap(
1521  rt_raster raster, int nband,
1522  rt_colormap colormap
1523 ) {
1524  _rti_colormap_arg arg = NULL;
1525  rt_raster rtnraster = NULL;
1526  rt_band band = NULL;
1527  int i = 0;
1528  int j = 0;
1529  int k = 0;
1530 
1531  assert(colormap != NULL);
1532 
1533  /* empty raster */
1534  if (rt_raster_is_empty(raster))
1535  return NULL;
1536 
1537  /* no colormap entries */
1538  if (colormap->nentry < 1) {
1539  rterror("rt_raster_colormap: colormap must have at least one entry");
1540  return NULL;
1541  }
1542 
1543  /* nband is valid */
1544  if (!rt_raster_has_band(raster, nband)) {
1545  rterror("rt_raster_colormap: raster has no band at index %d", nband);
1546  return NULL;
1547  }
1548 
1549  band = rt_raster_get_band(raster, nband);
1550  if (band == NULL) {
1551  rterror("rt_raster_colormap: Could not get band at index %d", nband);
1552  return NULL;
1553  }
1554 
1555  /* init internal variables */
1556  arg = _rti_colormap_arg_init(raster);
1557  if (arg == NULL) {
1558  rterror("rt_raster_colormap: Could not initialize internal variables");
1559  return NULL;
1560  }
1561 
1562  /* handle NODATA */
1563  if (rt_band_get_hasnodata_flag(band)) {
1564  arg->hasnodata = 1;
1565  rt_band_get_nodata(band, &(arg->nodataval));
1566  }
1567 
1568  /* # of colors */
1569  if (colormap->ncolor < 1) {
1570  rterror("rt_raster_colormap: At least one color must be provided");
1571  _rti_colormap_arg_destroy(arg);
1572  return NULL;
1573  }
1574  else if (colormap->ncolor > 4) {
1575  rtinfo("More than four colors indicated. Using only the first four colors");
1576  colormap->ncolor = 4;
1577  }
1578 
1579  /* find non-NODATA entries */
1580  arg->npos = 0;
1581  arg->pos = rtalloc(sizeof(uint16_t) * colormap->nentry);
1582  if (arg->pos == NULL) {
1583  rterror("rt_raster_colormap: Could not allocate memory for valid entries");
1584  _rti_colormap_arg_destroy(arg);
1585  return NULL;
1586  }
1587  for (i = 0, j = 0; i < colormap->nentry; i++) {
1588  /* special handling for NODATA entries */
1589  if (colormap->entry[i].isnodata) {
1590  /* first NODATA entry found, use it */
1591  if (arg->nodataentry == NULL)
1592  arg->nodataentry = &(colormap->entry[i]);
1593  else
1594  rtwarn("More than one colormap entry found for NODATA value. Only using first NOTDATA entry");
1595 
1596  continue;
1597  }
1598 
1599  (arg->npos)++;
1600  arg->pos[j++] = i;
1601  }
1602 
1603  /* INTERPOLATE and only one non-NODATA entry */
1604  if (colormap->method == CM_INTERPOLATE && arg->npos < 2) {
1605  rtwarn("Method INTERPOLATE requires at least two non-NODATA colormap entries. Using NEAREST instead");
1606  colormap->method = CM_NEAREST;
1607  }
1608 
1609  /* NODATA entry but band has no NODATA value */
1610  if (!arg->hasnodata && arg->nodataentry != NULL) {
1611  rtinfo("Band at index %d has no NODATA value. Ignoring NODATA entry", nband);
1612  arg->nodataentry = NULL;
1613  }
1614 
1615  /* allocate expr */
1616  arg->nexpr = arg->npos;
1617 
1618  /* INTERPOLATE needs one less than the number of entries */
1619  if (colormap->method == CM_INTERPOLATE)
1620  arg->nexpr -= 1;
1621  /* EXACT requires a no matching expression */
1622  else if (colormap->method == CM_EXACT)
1623  arg->nexpr += 1;
1624 
1625  /* NODATA entry exists, add expression */
1626  if (arg->nodataentry != NULL)
1627  arg->nexpr += 1;
1628  arg->expr = rtalloc(sizeof(rt_reclassexpr) * arg->nexpr);
1629  if (arg->expr == NULL) {
1630  rterror("rt_raster_colormap: Could not allocate memory for reclass expressions");
1631  _rti_colormap_arg_destroy(arg);
1632  return NULL;
1633  }
1634  RASTER_DEBUGF(4, "nexpr = %d", arg->nexpr);
1635  RASTER_DEBUGF(4, "expr @ %p", arg->expr);
1636 
1637  for (i = 0; i < arg->nexpr; i++) {
1638  arg->expr[i] = rtalloc(sizeof(struct rt_reclassexpr_t));
1639  if (arg->expr[i] == NULL) {
1640  rterror("rt_raster_colormap: Could not allocate memory for reclass expression");
1641  _rti_colormap_arg_destroy(arg);
1642  return NULL;
1643  }
1644  }
1645 
1646  /* reclassify bands */
1647  /* by # of colors */
1648  for (i = 0; i < colormap->ncolor; i++) {
1649  k = 0;
1650 
1651  /* handle NODATA entry first */
1652  if (arg->nodataentry != NULL) {
1653  arg->expr[k]->src.min = arg->nodataentry->value;
1654  arg->expr[k]->src.max = arg->nodataentry->value;
1655  arg->expr[k]->src.inc_min = 1;
1656  arg->expr[k]->src.inc_max = 1;
1657  arg->expr[k]->src.exc_min = 0;
1658  arg->expr[k]->src.exc_max = 0;
1659 
1660  arg->expr[k]->dst.min = arg->nodataentry->color[i];
1661  arg->expr[k]->dst.max = arg->nodataentry->color[i];
1662 
1663  arg->expr[k]->dst.inc_min = 1;
1664  arg->expr[k]->dst.inc_max = 1;
1665  arg->expr[k]->dst.exc_min = 0;
1666  arg->expr[k]->dst.exc_max = 0;
1667 
1668  RASTER_DEBUGF(4, "NODATA expr[%d]->src (min, max, in, ix, en, ex) = (%f, %f, %d, %d, %d, %d)",
1669  k,
1670  arg->expr[k]->src.min,
1671  arg->expr[k]->src.max,
1672  arg->expr[k]->src.inc_min,
1673  arg->expr[k]->src.inc_max,
1674  arg->expr[k]->src.exc_min,
1675  arg->expr[k]->src.exc_max
1676  );
1677  RASTER_DEBUGF(4, "NODATA expr[%d]->dst (min, max, in, ix, en, ex) = (%f, %f, %d, %d, %d, %d)",
1678  k,
1679  arg->expr[k]->dst.min,
1680  arg->expr[k]->dst.max,
1681  arg->expr[k]->dst.inc_min,
1682  arg->expr[k]->dst.inc_max,
1683  arg->expr[k]->dst.exc_min,
1684  arg->expr[k]->dst.exc_max
1685  );
1686 
1687  k++;
1688  }
1689 
1690  /* by non-NODATA entry */
1691  for (j = 0; j < arg->npos; j++) {
1692  if (colormap->method == CM_INTERPOLATE) {
1693  if (j == arg->npos - 1)
1694  continue;
1695 
1696  arg->expr[k]->src.min = colormap->entry[arg->pos[j + 1]].value;
1697  arg->expr[k]->src.inc_min = 1;
1698  arg->expr[k]->src.exc_min = 0;
1699 
1700  arg->expr[k]->src.max = colormap->entry[arg->pos[j]].value;
1701  arg->expr[k]->src.inc_max = 1;
1702  arg->expr[k]->src.exc_max = 0;
1703 
1704  arg->expr[k]->dst.min = colormap->entry[arg->pos[j + 1]].color[i];
1705  arg->expr[k]->dst.max = colormap->entry[arg->pos[j]].color[i];
1706 
1707  arg->expr[k]->dst.inc_min = 1;
1708  arg->expr[k]->dst.exc_min = 0;
1709 
1710  arg->expr[k]->dst.inc_max = 1;
1711  arg->expr[k]->dst.exc_max = 0;
1712  }
1713  else if (colormap->method == CM_NEAREST) {
1714 
1715  /* NOT last entry */
1716  if (j != arg->npos - 1) {
1717  arg->expr[k]->src.min = ((colormap->entry[arg->pos[j]].value - colormap->entry[arg->pos[j + 1]].value) / 2.) + colormap->entry[arg->pos[j + 1]].value;
1718  arg->expr[k]->src.inc_min = 0;
1719  arg->expr[k]->src.exc_min = 0;
1720  }
1721  /* last entry */
1722  else {
1723  arg->expr[k]->src.min = colormap->entry[arg->pos[j]].value;
1724  arg->expr[k]->src.inc_min = 1;
1725  arg->expr[k]->src.exc_min = 1;
1726  }
1727 
1728  /* NOT first entry */
1729  if (j > 0) {
1730  arg->expr[k]->src.max = arg->expr[k - 1]->src.min;
1731  arg->expr[k]->src.inc_max = 1;
1732  arg->expr[k]->src.exc_max = 0;
1733  }
1734  /* first entry */
1735  else {
1736  arg->expr[k]->src.max = colormap->entry[arg->pos[j]].value;
1737  arg->expr[k]->src.inc_max = 1;
1738  arg->expr[k]->src.exc_max = 1;
1739  }
1740 
1741  arg->expr[k]->dst.min = colormap->entry[arg->pos[j]].color[i];
1742  arg->expr[k]->dst.inc_min = 1;
1743  arg->expr[k]->dst.exc_min = 0;
1744 
1745  arg->expr[k]->dst.max = colormap->entry[arg->pos[j]].color[i];
1746  arg->expr[k]->dst.inc_max = 1;
1747  arg->expr[k]->dst.exc_max = 0;
1748  }
1749  else if (colormap->method == CM_EXACT) {
1750  arg->expr[k]->src.min = colormap->entry[arg->pos[j]].value;
1751  arg->expr[k]->src.inc_min = 1;
1752  arg->expr[k]->src.exc_min = 0;
1753 
1754  arg->expr[k]->src.max = colormap->entry[arg->pos[j]].value;
1755  arg->expr[k]->src.inc_max = 1;
1756  arg->expr[k]->src.exc_max = 0;
1757 
1758  arg->expr[k]->dst.min = colormap->entry[arg->pos[j]].color[i];
1759  arg->expr[k]->dst.inc_min = 1;
1760  arg->expr[k]->dst.exc_min = 0;
1761 
1762  arg->expr[k]->dst.max = colormap->entry[arg->pos[j]].color[i];
1763  arg->expr[k]->dst.inc_max = 1;
1764  arg->expr[k]->dst.exc_max = 0;
1765  }
1766 
1767  RASTER_DEBUGF(4, "expr[%d]->src (min, max, in, ix, en, ex) = (%f, %f, %d, %d, %d, %d)",
1768  k,
1769  arg->expr[k]->src.min,
1770  arg->expr[k]->src.max,
1771  arg->expr[k]->src.inc_min,
1772  arg->expr[k]->src.inc_max,
1773  arg->expr[k]->src.exc_min,
1774  arg->expr[k]->src.exc_max
1775  );
1776 
1777  RASTER_DEBUGF(4, "expr[%d]->dst (min, max, in, ix, en, ex) = (%f, %f, %d, %d, %d, %d)",
1778  k,
1779  arg->expr[k]->dst.min,
1780  arg->expr[k]->dst.max,
1781  arg->expr[k]->dst.inc_min,
1782  arg->expr[k]->dst.inc_max,
1783  arg->expr[k]->dst.exc_min,
1784  arg->expr[k]->dst.exc_max
1785  );
1786 
1787  k++;
1788  }
1789 
1790  /* EXACT has one last expression for catching all uncaught values */
1791  if (colormap->method == CM_EXACT) {
1792  arg->expr[k]->src.min = 0;
1793  arg->expr[k]->src.inc_min = 1;
1794  arg->expr[k]->src.exc_min = 1;
1795 
1796  arg->expr[k]->src.max = 0;
1797  arg->expr[k]->src.inc_max = 1;
1798  arg->expr[k]->src.exc_max = 1;
1799 
1800  arg->expr[k]->dst.min = 0;
1801  arg->expr[k]->dst.inc_min = 1;
1802  arg->expr[k]->dst.exc_min = 0;
1803 
1804  arg->expr[k]->dst.max = 0;
1805  arg->expr[k]->dst.inc_max = 1;
1806  arg->expr[k]->dst.exc_max = 0;
1807 
1808  RASTER_DEBUGF(4, "expr[%d]->src (min, max, in, ix, en, ex) = (%f, %f, %d, %d, %d, %d)",
1809  k,
1810  arg->expr[k]->src.min,
1811  arg->expr[k]->src.max,
1812  arg->expr[k]->src.inc_min,
1813  arg->expr[k]->src.inc_max,
1814  arg->expr[k]->src.exc_min,
1815  arg->expr[k]->src.exc_max
1816  );
1817 
1818  RASTER_DEBUGF(4, "expr[%d]->dst (min, max, in, ix, en, ex) = (%f, %f, %d, %d, %d, %d)",
1819  k,
1820  arg->expr[k]->dst.min,
1821  arg->expr[k]->dst.max,
1822  arg->expr[k]->dst.inc_min,
1823  arg->expr[k]->dst.inc_max,
1824  arg->expr[k]->dst.exc_min,
1825  arg->expr[k]->dst.exc_max
1826  );
1827 
1828  k++;
1829  }
1830 
1831  /* call rt_band_reclass */
1832  arg->band = rt_band_reclass(band, PT_8BUI, 0, 0, arg->expr, arg->nexpr);
1833  if (arg->band == NULL) {
1834  rterror("rt_raster_colormap: Could not reclassify band");
1835  _rti_colormap_arg_destroy(arg);
1836  return NULL;
1837  }
1838 
1839  /* add reclassified band to raster */
1840  if (rt_raster_add_band(arg->raster, arg->band, rt_raster_get_num_bands(arg->raster)) < 0) {
1841  rterror("rt_raster_colormap: Could not add reclassified band to output raster");
1842  _rti_colormap_arg_destroy(arg);
1843  return NULL;
1844  }
1845  }
1846 
1847  rtnraster = arg->raster;
1848  arg->raster = NULL;
1849  _rti_colormap_arg_destroy(arg);
1850 
1851  return rtnraster;
1852 }
rt_util_clamp_to_32BSI
int32_t rt_util_clamp_to_32BSI(double value)
Definition: rt_util.c:69
_rti_iterator_arg_t::rtband
rt_band * rtband
Definition: rt_mapalgebra.c:396
_rti_colormap_arg_t::expr
rt_reclassexpr * expr
Definition: rt_mapalgebra.c:1436
_rti_iterator_arg_t::distance
struct _rti_iterator_arg_t::@9 distance
_rti_colormap_arg_t::pos
uint16_t * pos
Definition: rt_mapalgebra.c:1439
rt_colormap_t::nentry
uint16_t nentry
Definition: librtcore.h:2444
rt_util_clamp_to_32BUI
uint32_t rt_util_clamp_to_32BUI(double value)
Definition: rt_util.c:74
rt_band_t
Definition: librtcore.h:2264
rt_pixel_t::y
int y
Definition: librtcore.h:2286
rt_raster_get_x_offset
double rt_raster_get_x_offset(rt_raster raster)
Get raster x offset, in projection units.
Definition: rt_raster.c:213
rt_band_reclass
rt_band rt_band_reclass(rt_band srcband, rt_pixtype pixtype, uint32_t hasnodata, double nodataval, rt_reclassexpr *exprset, int exprcount)
Returns new band with values reclassified.
Definition: rt_mapalgebra.c:50
rt_raster_get_num_bands
int rt_raster_get_num_bands(rt_raster raster)
Definition: rt_raster.c:372
rt_reclassexpr_t::rt_reclassrange::exc_min
int exc_min
Definition: librtcore.h:2387
rt_raster_get_y_skew
double rt_raster_get_y_skew(rt_raster raster)
Get skew about the Y axis.
Definition: rt_raster.c:190
rt_reclassexpr_t::rt_reclassrange::inc_max
int inc_max
Definition: librtcore.h:2386
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
_rti_iterator_arg_t::dimension
struct _rti_iterator_arg_t::@10 dimension
rt_raster_t::numBands
uint16_t numBands
Definition: librtcore.h:2241
rt_pixel
struct rt_pixel_t * rt_pixel
Definition: librtcore.h:147
rt_pixel_set_to_array
rt_errorstate rt_pixel_set_to_array(rt_pixel npixel, int count, rt_mask mask, int x, int y, uint16_t distancex, uint16_t distancey, double ***value, int ***nodata, int *dimx, int *dimy)
Definition: rt_pixel.c:286
_rti_colormap_arg_init
static _rti_colormap_arg _rti_colormap_arg_init(rt_raster raster)
Definition: rt_mapalgebra.c:1444
_rti_iterator_arg_t::isnodata
int * isnodata
Definition: rt_mapalgebra.c:398
_rti_iterator_arg_init
static _rti_iterator_arg _rti_iterator_arg_init()
Definition: rt_mapalgebra.c:422
rt_reclassexpr_t::dst
struct rt_reclassexpr_t::rt_reclassrange dst
rt_iterator_arg_t::src_pixel
int ** src_pixel
Definition: librtcore.h:2415
_rti_iterator_arg_t::raster
rt_raster * raster
Definition: rt_mapalgebra.c:389
_rti_iterator_arg_empty_init
static int _rti_iterator_arg_empty_init(_rti_iterator_arg _param)
Definition: rt_mapalgebra.c:670
ovdump.band
band
Definition: ovdump.py:57
librtcore_internal.h
_rti_colormap_arg_t::nexpr
int nexpr
Definition: rt_mapalgebra.c:1435
rt_raster_from_two_rasters
rt_errorstate rt_raster_from_two_rasters(rt_raster rast1, rt_raster rast2, rt_extenttype extenttype, rt_raster *rtnraster, double *offset)
Definition: rt_raster.c:3444
_rti_iterator_arg_t::width
int * width
Definition: rt_mapalgebra.c:392
_rti_iterator_arg_t::y
uint16_t y
Definition: rt_mapalgebra.c:405
rt_extenttype
rt_extenttype
Definition: librtcore.h:200
_rti_iterator_arg_t
Definition: rt_mapalgebra.c:386
rt_colormap_entry_t::color
uint8_t color[4]
Definition: librtcore.h:2433
rt_reclassexpr_t::rt_reclassrange::inc_min
int inc_min
Definition: librtcore.h:2385
FLT_EQ
#define FLT_EQ(x, y)
Definition: librtcore.h:2185
rterror
void rterror(const char *fmt,...)
Wrappers used for reporting errors and info.
Definition: rt_context.c:199
rtalloc
void * rtalloc(size_t size)
Wrappers used for managing memory.
Definition: rt_context.c:171
_rti_iterator_arg_t::arg
rt_iterator_arg arg
Definition: rt_mapalgebra.c:418
_rti_iterator_arg_t::isempty
int * isempty
Definition: rt_mapalgebra.c:390
rt_iterator_arg_t::rows
uint32_t rows
Definition: librtcore.h:2404
rtrealloc
void * rtrealloc(void *mem, size_t size)
Definition: rt_context.c:179
rt_pixel_t::value
double value
Definition: librtcore.h:2289
rt_errorstate
rt_errorstate
Enum definitions.
Definition: librtcore.h:179
rt_colormap_t::ncolor
int ncolor
Definition: librtcore.h:2443
rt_mask_t
Definition: librtcore.h:2294
_rti_colormap_arg_t::raster
rt_raster raster
Definition: rt_mapalgebra.c:1428
rt_pixtype
rt_pixtype
Definition: librtcore.h:185
rt_band_destroy
void rt_band_destroy(rt_band band)
Destroy a raster band.
Definition: rt_band.c:242
_rti_colormap_arg_destroy
static void _rti_colormap_arg_destroy(_rti_colormap_arg arg)
Definition: rt_mapalgebra.c:1479
rt_raster_is_empty
int rt_raster_is_empty(rt_raster raster)
Return TRUE if the raster is empty.
Definition: rt_raster.c:1338
_rti_colormap_arg_t::hasnodata
int hasnodata
Definition: rt_mapalgebra.c:1432
_rti_colormap_arg
struct _rti_colormap_arg_t * _rti_colormap_arg
Definition: rt_mapalgebra.c:1426
rt_band_get_nodata
rt_errorstate rt_band_get_nodata(rt_band band, double *nodata)
Get NODATA value.
Definition: rt_band.c:1597
uint32_t
unsigned int uint32_t
Definition: uthash.h:78
rt_util_clamp_to_16BUI
uint16_t rt_util_clamp_to_16BUI(double value)
Definition: rt_util.c:64
rt_band_set_ownsdata_flag
void rt_band_set_ownsdata_flag(rt_band band, int flag)
Definition: rt_band.c:534
pixval.nband
nband
Definition: pixval.py:52
rt_reclassexpr_t::src
struct rt_reclassexpr_t::rt_reclassrange src
rt_pixel_t
Definition: librtcore.h:2284
_rti_iterator_arg_t::values
double ** values
Definition: rt_mapalgebra.c:414
_rti_colormap_arg_t::nodataval
double nodataval
Definition: rt_mapalgebra.c:1433
rt_iterator_arg_t::columns
uint32_t columns
Definition: librtcore.h:2406
rt_reclassexpr_t::rt_reclassrange::exc_max
int exc_max
Definition: librtcore.h:2388
rt_raster_iterator
rt_errorstate rt_raster_iterator(rt_iterator itrset, uint16_t itrcount, rt_extenttype extenttype, rt_raster customextent, rt_pixtype pixtype, uint8_t hasnodata, double nodataval, uint16_t distancex, uint16_t distancey, rt_mask mask, void *userarg, int(*callback)(rt_iterator_arg arg, void *userarg, double *value, int *nodata), rt_raster *rtnraster)
n-raster iterator.
Definition: rt_mapalgebra.c:812
rtwarn
void rtwarn(const char *fmt,...)
Definition: rt_context.c:224
rt_iterator_arg_t::rasters
uint16_t rasters
Definition: librtcore.h:2402
_rti_iterator_arg_t::empty
struct _rti_iterator_arg_t::@11 empty
rt_band_get_pixel
rt_errorstate rt_band_get_pixel(rt_band band, int x, int y, double *value, int *nodata)
Get pixel value.
Definition: rt_band.c:1088
rt_colormap_t::method
enum rt_colormap_t::@7 method
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_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_iterator_arg_t
Definition: librtcore.h:2400
rt_raster_clone
rt_raster rt_raster_clone(rt_raster raster, uint8_t deep)
Clone an existing raster.
Definition: rt_raster.c:1544
_rti_iterator_arg_t::height
int * height
Definition: rt_mapalgebra.c:393
rt_util_clamp_to_8BSI
int8_t rt_util_clamp_to_8BSI(double value)
Definition: rt_util.c:49
rt_util_dbl_trunc_warning
int rt_util_dbl_trunc_warning(double initialvalue, int32_t checkvalint, uint32_t checkvaluint, float checkvalfloat, double checkvaldouble, rt_pixtype pixtype)
Definition: rt_util.c:625
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
_rti_iterator_arg_t::nodataval
double * nodataval
Definition: rt_mapalgebra.c:399
rtinfo
void rtinfo(const char *fmt,...)
Definition: rt_context.c:211
rt_band_get_hasnodata_flag
int rt_band_get_hasnodata_flag(rt_band band)
Get hasnodata flag value.
Definition: rt_band.c:541
rt_raster_get_x_scale
double rt_raster_get_x_scale(rt_raster raster)
Get scale X in projection units.
Definition: rt_raster.c:150
rt_band_get_width
uint16_t rt_band_get_width(rt_band band)
Return width of this band.
Definition: rt_band.c:507
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:1351
rt_band_get_height
uint16_t rt_band_get_height(rt_band band)
Return height of this band.
Definition: rt_band.c:516
rt_colormap_t::entry
rt_colormap_entry entry
Definition: librtcore.h:2445
RASTER_DEBUGF
#define RASTER_DEBUGF(level, msg,...)
Definition: librtcore.h:299
rt_raster_get_srid
int32_t rt_raster_get_srid(rt_raster raster)
Get raster&#39;s SRID.
Definition: rt_raster.c:356
rt_raster_destroy
void rt_raster_destroy(rt_raster raster)
Release memory associated to a raster.
Definition: rt_raster.c:82
rt_iterator_t::raster
rt_raster raster
Definition: librtcore.h:2394
_rti_iterator_arg_t::hasnodata
int * hasnodata
Definition: rt_mapalgebra.c:397
_rti_iterator_arg_t::count
int count
Definition: rt_mapalgebra.c:387
_rti_colormap_arg_t::npos
int npos
Definition: rt_mapalgebra.c:1438
_rti_iterator_arg_t::offset
double ** offset
Definition: rt_mapalgebra.c:391
rt_pixel_t::nodata
uint8_t nodata
Definition: librtcore.h:2288
rt_raster_get_y_scale
double rt_raster_get_y_scale(rt_raster raster)
Get scale Y in projection units.
Definition: rt_raster.c:159
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_get_x_skew
double rt_raster_get_x_skew(rt_raster raster)
Get skew about the X axis.
Definition: rt_raster.c:181
rt_band_new_inline
rt_band rt_band_new_inline(uint16_t width, uint16_t height, rt_pixtype pixtype, uint32_t hasnodata, double nodataval, uint8_t *data)
Create an in-db rt_band with no data.
Definition: rt_band.c:58
rt_util_clamp_to_4BUI
uint8_t rt_util_clamp_to_4BUI(double value)
Definition: rt_util.c:44
rt_util_clamp_to_8BUI
uint8_t rt_util_clamp_to_8BUI(double value)
Definition: rt_util.c:54
_rti_iterator_arg_t::nodata
int ** nodata
Definition: rt_mapalgebra.c:415
_rti_iterator_arg_t::columns
uint32_t columns
Definition: rt_mapalgebra.c:410
rt_pixtype_size
int rt_pixtype_size(rt_pixtype pixtype)
Return size in bytes of a value in the given pixtype.
Definition: rt_pixel.c:39
rt_util_clamp_to_1BB
uint8_t rt_util_clamp_to_1BB(double value)
Definition: rt_util.c:34
_rti_iterator_arg_callback_clean
static void _rti_iterator_arg_callback_clean(_rti_iterator_arg _param)
Definition: rt_mapalgebra.c:746
rt_raster_get_width
uint16_t rt_raster_get_width(rt_raster raster)
Definition: rt_raster.c:121
rt_iterator_arg_t::values
double *** values
Definition: librtcore.h:2410
_rti_iterator_arg_t::minval
double * minval
Definition: rt_mapalgebra.c:400
rtdealloc
void rtdealloc(void *mem)
Definition: rt_context.c:186
rt_colormap_entry_t
Definition: librtcore.h:2430
_rti_colormap_arg_t
Definition: rt_mapalgebra.c:1427
_rti_iterator_arg_t::rows
uint32_t rows
Definition: rt_mapalgebra.c:409
_rti_iterator_arg_callback_init
static int _rti_iterator_arg_callback_init(_rti_iterator_arg _param)
Definition: rt_mapalgebra.c:700
rt_iterator_t::nbnodata
uint8_t nbnodata
Definition: librtcore.h:2396
RASTER_DEBUG
#define RASTER_DEBUG(level, msg)
Definition: librtcore.h:295
rt_raster_serialized_t
Struct definitions.
Definition: librtcore.h:2201
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:1612
librtcore.h
This library is the generic raster handling section of PostGIS.
_rti_iterator_arg_t::x
uint16_t x
Definition: rt_mapalgebra.c:404
rt_raster_t::bands
rt_band * bands
Definition: librtcore.h:2254
rt_iterator_t
Definition: librtcore.h:2393
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
_rti_iterator_arg
struct _rti_iterator_arg_t * _rti_iterator_arg
Definition: rt_mapalgebra.c:385
rt_util_clamp_to_2BUI
uint8_t rt_util_clamp_to_2BUI(double value)
Definition: rt_util.c:39
rt_band_get_isnodata_flag
int rt_band_get_isnodata_flag(rt_band band)
Get isnodata flag value.
Definition: rt_band.c:581
_rti_colormap_arg_t::band
rt_band band
Definition: rt_mapalgebra.c:1429
rt_raster_t
Definition: librtcore.h:2235
rt_util_clamp_to_32F
float rt_util_clamp_to_32F(double value)
Definition: rt_util.c:79
rt_band_set_pixel
rt_errorstate rt_band_set_pixel(rt_band band, int x, int y, double val, int *converted)
Set single pixel&#39;s value.
Definition: rt_band.c:841
rt_reclassexpr_t::rt_reclassrange::max
double max
Definition: librtcore.h:2384
_rti_iterator_arg_populate
static int _rti_iterator_arg_populate(_rti_iterator_arg _param, rt_iterator itrset, uint16_t itrcount, uint16_t distancex, uint16_t distancey, int *allnull, int *allempty)
Definition: rt_mapalgebra.c:531
genraster.value
int value
Definition: genraster.py:61
rt_raster_get_height
uint16_t rt_raster_get_height(rt_raster raster)
Definition: rt_raster.c:129
_rti_iterator_arg_destroy
static void _rti_iterator_arg_destroy(_rti_iterator_arg _param)
Definition: rt_mapalgebra.c:460
rt_iterator_arg_t::nodata
int *** nodata
Definition: librtcore.h:2412
rt_colormap_t
Definition: librtcore.h:2436
rt_colormap_entry_t::isnodata
int isnodata
Definition: librtcore.h:2431
uint8_t
unsigned char uint8_t
Definition: uthash.h:79
rt_pixel_t::x
int x
Definition: librtcore.h:2285
rt_raster_colormap
rt_raster rt_raster_colormap(rt_raster raster, int nband, rt_colormap colormap)
Returns a new raster with up to four 8BUI bands (RGBA) from applying a colormap to the user-specified...
Definition: rt_mapalgebra.c:1520
rt_band_get_nearest_pixel
int rt_band_get_nearest_pixel(rt_band band, int x, int y, uint16_t distancex, uint16_t distancey, int exclude_nodata_value, rt_pixel *npixels)
Get nearest pixel(s) with value (not NODATA) to specified pixel.
Definition: rt_band.c:1241
_rti_iterator_arg_t::band
struct _rti_iterator_arg_t::@8 band
rt_iterator_arg_t::dst_pixel
int dst_pixel[2]
Definition: librtcore.h:2418
rt_raster_get_y_offset
double rt_raster_get_y_offset(rt_raster raster)
Get raster y offset, in projection units.
Definition: rt_raster.c:222
rt_colormap_entry_t::value
double value
Definition: librtcore.h:2432
rt_reclassexpr_t
Definition: librtcore.h:2381
if
if(!(yy_init))
Definition: lwin_wkt_lex.c:1080
rt_reclassexpr_t::rt_reclassrange::min
double min
Definition: librtcore.h:2383
_rti_colormap_arg_t::nodataentry
rt_colormap_entry nodataentry
Definition: rt_mapalgebra.c:1431
rt_util_clamp_to_16BSI
int16_t rt_util_clamp_to_16BSI(double value)
Definition: rt_util.c:59