GRASS GIS 7 Programmer's Manual  7.5.svn(2017)-r71933
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
gdal.c
Go to the documentation of this file.
1 /*!
2  \file lib/raster/gdal.c
3 
4  \brief Raster Library - Utilization of GDAL library.
5 
6  (C) 2010 by the GRASS Development Team
7 
8  This program is free software under the GNU General Public License
9  (>=v2). Read the file COPYING that comes with GRASS for details.
10 
11  \author Glynn Clements
12 */
13 
14 #include <stdlib.h>
15 #include <string.h>
16 #include <unistd.h>
17 
18 #include <grass/config.h>
19 #include <grass/gis.h>
20 #include <grass/raster.h>
21 #include <grass/glocale.h>
22 
23 #include "R.h"
24 
25 #ifndef HAVE_GDAL
26 #undef GDAL_LINK
27 #endif
28 
29 #ifdef GDAL_LINK
30 
31 #ifdef GDAL_DYNAMIC
32 # if defined(__unix) || defined(__unix__)
33 # include <dlfcn.h>
34 # endif
35 # ifdef _WIN32
36 # include <windows.h>
37 # endif
38 #endif
39 
40 static void CPL_STDCALL(*pGDALAllRegister) (void);
41 static void CPL_STDCALL(*pGDALClose) (GDALDatasetH);
42 static GDALRasterBandH CPL_STDCALL(*pGDALGetRasterBand) (GDALDatasetH, int);
43 static GDALDatasetH CPL_STDCALL(*pGDALOpen) (const char *pszFilename,
44  GDALAccess eAccess);
45 static CPLErr CPL_STDCALL(*pGDALRasterIO) (GDALRasterBandH hRBand,
46  GDALRWFlag eRWFlag, int nDSXOff,
47  int nDSYOff, int nDSXSize,
48  int nDSYSize, void *pBuffer,
49  int nBXSize, int nBYSize,
50  GDALDataType eBDataType,
51  int nPixelSpace, int nLineSpace);
52 static GDALDriverH CPL_STDCALL(*pGDALGetDriverByName) (const char *);
53 static const char *CPL_STDCALL(*pGDALGetMetadataItem) (GDALMajorObjectH,
54  const char *,
55  const char *);
56 static GDALDatasetH CPL_STDCALL(*pGDALCreate) (GDALDriverH hDriver,
57  const char *, int, int, int,
58  GDALDataType, char **);
59 static GDALDatasetH CPL_STDCALL(*pGDALCreateCopy) (GDALDriverH, const char *,
60  GDALDatasetH, int, char **,
61  GDALProgressFunc, void *);
62 static CPLErr CPL_STDCALL(*pGDALSetRasterNoDataValue) (GDALRasterBandH,
63  double);
64 static CPLErr CPL_STDCALL(*pGDALSetGeoTransform) (GDALDatasetH, double *);
65 static CPLErr CPL_STDCALL(*pGDALSetProjection) (GDALDatasetH, const char *);
66 static const char *CPL_STDCALL(*pGDALGetDriverShortName) (GDALDriverH);
67 static GDALDriverH CPL_STDCALL(*pGDALGetDatasetDriver) (GDALDatasetH);
68 
69 #if GDAL_DYNAMIC
70 # if defined(__unix) && !defined(__unix__)
71 # define __unix__ __unix
72 # endif
73 
74 static void *library_h;
75 
76 static void *get_symbol(const char *name)
77 {
78  void *sym;
79 
80 # ifdef __unix__
81  sym = dlsym(library_h, name);
82 # endif
83 # ifdef _WIN32
84  sym = GetProcAddress((HINSTANCE) library_h, name);
85 # endif
86 
87  if (!sym)
88  G_fatal_error(_("Unable to locate symbol <%s>"), name);
89 
90  return sym;
91 }
92 
93 static void try_load_library(const char *name)
94 {
95 # ifdef __unix__
96  library_h = dlopen(name, RTLD_NOW);
97 # endif
98 # ifdef _WIN32
99  library_h = LoadLibrary(name);
100 # endif
101 }
102 
103 static void load_library(void)
104 {
105  static const char *const candidates[] = {
106 # ifdef __unix__
107  "libgdal.so.20",
108  "libgdal.so.1",
109  "libgdal.1.1.so",
110  "gdal.1.0.so",
111  "gdal.so.1.0",
112  "libgdal.so",
113  "libgdal1.6.0.so",
114  "libgdal1.7.0.so",
115 # endif
116 # ifdef _WIN32
117  "gdal202.dll",
118  "gdal201.dll",
119  "gdal200.dll",
120  "gdal111.dll",
121  "gdal110.dll",
122  "gdal19.dll",
123  "gdal18.dll",
124  "gdal17.dll",
125  "gdal16.dll",
126  "gdal15.dll",
127  "gdal11.dll",
128  "gdal.1.0.dll",
129  "libgdal-1.dll",
130  "gdal.dll",
131 # endif
132  NULL
133  };
134  int i;
135 
136  for (i = 0; candidates[i]; i++) {
137  try_load_library(candidates[i]);
138  if (library_h) {
139  G_debug(3, "found %s", candidates[i]);
140  return;
141  }
142  }
143 
144  G_fatal_error(_("Unable to load GDAL library"));
145 }
146 
147 static void init_gdal(void)
148 {
149  load_library();
150 
151 # if defined(_WIN32) && !defined(_WIN64)
152  pGDALAllRegister = get_symbol("_GDALAllRegister@0");
153  pGDALOpen = get_symbol("_GDALOpen@8");
154  pGDALClose = get_symbol("_GDALClose@4");
155  pGDALGetRasterBand = get_symbol("_GDALGetRasterBand@8");
156  pGDALRasterIO = get_symbol("_GDALRasterIO@48");
157  pGDALGetDriverByName = get_symbol("_GDALGetDriverByName@4");
158  pGDALGetMetadataItem = get_symbol("_GDALGetMetadataItem@12");
159  pGDALCreate = get_symbol("_GDALCreate@28");
160  pGDALCreateCopy = get_symbol("_GDALCreateCopy@28");
161  pGDALSetRasterNoDataValue = get_symbol("_GDALSetRasterNoDataValue@12");
162  pGDALSetGeoTransform = get_symbol("_GDALSetGeoTransform@8");
163  pGDALSetProjection = get_symbol("_GDALSetProjection@8");
164  pGDALGetDriverShortName = get_symbol("_GDALGetDriverShortName@4");
165  pGDALGetDatasetDriver = get_symbol("_GDALGetDatasetDriver@4");
166 #else
167  pGDALAllRegister = get_symbol("GDALAllRegister");
168  pGDALOpen = get_symbol("GDALOpen");
169  pGDALClose = get_symbol("GDALClose");
170  pGDALGetRasterBand = get_symbol("GDALGetRasterBand");
171  pGDALRasterIO = get_symbol("GDALRasterIO");
172  pGDALGetDriverByName = get_symbol("GDALGetDriverByName");
173  pGDALGetMetadataItem = get_symbol("GDALGetMetadataItem");
174  pGDALCreate = get_symbol("GDALCreate");
175  pGDALCreateCopy = get_symbol("GDALCreateCopy");
176  pGDALSetRasterNoDataValue = get_symbol("GDALSetRasterNoDataValue");
177  pGDALSetGeoTransform = get_symbol("GDALSetGeoTransform");
178  pGDALSetProjection = get_symbol("GDALSetProjection");
179  pGDALGetDriverShortName = get_symbol("GDALGetDriverShortName");
180  pGDALGetDatasetDriver = get_symbol("GDALGetDatasetDriver");
181 #endif
182 }
183 
184 #else /* GDAL_DYNAMIC */
185 
186 static void init_gdal(void)
187 {
188  pGDALAllRegister = &GDALAllRegister;
189  pGDALOpen = &GDALOpen;
190  pGDALClose = &GDALClose;
191  pGDALGetRasterBand = &GDALGetRasterBand;
192  pGDALRasterIO = &GDALRasterIO;
193  pGDALGetDriverByName = &GDALGetDriverByName;
194  pGDALGetMetadataItem = &GDALGetMetadataItem;
195  pGDALCreate = &GDALCreate;
196  pGDALCreateCopy = &GDALCreateCopy;
197  pGDALSetRasterNoDataValue = &GDALSetRasterNoDataValue;
198  pGDALSetGeoTransform = &GDALSetGeoTransform;
199  pGDALSetProjection = &GDALSetProjection;
200  pGDALGetDriverShortName = &GDALGetDriverShortName;
201  pGDALGetDatasetDriver = &GDALGetDatasetDriver;
202 }
203 
204 #endif /* GDAL_DYNAMIC */
205 
206 #endif /* GDAL_LINK */
207 
208 /*!
209  \brief Initialization
210 
211  Register all GDAL drivers.
212 */
213 void Rast_init_gdal(void)
214 {
215 #ifdef GDAL_LINK
216  static int initialized;
217 
218  if (G_is_initialized(&initialized))
219  return;
220 
221  init_gdal();
222  (*pGDALAllRegister) ();
223  G_initialize_done(&initialized);
224 #endif
225 }
226 
227 /*!
228  \brief Get GDAL link settings for given raster map
229 
230  \param name map name
231  \param mapset name of mapset
232 
233  \return pointer to GDAL_link structure
234  \return NULL if link not found
235 */
236 struct GDAL_link *Rast_get_gdal_link(const char *name, const char *mapset)
237 {
238 #ifdef GDAL_LINK
239  GDALDatasetH data;
240  GDALRasterBandH band;
241  GDALDataType type;
242  RASTER_MAP_TYPE req_type;
243 #endif
244  const char *filename;
245  int band_num;
246  struct GDAL_link *gdal;
247  RASTER_MAP_TYPE map_type;
248  FILE *fp;
249  struct Key_Value *key_val;
250  const char *p;
251  DCELL null_val;
252  int hflip, vflip;
253 
254  if (!G_find_raster2(name, mapset))
255  return NULL;
256 
257  map_type = Rast_map_type(name, mapset);
258  if (map_type < 0)
259  return NULL;
260 
261  fp = G_fopen_old_misc("cell_misc", "gdal", name, mapset);
262  if (!fp)
263  return NULL;
264  key_val = G_fread_key_value(fp);
265  fclose(fp);
266 
267  if (!key_val)
268  return NULL;
269 
270  filename = G_find_key_value("file", key_val);
271  if (!filename)
272  return NULL;
273 
274  p = G_find_key_value("band", key_val);
275  if (!p)
276  return NULL;
277  band_num = atoi(p);
278  if (!band_num)
279  return NULL;
280 
281  p = G_find_key_value("null", key_val);
282  if (!p)
283  return NULL;
284  if (strcmp(p, "none") == 0)
285  Rast_set_d_null_value(&null_val, 1);
286  else
287  null_val = atof(p);
288 
289  hflip = G_find_key_value("hflip", key_val) ? 1 : 0;
290  vflip = G_find_key_value("vflip", key_val) ? 1 : 0;
291 
292 #ifdef GDAL_LINK
293  p = G_find_key_value("type", key_val);
294  if (!p)
295  return NULL;
296  type = atoi(p);
297 
298  switch (type) {
299  case GDT_Byte:
300  case GDT_Int16:
301  case GDT_UInt16:
302  case GDT_Int32:
303  case GDT_UInt32:
304  req_type = CELL_TYPE;
305  break;
306  case GDT_Float32:
307  req_type = FCELL_TYPE;
308  break;
309  case GDT_Float64:
310  req_type = DCELL_TYPE;
311  break;
312  default:
313  return NULL;
314  }
315 
316  if (req_type != map_type)
317  return NULL;
318 
319  Rast_init_gdal();
320 
321  data = (*pGDALOpen) (filename, GA_ReadOnly);
322  if (!data)
323  return NULL;
324 
325  band = (*pGDALGetRasterBand) (data, band_num);
326  if (!band) {
327  (*pGDALClose) (data);
328  return NULL;
329  }
330 #endif
331 
332  gdal = G_calloc(1, sizeof(struct GDAL_link));
333 
334  gdal->filename = G_store(filename);
335  gdal->band_num = band_num;
336  gdal->null_val = null_val;
337  gdal->hflip = hflip;
338  gdal->vflip = vflip;
339 #ifdef GDAL_LINK
340  gdal->data = data;
341  gdal->band = band;
342  gdal->type = type;
343 #endif
344 
345  return gdal;
346 }
347 
348 struct GDAL_Options
349 {
350  const char *dir;
351  const char *ext;
352  const char *format;
353  char **options;
354 };
355 
356 static struct state
357 {
358  int initialized;
359  struct GDAL_Options opts;
360  struct Key_Value *projinfo, *projunits;
361  char *srswkt;
362 } state;
363 
364 static struct state *st = &state;
365 
366 static void read_gdal_options(void)
367 {
368  FILE *fp;
369  struct Key_Value *key_val;
370  const char *p;
371 
372  fp = G_fopen_old("", "GDAL", G_mapset());
373  if (!fp)
374  G_fatal_error(_("Unable to open GDAL file"));
375  key_val = G_fread_key_value(fp);
376  fclose(fp);
377 
378  p = G_find_key_value("directory", key_val);
379  if (!p)
380  p = "gdal";
381  if (*p == '/') {
382  st->opts.dir = G_store(p);
383  }
384  else {
385  char path[GPATH_MAX];
386 
387  G_file_name(path, p, "", G_mapset());
388  st->opts.dir = G_store(path);
389  if (access(path, 0) != 0)
391  }
392 
393  p = G_find_key_value("extension", key_val);
394  st->opts.ext = G_store(p ? p : "");
395 
396  p = G_find_key_value("format", key_val);
397  st->opts.format = G_store(p ? p : "GTiff");
398 
399  p = G_find_key_value("options", key_val);
400  st->opts.options = p ? G_tokenize(p, ",") : NULL;
401 
402  G_free_key_value(key_val);
403 }
404 
405 /*!
406  \brief Create GDAL settings for given raster map
407 
408  \param name map name
409  \param map_type map type (CELL, FCELL, DCELL)
410 
411  \return pointer to allocated GDAL_link structure
412  \return NULL on error
413 */
414 struct GDAL_link *Rast_create_gdal_link(const char *name,
415  RASTER_MAP_TYPE map_type)
416 {
417 #ifdef GDAL_LINK
418  char path[GPATH_MAX];
419  GDALDriverH driver;
420  double transform[6];
421  struct GDAL_link *gdal;
422  FILE *fp;
423  struct Key_Value *key_val;
424  char buf[32];
425 
427 
428  Rast_init_gdal();
429 
430  if (!G_is_initialized(&st->initialized)) {
431  read_gdal_options();
432  st->projinfo = G_get_projinfo();
433  st->projunits = G_get_projunits();
434 #if 0
435  /* We cannot use GPJ_grass_to_wkt() here because that would create a
436  circular dependency between libgis and libgproj */
437  if (st->projinfo && st->projunits)
438  st->srswkt = GPJ_grass_to_wkt(st->projinfo, st->projunits);
439 #endif
440  G_initialize_done(&st->initialized);
441  }
442 
443  gdal = G_calloc(1, sizeof(struct GDAL_link));
444 
445  sprintf(path, "%s/%s%s", st->opts.dir, name, st->opts.ext);
446  gdal->filename = G_store(path);
447  gdal->band_num = 1;
448  gdal->hflip = 0;
449  gdal->vflip = 0;
450 
451  switch (map_type) {
452  case CELL_TYPE:
453  switch (R__.nbytes) {
454  case 1:
455  gdal->type = GDT_Byte;
456  gdal->null_val = (DCELL) 0xFF;
457  break;
458  case 2:
459  gdal->type = GDT_UInt16;
460  gdal->null_val = (DCELL) 0xFFFF;
461  break;
462  case 3:
463  case 4:
464  gdal->type = GDT_Int32;
465  gdal->null_val = (DCELL) 0x80000000U;
466  break;
467  }
468  break;
469  case FCELL_TYPE:
470  gdal->type = GDT_Float32;
471  Rast_set_d_null_value(&gdal->null_val, 1);
472  break;
473  case DCELL_TYPE:
474  gdal->type = GDT_Float64;
475  Rast_set_d_null_value(&gdal->null_val, 1);
476  break;
477  default:
478  G_fatal_error(_("Invalid map type <%d>"), map_type);
479  break;
480  }
481 
482  driver = (*pGDALGetDriverByName) (st->opts.format);
483  if (!driver)
484  G_fatal_error(_("Unable to get <%s> driver"), st->opts.format);
485 
486  /* Does driver support GDALCreate ? */
487  if ((*pGDALGetMetadataItem) (driver, GDAL_DCAP_CREATE, NULL)) {
488  gdal->data =
489  (*pGDALCreate)(driver, gdal->filename,
491  1, gdal->type, st->opts.options);
492  if (!gdal->data)
493  G_fatal_error(_("Unable to create <%s> dataset using <%s> driver"),
494  name, st->opts.format);
495  }
496  /* If not - create MEM driver for intermediate dataset.
497  * Check if raster can be created at all (with GDALCreateCopy) */
498  else if ((*pGDALGetMetadataItem) (driver, GDAL_DCAP_CREATECOPY, NULL)) {
499  GDALDriverH mem_driver;
500 
501  G_message(_("Driver <%s> does not support direct writing. "
502  "Using MEM driver for intermediate dataset."),
503  st->opts.format);
504 
505  mem_driver = (*pGDALGetDriverByName) ("MEM");
506  if (!mem_driver)
507  G_fatal_error(_("Unable to get in-memory raster driver"));
508 
509  gdal->data =
510  (*pGDALCreate)(mem_driver, "",
512  1, gdal->type, st->opts.options);
513  if (!gdal->data)
514  G_fatal_error(_("Unable to create <%s> dataset using memory driver"),
515  name);
516  }
517  else
518  G_fatal_error(_("Driver <%s> does not support creating rasters"),
519  st->opts.format);
520 
521  gdal->band = (*pGDALGetRasterBand) (gdal->data, gdal->band_num);
522 
523  (*pGDALSetRasterNoDataValue) (gdal->band, gdal->null_val);
524 
525  /* Set Geo Transform */
526  transform[0] = R__.wr_window.west;
527  transform[1] = R__.wr_window.ew_res;
528  transform[2] = 0.0;
529  transform[3] = R__.wr_window.north;
530  transform[4] = 0.0;
531  transform[5] = -R__.wr_window.ns_res;
532 
533  if ((*pGDALSetGeoTransform) (gdal->data, transform) >= CE_Failure)
534  G_warning(_("Unable to set geo transform"));
535 
536  if (st->srswkt)
537  if ((*pGDALSetProjection) (gdal->data, st->srswkt) == CE_Failure)
538  G_warning(_("Unable to set projection"));
539 
540  fp = G_fopen_new_misc("cell_misc", "gdal", name);
541  if (!fp)
542  G_fatal_error(_("Unable to create cell_misc/%s/gdal file"), name);
543 
544  key_val = G_create_key_value();
545 
546  G_set_key_value("file", gdal->filename, key_val);
547 
548  sprintf(buf, "%d", gdal->band_num);
549  G_set_key_value("band", buf, key_val);
550 
551  sprintf(buf, "%.22g", gdal->null_val);
552  G_set_key_value("null", buf, key_val);
553 
554  sprintf(buf, "%d", gdal->type);
555  G_set_key_value("type", buf, key_val);
556 
557  if (G_fwrite_key_value(fp, key_val) < 0)
558  G_fatal_error(_("Error writing cell_misc/%s/gdal file"), name);
559 
560  G_free_key_value(key_val);
561 
562  fclose(fp);
563 
564  return gdal;
565 #else
566  return NULL;
567 #endif
568 }
569 
570 /*!
571  \brief Close existing GDAL link
572 
573  \param gdal pointer to GDAL_link to be closed
574 */
575 void Rast_close_gdal_link(struct GDAL_link *gdal)
576 {
577 #ifdef GDAL_LINK
578  (*pGDALClose) (gdal->data);
579 #endif
580  G_free(gdal->filename);
581  G_free(gdal);
582 }
583 
584 /*!
585  \brief Close existing GDAL link and write out data
586 
587  \param gdal pointer to GDAL_link to be closed
588 
589  \return 1 on success
590  \return -1 on failure
591 */
593 {
594  int stat = 1;
595 
596 #ifdef GDAL_LINK
597  GDALDriverH src_drv = (*pGDALGetDatasetDriver) (gdal->data);
598 
599  if (G_strcasecmp((*pGDALGetDriverShortName) (src_drv), "MEM") == 0) {
600  GDALDriverH dst_drv = (*pGDALGetDriverByName) (st->opts.format);
601  GDALDatasetH dst =
602  (*pGDALCreateCopy) (dst_drv, gdal->filename, gdal->data, FALSE,
603  st->opts.options, NULL, NULL);
604 
605  if (!dst) {
606  G_warning(_("Unable to create output file <%s> using driver <%s>"),
607  gdal->filename, st->opts.format);
608  stat = -1;
609  }
610  (*pGDALClose) (dst);
611  }
612 
613  (*pGDALClose) (gdal->data);
614 
615 #endif
616  G_free(gdal->filename);
617  G_free(gdal);
618 
619  return stat;
620 }
621 
622 #ifdef GDAL_LINK
623 /*!
624  \brief Input/output function for GDAL links
625 
626  See GDAL's RasterIO for details.
627 */
628 CPLErr Rast_gdal_raster_IO(GDALRasterBandH band, GDALRWFlag rw_flag,
629  int x_off, int y_off, int x_size, int y_size,
630  void *buffer, int buf_x_size, int buf_y_size,
631  GDALDataType buf_type, int pixel_size,
632  int line_size)
633 {
634  return (*pGDALRasterIO) (band, rw_flag, x_off, y_off, x_size, y_size,
635  buffer, buf_x_size, buf_y_size, buf_type,
636  pixel_size, line_size);
637 }
638 #endif
int G_make_mapset_element(const char *p_element)
Create element in the current mapset.
Definition: mapset_msc.c:38
#define CELL_TYPE
Definition: raster.h:11
int G_strcasecmp(const char *x, const char *y)
String compare ignoring case (upper or lower)
Definition: strings.c:46
const char * G_find_key_value(const char *key, const struct Key_Value *kv)
Find given key (case sensitive)
Definition: key_value1.c:84
struct driver * driver
Definition: driver/init.c:25
void G_free(void *buf)
Free allocated memory.
Definition: gis/alloc.c:149
Definition: R.h:72
struct Key_Value * G_get_projunits(void)
Gets units information for location.
Definition: get_projinfo.c:29
void Rast__init_window(void)
struct Key_Value * G_fread_key_value(FILE *fd)
Read key/values pairs from file.
Definition: key_value2.c:49
int Rast_close_gdal_write_link(struct GDAL_link *gdal)
Close existing GDAL link and write out data.
Definition: gdal.c:592
double west
Extent coordinates (west)
Definition: gis.h:442
double DCELL
Definition: gis.h:581
char * GPJ_grass_to_wkt(const struct Key_Value *proj_info, const struct Key_Value *proj_units, int esri_style, int prettify)
Converts a GRASS co-ordinate system representation to WKT style.
Definition: convert.c:87
char ** G_tokenize(const char *buf, const char *delim)
Tokenize string.
Definition: gis/token.c:48
char * G_store(const char *s)
Copy string to allocated memory.
Definition: strings.c:86
char * dst
Definition: lz4.h:354
FILE * G_fopen_old_misc(const char *dir, const char *element, const char *name, const char *mapset)
open a database file for reading
Definition: open_misc.c:210
int G_is_initialized(int *p)
Definition: counter.c:59
int G_fwrite_key_value(FILE *fd, const struct Key_Value *kv)
Write key/value pairs to file.
Definition: key_value2.c:25
#define NULL
Definition: ccmath.h:32
void G_initialize_done(int *p)
Definition: counter.c:76
void G_free_key_value(struct Key_Value *kv)
Free allocated Key_Value structure.
Definition: key_value1.c:103
void G_fatal_error(const char *msg,...)
Print a fatal error message to stderr.
Definition: gis/error.c:159
struct Cell_head wr_window
Definition: R.h:84
struct state * st
Definition: parser.c:103
char * G_file_name(char *path, const char *element, const char *name, const char *mapset)
Builds full path names to GIS data files.
Definition: file_name.c:38
#define DCELL_TYPE
Definition: raster.h:13
RASTER_MAP_TYPE Rast_map_type(const char *name, const char *mapset)
Determine raster data type.
Definition: raster/open.c:869
double north
Extent coordinates (north)
Definition: gis.h:436
struct GDAL_link * Rast_create_gdal_link(const char *name, RASTER_MAP_TYPE map_type)
Create GDAL settings for given raster map.
Definition: gdal.c:414
#define FALSE
Definition: gis.h:53
int G_debug(int level, const char *msg,...)
Print debugging message.
Definition: debug.c:65
fclose(fd)
#define GPATH_MAX
Definition: gis.h:151
struct Key_Value * G_get_projinfo(void)
Gets projection information for location.
Definition: get_projinfo.c:58
Definition: gis.h:479
void Rast_set_d_null_value(DCELL *dcellVals, int numVals)
To set a number of DCELL raster values to NULL.
Definition: null_val.c:155
int cols
Number of columns for 2D data.
Definition: gis.h:409
double ns_res
Resolution - north to south cell size for 2D data.
Definition: gis.h:430
void Rast_init_gdal(void)
Initialization.
Definition: gdal.c:213
Definition: path.h:16
FILE * G_fopen_new_misc(const char *dir, const char *element, const char *name)
open a new database file
Definition: open_misc.c:182
#define _(str)
Definition: glocale.h:13
void G_message(const char *msg,...)
Print a message to stderr.
Definition: gis/error.c:89
int RASTER_MAP_TYPE
Definition: raster.h:25
const char * G_find_raster2(const char *name, const char *mapset)
Find a raster map (look but don&#39;t touch)
Definition: find_rast.c:76
#define FCELL_TYPE
Definition: raster.h:12
void Rast_close_gdal_link(struct GDAL_link *gdal)
Close existing GDAL link.
Definition: gdal.c:575
const char * G_mapset(void)
Get current mapset name.
Definition: gis/mapset.c:33
int
Reads the categories file for map name in mapset and stores the categories in the pcats structure...
void G_set_key_value(const char *key, const char *value, struct Key_Value *kv)
Set value for given key.
Definition: key_value1.c:38
int nbytes
Definition: R.h:78
FILE * G_fopen_old(const char *element, const char *name, const char *mapset)
Open a database file for reading.
Definition: gis/open.c:253
const char * name
Definition: named_colr.c:7
double ew_res
Resolution - east to west cell size for 2D data.
Definition: gis.h:426
int rows
Number of rows for 2D data.
Definition: gis.h:405
struct state state
Definition: parser.c:102
void G_warning(const char *msg,...)
Print a warning message to stderr.
Definition: gis/error.c:203
struct GDAL_link * Rast_get_gdal_link(const char *name, const char *mapset)
Get GDAL link settings for given raster map.
Definition: gdal.c:236
struct Key_Value * G_create_key_value(void)
Allocate and initialize Key_Value structure.
Definition: key_value1.c:23