GRASS GIS 7 Programmer's Manual  7.7.svn(2018)-r73373
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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,
580  RASTER_MAP_TYPE map_type)
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 
652  /* for writing fcb->data is allocated to be R__.wr_window.cols *
653  sizeof(CELL or DCELL or FCELL) */
654  fcb->data = (unsigned char *)G_calloc(R__.wr_window.cols,
655  Rast_cell_size(fcb->map_type));
656 
657  /*
658  * copy current window into cell header
659  * set format to cell/supercell
660  * for compressed writing
661  * allocate space to hold the row address array
662  */
663  fcb->cellhd = R__.wr_window;
664 
665  /* change open_mode to OPEN_NEW_UNCOMPRESSED if R__.compression_type == 0 ? */
666 
667  if (open_mode == OPEN_NEW_COMPRESSED && fcb->map_type == CELL_TYPE) {
668  fcb->row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
669  G_zero(fcb->row_ptr, (fcb->cellhd.rows + 1) * sizeof(off_t));
672 
673  fcb->nbytes = 1; /* to the minimum */
674  }
675  else {
676  fcb->nbytes = nbytes;
677  if (open_mode == OPEN_NEW_COMPRESSED) {
678  fcb->row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
679  G_zero(fcb->row_ptr, (fcb->cellhd.rows + 1) * sizeof(off_t));
682  }
683  else
684  fcb->cellhd.compressed = 0;
685 
686  if (fcb->map_type != CELL_TYPE) {
687  Rast_quant_init(&(fcb->quant));
688  }
689  }
690  if (open_mode == OPEN_NEW_COMPRESSED && fcb->map_type != CELL_TYPE &&
691  fcb->cellhd.compressed == 1) {
692  /* fp maps do not use RLE */
693  fcb->cellhd.compressed = 2;
694  }
695 
696  /* save name and mapset, and tempfile name */
697  fcb->name = map;
698  fcb->mapset = mapset;
699  fcb->temp_name = tempname;
700 
701  /* next row to be written (in order) is zero */
702  fcb->cur_row = 0;
703 
704  /* open a null tempfile name */
705  tempname = G_tempfile();
706  fcb->null_fd = creat(tempname, 0666);
707  if (fcb->null_fd < 0) {
708  int err = errno;
709  G_free(tempname);
710  G_free(fcb->name);
711  G_free(fcb->mapset);
712  G_free(fcb->temp_name);
713  close(cell_fd);
714  G_fatal_error(_("No temp files available: %s"), strerror(err));
715  }
716 
717  fcb->null_temp_name = tempname;
718 
719  fcb->null_row_ptr = NULL;
720  if (R__.compress_nulls) {
721  fcb->null_row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
722  G_zero(fcb->null_row_ptr, (fcb->cellhd.rows + 1) * sizeof(off_t));
724  }
725 
726  /* next row to be written (in order) is zero */
727  fcb->null_cur_row = 0;
728 
729  /* allocate null bitstream buffer for writing */
731 
732  /* init cell stats */
733  /* now works only for int maps */
734  if (fcb->map_type == CELL_TYPE)
735  if ((fcb->want_histogram = R__.want_histogram))
737 
738  /* init range and if map is double/float init d/f_range */
739  Rast_init_range(&fcb->range);
740 
741  if (fcb->map_type != CELL_TYPE)
743 
744  /* mark file as open for write */
745  fcb->open_mode = open_mode;
746  fcb->io_error = 0;
747 
748  return fd;
749 }
750 
751 int Rast__open_null_write(const char *name)
752 {
753  char xname[GNAME_MAX], xmapset[GMAPSET_MAX];
754  struct fileinfo *fcb;
755  int fd;
756  char *tempname;
757  char *map;
758  char *mapset;
759 
760  Rast__init();
761 
762  if (!G_find_raster2(name, G_mapset()))
763  G_fatal_error(_("Raster map <%s> does not exist in the current mapset (%s)"),
764  name, G_mapset());
765 
766  if (G_unqualified_name(name, G_mapset(), xname, xmapset) < 0)
767  G_fatal_error(_("Raster map <%s> is not in the current mapset (%s)"),
768  name, G_mapset());
769  map = G_store(xname);
770  mapset = G_store(xmapset);
771 
772  fd = new_fileinfo();
773  fcb = &R__.fileinfo[fd];
774 
775  G_zero(fcb, sizeof(*fcb));
776 
777  fcb->name = map;
778  fcb->mapset = mapset;
779 
780  Rast_get_cellhd(map, mapset, &fcb->cellhd);
781 
782  /* open a null tempfile name */
783  tempname = G_tempfile();
784  fcb->null_fd = creat(tempname, 0666);
785  if (fcb->null_fd < 0) {
786  int err = errno;
787  G_free(tempname);
788  G_free(fcb->name);
789  G_free(fcb->mapset);
790  G_fatal_error(_("No temp files available: %s"), strerror(err));
791  }
792  fcb->null_temp_name = tempname;
793 
794  if (R__.compress_nulls) {
795  fcb->null_row_ptr = G_calloc(fcb->cellhd.rows + 1, sizeof(off_t));
796  G_zero(fcb->null_row_ptr, (fcb->cellhd.rows + 1) * sizeof(off_t));
798  }
799 
800  /* allocate null bitstream buffer for writing */
802 
803  return fd;
804 }
805 
806 /*!
807  \brief Set raster map floating-point data format.
808 
809  This controls the storage type for floating-point maps. It affects
810  subsequent calls to G_open_fp_map_new(). The <i>type</i> must be
811  one of FCELL_TYPE (float) or DCELL_TYPE (double). The use of this
812  routine by applications is discouraged since its use would override
813  user preferences.
814 
815  \param type raster data type
816 
817  \return void
818  */
820 {
821  Rast__init();
822 
823  switch (map_type) {
824  case FCELL_TYPE:
825  case DCELL_TYPE:
826  R__.fp_type = map_type;
827  break;
828  default:
829  G_fatal_error(_("Rast_set_fp_type(): can only be called with FCELL_TYPE or DCELL_TYPE"));
830  break;
831  }
832 }
833 
834 /*!
835  \brief Check if raster map is floating-point
836 
837  Returns true (1) if raster map <i>name</i> in <i>mapset</i>
838  is a floating-point dataset; false(0) otherwise.
839 
840  \param name map name
841  \param mapset mapset name
842 
843  \return 1 floating-point
844  \return 0 int
845  */
846 int Rast_map_is_fp(const char *name, const char *mapset)
847 {
848  char path[GPATH_MAX];
849  const char *xmapset;
850 
851  xmapset = G_find_raster2(name, mapset);
852  if (!xmapset)
853  G_fatal_error(_("Raster map <%s> not found"),
854  G_fully_qualified_name(name, mapset));
855 
856  G_file_name(path, "fcell", name, xmapset);
857  if (access(path, 0) == 0)
858  return 1;
859 
860  G_file_name(path, "g3dcell", name, xmapset);
861  if (access(path, 0) == 0)
862  return 1;
863 
864  return 0;
865 }
866 
867 /*!
868  \brief Determine raster data type
869 
870  Determines if the raster map is floating point or integer. Returns
871  DCELL_TYPE for double maps, FCELL_TYPE for float maps, CELL_TYPE for
872  integer maps, -1 if error has occurred
873 
874  \param name map name
875  \param mapset mapset where map <i>name</i> lives
876 
877  \return raster data type
878  */
879 RASTER_MAP_TYPE Rast_map_type(const char *name, const char *mapset)
880 {
881  char path[GPATH_MAX];
882  const char *xmapset;
883 
884  xmapset = G_find_raster2(name, mapset);
885  if (!xmapset) {
886  if (mapset && *mapset)
887  G_fatal_error(_("Raster map <%s> not found in mapset <%s>"),
888  name, mapset);
889  else
890  G_fatal_error(_("Raster map <%s> not found"), name);
891  }
892 
893  G_file_name(path, "fcell", name, xmapset);
894 
895  if (access(path, 0) == 0)
896  return Rast__check_fp_type(name, xmapset);
897 
898  G_file_name(path, "g3dcell", name, xmapset);
899 
900  if (access(path, 0) == 0)
901  return DCELL_TYPE;
902 
903  return CELL_TYPE;
904 }
905 
906 /*!
907  \brief Determine raster type from descriptor
908 
909  Determines if the raster map is floating point or integer. Returns
910  DCELL_TYPE for double maps, FCELL_TYPE for float maps, CELL_TYPE for
911  integer maps, -1 if error has occurred
912 
913  \param fd file descriptor
914 
915  \return raster data type
916  */
918 {
919  struct fileinfo *fcb = &R__.fileinfo[fd];
920 
921  return fcb->map_type;
922 }
923 
924 /*!
925  \brief Determines whether the floating points cell file has double or float type
926 
927  \param name map name
928  \param mapset mapset where map <i>name</i> lives
929 
930  \return raster type (fcell, dcell)
931  */
932 RASTER_MAP_TYPE Rast__check_fp_type(const char *name, const char *mapset)
933 {
934  char path[GPATH_MAX];
935  struct Key_Value *format_keys;
936  const char *str, *str1;
937  RASTER_MAP_TYPE map_type;
938  const char *xmapset;
939 
940  xmapset = G_find_raster2(name, mapset);
941  if (!xmapset)
942  G_fatal_error(_("Raster map <%s> not found"),
943  G_fully_qualified_name(name, mapset));
944 
945  G_file_name_misc(path, "cell_misc", FORMAT_FILE, name, xmapset);
946 
947  if (access(path, 0) != 0)
948  G_fatal_error(_("Unable to find '%s'"), path);
949 
950  format_keys = G_read_key_value_file(path);
951 
952  if ((str = G_find_key_value("type", format_keys)) != NULL) {
953  if (strcmp(str, "double") == 0)
954  map_type = DCELL_TYPE;
955  else if (strcmp(str, "float") == 0)
956  map_type = FCELL_TYPE;
957  else {
958  G_free_key_value(format_keys);
959  G_fatal_error(_("Invalid type: field '%s' in file '%s'"), str, path);
960  }
961  }
962  else {
963  G_free_key_value(format_keys);
964  G_fatal_error(_("Missing type: field in file '%s'"), path);
965  }
966 
967  if ((str1 = G_find_key_value("byte_order", format_keys)) != NULL) {
968  if (strcmp(str1, "xdr") != 0)
969  G_warning(_("Raster map <%s> is not xdr: byte_order: %s"),
970  name, str);
971  /* here read and translate byte order if not using xdr */
972  }
973  G_free_key_value(format_keys);
974  return map_type;
975 }
976 
977 /*!
978  \brief Opens a new raster map
979 
980  Opens a new raster map of type <i>wr_type</i>
981 
982  See warnings and notes for Rast_open_new().
983 
984  Supported data types:
985  - CELL_TYPE
986  - FCELL_TYPE
987  - DCELL_TYPE
988 
989  On CELL_TYPE calls Rast_open_new() otherwise Rast_open_fp_new().
990 
991  \param name map name
992  \param wr_type raster data type
993 
994  \return nonnegative file descriptor (int)
995  */
996 int Rast_open_new(const char *name, RASTER_MAP_TYPE wr_type)
997 {
998  return open_raster_new(name, OPEN_NEW_COMPRESSED, wr_type);
999 }
1000 
1001 /*!
1002  \brief Opens a new raster map (uncompressed)
1003 
1004  See Rast_open_new().
1005 
1006  \param name map name
1007  \param wr_type raster data type
1008 
1009  \return nonnegative file descriptor (int)
1010  */
1011 int Rast_open_new_uncompressed(const char *name, RASTER_MAP_TYPE wr_type)
1012 {
1013  return open_raster_new(name, OPEN_NEW_UNCOMPRESSED, wr_type);
1014 }
1015 
1016 /*!
1017  \brief Sets quant translation rules for raster map opened for
1018  reading.
1019 
1020  Returned by Rast_open_old(). After calling this function,
1021  Rast_get_c_row() and Rast_get_c_row() will use rules defined by q
1022  (instead of using rules defined in map's quant file) to convert floats to
1023  ints.
1024 
1025  \param fd file descriptor (cell file)
1026  \param q pointer to Quant structure
1027 
1028  \return void
1029  */
1030 void Rast_set_quant_rules(int fd, struct Quant *q)
1031 {
1032  struct fileinfo *fcb = &R__.fileinfo[fd];
1033  CELL cell;
1034  DCELL dcell;
1035  struct Quant_table *p;
1036 
1037  if (fcb->open_mode != OPEN_OLD)
1038  G_fatal_error(_("Rast_set_quant_rules() can be called only for "
1039  "raster maps opened for reading"));
1040 
1041  /* copy all info from q to fcb->quant) */
1042  Rast_quant_init(&fcb->quant);
1043  if (q->truncate_only) {
1044  Rast_quant_truncate(&fcb->quant);
1045  return;
1046  }
1047 
1048  for (p = &(q->table[q->nofRules - 1]); p >= q->table; p--)
1049  Rast_quant_add_rule(&fcb->quant, p->dLow, p->dHigh, p->cLow,
1050  p->cHigh);
1051  if (Rast_quant_get_neg_infinite_rule(q, &dcell, &cell) > 0)
1052  Rast_quant_set_neg_infinite_rule(&fcb->quant, dcell, cell);
1053  if (Rast_quant_get_pos_infinite_rule(q, &dcell, &cell) > 0)
1054  Rast_quant_set_pos_infinite_rule(&fcb->quant, dcell, cell);
1055 }
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 &quot;dLeft&quot; and &quot;c&quot; 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:932
#define GMAPSET_MAX
Definition: gis.h:149
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:1030
Definition: R.h:88
char * name
Definition: R.h:78
2D/3D raster map header (used also for region)
Definition: gis.h:390
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:581
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 &quot;dRight&quot; and &quot;c&quot; 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:396
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
fd
Definition: d/range.c:69
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:403
#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:879
void Rast_get_cellhd(const char *name, const char *mapset, struct Cell_head *cellhd)
Read the raster header.
Definition: get_cellhd.c:44
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:415
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:424
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:819
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:151
#define XDR_DOUBLE_NBYTES
Definition: R.h:8
int proj
Projection code.
Definition: gis.h:422
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:479
#define OPEN_OLD
Definition: R.h:108
float FCELL
Definition: gis.h:582
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:751
int Rast_open_new_uncompressed(const char *name, RASTER_MAP_TYPE wr_type)
Opens a new raster map (uncompressed)
Definition: raster/open.c:1011
RASTER_MAP_TYPE Rast_get_map_type(int fd)
Determine raster type from descriptor.
Definition: raster/open.c:917
int cols
Number of columns for 2D data.
Definition: gis.h:409
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:148
int CELL
Definition: gis.h:580
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:996
int Rast_map_is_fp(const char *name, const char *mapset)
Check if raster map is floating-point.
Definition: raster/open.c:846
struct Quant quant
Definition: R.h:81
const char * G_mapset(void)
Get current mapset name.
Definition: gis/mapset.c:33
close(fd)
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 &quot;dRight&quot; 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:405
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 &quot;dLeft&quot; 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:237
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