GRASS Programmer's Manual  6.5.svn(2014)-r66266
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
adj_cellhd.c
Go to the documentation of this file.
1 
17 #include <grass/gis.h>
18 #include <grass/glocale.h>
19 
20 
43 char *G_adjust_Cell_head(struct Cell_head *cellhd, int row_flag, int col_flag)
44 {
45  if (!row_flag) {
46  if (cellhd->ns_res <= 0)
47  return (_("Illegal n-s resolution value"));
48  }
49  else {
50  if (cellhd->rows <= 0)
51  return (_("Illegal row value"));
52  }
53  if (!col_flag) {
54  if (cellhd->ew_res <= 0)
55  return (_("Illegal e-w resolution value"));
56  }
57  else {
58  if (cellhd->cols <= 0)
59  return (_("Illegal col value"));
60  }
61 
62  /* for lat/lon, check north,south. force east larger than west */
63  if (cellhd->proj == PROJECTION_LL) {
64  double epsilon_ns, epsilon_ew;
65 
66  /* TODO: find good thresholds */
67  epsilon_ns = 1. / cellhd->rows * 0.001;
68  epsilon_ew = .000001; /* epsilon_ew calculation doesn't work due to cellhd->cols update/global wraparound below */
69 
70  G_debug(3, "G_adjust_Cell_head: epsilon_ns: %g, epsilon_ew: %g",
71  epsilon_ns, epsilon_ew);
72 
73  /* TODO: once working, change below G_warning to G_debug */
74 
75  /* fix rounding problems if input map slightly exceeds the world definition -180 90 180 -90 */
76  if (cellhd->north > 90.0) {
77  if (((cellhd->north - 90.0) < epsilon_ns) &&
78  ((cellhd->north - 90.0) > GRASS_EPSILON)) {
79  G_warning(_("Fixing subtle input data rounding error of north boundary (%g>%g)"),
80  cellhd->north - 90.0, epsilon_ns);
81  cellhd->north = 90.0;
82  }
83  else
84  return (_("Illegal latitude for North"));
85  }
86 
87  if (cellhd->south < -90.0) {
88  if (((cellhd->south + 90.0) < epsilon_ns) &&
89  ((cellhd->south + 90.0) < GRASS_EPSILON)) {
90  G_warning(_("Fixing subtle input data rounding error of south boundary (%g>%g)"),
91  cellhd->south + 90.0, epsilon_ns);
92  cellhd->south = -90.0;
93  }
94  else
95  return (_("Illegal latitude for South"));
96  }
97 
98 #if 0
99  /* DISABLED: this breaks global wrap-around */
100 
101  G_debug(3,
102  "G_adjust_Cell_head() cellhd->west: %f, devi: %g, eps: %g",
103  cellhd->west, cellhd->west + 180.0, epsilon_ew);
104 
105  if ((cellhd->west < -180.0) && ((cellhd->west + 180.0) < epsilon_ew)
106  && ((cellhd->west + 180.0) < GRASS_EPSILON)) {
107  G_warning(_("Fixing subtle input data rounding error of west boundary (%g>%g)"),
108  cellhd->west + 180.0, epsilon_ew);
109  cellhd->west = -180.0;
110  }
111 
112  G_debug(3,
113  "G_adjust_Cell_head() cellhd->east: %f, devi: %g, eps: %g",
114  cellhd->east, cellhd->east - 180.0, epsilon_ew);
115 
116  if ((cellhd->east > 180.0) && ((cellhd->east - 180.0) > epsilon_ew)
117  && ((cellhd->east - 180.0) > GRASS_EPSILON)) {
118  G_warning(_("Fixing subtle input data rounding error of east boundary (%g>%g)"),
119  cellhd->east - 180.0, epsilon_ew);
120  cellhd->east = 180.0;
121  }
122 #endif
123 
124  while (cellhd->east <= cellhd->west)
125  cellhd->east += 360.0;
126  }
127 
128  /* check the edge values */
129  if (cellhd->north <= cellhd->south) {
130  if (cellhd->proj == PROJECTION_LL)
131  return (_("North must be north of South"));
132  else
133  return (_("North must be larger than South"));
134  }
135  if (cellhd->east <= cellhd->west)
136  return (_("East must be larger than West"));
137 
138  /* compute rows and columns, if not set */
139  if (!row_flag) {
140  cellhd->rows =
141  (cellhd->north - cellhd->south +
142  cellhd->ns_res / 2.0) / cellhd->ns_res;
143  if (cellhd->rows == 0)
144  cellhd->rows = 1;
145  }
146  if (!col_flag) {
147  cellhd->cols =
148  (cellhd->east - cellhd->west +
149  cellhd->ew_res / 2.0) / cellhd->ew_res;
150  if (cellhd->cols == 0)
151  cellhd->cols = 1;
152  }
153 
154  if (cellhd->cols < 0 || cellhd->rows < 0) {
155  return (_("Invalid coordinates"));
156  }
157 
158 
159  /* (re)compute the resolutions */
160  cellhd->ns_res = (cellhd->north - cellhd->south) / cellhd->rows;
161  cellhd->ew_res = (cellhd->east - cellhd->west) / cellhd->cols;
162 
163  return NULL;
164 }
165 
166 
194 char *G_adjust_Cell_head3(struct Cell_head *cellhd, int row_flag,
195  int col_flag, int depth_flag)
196 {
197  if (!row_flag) {
198  if (cellhd->ns_res <= 0)
199  return (_("Illegal n-s resolution value"));
200  if (cellhd->ns_res3 <= 0)
201  return (_("Illegal n-s3 resolution value"));
202  }
203  else {
204  if (cellhd->rows <= 0)
205  return (_("Illegal row value"));
206  if (cellhd->rows3 <= 0)
207  return (_("Illegal row3 value"));
208  }
209  if (!col_flag) {
210  if (cellhd->ew_res <= 0)
211  return (_("Illegal e-w resolution value"));
212  if (cellhd->ew_res3 <= 0)
213  return (_("Illegal e-w3 resolution value"));
214  }
215  else {
216  if (cellhd->cols <= 0)
217  return (_("Illegal col value"));
218  if (cellhd->cols3 <= 0)
219  return (_("Illegal col3 value"));
220  }
221  if (!depth_flag) {
222  if (cellhd->tb_res <= 0)
223  return (_("Illegal t-b3 resolution value"));
224  }
225  else {
226  if (cellhd->depths <= 0)
227  return (_("Illegal depths value"));
228  }
229 
230  /* for lat/lon, check north,south. force east larger than west */
231  if (cellhd->proj == PROJECTION_LL) {
232  double epsilon_ns, epsilon_ew;
233 
234  /* TODO: find good thresholds */
235  epsilon_ns = 1. / cellhd->rows * 0.001;
236  epsilon_ew = .000001; /* epsilon_ew calculation doesn't work due to cellhd->cols update/global wraparound below */
237 
238  G_debug(3, "G_adjust_Cell_head: epsilon_ns: %g, epsilon_ew: %g",
239  epsilon_ns, epsilon_ew);
240 
241  /* TODO: once working, change below G_warning to G_debug */
242 
243  /* fix rounding problems if input map slightly exceeds the world definition -180 90 180 -90 */
244  if (cellhd->north > 90.0) {
245  if (((cellhd->north - 90.0) < epsilon_ns) &&
246  ((cellhd->north - 90.0) > GRASS_EPSILON)) {
247  G_warning(_("Fixing subtle input data rounding error of north boundary (%g>%g)"),
248  cellhd->north - 90.0, epsilon_ns);
249  cellhd->north = 90.0;
250  }
251  else
252  return (_("Illegal latitude for North"));
253  }
254 
255  if (cellhd->south < -90.0) {
256  if (((cellhd->south + 90.0) < epsilon_ns) &&
257  ((cellhd->south + 90.0) < GRASS_EPSILON)) {
258  G_warning(_("Fixing subtle input data rounding error of south boundary (%g>%g)"),
259  cellhd->south + 90.0, epsilon_ns);
260  cellhd->south = -90.0;
261  }
262  else
263  return (_("Illegal latitude for South"));
264  }
265 
266 #if 0
267  /* DISABLED: this breaks global wrap-around */
268 
269  G_debug(3,
270  "G_adjust_Cell_head3() cellhd->west: %f, devi: %g, eps: %g",
271  cellhd->west, cellhd->west + 180.0, epsilon_ew);
272 
273  if ((cellhd->west < -180.0) && ((cellhd->west + 180.0) < epsilon_ew)
274  && ((cellhd->west + 180.0) < GRASS_EPSILON)) {
275  G_warning(_("Fixing subtle input data rounding error of west boundary (%g>%g)"),
276  cellhd->west + 180.0, epsilon_ew);
277  cellhd->west = -180.0;
278  }
279 
280  G_debug(3,
281  "G_adjust_Cell_head3() cellhd->east: %f, devi: %g, eps: %g",
282  cellhd->east, cellhd->east - 180.0, epsilon_ew);
283 
284  if ((cellhd->east > 180.0) && ((cellhd->east - 180.0) > epsilon_ew)
285  && ((cellhd->east - 180.0) > GRASS_EPSILON)) {
286  G_warning(_("Fixing subtle input data rounding error of east boundary (%g>%g)"),
287  cellhd->east - 180.0, epsilon_ew);
288  cellhd->east = 180.0;
289  }
290 #endif
291 
292  while (cellhd->east <= cellhd->west)
293  cellhd->east += 360.0;
294  }
295 
296  /* check the edge values */
297  if (cellhd->north <= cellhd->south) {
298  if (cellhd->proj == PROJECTION_LL)
299  return (_("North must be north of South"));
300  else
301  return (_("North must be larger than South"));
302  }
303  if (cellhd->east <= cellhd->west)
304  return (_("East must be larger than West"));
305  if (cellhd->top <= cellhd->bottom)
306  return (_("Top must be larger than Bottom"));
307 
308 
309  /* compute rows and columns, if not set */
310  if (!row_flag) {
311  cellhd->rows =
312  (cellhd->north - cellhd->south +
313  cellhd->ns_res / 2.0) / cellhd->ns_res;
314  if (cellhd->rows == 0)
315  cellhd->rows = 1;
316 
317  cellhd->rows3 =
318  (cellhd->north - cellhd->south +
319  cellhd->ns_res3 / 2.0) / cellhd->ns_res3;
320  if (cellhd->rows3 == 0)
321  cellhd->rows3 = 1;
322  }
323  if (!col_flag) {
324  cellhd->cols =
325  (cellhd->east - cellhd->west +
326  cellhd->ew_res / 2.0) / cellhd->ew_res;
327  if (cellhd->cols == 0)
328  cellhd->cols = 1;
329 
330  cellhd->cols3 =
331  (cellhd->east - cellhd->west +
332  cellhd->ew_res3 / 2.0) / cellhd->ew_res3;
333  if (cellhd->cols3 == 0)
334  cellhd->cols3 = 1;
335  }
336 
337  if (!depth_flag) {
338  cellhd->depths =
339  (cellhd->top - cellhd->bottom +
340  cellhd->tb_res / 2.0) / cellhd->tb_res;
341  if (cellhd->depths == 0)
342  cellhd->depths = 1;
343 
344  }
345 
346  if (cellhd->cols < 0 || cellhd->rows < 0 || cellhd->cols3 < 0 ||
347  cellhd->rows3 < 0 || cellhd->depths < 0) {
348  return (_("Invalid coordinates"));
349  }
350 
351  /* (re)compute the resolutions */
352  cellhd->ns_res = (cellhd->north - cellhd->south) / cellhd->rows;
353  cellhd->ns_res3 = (cellhd->north - cellhd->south) / cellhd->rows3;
354  cellhd->ew_res = (cellhd->east - cellhd->west) / cellhd->cols;
355  cellhd->ew_res3 = (cellhd->east - cellhd->west) / cellhd->cols3;
356  cellhd->tb_res = (cellhd->top - cellhd->bottom) / cellhd->depths;
357 
358  return NULL;
359 }
G_adjust_Cell_head3
char * G_adjust_Cell_head3(struct Cell_head *cellhd, int row_flag, int col_flag, int depth_flag)
Adjust cell header for 3D values.
Definition: adj_cellhd.c:194
G_adjust_Cell_head
char * G_adjust_Cell_head(struct Cell_head *cellhd, int row_flag, int col_flag)
Adjust cell header.
Definition: adj_cellhd.c:43
NULL
return NULL
Definition: dbfopen.c:1394
G_warning
G_warning("category support for [%s] in mapset [%s] %s", name, mapset, type)
G_debug
int G_debug(int level, const char *msg,...)
Print debugging message.
Definition: gis/debug.c:51