GRASS GIS 7 Programmer's Manual  7.9.dev(2020)-ef56e0001
raster/open.c
Go to the documentation of this file.
1 /*!
2  * \file lib/raster/open.c
3  *
4  * \brief Raster Library - Open raster file
5  *
6  * (C) 1999-2009 by the GRASS Development Team
7  *
8  * This program is free software under the GNU General Public
9  * License (>=v2). Read the file COPYING that comes with GRASS
10  * for details.
11  *
12  * \author USACERL and many others
13  */
14 
15 #include <unistd.h>
16 #include <string.h>
17 #include <sys/types.h>
18 #include <sys/stat.h>
19 #include <fcntl.h>
20 #include <errno.h>
21 
22 #include <grass/config.h>
23 #include <grass/gis.h>
24 #include <grass/raster.h>
25 #include <grass/glocale.h>
26 
27 #include "R.h"
28 #define FORMAT_FILE "f_format"
29 #define NULL_FILE "null"
30 /* cmpressed null file */
31 #define NULLC_FILE "nullcmpr"
32 
33 static int new_fileinfo(void)
34 {
35  int oldsize = R__.fileinfo_count;
36  int newsize = oldsize;
37  int i;
38 
39  for (i = 0; i < oldsize; i++)
40  if (R__.fileinfo[i].open_mode <= 0) {
41  memset(&R__.fileinfo[i], 0, sizeof(struct fileinfo));
42  R__.fileinfo[i].open_mode = -1;
43  return i;
44  }
45 
46  if (newsize < 20)
47  newsize += 20;
48  else
49  newsize *= 2;
50 
51  R__.fileinfo = G_realloc(R__.fileinfo, newsize * sizeof(struct fileinfo));
52 
53  /* Mark all cell files as closed */
54  for (i = oldsize; i < newsize; i++) {
55  memset(&R__.fileinfo[i], 0, sizeof(struct fileinfo));
56  R__.fileinfo[i].open_mode = -1;
57  }
58 
59  R__.fileinfo_count = newsize;
60 
61  return oldsize;
62 }
63 
64 /*!
65  * \brief Open raster file
66  *
67  * Arrange for the NULL-value bitmap to be read as well as the raster
68  * map. If no NULL-value bitmap exists, arrange for the production of
69  * NULL-values based on zeros in the raster map. If the map is
70  * floating-point, arrange for quantization to integer for
71  * Rast_get_c_row(), et. al., by reading the quantization rules
72  * for the map using Rast_read_quant(). If the programmer wants to read
73  * the floating point map using uing quant rules other than the ones
74  * stored in map's quant file, he/she should call Rast_set_quant_rules()
75  * after the call to Rast_open_old().
76  *
77  * \param name map name
78  * \param open_mode mode
79  * \param map_type map type (CELL, FCELL, DCELL)
80  *
81  * \return open file descriptor ( >= 0) if successful
82  */
83 
84 static int open_raster_new(const char *name, int open_mode,
85  RASTER_MAP_TYPE map_type);
86 
87 /*!
88  \brief Open an existing integer raster map (cell)
89 
90  Opens the existing cell file <i>name</i> in the <i>mapset</i> for
91  reading by Rast_get_row() with mapping into the current window.
92 
93  This routine opens the raster map <i>name</i> in <i>mapset</i> for
94  reading. A nonnegative file descriptor is returned if the open is
95  successful. Otherwise a diagnostic message is printed and a negative
96  value is returned. This routine does quite a bit of work. Since
97  GRASS users expect that all raster maps will be resampled into the
98  current region, the resampling index for the raster map is prepared
99  by this routine after the file is opened. The resampling is based on
100  the active module region (see also \ref The_Region}. Preparation
101  required for reading the various raster file formats (see \ref
102  Raster_File_Format for an explanation of the various raster file
103  formats) is also done.
104 
105  Diagnostics: warning message printed if open fails.
106 
107  \param name map name
108  \param mapset mapset name where raster map <i>name</i> lives
109 
110  \return nonnegative file descriptor (int)
111  */
112 int Rast_open_old(const char *name, const char *mapset)
113 {
114  int fd = Rast__open_old(name, mapset);
115 
116  /* turn on auto masking, if not already on */
118  /*
119  if(R__.auto_mask <= 0)
120  R__.mask_buf = Rast_allocate_c_buf();
121  now we don't ever free it!, so no need to allocate it (Olga)
122  */
123  /* mask_buf is used for reading MASK file when mask is set and
124  for reading map rows when the null file doesn't exist */
125 
126  return fd;
127 }
128 
129 /*! \brief Lower level function, open cell files, supercell files,
130  and the MASK file.
131 
132  Actions:
133  - opens the named cell file, following reclass reference if
134  named layer is a reclass layer.
135  - creates the required mapping between the data and the window
136  for use by the get_map_row family of routines.
137 
138  Diagnostics: Errors other than actual open failure will cause a
139  diagnostic to be delivered through G_warning() open failure messages
140  are left to the calling routine since the masking logic will want to
141  issue a different warning.
142 
143  Note: This routine does NOT open the MASK layer. If it did we would
144  get infinite recursion. This routine is called to open the mask by
145  Rast__check_for_auto_masking() which is called by Rast_open_old().
146 
147  \param name map name
148  \param mapset mapset of cell file to be opened
149 
150  \return open file descriptor
151  */
152 int Rast__open_old(const char *name, const char *mapset)
153 {
154  struct fileinfo *fcb;
155  int cell_fd, fd;
156  char *cell_dir;
157  const char *r_name;
158  const char *r_mapset;
159  struct Cell_head cellhd;
160  int CELL_nbytes = 0; /* bytes per cell in CELL map */
161  int INTERN_SIZE;
162  int reclass_flag;
163  int MAP_NBYTES;
164  RASTER_MAP_TYPE MAP_TYPE;
165  struct Reclass reclass;
166  char xname[GNAME_MAX], xmapset[GMAPSET_MAX];
167  struct GDAL_link *gdal;
168  struct R_vrt *vrt;
169 
170  Rast__init();
171 
172  G_unqualified_name(name, mapset, xname, xmapset);
173  name = xname;
174  mapset = xmapset;
175 
176  if (!G_find_raster2(name, mapset))
177  G_fatal_error(_("Raster map <%s> not found"),
178  G_fully_qualified_name(name, mapset));
179 
180  /* Check for reclassification */
181  reclass_flag = Rast_get_reclass(name, mapset, &reclass);
182 
183  switch (reclass_flag) {
184  case 0:
185  r_name = name;
186  r_mapset = mapset;
187  break;
188  case 1:
189  r_name = reclass.name;
190  r_mapset = reclass.mapset;
191  if (!G_find_raster2(r_name, r_mapset))
192  G_fatal_error(_("Unable to open raster map <%s@%s> since it is a reclass "
193  "of raster map <%s@%s> which does not exist"),
194  name, mapset, r_name, r_mapset);
195  break;
196  default: /* Error reading cellhd/reclass file */
197  G_fatal_error(_("Error reading reclass file for raster map <%s>"),
198  G_fully_qualified_name(name, mapset));
199  break;
200  }
201 
202  /* read the cell header */
203  Rast_get_cellhd(r_name, r_mapset, &cellhd);
204 
205  /* now check the type */
206  MAP_TYPE = Rast_map_type(r_name, r_mapset);
207  if (MAP_TYPE < 0)
208  G_fatal_error(_("Error reading map type for raster map <%s>"),
209  G_fully_qualified_name(name, mapset));
210 
211  if (MAP_TYPE == CELL_TYPE)
212  /* set the number of bytes for CELL map */
213  {
214  CELL_nbytes = cellhd.format + 1;
215  if (CELL_nbytes < 1)
216  G_fatal_error(_("Raster map <%s@%s>: format field in header file invalid"),
217  r_name, r_mapset);
218  }
219 
220  /* compressor */
221  if (MAP_TYPE != CELL_TYPE) {
222  /* fp maps do not use RLE */
223  /* previously, compressed simply meant yes (ZLIB) or no
224  * now compressed encodes compressor type
225  * 0: not compressed
226  * 1, 2: ZLIB
227  * 3: LZ4
228  * 4: BZIP2
229  * etc */
230  if (cellhd.compressed == 1)
231  cellhd.compressed = 2;
232  }
233  /* test if compressor type is supported */
234  if (!G_check_compressor(cellhd.compressed)) {
235  G_fatal_error(_("Compression with %s is not supported in this GRASS GIS installation"), G_compressor_name(cellhd.compressed));
236  }
237 
238  if (cellhd.proj != R__.rd_window.proj)
239  G_fatal_error(_("Raster map <%s> is in different projection than current region. "
240  "Found <%s>, should be <%s>."),
241  G_fully_qualified_name(name, mapset),
242  G_projection_name(cellhd.proj),
244 
245  if (cellhd.zone != R__.rd_window.zone)
246  G_fatal_error(_("Raster map <%s> is in different zone (%d) than current region (%d)"),
247  G_fully_qualified_name(name, mapset), cellhd.zone, R__.rd_window.zone);
248 
249  /* when map is int warn if too large cell size */
250  if (MAP_TYPE == CELL_TYPE && (unsigned int)CELL_nbytes > sizeof(CELL))
251  G_fatal_error(_("Raster map <%s>: bytes per cell (%d) too large"),
252  G_fully_qualified_name(name, mapset), CELL_nbytes);
253 
254  /* record number of bytes per cell */
255  if (MAP_TYPE == FCELL_TYPE) {
256  cell_dir = "fcell";
257  INTERN_SIZE = sizeof(FCELL);
258  MAP_NBYTES = XDR_FLOAT_NBYTES;
259  }
260  else if (MAP_TYPE == DCELL_TYPE) {
261  cell_dir = "fcell";
262  INTERN_SIZE = sizeof(DCELL);
263  MAP_NBYTES = XDR_DOUBLE_NBYTES;
264  }
265  else { /* integer */
266  cell_dir = "cell";
267  INTERN_SIZE = sizeof(CELL);
268  MAP_NBYTES = CELL_nbytes;
269  }
270 
271  gdal = Rast_get_gdal_link(r_name, r_mapset);
272  vrt = Rast_get_vrt(r_name, r_mapset);
273  cell_fd = -1;
274  if (gdal) {
275 #ifdef HAVE_GDAL
276  cell_fd = -1;
277 #else
278  G_fatal_error(_("Raster map <%s@%s> is a GDAL link but GRASS is compiled without GDAL support"),
279  r_name, r_mapset);
280 #endif
281  }
282  else if (vrt) {
283  cell_fd = -1;
284  }
285  else {
286  /* now actually open file for reading */
287  cell_fd = G_open_old(cell_dir, r_name, r_mapset);
288  if (cell_fd < 0)
289  G_fatal_error(_("Unable to open %s file for raster map <%s@%s>"),
290  cell_dir, r_name, r_mapset);
291  }
292 
293  fd = new_fileinfo();
294  fcb = &R__.fileinfo[fd];
295  fcb->data_fd = cell_fd;
296 
297  fcb->map_type = MAP_TYPE;
298 
299  /* Save cell header */
300  fcb->cellhd = cellhd;
301 
302  /* allocate null bitstream buffers for reading null rows */
303  fcb->null_fd = -1;
304  fcb->null_cur_row = -1;
305  fcb->null_bits = Rast__allocate_null_bits(cellhd.cols);
306 
307  /* mark closed */
308  fcb->open_mode = -1;
309 
310  /* save name and mapset */
311  fcb->name = G_store(name);
312  fcb->mapset = G_store(mapset);
313 
314  /* mark no data row in memory */
315  fcb->cur_row = -1;
316 
317  /* if reclass, copy reclass structure */
318  if ((fcb->reclass_flag = reclass_flag))
319  fcb->reclass = reclass;
320 
321  fcb->gdal = gdal;
322  fcb->vrt = vrt;
323  if (!gdal && !vrt) {
324  /* check for compressed data format, making initial reads if necessary */
325  if (Rast__check_format(fd) < 0) {
326  close(cell_fd); /* warning issued by check_format() */
327  G_fatal_error(_("Error reading format for <%s@%s>"),
328  r_name, r_mapset);
329  }
330  }
331 
332  if (!vrt) {
333  /* create the mapping from cell file to window */
335  }
336 
337  /*
338  * allocate the data buffer
339  * number of bytes per cell is cellhd.format+1
340  */
341 
342  /* for reading fcb->data is allocated to be fcb->cellhd.cols * fcb->nbytes
343  (= XDR_FLOAT/DOUBLE_NBYTES) */
344  fcb->data = (unsigned char *)G_calloc(fcb->cellhd.cols, MAP_NBYTES);
345 
346  /* initialize/read in quant rules for float point maps */
347  if (fcb->map_type != CELL_TYPE) {
348  if (fcb->reclass_flag)
350  &(fcb->quant));
351  else
352  Rast_read_quant(fcb->name, fcb->mapset, &(fcb->quant));
353  }
354 
355  /* now mark open for read: this must follow create_window_mapping() */
356  fcb->open_mode = OPEN_OLD;
357  fcb->io_error = 0;
358  fcb->map_type = MAP_TYPE;
359  fcb->nbytes = MAP_NBYTES;
360  fcb->null_row_ptr = NULL;
361 
362  if (!gdal && !vrt) {
363  /* First, check for compressed null file */
364  fcb->null_fd = G_open_old_misc("cell_misc", NULL_FILE, r_name, r_mapset);
365  if (fcb->null_fd < 0) {
366  fcb->null_fd = G_open_old_misc("cell_misc", NULLC_FILE, r_name, r_mapset);
367  if (fcb->null_fd >= 0) {
368  fcb->null_row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
369  if (Rast__read_null_row_ptrs(fd, fcb->null_fd) < 0) {
370  close(fcb->null_fd);
371  fcb->null_fd = -1;
372  G_free(fcb->null_row_ptr);
373  fcb->null_row_ptr = NULL;
374  }
375  }
376  }
377  fcb->null_file_exists = fcb->null_fd >= 0;
378  }
379 
380  return fd;
381 }
382 
383 /*!
384  \brief Opens a new cell file in a database (compressed)
385 
386  Opens a new cell file <i>name</i> in the current mapset for writing
387  by Rast_put_row().
388 
389  The file is created and filled with no data it is assumed that the
390  new cell file is to conform to the current window.
391 
392  The file must be written sequentially. Use Rast_open_new_random()
393  for non sequential writes.
394 
395  Note: the open actually creates a temporary file Rast_close() will
396  move the temporary file to the cell file and write out the necessary
397  support files (cellhd, cats, hist, etc.).
398 
399  Diagnostics: warning message printed if open fails
400 
401  Warning: calls to Rast_set_window() made after opening a new cell file
402  may create confusion and should be avoided the new cell file will be
403  created to conform to the window at the time of the open.
404 
405  \param name map name
406 
407  \return open file descriptor ( >= 0) if successful
408  \return negative integer if error
409  */
410 int Rast_open_c_new(const char *name)
411 {
412  return open_raster_new(name, OPEN_NEW_COMPRESSED, CELL_TYPE);
413 }
414 
415 /*!
416  \brief Opens a new cell file in a database (uncompressed)
417 
418  See also Rast_open_new().
419 
420  \param name map name
421 
422  \return open file descriptor ( >= 0) if successful
423  \return negative integer if error
424  */
426 {
427  return open_raster_new(name, OPEN_NEW_UNCOMPRESSED, CELL_TYPE);
428 }
429 
430 /*!
431  \brief Save histogram for newly create raster map (cell)
432 
433  If newly created cell files should have histograms, set flag=1
434  otherwise set flag=0. Applies to subsequent opens.
435 
436  \param flag flag indicator
437  */
438 void Rast_want_histogram(int flag)
439 {
440  R__.want_histogram = flag;
441 }
442 
443 /*!
444  \brief Sets the format for subsequent opens on new integer cell files
445  (uncompressed and random only).
446 
447  Warning: subsequent put_row calls will only write n+1 bytes per
448  cell. If the data requires more, the cell file will be written
449  incorrectly (but with n+1 bytes per cell)
450 
451  When writing float map: format is -1
452 
453  \param n format
454  */
456 /* sets the format for integer raster map */
457 {
458  R__.nbytes = n + 1;
459  if (R__.nbytes <= 0)
460  R__.nbytes = 1;
461  if (R__.nbytes > sizeof(CELL))
462  R__.nbytes = sizeof(CELL);
463 }
464 
465 /*!
466  \brief Get cell value format
467 
468  \param v cell
469 
470  \return cell format
471  */
473 {
474  unsigned int i;
475 
476  if (v >= 0)
477  for (i = 0; i < sizeof(CELL); i++)
478  if (!(v /= 256))
479  return i;
480  return sizeof(CELL) - 1;
481 }
482 
483 /*!
484  \brief Opens new fcell file in a database
485 
486  Opens a new floating-point map <i>name</i> in the current mapset for
487  writing. The type of the file (i.e. either double or float) is
488  determined and fixed at this point. The default is FCELL_TYPE. In
489  order to change this default
490 
491  Use Rast_set_fp_type() where type is one of DCELL_TYPE or FCELL_TYPE.
492 
493  See warnings and notes for Rast_open_new().
494 
495  \param name map name
496 
497  \return nonnegative file descriptor (int)
498  */
499 int Rast_open_fp_new(const char *name)
500 {
501  return open_raster_new(name, OPEN_NEW_COMPRESSED, R__.fp_type);
502 }
503 
504 /*!
505  \brief Opens new fcell file in a database (uncompressed)
506 
507  See Rast_open_fp_new() for details.
508 
509  \param name map name
510 
511  \return nonnegative file descriptor (int)
512  */
514 {
515  return open_raster_new(name, OPEN_NEW_UNCOMPRESSED, R__.fp_type);
516 }
517 
518 #ifdef HAVE_GDAL
519 static int open_raster_new_gdal(char *map, char *mapset,
520  RASTER_MAP_TYPE map_type)
521 {
522  int fd;
523  struct fileinfo *fcb;
524 
525  fd = new_fileinfo();
526  fcb = &R__.fileinfo[fd];
527  fcb->data_fd = -1;
528 
529  /* mark closed */
530  fcb->map_type = map_type;
531  fcb->open_mode = -1;
532 
533  fcb->gdal = Rast_create_gdal_link(map, map_type);
534  if (!fcb->gdal)
535  G_fatal_error(_("Unable to create GDAL link"));
536 
537  fcb->cellhd = R__.wr_window;
538  fcb->cellhd.compressed = 0;
539  fcb->nbytes = Rast_cell_size(fcb->map_type);
540  /* for writing fcb->data is allocated to be R__.wr_window.cols *
541  sizeof(CELL or DCELL or FCELL) */
542  fcb->data = G_calloc(R__.wr_window.cols, fcb->nbytes);
543 
544  fcb->name = map;
545  fcb->mapset = mapset;
546  fcb->cur_row = 0;
547 
548  fcb->row_ptr = NULL;
549  fcb->temp_name = NULL;
550  fcb->null_temp_name = NULL;
551  fcb->null_cur_row = 0;
552  fcb->null_bits = NULL;
553  fcb->null_fd = -1;
554  fcb->null_row_ptr = NULL;
555 
556  if (fcb->map_type != CELL_TYPE)
557  Rast_quant_init(&(fcb->quant));
558 
559  /* init cell stats */
560  /* now works only for int maps */
561  if (fcb->map_type == CELL_TYPE)
562  if ((fcb->want_histogram = R__.want_histogram))
564 
565  /* init range and if map is double/float init d/f_range */
566  Rast_init_range(&fcb->range);
567 
568  if (fcb->map_type != CELL_TYPE)
570 
571  /* mark file as open for write */
573  fcb->io_error = 0;
574 
575  return fd;
576 }
577 #endif /* HAVE_GDAL */
578 
579 static int open_raster_new(const char *name, int open_mode,
581 {
582  char xname[GNAME_MAX], xmapset[GMAPSET_MAX];
583  struct fileinfo *fcb;
584  int fd, cell_fd;
585  char *tempname;
586  char *map;
587  char *mapset;
588  const char *cell_dir;
589  int nbytes;
590 
591  Rast__init();
592 
593  switch (map_type) {
594  case CELL_TYPE:
595  cell_dir = "cell";
596  nbytes = R__.nbytes;
597  break;
598  case FCELL_TYPE:
599  nbytes = XDR_FLOAT_NBYTES;
600  cell_dir = "fcell";
601  break;
602  case DCELL_TYPE:
603  nbytes = XDR_DOUBLE_NBYTES;
604  cell_dir = "fcell";
605  break;
606  default:
607  G_fatal_error(_("Invalid map type <%d>"), map_type);
608  break;
609  }
610 
611  if (G_unqualified_name(name, G_mapset(), xname, xmapset) < 0)
612  G_fatal_error(_("Raster map <%s> is not in the current mapset (%s)"),
613  name, G_mapset());
614  map = G_store(xname);
615  mapset = G_store(xmapset);
616 
617  /* check for legal grass name */
618  if (G_legal_filename(map) < 0)
619  G_fatal_error(_("<%s> is an illegal file name"), map);
620 
621 #ifdef HAVE_GDAL
622  if (G_find_file2("", "GDAL", G_mapset()))
623  return open_raster_new_gdal(map, mapset, map_type);
624 #endif
625 
626  /* open a tempfile name */
627  tempname = G_tempfile();
628  cell_fd = creat(tempname, 0666);
629  if (cell_fd < 0) {
630  int err = errno;
631  G_free(mapset);
632  G_free(tempname);
633  G_free(map);
634  G_fatal_error(_("No temp files available: %s"), strerror(err));
635  }
636 
637  fd = new_fileinfo();
638  fcb = &R__.fileinfo[fd];
639  fcb->data_fd = cell_fd;
640 
641  /*
642  * since we are bypassing the normal open logic
643  * must create the cell element
644  */
645  G_make_mapset_element(cell_dir);
646 
647  /* mark closed */
648  fcb->map_type = map_type;
649  fcb->open_mode = -1;
650  fcb->gdal = NULL;
651  fcb->vrt = NULL;
652 
653  /* for writing fcb->data is allocated to be R__.wr_window.cols *
654  sizeof(CELL or DCELL or FCELL) */
655  fcb->data = (unsigned char *)G_calloc(R__.wr_window.cols,
656  Rast_cell_size(fcb->map_type));
657 
658  /*
659  * copy current window into cell header
660  * set format to cell/supercell
661  * for compressed writing
662  * allocate space to hold the row address array
663  */
664  fcb->cellhd = R__.wr_window;
665 
666  /* change open_mode to OPEN_NEW_UNCOMPRESSED if R__.compression_type == 0 ? */
667 
668  if (open_mode == OPEN_NEW_COMPRESSED && fcb->map_type == CELL_TYPE) {
669  fcb->row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
670  G_zero(fcb->row_ptr, (fcb->cellhd.rows + 1) * sizeof(off_t));
673 
674  fcb->nbytes = 1; /* to the minimum */
675  }
676  else {
677  fcb->nbytes = nbytes;
679  fcb->row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
680  G_zero(fcb->row_ptr, (fcb->cellhd.rows + 1) * sizeof(off_t));
683  }
684  else
685  fcb->cellhd.compressed = 0;
686 
687  if (fcb->map_type != CELL_TYPE) {
688  Rast_quant_init(&(fcb->quant));
689  }
690  }
691  if (open_mode == OPEN_NEW_COMPRESSED && fcb->map_type != CELL_TYPE &&
692  fcb->cellhd.compressed == 1) {
693  /* fp maps do not use RLE */
694  fcb->cellhd.compressed = 2;
695  }
696 
697  /* save name and mapset, and tempfile name */
698  fcb->name = map;
699  fcb->mapset = mapset;
700  fcb->temp_name = tempname;
701 
702  /* next row to be written (in order) is zero */
703  fcb->cur_row = 0;
704 
705  /* open a null tempfile name */
706  tempname = G_tempfile();
707  fcb->null_fd = creat(tempname, 0666);
708  if (fcb->null_fd < 0) {
709  int err = errno;
710  G_free(tempname);
711  G_free(fcb->name);
712  G_free(fcb->mapset);
713  G_free(fcb->temp_name);
714  close(cell_fd);
715  G_fatal_error(_("No temp files available: %s"), strerror(err));
716  }
717 
718  fcb->null_temp_name = tempname;
719 
720  fcb->null_row_ptr = NULL;
721  if (R__.compress_nulls) {
722  fcb->null_row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
723  G_zero(fcb->null_row_ptr, (fcb->cellhd.rows + 1) * sizeof(off_t));
725  }
726 
727  /* next row to be written (in order) is zero */
728  fcb->null_cur_row = 0;
729 
730  /* allocate null bitstream buffer for writing */
732 
733  /* init cell stats */
734  /* now works only for int maps */
735  if (fcb->map_type == CELL_TYPE)
736  if ((fcb->want_histogram = R__.want_histogram))
738 
739  /* init range and if map is double/float init d/f_range */
740  Rast_init_range(&fcb->range);
741 
742  if (fcb->map_type != CELL_TYPE)
744 
745  /* mark file as open for write */
746  fcb->open_mode = open_mode;
747  fcb->io_error = 0;
748 
749  return fd;
750 }
751 
752 int Rast__open_null_write(const char *name)
753 {
754  char xname[GNAME_MAX], xmapset[GMAPSET_MAX];
755  struct fileinfo *fcb;
756  int fd;
757  char *tempname;
758  char *map;
759  char *mapset;
760 
761  Rast__init();
762 
763  if (!G_find_raster2(name, G_mapset()))
764  G_fatal_error(_("Raster map <%s> does not exist in the current mapset (%s)"),
765  name, G_mapset());
766 
767  if (G_unqualified_name(name, G_mapset(), xname, xmapset) < 0)
768  G_fatal_error(_("Raster map <%s> is not in the current mapset (%s)"),
769  name, G_mapset());
770  map = G_store(xname);
771  mapset = G_store(xmapset);
772 
773  fd = new_fileinfo();
774  fcb = &R__.fileinfo[fd];
775 
776  G_zero(fcb, sizeof(*fcb));
777 
778  fcb->name = map;
779  fcb->mapset = mapset;
780 
781  Rast_get_cellhd(map, mapset, &fcb->cellhd);
782 
783  /* open a null tempfile name */
784  tempname = G_tempfile();
785  fcb->null_fd = creat(tempname, 0666);
786  if (fcb->null_fd < 0) {
787  int err = errno;
788  G_free(tempname);
789  G_free(fcb->name);
790  G_free(fcb->mapset);
791  G_fatal_error(_("No temp files available: %s"), strerror(err));
792  }
793  fcb->null_temp_name = tempname;
794 
795  if (R__.compress_nulls) {
796  fcb->null_row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
797  G_zero(fcb->null_row_ptr, (fcb->cellhd.rows + 1) * sizeof(off_t));
799  }
800 
801  /* allocate null bitstream buffer for writing */
803 
804  return fd;
805 }
806 
807 /*!
808  \brief Set raster map floating-point data format.
809 
810  This controls the storage type for floating-point maps. It affects
811  subsequent calls to G_open_fp_map_new(). The <i>type</i> must be
812  one of FCELL_TYPE (float) or DCELL_TYPE (double). The use of this
813  routine by applications is discouraged since its use would override
814  user preferences.
815 
816  \param type raster data type
817 
818  \return void
819  */
821 {
822  Rast__init();
823 
824  switch (map_type) {
825  case FCELL_TYPE:
826  case DCELL_TYPE:
827  R__.fp_type = map_type;
828  break;
829  default:
830  G_fatal_error(_("Rast_set_fp_type(): can only be called with FCELL_TYPE or DCELL_TYPE"));
831  break;
832  }
833 }
834 
835 /*!
836  \brief Check if raster map is floating-point
837 
838  Returns true (1) if raster map <i>name</i> in <i>mapset</i>
839  is a floating-point dataset; false(0) otherwise.
840 
841  \param name map name
842  \param mapset mapset name
843 
844  \return 1 floating-point
845  \return 0 int
846  */
847 int Rast_map_is_fp(const char *name, const char *mapset)
848 {
849  char path[GPATH_MAX];
850  const char *xmapset;
851 
852  xmapset = G_find_raster2(name, mapset);
853  if (!xmapset)
854  G_fatal_error(_("Raster map <%s> not found"),
855  G_fully_qualified_name(name, mapset));
856 
857  G_file_name(path, "fcell", name, xmapset);
858  if (access(path, 0) == 0)
859  return 1;
860 
861  G_file_name(path, "g3dcell", name, xmapset);
862  if (access(path, 0) == 0)
863  return 1;
864 
865  return 0;
866 }
867 
868 /*!
869  \brief Determine raster data type
870 
871  Determines if the raster map is floating point or integer. Returns
872  DCELL_TYPE for double maps, FCELL_TYPE for float maps, CELL_TYPE for
873  integer maps, -1 if error has occurred
874 
875  \param name map name
876  \param mapset mapset where map <i>name</i> lives
877 
878  \return raster data type
879  */
880 RASTER_MAP_TYPE Rast_map_type(const char *name, const char *mapset)
881 {
882  char path[GPATH_MAX];
883  const char *xmapset;
884 
885  xmapset = G_find_raster2(name, mapset);
886  if (!xmapset) {
887  if (mapset && *mapset)
888  G_fatal_error(_("Raster map <%s> not found in mapset <%s>"),
889  name, mapset);
890  else
891  G_fatal_error(_("Raster map <%s> not found"), name);
892  }
893 
894  G_file_name(path, "fcell", name, xmapset);
895 
896  if (access(path, 0) == 0)
897  return Rast__check_fp_type(name, xmapset);
898 
899  G_file_name(path, "g3dcell", name, xmapset);
900 
901  if (access(path, 0) == 0)
902  return DCELL_TYPE;
903 
904  return CELL_TYPE;
905 }
906 
907 /*!
908  \brief Determine raster type from descriptor
909 
910  Determines if the raster map is floating point or integer. Returns
911  DCELL_TYPE for double maps, FCELL_TYPE for float maps, CELL_TYPE for
912  integer maps, -1 if error has occurred
913 
914  \param fd file descriptor
915 
916  \return raster data type
917  */
919 {
920  struct fileinfo *fcb = &R__.fileinfo[fd];
921 
922  return fcb->map_type;
923 }
924 
925 /*!
926  \brief Determines whether the floating points cell file has double or float type
927 
928  \param name map name
929  \param mapset mapset where map <i>name</i> lives
930 
931  \return raster type (fcell, dcell)
932  */
934 {
935  char path[GPATH_MAX];
936  struct Key_Value *format_keys;
937  const char *str, *str1;
938  RASTER_MAP_TYPE map_type;
939  const char *xmapset;
940 
941  xmapset = G_find_raster2(name, mapset);
942  if (!xmapset)
943  G_fatal_error(_("Raster map <%s> not found"),
944  G_fully_qualified_name(name, mapset));
945 
946  G_file_name_misc(path, "cell_misc", FORMAT_FILE, name, xmapset);
947 
948  if (access(path, 0) != 0)
949  G_fatal_error(_("Unable to find '%s'"), path);
950 
951  format_keys = G_read_key_value_file(path);
952 
953  if ((str = G_find_key_value("type", format_keys)) != NULL) {
954  if (strcmp(str, "double") == 0)
955  map_type = DCELL_TYPE;
956  else if (strcmp(str, "float") == 0)
957  map_type = FCELL_TYPE;
958  else {
959  G_free_key_value(format_keys);
960  G_fatal_error(_("Invalid type: field '%s' in file '%s'"), str, path);
961  }
962  }
963  else {
964  G_free_key_value(format_keys);
965  G_fatal_error(_("Missing type: field in file '%s'"), path);
966  }
967 
968  if ((str1 = G_find_key_value("byte_order", format_keys)) != NULL) {
969  if (strcmp(str1, "xdr") != 0)
970  G_warning(_("Raster map <%s> is not xdr: byte_order: %s"),
971  name, str);
972  /* here read and translate byte order if not using xdr */
973  }
974  G_free_key_value(format_keys);
975  return map_type;
976 }
977 
978 /*!
979  \brief Opens a new raster map
980 
981  Opens a new raster map of type <i>wr_type</i>
982 
983  See warnings and notes for Rast_open_new().
984 
985  Supported data types:
986  - CELL_TYPE
987  - FCELL_TYPE
988  - DCELL_TYPE
989 
990  On CELL_TYPE calls Rast_open_new() otherwise Rast_open_fp_new().
991 
992  \param name map name
993  \param wr_type raster data type
994 
995  \return nonnegative file descriptor (int)
996  */
997 int Rast_open_new(const char *name, RASTER_MAP_TYPE wr_type)
998 {
999  return open_raster_new(name, OPEN_NEW_COMPRESSED, wr_type);
1000 }
1001 
1002 /*!
1003  \brief Opens a new raster map (uncompressed)
1004 
1005  See Rast_open_new().
1006 
1007  \param name map name
1008  \param wr_type raster data type
1009 
1010  \return nonnegative file descriptor (int)
1011  */
1013 {
1014  return open_raster_new(name, OPEN_NEW_UNCOMPRESSED, wr_type);
1015 }
1016 
1017 /*!
1018  \brief Sets quant translation rules for raster map opened for
1019  reading.
1020 
1021  Returned by Rast_open_old(). After calling this function,
1022  Rast_get_c_row() and Rast_get_c_row() will use rules defined by q
1023  (instead of using rules defined in map's quant file) to convert floats to
1024  ints.
1025 
1026  \param fd file descriptor (cell file)
1027  \param q pointer to Quant structure
1028 
1029  \return void
1030  */
1031 void Rast_set_quant_rules(int fd, struct Quant *q)
1032 {
1033  struct fileinfo *fcb = &R__.fileinfo[fd];
1034  CELL cell;
1035  DCELL dcell;
1036  struct Quant_table *p;
1037 
1038  if (fcb->open_mode != OPEN_OLD)
1039  G_fatal_error(_("Rast_set_quant_rules() can be called only for "
1040  "raster maps opened for reading"));
1041 
1042  /* copy all info from q to fcb->quant) */
1043  Rast_quant_init(&fcb->quant);
1044  if (q->truncate_only) {
1045  Rast_quant_truncate(&fcb->quant);
1046  return;
1047  }
1048 
1049  for (p = &(q->table[q->nofRules - 1]); p >= q->table; p--)
1050  Rast_quant_add_rule(&fcb->quant, p->dLow, p->dHigh, p->cLow,
1051  p->cHigh);
1052  if (Rast_quant_get_neg_infinite_rule(q, &dcell, &cell) > 0)
1053  Rast_quant_set_neg_infinite_rule(&fcb->quant, dcell, cell);
1054  if (Rast_quant_get_pos_infinite_rule(q, &dcell, &cell) > 0)
1055  Rast_quant_set_pos_infinite_rule(&fcb->quant, dcell, cell);
1056 }
int G_make_mapset_element(const char *p_element)
Create element in the current mapset.
Definition: mapset_msc.c:38
int Rast__check_for_auto_masking(void)
Checks for auto masking.
Definition: auto_mask.c:37
#define CELL_TYPE
Definition: raster.h:11
int Rast__write_null_row_ptrs(int fd, int null_fd)
#define OPEN_NEW_COMPRESSED
Definition: R.h:109
struct Cell_head rd_window
Definition: R.h:99
#define NULL_FILE
Definition: raster/open.c:29
int nbytes
Definition: R.h:73
const char * G_find_key_value(const char *key, const struct Key_Value *kv)
Find given key (case sensitive)
Definition: key_value1.c:84
int Rast_quant_get_neg_infinite_rule(const struct Quant *q, DCELL *dLeft, CELL *c)
Returns in "dLeft" and "c" the rule values.
Definition: quant.c:399
int Rast_get_reclass(const char *name, const char *mapset, struct Reclass *reclass)
Get reclass.
Definition: reclass.c:140
void G_free(void *buf)
Free allocated memory.
Definition: gis/alloc.c:149
void Rast__create_window_mapping(int fd)
Create window mapping.
char * name
Definition: raster.h:33
int want_histogram
Definition: R.h:62
int Rast_open_fp_new(const char *name)
Opens new fcell file in a database.
Definition: raster/open.c:499
RASTER_MAP_TYPE Rast__check_fp_type(const char *name, const char *mapset)
Determines whether the floating points cell file has double or float type.
Definition: raster/open.c:933
#define GMAPSET_MAX
Definition: gis.h:165
void Rast_quant_truncate(struct Quant *quant)
Sets the quant rules to perform simple truncation on floats.
Definition: quant.c:221
RASTER_MAP_TYPE map_type
Definition: R.h:74
void Rast_set_quant_rules(int fd, struct Quant *q)
Sets quant translation rules for raster map opened for reading.
Definition: raster/open.c:1031
Definition: R.h:88
char * name
Definition: R.h:78
2D/3D raster map header (used also for region)
Definition: gis.h:409
int data_fd
Definition: R.h:83
int Rast__check_format(int fd)
Definition: raster/format.c:63
void Rast_init_cell_stats(struct Cell_stats *s)
Initialize cell stats.
Definition: cell_stats.c:39
off_t * row_ptr
Definition: R.h:64
const char * G_projection_name(int n)
Get projection name.
Definition: proj2.c:55
double DCELL
Definition: gis.h:600
char * G_compressor_name(int number)
Definition: compress.c:118
int io_error
Definition: R.h:80
RASTER_MAP_TYPE fp_type
Definition: R.h:90
void Rast_want_histogram(int flag)
Save histogram for newly create raster map (cell)
Definition: raster/open.c:438
void Rast__init(void)
Definition: raster/init.c:65
int Rast__read_null_row_ptrs(int fd, int null_fd)
char * temp_name
Definition: R.h:75
char * G_store(const char *s)
Copy string to allocated memory.
Definition: strings.c:86
unsigned char * data
Definition: R.h:70
int G_unqualified_name(const char *name, const char *mapset, char *xname, char *xmapset)
Returns unqualified map name (without @ mapset)
Definition: nme_in_mps.c:134
DCELL dHigh
Definition: raster.h:79
char * G_tempfile(void)
Returns a temporary file name.
Definition: tempfile.c:62
struct R_vrt * vrt
Definition: R.h:85
int Rast_open_fp_new_uncompressed(const char *name)
Opens new fcell file in a database (uncompressed)
Definition: raster/open.c:513
#define NULLC_FILE
Definition: raster/open.c:31
int Rast_quant_get_pos_infinite_rule(const struct Quant *q, DCELL *dRight, CELL *c)
Returns in "dRight" and "c" the rule values.
Definition: quant.c:447
#define NULL
Definition: ccmath.h:32
unsigned char * null_bits
Definition: R.h:72
struct GDAL_link * gdal
Definition: R.h:82
int format
Max number of bytes per raster data value minus 1 (raster header only)
Definition: gis.h:415
struct R_vrt * Rast_get_vrt(const char *vname, const char *vmapset)
Definition: vrt.c:47
void Rast_quant_add_rule(struct Quant *q, DCELL dLow, DCELL dHigh, CELL cLow, CELL cHigh)
Adds a new rule to the set of quantization rules.
Definition: quant.c:478
struct Range range
Definition: R.h:60
#define OPEN_NEW_UNCOMPRESSED
Definition: R.h:110
void G_free_key_value(struct Key_Value *kv)
Free allocated Key_Value structure.
Definition: key_value1.c:103
struct Reclass reclass
Definition: R.h:58
int truncate_only
Definition: raster.h:86
void G_fatal_error(const char *msg,...)
Print a fatal error message to stderr.
Definition: gis/error.c:160
int G_open_old_misc(const char *dir, const char *element, const char *name, const char *mapset)
open a database file for reading
Definition: open_misc.c:134
struct Cell_head wr_window
Definition: R.h:100
#define FORMAT_FILE
Definition: raster/open.c:28
SYMBOL * err(FILE *fp, SYMBOL *s, char *msg)
Definition: symbol/read.c:220
Definition: raster.h:31
struct Key_Value * G_read_key_value_file(const char *file)
Read key/values pairs from file.
Definition: key_value3.c:53
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
CELL cLow
Definition: raster.h:80
int compressed
Compression mode (raster header only)
Definition: gis.h:422
#define DCELL_TYPE
Definition: raster.h:13
size_t Rast_cell_size(RASTER_MAP_TYPE data_type)
Returns size of a raster cell in bytes.
Definition: alloc_cell.c:39
RASTER_MAP_TYPE Rast_map_type(const char *name, const char *mapset)
Determine raster data type.
Definition: raster/open.c:880
void Rast_get_cellhd(const char *name, const char *mapset, struct Cell_head *cellhd)
Read the raster header.
Definition: get_cellhd.c:41
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:422
struct Cell_stats statf
Definition: R.h:59
void Rast_init_range(struct Range *range)
Initialize range structure.
Definition: range.c:678
struct Quant_table * table
Definition: raster.h:107
Definition: raster.h:84
DCELL dLow
Definition: raster.h:78
int Rast_open_c_new_uncompressed(const char *name)
Opens a new cell file in a database (uncompressed)
Definition: raster/open.c:425
int null_fd
Definition: R.h:71
void Rast_set_cell_format(int n)
Sets the format for subsequent opens on new integer cell files (uncompressed and random only)...
Definition: raster/open.c:455
char * mapset
Definition: R.h:79
int open_mode
Definition: R.h:56
int zone
Projection zone (UTM)
Definition: gis.h:443
int cur_row
Definition: R.h:67
void Rast_set_fp_type(RASTER_MAP_TYPE map_type)
Set raster map floating-point data format.
Definition: raster/open.c:820
int Rast_open_c_new(const char *name)
Opens a new cell file in a database (compressed)
Definition: raster/open.c:410
#define GPATH_MAX
Definition: gis.h:167
#define XDR_DOUBLE_NBYTES
Definition: R.h:8
int proj
Projection code.
Definition: gis.h:441
int fileinfo_count
Definition: R.h:102
struct fileinfo * fileinfo
Definition: R.h:103
char * null_temp_name
Definition: R.h:76
int Rast_read_quant(const char *name, const char *mapset, struct Quant *quant)
Reads quantization rules for name in mapset and stores them in the quantization structure. If the map is in another mapset, first checks for quant2 table for this map in current mapset.
Definition: quant_rw.c:186
Definition: gis.h:498
#define OPEN_OLD
Definition: R.h:108
float FCELL
Definition: gis.h:601
int G_open_old(const char *element, const char *name, const char *mapset)
Open a database file for reading.
Definition: gis/open.c:170
int Rast__open_null_write(const char *name)
Definition: raster/open.c:752
int Rast_open_new_uncompressed(const char *name, RASTER_MAP_TYPE wr_type)
Opens a new raster map (uncompressed)
Definition: raster/open.c:1012
RASTER_MAP_TYPE Rast_get_map_type(int fd)
Determine raster type from descriptor.
Definition: raster/open.c:918
int cols
Number of columns for 2D data.
Definition: gis.h:428
Definition: R.h:47
int null_file_exists
Definition: R.h:77
int G_check_compressor(int number)
Definition: compress.c:140
int compress_nulls
Definition: R.h:96
off_t * null_row_ptr
Definition: R.h:84
struct Cell_head cellhd
Definition: R.h:57
#define GNAME_MAX
Definition: gis.h:164
int CELL
Definition: gis.h:599
Definition: R.h:54
Definition: path.h:16
#define XDR_FLOAT_NBYTES
Definition: R.h:7
void Rast_quant_init(struct Quant *quant)
Initialize the structure.
Definition: quant.c:179
#define _(str)
Definition: glocale.h:13
int RASTER_MAP_TYPE
Definition: raster.h:25
int null_cur_row
Definition: R.h:68
char * G_fully_qualified_name(const char *name, const char *mapset)
Get fully qualified element name.
Definition: nme_in_mps.c:101
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
int Rast_open_new(const char *name, RASTER_MAP_TYPE wr_type)
Opens a new raster map.
Definition: raster/open.c:997
int Rast_map_is_fp(const char *name, const char *mapset)
Check if raster map is floating-point.
Definition: raster/open.c:847
struct Quant quant
Definition: R.h:81
const char * G_mapset(void)
Get current mapset name.
Definition: gis/mapset.c:33
char * G_file_name_misc(char *path, const char *dir, const char *element, const char *name, const char *mapset)
Builds full path names to GIS misc data files.
Definition: file_name.c:55
int want_histogram
Definition: R.h:93
void G_zero(void *buf, int i)
Zero out a buffer, buf, of length i.
Definition: gis/zero.c:23
int Rast_open_old(const char *name, const char *mapset)
Open an existing integer raster map (cell)
Definition: raster/open.c:112
int nbytes
Definition: R.h:94
const char * name
Definition: named_colr.c:7
unsigned char * Rast__allocate_null_bits(int cols)
Allocates memory for null bits.
Definition: alloc_cell.c:135
char * mapset
Definition: raster.h:34
void Rast_quant_set_pos_infinite_rule(struct Quant *q, DCELL dRight, CELL c)
Defines a rule for values "dRight" and larger.
Definition: quant.c:421
CELL cHigh
Definition: raster.h:81
int nofRules
Definition: raster.h:94
int rows
Number of rows for 2D data.
Definition: gis.h:424
struct FPRange fp_range
Definition: R.h:61
int Rast__write_row_ptrs(int fd)
void G_warning(const char *msg,...)
Print a warning message to stderr.
Definition: gis/error.c:204
void Rast_quant_set_neg_infinite_rule(struct Quant *q, DCELL dLeft, CELL c)
Defines a rule for values "dLeft" and smaller.
Definition: quant.c:373
const char * G_find_file2(const char *element, const char *name, const char *mapset)
Searches for a file from the mapset search list or in a specified mapset. (look but don&#39;t touch) ...
Definition: find_file.c:247
int Rast__open_old(const char *name, const char *mapset)
Lower level function, open cell files, supercell files, and the MASK file.
Definition: raster/open.c:152
void Rast_init_fp_range(struct FPRange *range)
Initialize fp range.
Definition: range.c:730
struct GDAL_link * Rast_get_gdal_link(const char *name, const char *mapset)
Get GDAL link settings for given raster map.
Definition: gdal.c:241
int Rast_get_cell_format(CELL v)
Get cell value format.
Definition: raster/open.c:472
int compression_type
Definition: R.h:95
int reclass_flag
Definition: R.h:63