GRASS GIS 7 Programmer's Manual  7.7.svn(2018)-r73679
 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/gprojects.h>
22 #include <grass/glocale.h>
23 
24 #include "R.h"
25 
26 #ifndef HAVE_GDAL
27 #undef GDAL_LINK
28 #endif
29 
30 #ifdef GDAL_LINK
31 
32 #ifdef GDAL_DYNAMIC
33 # if defined(__unix) || defined(__unix__)
34 # include <dlfcn.h>
35 # endif
36 # ifdef _WIN32
37 # include <windows.h>
38 # endif
39 #endif
40 
41 static void CPL_STDCALL(*pGDALAllRegister) (void);
42 static void CPL_STDCALL(*pGDALClose) (GDALDatasetH);
43 static GDALRasterBandH CPL_STDCALL(*pGDALGetRasterBand) (GDALDatasetH, int);
44 static GDALDatasetH CPL_STDCALL(*pGDALOpen) (const char *pszFilename,
45  GDALAccess eAccess);
46 static CPLErr CPL_STDCALL(*pGDALRasterIO) (GDALRasterBandH hRBand,
47  GDALRWFlag eRWFlag, int nDSXOff,
48  int nDSYOff, int nDSXSize,
49  int nDSYSize, void *pBuffer,
50  int nBXSize, int nBYSize,
51  GDALDataType eBDataType,
52  int nPixelSpace, int nLineSpace);
53 static GDALDriverH CPL_STDCALL(*pGDALGetDriverByName) (const char *);
54 static const char *CPL_STDCALL(*pGDALGetMetadataItem) (GDALMajorObjectH,
55  const char *,
56  const char *);
57 static GDALDatasetH CPL_STDCALL(*pGDALCreate) (GDALDriverH hDriver,
58  const char *, int, int, int,
59  GDALDataType, char **);
60 static GDALDatasetH CPL_STDCALL(*pGDALCreateCopy) (GDALDriverH, const char *,
61  GDALDatasetH, int, char **,
62  GDALProgressFunc, void *);
63 static CPLErr CPL_STDCALL(*pGDALSetRasterNoDataValue) (GDALRasterBandH,
64  double);
65 static CPLErr CPL_STDCALL(*pGDALSetGeoTransform) (GDALDatasetH, double *);
66 static CPLErr CPL_STDCALL(*pGDALSetProjection) (GDALDatasetH, const char *);
67 static const char *CPL_STDCALL(*pGDALGetDriverShortName) (GDALDriverH);
68 static GDALDriverH CPL_STDCALL(*pGDALGetDatasetDriver) (GDALDatasetH);
69 
70 #if GDAL_DYNAMIC
71 # if defined(__unix) && !defined(__unix__)
72 # define __unix__ __unix
73 # endif
74 
75 static void *library_h;
76 
77 static void *get_symbol(const char *name)
78 {
79  void *sym;
80 
81 # ifdef __unix__
82  sym = dlsym(library_h, name);
83 # endif
84 # ifdef _WIN32
85  sym = GetProcAddress((HINSTANCE) library_h, name);
86 # endif
87 
88  if (!sym)
89  G_fatal_error(_("Unable to locate symbol <%s>"), name);
90 
91  return sym;
92 }
93 
94 static void try_load_library(const char *name)
95 {
96 # ifdef __unix__
97  library_h = dlopen(name, RTLD_NOW);
98 # endif
99 # ifdef _WIN32
100  library_h = LoadLibrary(name);
101 # endif
102 }
103 
104 static void load_library(void)
105 {
106  static const char *const candidates[] = {
107 # ifdef __unix__
108  "libgdal.so.20",
109  "libgdal.so.1",
110  "libgdal.1.1.so",
111  "gdal.1.0.so",
112  "gdal.so.1.0",
113  "libgdal.so",
114  "libgdal1.6.0.so",
115  "libgdal1.7.0.so",
116 # endif
117 # ifdef _WIN32
118  "gdal203.dll",
119  "gdal202.dll",
120  "gdal201.dll",
121  "gdal200.dll",
122  "gdal111.dll",
123  "gdal110.dll",
124  "gdal19.dll",
125  "gdal18.dll",
126  "gdal17.dll",
127  "gdal16.dll",
128  "gdal15.dll",
129  "gdal11.dll",
130  "gdal.1.0.dll",
131  "libgdal-1.dll",
132  "gdal.dll",
133 # endif
134  NULL
135  };
136  int i;
137 
138  for (i = 0; candidates[i]; i++) {
139  try_load_library(candidates[i]);
140  if (library_h) {
141  G_debug(3, "found %s", candidates[i]);
142  return;
143  }
144  }
145 
146  G_fatal_error(_("Unable to load GDAL library"));
147 }
148 
149 static void init_gdal(void)
150 {
151  load_library();
152 
153 # if defined(_WIN32) && !defined(_WIN64)
154  pGDALAllRegister = get_symbol("_GDALAllRegister@0");
155  pGDALOpen = get_symbol("_GDALOpen@8");
156  pGDALClose = get_symbol("_GDALClose@4");
157  pGDALGetRasterBand = get_symbol("_GDALGetRasterBand@8");
158  pGDALRasterIO = get_symbol("_GDALRasterIO@48");
159  pGDALGetDriverByName = get_symbol("_GDALGetDriverByName@4");
160  pGDALGetMetadataItem = get_symbol("_GDALGetMetadataItem@12");
161  pGDALCreate = get_symbol("_GDALCreate@28");
162  pGDALCreateCopy = get_symbol("_GDALCreateCopy@28");
163  pGDALSetRasterNoDataValue = get_symbol("_GDALSetRasterNoDataValue@12");
164  pGDALSetGeoTransform = get_symbol("_GDALSetGeoTransform@8");
165  pGDALSetProjection = get_symbol("_GDALSetProjection@8");
166  pGDALGetDriverShortName = get_symbol("_GDALGetDriverShortName@4");
167  pGDALGetDatasetDriver = get_symbol("_GDALGetDatasetDriver@4");
168 #else
169  pGDALAllRegister = get_symbol("GDALAllRegister");
170  pGDALOpen = get_symbol("GDALOpen");
171  pGDALClose = get_symbol("GDALClose");
172  pGDALGetRasterBand = get_symbol("GDALGetRasterBand");
173  pGDALRasterIO = get_symbol("GDALRasterIO");
174  pGDALGetDriverByName = get_symbol("GDALGetDriverByName");
175  pGDALGetMetadataItem = get_symbol("GDALGetMetadataItem");
176  pGDALCreate = get_symbol("GDALCreate");
177  pGDALCreateCopy = get_symbol("GDALCreateCopy");
178  pGDALSetRasterNoDataValue = get_symbol("GDALSetRasterNoDataValue");
179  pGDALSetGeoTransform = get_symbol("GDALSetGeoTransform");
180  pGDALSetProjection = get_symbol("GDALSetProjection");
181  pGDALGetDriverShortName = get_symbol("GDALGetDriverShortName");
182  pGDALGetDatasetDriver = get_symbol("GDALGetDatasetDriver");
183 #endif
184 }
185 
186 #else /* GDAL_DYNAMIC */
187 
188 static void init_gdal(void)
189 {
190  pGDALAllRegister = &GDALAllRegister;
191  pGDALOpen = &GDALOpen;
192  pGDALClose = &GDALClose;
193  pGDALGetRasterBand = &GDALGetRasterBand;
194  pGDALRasterIO = &GDALRasterIO;
195  pGDALGetDriverByName = &GDALGetDriverByName;
196  pGDALGetMetadataItem = &GDALGetMetadataItem;
197  pGDALCreate = &GDALCreate;
198  pGDALCreateCopy = &GDALCreateCopy;
199  pGDALSetRasterNoDataValue = &GDALSetRasterNoDataValue;
200  pGDALSetGeoTransform = &GDALSetGeoTransform;
201  pGDALSetProjection = &GDALSetProjection;
202  pGDALGetDriverShortName = &GDALGetDriverShortName;
203  pGDALGetDatasetDriver = &GDALGetDatasetDriver;
204 }
205 
206 #endif /* GDAL_DYNAMIC */
207 
208 #endif /* GDAL_LINK */
209 
210 /*!
211  \brief Initialization
212 
213  Register all GDAL drivers.
214 */
215 void Rast_init_gdal(void)
216 {
217 #ifdef GDAL_LINK
218  static int initialized;
219 
220  if (G_is_initialized(&initialized))
221  return;
222 
223  init_gdal();
224  (*pGDALAllRegister) ();
225  G_initialize_done(&initialized);
226 #endif
227 }
228 
229 /*!
230  \brief Get GDAL link settings for given raster map
231 
232  \param name map name
233  \param mapset name of mapset
234 
235  \return pointer to GDAL_link structure
236  \return NULL if link not found
237 */
238 struct GDAL_link *Rast_get_gdal_link(const char *name, const char *mapset)
239 {
240 #ifdef GDAL_LINK
241  GDALDatasetH data;
242  GDALRasterBandH band;
243  GDALDataType type;
244  RASTER_MAP_TYPE req_type;
245 #endif
246  const char *filename;
247  int band_num;
248  struct GDAL_link *gdal;
249  RASTER_MAP_TYPE map_type;
250  FILE *fp;
251  struct Key_Value *key_val;
252  const char *p;
253  DCELL null_val;
254  int hflip, vflip;
255 
256  if (!G_find_raster2(name, mapset))
257  return NULL;
258 
259  map_type = Rast_map_type(name, mapset);
260  if (map_type < 0)
261  return NULL;
262 
263  fp = G_fopen_old_misc("cell_misc", "gdal", name, mapset);
264  if (!fp)
265  return NULL;
266  key_val = G_fread_key_value(fp);
267  fclose(fp);
268 
269  if (!key_val)
270  return NULL;
271 
272  filename = G_find_key_value("file", key_val);
273  if (!filename)
274  return NULL;
275 
276  p = G_find_key_value("band", key_val);
277  if (!p)
278  return NULL;
279  band_num = atoi(p);
280  if (!band_num)
281  return NULL;
282 
283  p = G_find_key_value("null", key_val);
284  if (!p)
285  return NULL;
286  if (strcmp(p, "none") == 0)
287  Rast_set_d_null_value(&null_val, 1);
288  else
289  null_val = atof(p);
290 
291  hflip = G_find_key_value("hflip", key_val) ? 1 : 0;
292  vflip = G_find_key_value("vflip", key_val) ? 1 : 0;
293 
294 #ifdef GDAL_LINK
295  p = G_find_key_value("type", key_val);
296  if (!p)
297  return NULL;
298  type = atoi(p);
299 
300  switch (type) {
301  case GDT_Byte:
302  case GDT_Int16:
303  case GDT_UInt16:
304  case GDT_Int32:
305  case GDT_UInt32:
306  req_type = CELL_TYPE;
307  break;
308  case GDT_Float32:
309  req_type = FCELL_TYPE;
310  break;
311  case GDT_Float64:
312  req_type = DCELL_TYPE;
313  break;
314  default:
315  return NULL;
316  }
317 
318  if (req_type != map_type)
319  return NULL;
320 
321  Rast_init_gdal();
322 
323  data = (*pGDALOpen) (filename, GA_ReadOnly);
324  if (!data)
325  return NULL;
326 
327  band = (*pGDALGetRasterBand) (data, band_num);
328  if (!band) {
329  (*pGDALClose) (data);
330  return NULL;
331  }
332 #endif
333 
334  gdal = G_calloc(1, sizeof(struct GDAL_link));
335 
336  gdal->filename = G_store(filename);
337  gdal->band_num = band_num;
338  gdal->null_val = null_val;
339  gdal->hflip = hflip;
340  gdal->vflip = vflip;
341 #ifdef GDAL_LINK
342  gdal->data = data;
343  gdal->band = band;
344  gdal->type = type;
345 #endif
346 
347  return gdal;
348 }
349 
350 struct GDAL_Options
351 {
352  const char *dir;
353  const char *ext;
354  const char *format;
355  char **options;
356 };
357 
358 static struct state
359 {
360  int initialized;
361  struct GDAL_Options opts;
362  struct Key_Value *projinfo, *projunits, *projepsg;
363  char *srswkt;
364 } state;
365 
366 static struct state *st = &state;
367 
368 static void read_gdal_options(void)
369 {
370  FILE *fp;
371  struct Key_Value *key_val;
372  const char *p;
373 
374  fp = G_fopen_old("", "GDAL", G_mapset());
375  if (!fp)
376  G_fatal_error(_("Unable to open GDAL file"));
377  key_val = G_fread_key_value(fp);
378  fclose(fp);
379 
380  p = G_find_key_value("directory", key_val);
381  if (!p)
382  p = "gdal";
383  if (*p == '/') {
384  st->opts.dir = G_store(p);
385  }
386  else {
387  char path[GPATH_MAX];
388 
389  G_file_name(path, p, "", G_mapset());
390  st->opts.dir = G_store(path);
391  if (access(path, 0) != 0)
393  }
394 
395  p = G_find_key_value("extension", key_val);
396  st->opts.ext = G_store(p ? p : "");
397 
398  p = G_find_key_value("format", key_val);
399  st->opts.format = G_store(p ? p : "GTiff");
400 
401  p = G_find_key_value("options", key_val);
402  st->opts.options = p ? G_tokenize(p, ",") : NULL;
403 
404  G_free_key_value(key_val);
405 }
406 
407 /*!
408  \brief Create GDAL settings for given raster map
409 
410  \param name map name
411  \param map_type map type (CELL, FCELL, DCELL)
412 
413  \return pointer to allocated GDAL_link structure
414  \return NULL on error
415 */
416 struct GDAL_link *Rast_create_gdal_link(const char *name,
417  RASTER_MAP_TYPE map_type)
418 {
419 #ifdef GDAL_LINK
420  char path[GPATH_MAX];
421  GDALDriverH driver;
422  double transform[6];
423  struct GDAL_link *gdal;
424  FILE *fp;
425  struct Key_Value *key_val;
426  char buf[32];
427 
429 
430  Rast_init_gdal();
431 
432  if (!G_is_initialized(&st->initialized)) {
433  read_gdal_options();
434  st->projinfo = G_get_projinfo();
435  st->projunits = G_get_projunits();
436  st->projepsg = G_get_projepsg();
437  if (st->projinfo && st->projunits)
438  st->srswkt = GPJ_grass_to_wkt2(st->projinfo, st->projunits,
439  st->projepsg, 0, 0);
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:88
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 ** 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:599
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
char * GPJ_grass_to_wkt2(const struct Key_Value *proj_info, const struct Key_Value *proj_units, const struct Key_Value *proj_epsg, int esri_style, int prettify)
Converts a GRASS co-ordinate system representation to WKT style. EPSG code is preferred if available...
Definition: convert.c:152
#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:160
struct Cell_head wr_window
Definition: R.h:100
struct state * st
Definition: parser.c:102
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:879
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:416
#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
struct Key_Value * G_get_projepsg(void)
Gets EPSG information for the current location.
Definition: get_projinfo.c:85
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:215
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:90
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:94
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:101
void G_warning(const char *msg,...)
Print a warning message to stderr.
Definition: gis/error.c:204
struct GDAL_link * Rast_get_gdal_link(const char *name, const char *mapset)
Get GDAL link settings for given raster map.
Definition: gdal.c:238
struct Key_Value * G_create_key_value(void)
Allocate and initialize Key_Value structure.
Definition: key_value1.c:23