PostGIS  2.5.7dev-r@@SVN_REVISION@@
rt_pixel.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  * Copyright (C) 2013 Nathaniel Hunter Clay <clay.nathaniel@gmail.com>
14  *
15  * This program is free software; you can redistribute it and/or
16  * modify it under the terms of the GNU General Public License
17  * as published by the Free Software Foundation; either version 2
18  * of the License, or (at your option) any later version.
19  *
20  * This program is distributed in the hope that it will be useful,
21  * but WITHOUT ANY WARRANTY; without even the implied warranty of
22  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
23  * GNU General Public License for more details.
24  *
25  * You should have received a copy of the GNU General Public License
26  * along with this program; if not, write to the Free Software Foundation,
27  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
28  *
29  */
30 
31 #include "librtcore.h"
32 #include "librtcore_internal.h"
33 
34 /******************************************************************************
35 * rt_pixeltype
36 ******************************************************************************/
37 
38 int
40  int pixbytes = -1;
41 
42  switch (pixtype) {
43  case PT_1BB:
44  case PT_2BUI:
45  case PT_4BUI:
46  case PT_8BSI:
47  case PT_8BUI:
48  pixbytes = 1;
49  break;
50  case PT_16BSI:
51  case PT_16BUI:
52  pixbytes = 2;
53  break;
54  case PT_32BSI:
55  case PT_32BUI:
56  case PT_32BF:
57  pixbytes = 4;
58  break;
59  case PT_64BF:
60  pixbytes = 8;
61  break;
62  default:
63  rterror("rt_pixtype_size: Unknown pixeltype %d", pixtype);
64  pixbytes = -1;
65  break;
66  }
67 
68  RASTER_DEBUGF(3, "Pixel type = %s and size = %d bytes",
69  rt_pixtype_name(pixtype), pixbytes);
70 
71  return pixbytes;
72 }
73 
74 int
76  return rt_pixtype_size(pixtype);
77 }
78 
80 rt_pixtype_index_from_name(const char* pixname) {
81  assert(pixname && strlen(pixname) > 0);
82 
83  if (strcmp(pixname, "1BB") == 0)
84  return PT_1BB;
85  else if (strcmp(pixname, "2BUI") == 0)
86  return PT_2BUI;
87  else if (strcmp(pixname, "4BUI") == 0)
88  return PT_4BUI;
89  else if (strcmp(pixname, "8BSI") == 0)
90  return PT_8BSI;
91  else if (strcmp(pixname, "8BUI") == 0)
92  return PT_8BUI;
93  else if (strcmp(pixname, "16BSI") == 0)
94  return PT_16BSI;
95  else if (strcmp(pixname, "16BUI") == 0)
96  return PT_16BUI;
97  else if (strcmp(pixname, "32BSI") == 0)
98  return PT_32BSI;
99  else if (strcmp(pixname, "32BUI") == 0)
100  return PT_32BUI;
101  else if (strcmp(pixname, "32BF") == 0)
102  return PT_32BF;
103  else if (strcmp(pixname, "64BF") == 0)
104  return PT_64BF;
105 
106  return PT_END;
107 }
108 
109 const char*
111  switch (pixtype) {
112  case PT_1BB:
113  return "1BB";
114  case PT_2BUI:
115  return "2BUI";
116  case PT_4BUI:
117  return "4BUI";
118  case PT_8BSI:
119  return "8BSI";
120  case PT_8BUI:
121  return "8BUI";
122  case PT_16BSI:
123  return "16BSI";
124  case PT_16BUI:
125  return "16BUI";
126  case PT_32BSI:
127  return "32BSI";
128  case PT_32BUI:
129  return "32BUI";
130  case PT_32BF:
131  return "32BF";
132  case PT_64BF:
133  return "64BF";
134  default:
135  rterror("rt_pixtype_name: Unknown pixeltype %d", pixtype);
136  return "Unknown";
137  }
138 }
139 
147 double
149  switch (pixtype) {
150  case PT_1BB: {
151  return (double) rt_util_clamp_to_1BB((double) CHAR_MIN);
152  }
153  case PT_2BUI: {
154  return 0;
155  }
156  case PT_4BUI: {
157  return 0;
158  }
159  case PT_8BUI: {
160  return 0;
161  }
162  case PT_8BSI: {
163  return (double) rt_util_clamp_to_8BSI((double) SCHAR_MIN);
164  }
165  case PT_16BSI: {
166  return (double) rt_util_clamp_to_16BSI((double) SHRT_MIN);
167  }
168  case PT_16BUI: {
169  return 0;
170  }
171  case PT_32BSI: {
172  return (double) rt_util_clamp_to_32BSI((double) INT_MIN);
173  }
174  case PT_32BUI: {
175  return 0;
176  }
177  case PT_32BF: {
178  return (double) -FLT_MAX;
179  }
180  case PT_64BF: {
181  return (double) -DBL_MAX;
182  }
183  default: {
184  rterror("rt_pixtype_get_min_value: Unknown pixeltype %d", pixtype);
185  return (double) rt_util_clamp_to_8BUI((double) CHAR_MIN);
186  }
187  }
188 }
189 
201  rt_pixtype pixtype,
202  double val, double refval,
203  int *isequal
204 ) {
205  assert(isequal != NULL);
206  *isequal = 0;
207 
208  switch (pixtype) {
209  case PT_1BB:
210  if (rt_util_clamp_to_1BB(val) == rt_util_clamp_to_1BB(refval))
211  *isequal = 1;
212  break;
213  case PT_2BUI:
214  if (rt_util_clamp_to_2BUI(val) == rt_util_clamp_to_2BUI(refval))
215  *isequal = 1;
216  break;
217  case PT_4BUI:
218  if (rt_util_clamp_to_4BUI(val) == rt_util_clamp_to_4BUI(refval))
219  *isequal = 1;
220  break;
221  case PT_8BSI:
222  if (rt_util_clamp_to_8BSI(val) == rt_util_clamp_to_8BSI(refval))
223  *isequal = 1;
224  break;
225  case PT_8BUI:
226  if (rt_util_clamp_to_8BUI(val) == rt_util_clamp_to_8BUI(refval))
227  *isequal = 1;
228  break;
229  case PT_16BSI:
231  *isequal = 1;
232  break;
233  case PT_16BUI:
235  *isequal = 1;
236  break;
237  case PT_32BSI:
239  *isequal = 1;
240  break;
241  case PT_32BUI:
243  *isequal = 1;
244  break;
245  case PT_32BF:
247  *isequal = 1;
248  break;
249  case PT_64BF:
250  if (FLT_EQ(val, refval))
251  *isequal = 1;
252  break;
253  default:
254  rterror("rt_pixtype_compare_clamped_values: Unknown pixeltype %d", pixtype);
255  return ES_ERROR;
256  }
257 
258  return ES_NONE;
259 }
260 
261 /******************************************************************************
262 * rt_pixel
263 ******************************************************************************/
264 
265 /*
266  * Convert an array of rt_pixel objects to two 2D arrays of value and NODATA.
267  * The dimensions of the returned 2D array are [Y][X], going by row Y and
268  * then column X.
269  *
270  * @param npixel : array of rt_pixel objects
271  * @param count : number of elements in npixel
272  * @param mask : mask to be respected when retruning array
273  * @param x : the column of the center pixel (0-based)
274  * @param y : the line of the center pixel (0-based)
275  * @param distancex : the number of pixels around the specified pixel
276  * along the X axis
277  * @param distancey : the number of pixels around the specified pixel
278  * along the Y axis
279  * @param value : pointer to pointer for 2D value array
280  * @param nodata : pointer to pointer for 2D NODATA array
281  * @param dimx : size of value and nodata along the X axis
282  * @param dimy : size of value and nodata along the Y axis
283  *
284  * @return ES_NONE on success, ES_ERROR on error
285  */
287  rt_pixel npixel, uint32_t count, rt_mask mask,
288  int x, int y,
289  uint16_t distancex, uint16_t distancey,
290  double ***value,
291  int ***nodata,
292  int *dimx, int *dimy
293 ) {
294  uint32_t i;
295  uint32_t j;
296  uint32_t dim[2] = {0};
297  double **values = NULL;
298  int **nodatas = NULL;
299  int zero[2] = {0};
300  int _x;
301  int _y;
302 
303  assert(npixel != NULL && count > 0);
304  assert(value != NULL);
305  assert(nodata != NULL);
306 
307  /* dimensions */
308  dim[0] = distancex * 2 + 1;
309  dim[1] = distancey * 2 + 1;
310  RASTER_DEBUGF(4, "dimensions = %d x %d", dim[0], dim[1]);
311 
312  /* make sure that the dimx and dimy match mask */
313  if (mask != NULL) {
314  if (dim[0] != mask->dimx || dim[1] != mask->dimy) {
315  rterror("rt_pixel_set_array: mask dimensions %d x %d do not match given dims %d x %d", mask->dimx, mask->dimy, dim[0], dim[1]);
316  return ES_ERROR;
317  }
318 
319  if (mask->values == NULL || mask->nodata == NULL) {
320  rterror("rt_pixel_set_array: Invalid mask");
321  return ES_ERROR;
322  }
323 
324  }
325 
326  /* establish 2D arrays (Y axis) */
327  values = rtalloc(sizeof(double *) * dim[1]);
328  nodatas = rtalloc(sizeof(int *) * dim[1]);
329 
330  if (values == NULL || nodatas == NULL) {
331  rterror("rt_pixel_set_to_array: Could not allocate memory for 2D array");
332  return ES_ERROR;
333  }
334 
335  /* initialize X axis */
336  for (i = 0; i < dim[1]; i++) {
337  values[i] = rtalloc(sizeof(double) * dim[0]);
338  nodatas[i] = rtalloc(sizeof(int) * dim[0]);
339 
340  if (values[i] == NULL || nodatas[i] == NULL) {
341  rterror("rt_pixel_set_to_array: Could not allocate memory for dimension of 2D array");
342 
343  if (values[i] == NULL) {
344  for (j = 0; j < i; j++) {
345  rtdealloc(values[j]);
346  rtdealloc(nodatas[j]);
347  }
348  }
349  else {
350  for (j = 0; j <= i; j++) {
351  rtdealloc(values[j]);
352  if (j < i)
353  rtdealloc(nodatas[j]);
354  }
355  }
356 
357  rtdealloc(values);
358  rtdealloc(nodatas);
359 
360  return ES_ERROR;
361  }
362 
363  /* set values to 0 */
364  memset(values[i], 0, sizeof(double) * dim[0]);
365 
366  /* set nodatas to 1 */
367  for (j = 0; j < dim[0]; j++)
368  nodatas[i][j] = 1;
369  }
370 
371  /* get 0,0 of grid */
372  zero[0] = x - distancex;
373  zero[1] = y - distancey;
374 
375  /* populate 2D arrays */
376  for (i = 0; i < count; i++) {
377  if (npixel[i].nodata)
378  continue;
379 
380  _x = npixel[i].x - zero[0];
381  _y = npixel[i].y - zero[1];
382 
383  RASTER_DEBUGF(4, "absolute x,y: %d x %d", npixel[i].x, npixel[i].y);
384  RASTER_DEBUGF(4, "relative x,y: %d x %d", _x, _y);
385 
386  /* no mask */
387  if (mask == NULL) {
388  values[_y][_x] = npixel[i].value;
389  nodatas[_y][_x] = 0;
390  }
391  /* mask */
392  else {
393  /* unweighted (boolean) mask */
394  if (mask->weighted == 0) {
395  /* pixel is set to zero or nodata */
396  if (FLT_EQ(mask->values[_y][_x],0) || mask->nodata[_y][_x] == 1) {
397  values[_y][_x] = 0;
398  nodatas[_y][_x] = 1;
399  }
400  /* use pixel */
401  else {
402  values[_y][_x] = npixel[i].value;
403  nodatas[_y][_x] = 0;
404  }
405  }
406  /* weighted mask */
407  else {
408  /* nodata */
409  if(mask->nodata[_y][_x] == 1) {
410  values[_y][_x] = 0;
411  nodatas[_y][_x] = 1;
412  }
413  /* apply weight to pixel value */
414  else {
415  values[_y][_x] = npixel[i].value * mask->values[_y][_x];
416  nodatas[_y][_x] = 0;
417  }
418  }
419  }
420 
421  RASTER_DEBUGF(4, "(x, y, nodata, value) = (%d, %d, %d, %f)", _x, _y, nodatas[_y][_x], values[_y][_x]);
422  }
423 
424  *value = &(*values);
425  *nodata = &(*nodatas);
426  if (dimx != NULL)
427  *dimx = dim[0];
428  if (dimy != NULL)
429  *dimy = dim[1];
430 
431  return ES_NONE;
432 }
void rterror(const char *fmt,...)
Wrappers used for reporting errors and info.
Definition: rt_context.c:199
void * rtalloc(size_t size)
Wrappers used for managing memory.
Definition: rt_context.c:171
#define RASTER_DEBUGF(level, msg,...)
Definition: librtcore.h:299
int8_t rt_util_clamp_to_8BSI(double value)
Definition: rt_util.c:49
uint8_t rt_util_clamp_to_1BB(double value)
Definition: rt_util.c:34
int32_t rt_util_clamp_to_32BSI(double value)
Definition: rt_util.c:69
rt_pixtype
Definition: librtcore.h:185
@ PT_32BUI
Definition: librtcore.h:194
@ PT_2BUI
Definition: librtcore.h:187
@ PT_32BSI
Definition: librtcore.h:193
@ PT_END
Definition: librtcore.h:197
@ PT_4BUI
Definition: librtcore.h:188
@ PT_32BF
Definition: librtcore.h:195
@ PT_1BB
Definition: librtcore.h:186
@ PT_16BUI
Definition: librtcore.h:192
@ PT_8BSI
Definition: librtcore.h:189
@ PT_16BSI
Definition: librtcore.h:191
@ PT_64BF
Definition: librtcore.h:196
@ PT_8BUI
Definition: librtcore.h:190
#define FLT_EQ(x, y)
Definition: librtcore.h:2234
uint8_t rt_util_clamp_to_2BUI(double value)
Definition: rt_util.c:39
uint8_t rt_util_clamp_to_8BUI(double value)
Definition: rt_util.c:54
rt_errorstate
Enum definitions.
Definition: librtcore.h:179
@ ES_NONE
Definition: librtcore.h:180
@ ES_ERROR
Definition: librtcore.h:181
int16_t rt_util_clamp_to_16BSI(double value)
Definition: rt_util.c:59
uint8_t rt_util_clamp_to_4BUI(double value)
Definition: rt_util.c:44
void rtdealloc(void *mem)
Definition: rt_context.c:186
uint16_t rt_util_clamp_to_16BUI(double value)
Definition: rt_util.c:64
uint32_t rt_util_clamp_to_32BUI(double value)
Definition: rt_util.c:74
float rt_util_clamp_to_32F(double value)
Definition: rt_util.c:79
This library is the generic raster handling section of PostGIS.
int value
Definition: genraster.py:61
int count
Definition: genraster.py:56
rt_pixtype rt_pixtype_index_from_name(const char *pixname)
Definition: rt_pixel.c:80
double rt_pixtype_get_min_value(rt_pixtype pixtype)
Return minimum value possible for pixel type.
Definition: rt_pixel.c:148
int rt_pixtype_alignment(rt_pixtype pixtype)
Return alignment requirements for data in the given pixel type.
Definition: rt_pixel.c:75
rt_errorstate rt_pixtype_compare_clamped_values(rt_pixtype pixtype, double val, double refval, int *isequal)
Returns 1 if clamped values are equal, 0 if not equal, -1 if error.
Definition: rt_pixel.c:200
rt_errorstate rt_pixel_set_to_array(rt_pixel npixel, uint32_t 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
const char * rt_pixtype_name(rt_pixtype pixtype)
Definition: rt_pixel.c:110
int rt_pixtype_size(rt_pixtype pixtype)
Return size in bytes of a value in the given pixtype.
Definition: rt_pixel.c:39
double ** values
Definition: librtcore.h:2346
uint16_t dimy
Definition: librtcore.h:2345
uint16_t dimx
Definition: librtcore.h:2344
int ** nodata
Definition: librtcore.h:2347
int weighted
Definition: librtcore.h:2348
double value
Definition: librtcore.h:2338
unsigned int uint32_t
Definition: uthash.h:78