~ubuntu-branches/ubuntu/vivid/grass/vivid-proposed

« back to all changes in this revision

Viewing changes to lib/raster3d/region.c

  • Committer: Package Import Robot
  • Author(s): Bas Couwenberg
  • Date: 2015-02-20 23:12:08 UTC
  • mfrom: (8.2.6 experimental)
  • Revision ID: package-import@ubuntu.com-20150220231208-1u6qvqm84v430b10
Tags: 7.0.0-1~exp1
* New upstream release.
* Update python-ctypes-ternary.patch to use if/else instead of and/or.
* Drop check4dev patch, rely on upstream check.
* Add build dependency on libpq-dev to grass-dev for libpq-fe.h.
* Drop patches applied upstream, refresh remaining patches.
* Update symlinks for images switched from jpg to png.

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#include <stdio.h>
 
2
#include <math.h>
 
3
 
 
4
#include <grass/raster.h>
 
5
#include "raster3d_intern.h"
 
6
 
 
7
 
 
8
/*---------------------------------------------------------------------------*/
 
9
 
 
10
 
 
11
/*!
 
12
 * \brief
 
13
 *
 
14
 *  Returns in <em>region2d</em> the <em>2d</em> portion of <em>region3d</em>.
 
15
 *
 
16
 *  \param region3d
 
17
 *  \param region2d
 
18
 *  \return void
 
19
 */
 
20
 
 
21
void Rast3d_extract2d_region(RASTER3D_Region * region3d, struct Cell_head *region2d)
 
22
{
 
23
    region2d->proj = region3d->proj;
 
24
    region2d->zone = region3d->zone;
 
25
 
 
26
    region2d->north = region3d->north;
 
27
    region2d->south = region3d->south;
 
28
    region2d->east = region3d->east;
 
29
    region2d->west = region3d->west;
 
30
 
 
31
    region2d->rows = region3d->rows;
 
32
    region2d->cols = region3d->cols;
 
33
 
 
34
    region2d->ns_res = region3d->ns_res;
 
35
    region2d->ew_res = region3d->ew_res;
 
36
}
 
37
 
 
38
/*!
 
39
 * \brief
 
40
 *
 
41
 *  Returns in <em>region2d</em> the <em>2d</em> portion of <em>region3d</em>.
 
42
 *
 
43
 *  \param region3d
 
44
 *  \param region2d
 
45
 *  \return void
 
46
 */
 
47
 
 
48
void Rast3d_region_to_cell_head(RASTER3D_Region * region3d, struct Cell_head *region2d)
 
49
{
 
50
    region2d->proj = region3d->proj;
 
51
    region2d->zone = region3d->zone;
 
52
 
 
53
    region2d->north = region3d->north;
 
54
    region2d->south = region3d->south;
 
55
    region2d->east = region3d->east;
 
56
    region2d->west = region3d->west;
 
57
    region2d->top = region3d->top;
 
58
    region2d->bottom = region3d->bottom;
 
59
 
 
60
    region2d->rows = region3d->rows;
 
61
    region2d->rows3 = region3d->rows;
 
62
    region2d->cols = region3d->cols;
 
63
    region2d->cols3 = region3d->cols;
 
64
    region2d->depths = region3d->depths;
 
65
 
 
66
    region2d->ns_res = region3d->ns_res;
 
67
    region2d->ns_res3 = region3d->ns_res;
 
68
    region2d->ew_res = region3d->ew_res;
 
69
    region2d->ew_res3 = region3d->ew_res;
 
70
    region2d->tb_res = region3d->tb_res;
 
71
}
 
72
 
 
73
/*---------------------------------------------------------------------------*/
 
74
 
 
75
 
 
76
/*!
 
77
 * \brief
 
78
 *
 
79
 * Replaces the <em>2d</em> portion of <em>region3d</em> with the
 
80
 * values stored in <em>region2d</em>.
 
81
 *
 
82
 *  \param region2d
 
83
 *  \param region3d
 
84
 *  \return void
 
85
 */
 
86
 
 
87
void
 
88
Rast3d_incorporate2d_region(struct Cell_head *region2d, RASTER3D_Region * region3d)
 
89
{
 
90
    region3d->proj = region2d->proj;
 
91
    region3d->zone = region2d->zone;
 
92
 
 
93
    region3d->north = region2d->north;
 
94
    region3d->south = region2d->south;
 
95
    region3d->east = region2d->east;
 
96
    region3d->west = region2d->west;
 
97
 
 
98
    region3d->rows = region2d->rows;
 
99
    region3d->cols = region2d->cols;
 
100
 
 
101
    region3d->ns_res = region2d->ns_res;
 
102
    region3d->ew_res = region2d->ew_res;
 
103
}
 
104
 
 
105
/*!
 
106
 * \brief
 
107
 *
 
108
 * Replaces the <em>2d</em> portion of <em>region3d</em> with the
 
109
 * values stored in <em>region2d</em>.
 
110
 *
 
111
 *  \param region2d
 
112
 *  \param region3d
 
113
 *  \return void
 
114
 */
 
115
 
 
116
void
 
117
Rast3d_region_from_to_cell_head(struct Cell_head *region2d, RASTER3D_Region * region3d)
 
118
{
 
119
    region3d->proj = region2d->proj;
 
120
    region3d->zone = region2d->zone;
 
121
 
 
122
    region3d->north = region2d->north;
 
123
    region3d->south = region2d->south;
 
124
    region3d->east = region2d->east;
 
125
    region3d->west = region2d->west;
 
126
    region3d->top = region2d->top;
 
127
    region3d->bottom = region2d->bottom;
 
128
 
 
129
    region3d->rows = region2d->rows3;
 
130
    region3d->cols = region2d->cols3;
 
131
    region3d->depths = region2d->depths;
 
132
 
 
133
    region3d->ns_res = region2d->ns_res3;
 
134
    region3d->ew_res = region2d->ew_res3;
 
135
    region3d->tb_res = region2d->tb_res;
 
136
}
 
137
 
 
138
/*---------------------------------------------------------------------------*/
 
139
 
 
140
 
 
141
/*!
 
142
 * \brief
 
143
 *
 
144
 * Computes an adjusts the resolutions in the region structure from the region
 
145
 * boundaries and number of cells per dimension.
 
146
 *
 
147
 *  \param region
 
148
 *  \return void
 
149
 */
 
150
 
 
151
void Rast3d_adjust_region(RASTER3D_Region * region)
 
152
{
 
153
    struct Cell_head region2d;
 
154
 
 
155
    Rast3d_region_to_cell_head(region, &region2d);
 
156
    G_adjust_Cell_head3(&region2d, 1, 1, 1);
 
157
    Rast3d_region_from_to_cell_head(&region2d, region);
 
158
 
 
159
    if (region->depths <= 0)
 
160
    Rast3d_fatal_error("Rast3d_adjust_region: depths <= 0");
 
161
    region->tb_res = (region->top - region->bottom) / region->depths;
 
162
}
 
163
 
 
164
/*---------------------------------------------------------------------------*/
 
165
 
 
166
 
 
167
/*!
 
168
 * \brief
 
169
 *
 
170
 * Computes an adjusts the number of cells per dimension in the region
 
171
 * structure from the region boundaries and resolutions.
 
172
 *
 
173
 *  \param region
 
174
 *  \return void
 
175
 */
 
176
 
 
177
void Rast3d_adjust_region_res(RASTER3D_Region * region)
 
178
{
 
179
    struct Cell_head region2d;
 
180
 
 
181
    Rast3d_region_to_cell_head(region, &region2d);
 
182
    G_adjust_Cell_head3(&region2d, 1, 1, 1);
 
183
    Rast3d_region_from_to_cell_head(&region2d, region);
 
184
 
 
185
    if (region->tb_res <= 0)
 
186
    Rast3d_fatal_error("Rast3d_adjust_region_res: tb_res <= 0");
 
187
 
 
188
    region->depths = (region->top - region->bottom + region->tb_res / 2.0) /
 
189
    region->tb_res;
 
190
    if (region->depths == 0)
 
191
    region->depths = 1;
 
192
}
 
193
 
 
194
/*---------------------------------------------------------------------------*/
 
195
 
 
196
 
 
197
/*!
 
198
 * \brief
 
199
 *
 
200
 * Copies the values of <em>regionSrc</em> into <em>regionDst</em>.
 
201
 *
 
202
 *  \param regionDest
 
203
 *  \param regionSrc
 
204
 *  \return void
 
205
 */
 
206
 
 
207
void Rast3d_region_copy(RASTER3D_Region * regionDest, RASTER3D_Region * regionSrc)
 
208
{
 
209
    regionDest->proj = regionSrc->proj;
 
210
    regionDest->zone = regionSrc->zone;
 
211
 
 
212
    regionDest->north = regionSrc->north;
 
213
    regionDest->south = regionSrc->south;
 
214
    regionDest->east = regionSrc->east;
 
215
    regionDest->west = regionSrc->west;
 
216
    regionDest->top = regionSrc->top;
 
217
    regionDest->bottom = regionSrc->bottom;
 
218
 
 
219
    regionDest->rows = regionSrc->rows;
 
220
    regionDest->cols = regionSrc->cols;
 
221
    regionDest->depths = regionSrc->depths;
 
222
 
 
223
    regionDest->ns_res = regionSrc->ns_res;
 
224
    regionDest->ew_res = regionSrc->ew_res;
 
225
    regionDest->tb_res = regionSrc->tb_res;
 
226
}
 
227
 
 
228
 
 
229
/*---------------------------------------------------------------------------*/
 
230
 
 
231
int
 
232
Rast3d_read_region_map(const char *name, const char *mapset, RASTER3D_Region * region)
 
233
{
 
234
    char fullName[GPATH_MAX];
 
235
    char xname[GNAME_MAX], xmapset[GMAPSET_MAX];
 
236
 
 
237
    if (G_name_is_fully_qualified(name, xname, xmapset))
 
238
    Rast3d_filename(fullName, RASTER3D_HEADER_ELEMENT, xname, xmapset);
 
239
    else {
 
240
    if (!mapset || !*mapset)
 
241
        mapset = G_find_raster3d(name, "");
 
242
    Rast3d_filename(fullName, RASTER3D_HEADER_ELEMENT, name, mapset);
 
243
    }
 
244
    return Rast3d_read_window(region, fullName);
 
245
}
 
246
 
 
247
/*---------------------------------------------------------------------------*/
 
248
 
 
249
 
 
250
/*!
 
251
 * \brief
 
252
 *
 
253
 *  Returns 1 if region-coordinates <em>(north, east, top)</em> are
 
254
 * inside the region of <em>map</em>. Returns 0 otherwise.
 
255
 *
 
256
 *  \param region
 
257
 *  \param north
 
258
 *  \param east
 
259
 *  \param top
 
260
 *  \return int
 
261
 */
 
262
 
 
263
int Rast3d_is_valid_location(RASTER3D_Region *region, double north, double east, double top)
 
264
{
 
265
    return ((north >= region->south) && (north <= region->north) &&
 
266
        (east >= region->west) && (east <= region->east) &&
 
267
        (((top >= region->bottom) && (top <= region->top)) ||
 
268
         ((top <= region->bottom) && (top >= region->top))));
 
269
}
 
270
 
 
271
/*---------------------------------------------------------------------------*/
 
272
 
 
273
/*!
 
274
 * \brief
 
275
 *
 
276
 *  Converts region-coordinates <em>(north, east,
 
277
 *  top)</em> into cell-coordinates <em>(x, y, z)</em>.
 
278
 *
 
279
 *  \param region
 
280
 *  \param north
 
281
 *  \param east
 
282
 *  \param top
 
283
 *  \param x
 
284
 *  \param y
 
285
 *  \param z
 
286
 *  \return void
 
287
 */
 
288
 
 
289
void
 
290
Rast3d_location2coord(RASTER3D_Region *region, double north, double east, double top,
 
291
           int *x, int *y, int *z)
 
292
{
 
293
    double col, row, depth;
 
294
    LOCATION_TO_COORD(region, north, east, top, &col, &row, &depth);
 
295
 
 
296
    *x = (int)floor(col);
 
297
    *y = (int)floor(row);
 
298
    *z = (int)floor(depth);
 
299
}
 
300
 
 
301
/*!
 
302
 * \brief
 
303
 *
 
304
 *  Converts region-coordinates <em>(north, east,
 
305
 *  top)</em> into cell-coordinates <em>(x, y, z)</em>.
 
306
 *
 
307
 * <b>Note:</b> The results are <i>double</i> numbers. Casting them to
 
308
 * <i>int</i> will give the column, row and depth number.
 
309
 *
 
310
 *  \param region
 
311
 *  \param north
 
312
 *  \param east
 
313
 *  \param top
 
314
 *  \param x
 
315
 *  \param y
 
316
 *  \param z
 
317
 *  \return void
 
318
 */
 
319
 
 
320
void
 
321
Rast3d_location2coord_double(RASTER3D_Region *region, double north, double east, double top,
 
322
           double *x, double *y, double *z)
 
323
{
 
324
    LOCATION_TO_COORD(region, north, east, top, x, y, z);
 
325
 
 
326
    G_debug(4, "Rast3d_location2coord_double x %f y %f z %f\n", *x, *y, *z);
 
327
}
 
328
 
 
329
/*!
 
330
 * \brief
 
331
 *
 
332
 *  Converts region-coordinates <em>(north, east,
 
333
 *  top)</em> into cell-coordinates <em>(x, y, z)</em>.
 
334
 *  This function calls Rast3d_fatal_error in case location is not in window.
 
335
 *
 
336
 *  \param region
 
337
 *  \param north
 
338
 *  \param east
 
339
 *  \param top
 
340
 *  \param x
 
341
 *  \param y
 
342
 *  \param z
 
343
 *  \return void
 
344
 */
 
345
 
 
346
void
 
347
Rast3d_location2coord2(RASTER3D_Region *region, double north, double east, double top,
 
348
           int *x, int *y, int *z)
 
349
{
 
350
    if (!Rast3d_is_valid_location(region, north, east, top))
 
351
    Rast3d_fatal_error("Rast3d_location2coord2: location not in region");
 
352
 
 
353
    double col, row, depth;
 
354
    LOCATION_TO_COORD(region, north, east, top, &col, &row, &depth);
 
355
 
 
356
    *x = (int)floor(col);
 
357
    *y = (int)floor(row);
 
358
    *z = (int)floor(depth);
 
359
}
 
360
 
 
361
/*!
 
362
 * \brief
 
363
 *
 
364
 *  Converts cell-coordinates <em>(x, y, z)</em> into region-coordinates
 
365
 * <em>(north, east, top)</em>.
 
366
 *
 
367
 *  * <b>Note:</b> x, y and z is a double:
 
368
 *  - x+0.0 will return the easting for the western edge of the column.
 
369
 *  - x+0.5 will return the easting for the center of the column.
 
370
 *  - x+1.0 will return the easting for the eastern edge of the column.
 
371
 *
 
372
 *  - y+0.0 will return the northing for the northern edge of the row.
 
373
 *  - y+0.5 will return the northing for the center of the row.
 
374
 *  - y+1.0 will return the northing for the southern edge of the row.
 
375
 *
 
376
 *  - z+0.0 will return the top for the lower edge of the depth.
 
377
 *  - z+0.5 will return the top for the center of the depth.
 
378
 *  - z+1.0 will return the top for the upper edge of the column.
 
379
 *
 
380
*
 
381
 *  \param region
 
382
 *  \param x
 
383
 *  \param y
 
384
 *  \param z
 
385
 *  \param north
 
386
 *  \param east
 
387
 *  \param top
 
388
 *  \return void
 
389
 */
 
390
 
 
391
void
 
392
Rast3d_coord2location(RASTER3D_Region * region, double x, double y, double z, double *north, double *east, double *top)
 
393
{
 
394
    COORD_TO_LOCATION(region, x, y, z, north, east, top);
 
395
 
 
396
    G_debug(4, "Rast3d_coord2location north %g east %g top %g\n", *north, *east, *top);
 
397
}