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