GRASS GIS 7 Programmer's Manual  7.9.dev(2021)-0140b987a
vrt.c
Go to the documentation of this file.
1 /*!
2  \file lib/raster/vrt.c
3 
4  \brief Raster Library - virtual GRASS raster maps.
5 
6  (C) 2010 by the GRASS Development Team
7 
8  This program is free software under the GNU General Public License
9  (>=v2). Read the file COPYING that comes with GRASS for details.
10 
11  \author Markus Metz
12 */
13 
14 #include <grass/gis.h>
15 #include <grass/raster.h>
16 #include <grass/gprojects.h>
17 #include <grass/glocale.h>
18 
19 #include "R.h"
20 
21 int cmp_wnd(const void *a, const void *b)
22 {
23  struct Cell_head *cellhda = &((struct tileinfo *) a)->cellhd;
24  struct Cell_head *cellhdb = &((struct tileinfo *) b)->cellhd;
25 
26  /* sort from descending N to S, then ascending from W to E */
27  if (cellhda->south > cellhdb->south)
28  return -1;
29  if (cellhda->south < cellhdb->south)
30  return 1;
31  if (cellhda->north > cellhdb->north)
32  return -1;
33  if (cellhda->north < cellhdb->north)
34  return 1;
35  if (cellhda->west < cellhdb->west)
36  return -1;
37  if (cellhda->west > cellhdb->west)
38  return 1;
39  if (cellhda->east < cellhdb->east)
40  return -1;
41  if (cellhda->east > cellhdb->east)
42  return 1;
43 
44  return 0;
45 }
46 
47 struct R_vrt *Rast_get_vrt(const char *vname, const char *vmapset)
48 {
49  FILE *fp;
50  int talloc, tilecount;
51  struct tileinfo *ti;
52  struct R_vrt *vrt;
53  struct Cell_head *rd_window = &R__.rd_window;
54  struct ilist *tlist;
55 
56  tilecount = 0;
57  ti = NULL;
58 
59  if (!G_find_raster2(vname, vmapset))
60  return NULL;
61 
62  fp = G_fopen_old_misc("cell_misc", "vrt", vname, vmapset);
63  if (!fp)
64  return NULL;
65 
66  tlist = G_new_ilist();
67  talloc = 0;
68  while (1) {
69  char buf[GNAME_MAX];
70  char *name;
71  const char *mapset;
72  struct tileinfo *p;
73 
74  if (!G_getl2(buf, sizeof(buf), fp))
75  break;
76 
77  /* Ignore empty lines */
78  if (!*buf)
79  continue;
80 
81  name = buf;
82  if ((mapset = G_find_raster(name, "")) == NULL)
83  G_fatal_error(_("Tile raster map <%s> not found"), name);
84 
85  if (strcmp(name, vname) == 0)
86  G_fatal_error(_("A virtual raster can not contain itself"));
87 
88  if (tilecount >= talloc) {
89  talloc += 100;
90  ti = G_realloc(ti, talloc * sizeof(struct tileinfo));
91  }
92  p = &ti[tilecount];
93 
94  p->name = G_store(name);
95  p->mapset = G_store(mapset);
96  Rast_get_cellhd(p->name, p->mapset, &(p->cellhd));
97  p->clist = NULL;
98 
99  if (rd_window->proj == PROJECTION_LL) {
100  while (p->cellhd.west >= rd_window->east) {
101  p->cellhd.west -= 360.0;
102  p->cellhd.east -= 360.0;
103  }
104  while (p->cellhd.east <= rd_window->west) {
105  p->cellhd.west += 360.0;
106  p->cellhd.east += 360.0;
107  }
108  }
109 
110  if (p->cellhd.north > rd_window->south &&
111  p->cellhd.south <= rd_window->north &&
112  p->cellhd.west < rd_window->east &&
113  p->cellhd.east >= rd_window->west) {
114 
115  int col;
116  double east;
117 
118  G_ilist_add(tlist, tilecount);
119 
120  p->clist = G_new_ilist();
121  for (col = 0; col < rd_window->cols; col++) {
122  east = rd_window->west + rd_window->ew_res * (col + 0.5);
123 
124  if (rd_window->proj == PROJECTION_LL) {
125  while (east > p->cellhd.east)
126  east -= 360;
127  while (east < p->cellhd.west)
128  east += 360;
129  }
130  if (east >= p->cellhd.west && east < p->cellhd.east)
131  G_ilist_add(p->clist, col);
132  }
133  }
134  tilecount++;
135  }
136 
137  if (tilecount > 1)
138  qsort(ti, tilecount, sizeof(struct tileinfo), cmp_wnd);
139 
140  fclose(fp);
141 
142  vrt = G_calloc(1, sizeof(struct R_vrt));
143  vrt->tilecount = tilecount;
144  vrt->tileinfo = ti;
145  vrt->tlist = tlist;
146 
147  return vrt;
148 }
149 
150 void Rast_close_vrt(struct R_vrt *vrt)
151 {
152  int i;
153 
154  for (i = 0; i < vrt->tilecount; i++) {
155  struct tileinfo *p;
156 
157  p = &(vrt->tileinfo[i]);
158 
159  G_free(p->name);
160  G_free(p->mapset);
161  if (p->clist)
162  G_free_ilist(p->clist);
163  }
164  G_free(vrt->tileinfo);
165  G_free_ilist(vrt->tlist);
166  G_free(vrt);
167 }
168 
169 /* must only be called by get_map_row_nomask()
170  * move to get_row.c as read_data_vrt() ? */
171 int Rast_get_vrt_row(int fd, void *buf, int row, RASTER_MAP_TYPE data_type)
172 {
173  struct fileinfo *fcb = &R__.fileinfo[fd];
174  struct R_vrt *vrt = fcb->vrt;
175  struct tileinfo *ti = vrt->tileinfo;
176  struct Cell_head *rd_window = &R__.rd_window;
177  double rown, rows;
178  int i, j;
179  int have_tile;
180  void *tmpbuf;
181  size_t size = Rast_cell_size(data_type);
182 
183  rown = rd_window->north - rd_window->ns_res * row;
184  rows = rd_window->north - rd_window->ns_res * (row + 1);
185 
186  Rast_set_null_value(buf, rd_window->cols, data_type);
187  tmpbuf = Rast_allocate_input_buf(data_type);
188  have_tile = 0;
189 
190  for (i = 0; i < vrt->tlist->n_values; i++) {
191  struct tileinfo *p = &ti[vrt->tlist->value[i]];
192 
193  if (p->cellhd.north > rows && p->cellhd.south <= rown) {
194  int tfd;
195  void *p1, *p2;
196 
197  /* recurse into get_map_row(), collect data for all tiles
198  * a mask is applied to the collected data
199  * after this function returns */
200  Rast_set_null_value(tmpbuf, rd_window->cols, data_type);
201  /* avoid Rast__check_for_auto_masking() */
202  tfd = Rast__open_old(p->name, p->mapset);
203  Rast_get_row_nomask(tfd, tmpbuf, row, data_type);
204  Rast_unopen(tfd);
205 
206  p1 = buf;
207  p2 = tmpbuf;
208  /* restrict to start and end col ? */
209  for (j = 0; j < p->clist->n_values; j++) {
210  p1 = (unsigned char *)buf + size * p->clist->value[j];
211  p2 = (unsigned char *)tmpbuf + size * p->clist->value[j];
212 
213  if (!Rast_is_null_value(p2, data_type)) {
214  switch (data_type) {
215  case CELL_TYPE:
216  *(CELL *) p1 = *(CELL *) p2;
217  break;
218  case FCELL_TYPE:
219  *(FCELL *) p1 = *(FCELL *) p2;
220  break;
221  case DCELL_TYPE:
222  *(DCELL *) p1 = *(DCELL *) p2;
223  break;
224  default:
225  break;
226  }
227  }
228  }
229  have_tile = 1;
230  }
231  }
232  G_free(tmpbuf);
233 
234  return have_tile;
235 }
int G_getl2(char *, int, FILE *)
Gets a line of text from a file of any pedigree.
Definition: getl.c:64
int tilecount
Definition: R.h:49
#define CELL_TYPE
Definition: raster.h:11
struct Cell_head rd_window
Definition: R.h:99
struct ilist * clist
Definition: R.h:44
void Rast_set_null_value(void *, int, RASTER_MAP_TYPE)
To set one or more raster values to null.
Definition: null_val.c:98
void void void void G_fatal_error(const char *,...) __attribute__((format(printf
int Rast_get_vrt_row(int fd, void *buf, int row, RASTER_MAP_TYPE data_type)
Definition: vrt.c:171
struct tileinfo * tileinfo
Definition: R.h:50
Definition: R.h:88
2D/3D raster map header (used also for region)
Definition: gis.h:409
double west
Extent coordinates (west)
Definition: gis.h:461
double DCELL
Definition: gis.h:600
void Rast_get_row_nomask(int, void *, int, RASTER_MAP_TYPE)
Read raster row without masking.
void Rast_close_vrt(struct R_vrt *vrt)
Definition: vrt.c:150
int n_values
Number of values in the list.
Definition: gis.h:695
void G_free(void *)
Free allocated memory.
Definition: gis/alloc.c:149
int Rast__open_old(const char *, const char *)
Lower level function, open cell files, supercell files, and the MASK file.
Definition: raster/open.c:152
struct R_vrt * vrt
Definition: R.h:85
#define NULL
Definition: ccmath.h:32
struct R_vrt * Rast_get_vrt(const char *vname, const char *vmapset)
Definition: vrt.c:47
void * Rast_allocate_input_buf(RASTER_MAP_TYPE)
Definition: alloc_cell.c:157
struct Cell_head cellhd
Definition: R.h:43
#define G_calloc(m, n)
Definition: defs/gis.h:113
void G_ilist_add(struct ilist *, int)
Add item to ilist.
Definition: ilist.c:77
struct ilist * tlist
Definition: R.h:51
#define DCELL_TYPE
Definition: raster.h:13
double north
Extent coordinates (north)
Definition: gis.h:455
double b
Definition: r_raster.c:39
void Rast_unopen(int)
Unopen a raster map.
Definition: raster/close.c:132
size_t Rast_cell_size(RASTER_MAP_TYPE)
Returns size of a raster cell in bytes.
Definition: alloc_cell.c:39
double south
Extent coordinates (south)
Definition: gis.h:457
char * name
Definition: R.h:41
FILE * G_fopen_old_misc(const char *, const char *, const char *, const char *)
open a database file for reading
Definition: open_misc.c:210
int cmp_wnd(const void *a, const void *b)
Definition: vrt.c:21
#define PROJECTION_LL
Projection code - Latitude-Longitude.
Definition: gis.h:113
int proj
Projection code.
Definition: gis.h:441
struct fileinfo * fileinfo
Definition: R.h:103
float FCELL
Definition: gis.h:601
const char * G_find_raster(char *, const char *)
Find a raster map.
Definition: find_rast.c:55
int cols
Number of columns for 2D data.
Definition: gis.h:428
Definition: R.h:47
void Rast_get_cellhd(const char *, const char *, struct Cell_head *)
Read the raster header.
Definition: get_cellhd.c:41
double ns_res
Resolution - north to south cell size for 2D data.
Definition: gis.h:449
#define GNAME_MAX
Definition: gis.h:164
int CELL
Definition: gis.h:599
Definition: R.h:54
double east
Extent coordinates (east)
Definition: gis.h:459
#define G_realloc(p, n)
Definition: defs/gis.h:114
#define _(str)
Definition: glocale.h:10
List of integers.
Definition: gis.h:686
int RASTER_MAP_TYPE
Definition: raster.h:25
#define FCELL_TYPE
Definition: raster.h:12
const char * G_find_raster2(const char *, const char *)
Find a raster map (look but don&#39;t touch)
Definition: find_rast.c:76
char * G_store(const char *)
Copy string to allocated memory.
Definition: strings.c:87
int * value
Array of values.
Definition: gis.h:691
Definition: R.h:39
const char * name
Definition: named_colr.c:7
char * mapset
Definition: R.h:42
double ew_res
Resolution - east to west cell size for 2D data.
Definition: gis.h:445
struct ilist * G_new_ilist()
Return a new integer list.
Definition: ilist.c:43
int rows
Number of rows for 2D data.
Definition: gis.h:424
void G_free_ilist(struct ilist *)
Free allocated memory of an integer list.
Definition: ilist.c:27
int Rast_is_null_value(const void *, RASTER_MAP_TYPE)
To check if a raster value is set to NULL.
Definition: null_val.c:179