GRASS GIS 8 Programmer's Manual  8.5.0dev(2024)-36359e2344
adj_cellhd.c
Go to the documentation of this file.
1 /*!
2  * \file lib/gis/adj_cellhd.c
3  *
4  * \brief GIS Library - CELL header adjustment.
5  *
6  * (C) 2001-2009 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 Original author CERL
12  */
13 
14 #include <math.h>
15 #include <string.h>
16 #include <grass/gis.h>
17 #include <grass/glocale.h>
18 
19 #define LL_TOLERANCE 10
20 
21 /* TODO: find good thresholds */
22 /* deviation measured in cells */
23 static double llepsilon = 0.01;
24 static double fpepsilon = 1.0e-9;
25 
26 static int ll_wrap(struct Cell_head *cellhd);
27 static int ll_check_ns(struct Cell_head *cellhd);
28 static int ll_check_ew(struct Cell_head *cellhd);
29 
30 /*!
31  * \brief Adjust cell header.
32  *
33  * This function fills in missing parts of the input cell header (or
34  * region). It also makes projection-specific adjustments. The
35  * <i>cellhd</i> structure must have its <i>north, south, east,
36  * west</i>, and <i>proj</i> fields set.
37  *
38  * If <i>row_flag</i> is true, then the north-south resolution is
39  * computed from the number of <i>rows</i> in the <i>cellhd</i>
40  * structure. Otherwise the number of <i>rows</i> is computed from the
41  * north-south resolution in the structure, similarly for
42  * <i>col_flag</i> and the number of columns and the east-west
43  * resolution.
44  *
45  * <b>Note:</b> 3D values are not adjusted.
46  *
47  * \param[in,out] cellhd pointer to Cell_head structure
48  * \param row_flag compute n-s resolution
49  * \param col_flag compute e-w resolution
50  */
51 void G_adjust_Cell_head(struct Cell_head *cellhd, int row_flag, int col_flag)
52 {
53  double old_res;
54 
55  if (!row_flag) {
56  if (cellhd->ns_res <= 0)
57  G_fatal_error(_("Illegal n-s resolution value: %g"),
58  cellhd->ns_res);
59  }
60  else {
61  if (cellhd->rows <= 0)
62  G_fatal_error(_("Illegal number of rows: %d"
63  " (resolution is %g)"),
64  cellhd->rows, cellhd->ns_res);
65  }
66  if (!col_flag) {
67  if (cellhd->ew_res <= 0)
68  G_fatal_error(_("Illegal e-w resolution value: %g"),
69  cellhd->ew_res);
70  }
71  else {
72  if (cellhd->cols <= 0)
73  G_fatal_error(_("Illegal number of columns: %d"
74  " (resolution is %g)"),
75  cellhd->cols, cellhd->ew_res);
76  }
77 
78  /* check the edge values */
79  if (cellhd->north <= cellhd->south) {
80  if (cellhd->proj == PROJECTION_LL)
81  G_fatal_error(_("North must be north of South,"
82  " but %g (north) <= %g (south"),
83  cellhd->north, cellhd->south);
84  else
85  G_fatal_error(_("North must be larger than South,"
86  " but %g (north) <= %g (south"),
87  cellhd->north, cellhd->south);
88  }
89 
90  ll_wrap(cellhd);
91 
92  if (cellhd->east <= cellhd->west)
93  G_fatal_error(_("East must be larger than West,"
94  " but %g (east) <= %g (west)"),
95  cellhd->east, cellhd->west);
96 
97  /* compute rows and columns, if not set */
98  if (!row_flag) {
99  cellhd->rows = (cellhd->north - cellhd->south + cellhd->ns_res / 2.0) /
100  cellhd->ns_res;
101  if (cellhd->rows == 0)
102  cellhd->rows = 1;
103  }
104  if (!col_flag) {
105  cellhd->cols = (cellhd->east - cellhd->west + cellhd->ew_res / 2.0) /
106  cellhd->ew_res;
107  if (cellhd->cols == 0)
108  cellhd->cols = 1;
109  }
110 
111  if (cellhd->cols < 0) {
112  G_fatal_error(_("Invalid coordinates: negative number of columns"));
113  }
114  if (cellhd->rows < 0) {
115  G_fatal_error(_("Invalid coordinates: negative number of rows"));
116  }
117 
118  /* (re)compute the resolutions */
119  old_res = cellhd->ns_res;
120  cellhd->ns_res = (cellhd->north - cellhd->south) / cellhd->rows;
121  if (old_res > 0 && fabs(old_res - cellhd->ns_res) / old_res > 0.01)
122  G_verbose_message(_("NS resolution has been changed"));
123 
124  old_res = cellhd->ew_res;
125  cellhd->ew_res = (cellhd->east - cellhd->west) / cellhd->cols;
126  if (old_res > 0 && fabs(old_res - cellhd->ew_res) / old_res > 0.01)
127  G_verbose_message(_("EW resolution has been changed"));
128 
129  if (fabs(cellhd->ns_res - cellhd->ew_res) / cellhd->ns_res > 0.01)
130  G_verbose_message(_("NS and EW resolutions are different"));
131 
132  ll_check_ns(cellhd);
133  ll_check_ew(cellhd);
134 }
135 
136 /*!
137  * \brief Adjust cell header for 3D values.
138  *
139  * This function fills in missing parts of the input cell header (or
140  * region). It also makes projection-specific adjustments. The
141  * <i>cellhd</i> structure must have its <i>north, south, east,
142  * west</i>, and <i>proj</i> fields set.
143  *
144  * If <i>row_flag</i> is true, then the north-south resolution is computed
145  * from the number of <i>rows</i> in the <i>cellhd</i> structure.
146  * Otherwise the number of <i>rows</i> is computed from the north-south
147  * resolution in the structure, similarly for <i>col_flag</i> and the
148  * number of columns and the east-west resolution.
149  *
150  * If <i>depth_flag</i> is true, top-bottom resolution is calculated
151  * from depths.
152  * If <i>depth_flag</i> are false, number of depths is calculated from
153  * top-bottom resolution.
154  *
155  * \warning This function can cause segmentation fault without any warning
156  * when it is called with Cell_head top and bottom set to zero.
157  *
158  * \param[in,out] cellhd pointer to Cell_head structure
159  * \param row_flag compute n-s resolution
160  * \param col_flag compute e-w resolution
161  * \param depth_flag compute t-b resolution
162  */
163 void G_adjust_Cell_head3(struct Cell_head *cellhd, int row_flag, int col_flag,
164  int depth_flag)
165 {
166  double old_res;
167 
168  if (!row_flag) {
169  if (cellhd->ns_res <= 0)
170  G_fatal_error(_("Illegal n-s resolution value: %g"),
171  cellhd->ns_res);
172  if (cellhd->ns_res3 <= 0)
173  G_fatal_error(_("Illegal n-s resolution value for 3D: %g"),
174  cellhd->ns_res3);
175  }
176  else {
177  if (cellhd->rows <= 0)
178  G_fatal_error(_("Illegal number of rows: %d"
179  " (resolution is %g)"),
180  cellhd->rows, cellhd->ns_res);
181  if (cellhd->rows3 <= 0)
182  G_fatal_error(_("Illegal number of rows for 3D: %d"
183  " (resolution is %g)"),
184  cellhd->rows3, cellhd->ns_res3);
185  }
186  if (!col_flag) {
187  if (cellhd->ew_res <= 0)
188  G_fatal_error(_("Illegal e-w resolution value: %g"),
189  cellhd->ew_res);
190  if (cellhd->ew_res3 <= 0)
191  G_fatal_error(_("Illegal e-w resolution value for 3D: %g"),
192  cellhd->ew_res3);
193  }
194  else {
195  if (cellhd->cols <= 0)
196  G_fatal_error(_("Illegal number of columns: %d"
197  " (resolution is %g)"),
198  cellhd->cols, cellhd->ew_res);
199  if (cellhd->cols3 <= 0)
200  G_fatal_error(_("Illegal number of columns for 3D: %d"
201  " (resolution is %g)"),
202  cellhd->cols3, cellhd->ew_res3);
203  }
204  if (!depth_flag) {
205  if (cellhd->tb_res <= 0)
206  G_fatal_error(_("Illegal t-b resolution value: %g"),
207  cellhd->tb_res);
208  }
209  else {
210  if (cellhd->depths <= 0)
211  G_fatal_error(_("Illegal depths value: %d"), cellhd->depths);
212  }
213 
214  /* check the edge values */
215  if (cellhd->north <= cellhd->south) {
216  if (cellhd->proj == PROJECTION_LL)
217  G_fatal_error(_("North must be north of South,"
218  " but %g (north) <= %g (south"),
219  cellhd->north, cellhd->south);
220  else
221  G_fatal_error(_("North must be larger than South,"
222  " but %g (north) <= %g (south"),
223  cellhd->north, cellhd->south);
224  }
225 
226  ll_wrap(cellhd);
227 
228  if (cellhd->east <= cellhd->west)
229  G_fatal_error(_("East must be larger than West,"
230  " but %g (east) <= %g (west)"),
231  cellhd->east, cellhd->west);
232 
233  if (cellhd->top <= cellhd->bottom)
234  G_fatal_error(_("Top must be larger than Bottom,"
235  " but %g (top) <= %g (bottom)"),
236  cellhd->top, cellhd->bottom);
237 
238  /* compute rows and columns, if not set */
239  if (!row_flag) {
240  cellhd->rows = (cellhd->north - cellhd->south + cellhd->ns_res / 2.0) /
241  cellhd->ns_res;
242  if (cellhd->rows == 0)
243  cellhd->rows = 1;
244 
245  cellhd->rows3 =
246  (cellhd->north - cellhd->south + cellhd->ns_res3 / 2.0) /
247  cellhd->ns_res3;
248  if (cellhd->rows3 == 0)
249  cellhd->rows3 = 1;
250  }
251  if (!col_flag) {
252  cellhd->cols = (cellhd->east - cellhd->west + cellhd->ew_res / 2.0) /
253  cellhd->ew_res;
254  if (cellhd->cols == 0)
255  cellhd->cols = 1;
256 
257  cellhd->cols3 = (cellhd->east - cellhd->west + cellhd->ew_res3 / 2.0) /
258  cellhd->ew_res3;
259  if (cellhd->cols3 == 0)
260  cellhd->cols3 = 1;
261  }
262 
263  if (!depth_flag) {
264  cellhd->depths = (cellhd->top - cellhd->bottom + cellhd->tb_res / 2.0) /
265  cellhd->tb_res;
266  if (cellhd->depths == 0)
267  cellhd->depths = 1;
268  }
269 
270  if (cellhd->cols < 0 || cellhd->cols3 < 0) {
271  G_fatal_error(_("Invalid coordinates: negative number of columns"));
272  }
273  if (cellhd->rows < 0 || cellhd->rows3 < 0) {
274  G_fatal_error(_("Invalid coordinates: negative number of rows"));
275  }
276  if (cellhd->depths < 0) {
277  G_fatal_error(_("Invalid coordinates: negative number of depths"));
278  }
279 
280  /* (re)compute the resolutions */
281  old_res = cellhd->ns_res;
282  cellhd->ns_res = (cellhd->north - cellhd->south) / cellhd->rows;
283  if (old_res > 0 && fabs(old_res - cellhd->ns_res) / old_res > 0.01)
284  G_verbose_message(_("NS resolution has been changed"));
285 
286  old_res = cellhd->ew_res;
287  cellhd->ew_res = (cellhd->east - cellhd->west) / cellhd->cols;
288  if (old_res > 0 && fabs(old_res - cellhd->ew_res) / old_res > 0.01)
289  G_verbose_message(_("EW resolution has been changed"));
290 
291  if (fabs(cellhd->ns_res - cellhd->ew_res) / cellhd->ns_res > 0.01)
292  G_verbose_message(_("NS and EW resolutions are different"));
293 
294  ll_check_ns(cellhd);
295  ll_check_ew(cellhd);
296 
297  cellhd->ns_res3 = (cellhd->north - cellhd->south) / cellhd->rows3;
298  cellhd->ew_res3 = (cellhd->east - cellhd->west) / cellhd->cols3;
299  cellhd->tb_res = (cellhd->top - cellhd->bottom) / cellhd->depths;
300 }
301 
302 static int ll_wrap(struct Cell_head *cellhd)
303 {
304  double shift;
305 
306  /* for lat/lon, force east larger than west, try to wrap to -180, 180 */
307  if (cellhd->proj != PROJECTION_LL)
308  return 0;
309 
310  if (cellhd->east <= cellhd->west) {
311  G_warning(_("East (%.15g) is not larger than West (%.15g)"),
312  cellhd->east, cellhd->west);
313 
314  while (cellhd->east <= cellhd->west)
315  cellhd->east += 360.0;
316  }
317 
318  /* with east larger than west,
319  * any 360 degree W-E extent can be represented within -360, 360
320  * but not within -180, 180 */
321 
322  /* try to shift to within -180, 180 */
323  shift = 0;
324  while (cellhd->west + shift >= 180) {
325  shift -= 360.0;
326  }
327  while (cellhd->east + shift <= -180) {
328  shift += 360.0;
329  }
330 
331  /* try to shift to within -360, 360 */
332  while (cellhd->east + shift > 360) {
333  shift -= 360.0;
334  }
335  while (cellhd->west + shift <= -360) {
336  shift += 360.0;
337  }
338 
339  if (shift) {
340  cellhd->west += shift;
341  cellhd->east += shift;
342  }
343 
344  /* very liberal thresholds */
345  if (cellhd->north > 90.0 + LL_TOLERANCE)
346  G_fatal_error(_("Illegal latitude for North: %g"), cellhd->north);
347  if (cellhd->south < -90.0 - LL_TOLERANCE)
348  G_fatal_error(_("Illegal latitude for South: %g"), cellhd->south);
349 
350 #if 0
351  /* disabled: allow W-E extents larger than 360 degree e.g. for display */
352  if (cellhd->west < -360.0 - LL_TOLERANCE) {
353  G_debug(1, "East: %g", cellhd->east);
354  G_fatal_error(_("Illegal longitude for West: %g"), cellhd->west);
355  }
356  if (cellhd->east > 360.0 + LL_TOLERANCE) {
357  G_debug(1, "West: %g", cellhd->west);
358  G_fatal_error(_("Illegal longitude for East: %g"), cellhd->east);
359  }
360 #endif
361 
362  return 1;
363 }
364 
365 static int ll_check_ns(struct Cell_head *cellhd)
366 {
367  int lladjust;
368  double diff;
369  int ncells;
370 
371  /* lat/lon checks */
372  if (cellhd->proj != PROJECTION_LL)
373  return 0;
374 
375  lladjust = 0;
376 
377  G_debug(3, "ll_check_ns: epsilon: %g", llepsilon);
378 
379  /* North, South: allow a half cell spill-over */
380 
381  diff = (cellhd->north - cellhd->south) / cellhd->ns_res;
382  ncells = (int)(diff + 0.5);
383  diff -= ncells;
384  if ((diff < 0 && diff < -fpepsilon) || (diff > 0 && diff > fpepsilon)) {
386  _("NS extent does not match NS resolution: %g cells difference"),
387  diff);
388  }
389 
390  /* north */
391  diff = (cellhd->north - 90) / cellhd->ns_res;
392  if (diff < 0)
393  diff = -diff;
394  if (cellhd->north < 90.0 && diff < 1.0) {
395  G_verbose_message(_("%g cells missing to reach 90 degree north"), diff);
396  if (diff < llepsilon && diff > fpepsilon) {
398  _("Subtle input data rounding error of north boundary (%g)"),
399  cellhd->north - 90.0);
400  /* check only, do not modify
401  cellhd->north = 90.0;
402  lladjust = 1;
403  */
404  }
405  }
406  if (cellhd->north > 90.0) {
407  if (diff <= 0.5 + llepsilon) {
408  G_important_message(_("90 degree north is exceeded by %g cells"),
409  diff);
410 
411  if (diff < llepsilon && diff > fpepsilon) {
412  G_verbose_message(_("Subtle input data rounding error of north "
413  "boundary (%g)"),
414  cellhd->north - 90.0);
415  G_debug(1, "North of north in seconds: %g",
416  (cellhd->north - 90.0) * 3600);
417  /* check only, do not modify
418  cellhd->north = 90.0;
419  lladjust = 1;
420  */
421  }
422 
423  diff = diff - 0.5;
424  if (diff < 0)
425  diff = -diff;
426  if (diff < llepsilon && diff > fpepsilon) {
427  G_verbose_message(_("Subtle input data rounding error of north "
428  "boundary (%g)"),
429  cellhd->north - 90.0 - cellhd->ns_res / 2.0);
430  G_debug(1, "North of north + 0.5 cells in seconds: %g",
431  (cellhd->north - 90.0 - cellhd->ns_res / 2.0) * 3600);
432  /* check only, do not modify
433  cellhd->north = 90.0 + cellhd->ns_res / 2.0;
434  lladjust = 1;
435  */
436  }
437  }
438  else
439  G_fatal_error(_("Illegal latitude for North: %g"), cellhd->north);
440  }
441 
442  /* south */
443  diff = (cellhd->south + 90) / cellhd->ns_res;
444  if (diff < 0)
445  diff = -diff;
446  if (cellhd->south > -90.0 && diff < 1.0) {
447  G_verbose_message(_("%g cells missing to reach 90 degree south"), diff);
448  if (diff < llepsilon && diff > fpepsilon) {
450  _("Subtle input data rounding error of south boundary (%g)"),
451  cellhd->south + 90.0);
452  /* check only, do not modify
453  cellhd->south = -90.0;
454  lladjust = 1;
455  */
456  }
457  }
458  if (cellhd->south < -90.0) {
459  if (diff <= 0.5 + llepsilon) {
460  G_important_message(_("90 degree south is exceeded by %g cells"),
461  diff);
462 
463  if (diff < llepsilon && diff > fpepsilon) {
464  G_verbose_message(_("Subtle input data rounding error of south "
465  "boundary (%g)"),
466  cellhd->south + 90);
467  G_debug(1, "South of south in seconds: %g",
468  (-cellhd->south - 90) * 3600);
469  /* check only, do not modify
470  cellhd->south = -90.0;
471  lladjust = 1;
472  */
473  }
474 
475  diff = diff - 0.5;
476  if (diff < 0)
477  diff = -diff;
478  if (diff < llepsilon && diff > fpepsilon) {
479  G_verbose_message(_("Subtle input data rounding error of south "
480  "boundary (%g)"),
481  cellhd->south + 90 + cellhd->ns_res / 2.0);
482  G_debug(1, "South of south + 0.5 cells in seconds: %g",
483  (-cellhd->south - 90 - cellhd->ns_res / 2.0) * 3600);
484  /* check only, do not modify
485  cellhd->south = -90.0 - cellhd->ns_res / 2.0;
486  lladjust = 1;
487  */
488  }
489  }
490  else
491  G_fatal_error(_("Illegal latitude for South: %g"), cellhd->south);
492  }
493 
494  if (lladjust)
495  cellhd->ns_res = (cellhd->north - cellhd->south) / cellhd->rows;
496 
497  return lladjust;
498 }
499 
500 static int ll_check_ew(struct Cell_head *cellhd)
501 {
502  int lladjust;
503  double diff;
504  int ncells;
505 
506  /* lat/lon checks */
507  if (cellhd->proj != PROJECTION_LL)
508  return 0;
509 
510  lladjust = 0;
511 
512  G_debug(3, "ll_check_ew: epsilon: %g", llepsilon);
513 
514  /* west - east, no adjustment */
515  diff = (cellhd->east - cellhd->west) / cellhd->ew_res;
516  ncells = (int)(diff + 0.5);
517  diff -= ncells;
518  if ((diff < 0 && diff < -fpepsilon) || (diff > 0 && diff > fpepsilon)) {
520  _("EW extent does not match EW resolution: %g cells difference"),
521  diff);
522  }
523  if (cellhd->east - cellhd->west > 360.0) {
524  diff = (cellhd->east - cellhd->west - 360.0) / cellhd->ew_res;
525  if (diff > fpepsilon)
526  G_important_message(_("360 degree EW extent is exceeded by %g cells"
527  " (East: %g, West: %g)"),
528  diff, cellhd->east, cellhd->west);
529  }
530  else if (cellhd->east - cellhd->west < 360.0) {
531  diff = (360.0 - (cellhd->east - cellhd->west)) / cellhd->ew_res;
532  if (diff < 1.0 && diff > fpepsilon)
534  _("%g cells missing to cover 360 degree EW extent"), diff);
535  }
536 
537  return lladjust;
538 }
539 
540 /*!
541  * \brief Adjust window for lat/lon.
542  *
543  * This function tries to automatically fix fp precision issues and
544  * adjust rounding errors for lat/lon.
545  *
546  * <b>Note:</b> 3D values are not adjusted.
547  *
548  * \param[in,out] cellhd pointer to Cell_head structure
549  * \return 1 if window was adjusted
550  * \return 0 if window was not adjusted
551  */
552 int G_adjust_window_ll(struct Cell_head *cellhd)
553 {
554  int ll_adjust, res_adj;
555  double dsec, dsec2;
556  char buf[100], buf2[100];
557  double diff;
558  double old, new;
559  struct Cell_head cellhds; /* everything in seconds, not degrees */
560 
561  /* lat/lon checks */
562  if (cellhd->proj != PROJECTION_LL)
563  return 0;
564 
565  /* put everything through ll_format + ll_scan */
566  G_llres_format(cellhd->ns_res, buf);
567  if (G_llres_scan(buf, &new) != 1)
568  G_fatal_error(_("Invalid NS resolution"));
569  cellhd->ns_res = new;
570 
571  G_llres_format(cellhd->ew_res, buf);
572  if (G_llres_scan(buf, &new) != 1)
573  G_fatal_error(_("Invalid EW resolution"));
574  cellhd->ew_res = new;
575 
576  G_lat_format(cellhd->north, buf);
577  if (G_lat_scan(buf, &new) != 1)
578  G_fatal_error(_("Invalid North"));
579  cellhd->north = new;
580 
581  G_lat_format(cellhd->south, buf);
582  if (G_lat_scan(buf, &new) != 1)
583  G_fatal_error(_("Invalid South"));
584  cellhd->south = new;
585 
586  G_lon_format(cellhd->west, buf);
587  if (G_lon_scan(buf, &new) != 1)
588  G_fatal_error(_("Invalid West"));
589  cellhd->west = new;
590 
591  G_lon_format(cellhd->east, buf);
592  if (G_lon_scan(buf, &new) != 1)
593  G_fatal_error(_("Invalid East"));
594  cellhd->east = new;
595 
596  /* convert to seconds */
597  cellhds = *cellhd;
598 
599  old = cellhds.ns_res * 3600;
600  sprintf(buf, "%f", old);
601  sscanf(buf, "%lf", &new);
602  cellhds.ns_res = new;
603 
604  old = cellhds.ew_res * 3600;
605  sprintf(buf, "%f", old);
606  sscanf(buf, "%lf", &new);
607  cellhds.ew_res = new;
608 
609  old = cellhds.north * 3600;
610  sprintf(buf, "%f", old);
611  sscanf(buf, "%lf", &new);
612  cellhds.north = new;
613 
614  old = cellhds.south * 3600;
615  sprintf(buf, "%f", old);
616  sscanf(buf, "%lf", &new);
617  cellhds.south = new;
618 
619  old = cellhds.west * 3600;
620  sprintf(buf, "%f", old);
621  sscanf(buf, "%lf", &new);
622  cellhds.west = new;
623 
624  old = cellhds.east * 3600;
625  sprintf(buf, "%f", old);
626  sscanf(buf, "%lf", &new);
627  cellhds.east = new;
628 
629  ll_adjust = 0;
630 
631  /* N - S */
632  /* resolution */
633  res_adj = 0;
634  old = cellhds.ns_res;
635 
636  if (old > 0.4) {
637  /* round to nearest 0.1 sec */
638  dsec = old * 10;
639  dsec2 = floor(dsec + 0.5);
640  new = dsec2 / 10;
641  diff = fabs(dsec2 - dsec) / dsec;
642  if (diff > 0 && diff < llepsilon) {
643  G_llres_format(old / 3600, buf);
644  G_llres_format(new / 3600, buf2);
645  if (strcmp(buf, buf2))
646  G_verbose_message(_("NS resolution rounded from %s to %s"), buf,
647  buf2);
648  ll_adjust = 1;
649  res_adj = 1;
650  cellhds.ns_res = new;
651  }
652  }
653 
654  if (res_adj) {
655  double n_off, s_off;
656 
657  old = cellhds.north;
658  dsec = old * 10;
659  dsec2 = floor(dsec + 0.5);
660  diff = fabs(dsec2 - dsec) / (cellhds.ns_res * 10);
661  n_off = diff;
662 
663  old = cellhds.south;
664  dsec = old * 10;
665  dsec2 = floor(dsec + 0.5);
666  diff = fabs(dsec2 - dsec) / (cellhds.ns_res * 10);
667  s_off = diff;
668 
669  if (n_off < llepsilon || n_off <= s_off) {
670  old = cellhds.north;
671  dsec = old * 10;
672  dsec2 = floor(dsec + 0.5);
673  new = dsec2 / 10;
674  diff = n_off;
675  if (diff > 0 && diff < llepsilon) {
676  G_lat_format(old / 3600, buf);
677  G_lat_format(new / 3600, buf2);
678  if (strcmp(buf, buf2))
679  G_verbose_message(_("North rounded from %s to %s"), buf,
680  buf2);
681  cellhds.north = new;
682  }
683 
684  old = cellhds.south;
685  new = cellhds.north - cellhds.ns_res *cellhds.rows;
686  diff = fabs(new - old) / cellhds.ns_res;
687  if (diff > 0) {
688  G_lat_format(old / 3600, buf);
689  G_lat_format(new / 3600, buf2);
690  if (strcmp(buf, buf2))
691  G_verbose_message(_("South adjusted from %s to %s"), buf,
692  buf2);
693  }
694  cellhds.south = new;
695  }
696  else {
697  old = cellhds.south;
698  dsec = old * 10;
699  dsec2 = floor(dsec + 0.5);
700  new = dsec2 / 10;
701  diff = s_off;
702  if (diff > 0 && diff < llepsilon) {
703  G_lat_format(old / 3600, buf);
704  G_lat_format(new / 3600, buf2);
705  if (strcmp(buf, buf2))
706  G_verbose_message(_("South rounded from %s to %s"), buf,
707  buf2);
708  cellhds.south = new;
709  }
710 
711  old = cellhds.north;
712  new = cellhds.south + cellhds.ns_res *cellhds.rows;
713  diff = fabs(new - old) / cellhds.ns_res;
714  if (diff > 0) {
715  G_lat_format(old / 3600, buf);
716  G_lat_format(new / 3600, buf2);
717  if (strcmp(buf, buf2))
718  G_verbose_message(_("North adjusted from %s to %s"), buf,
719  buf2);
720  }
721  cellhds.north = new;
722  }
723  }
724  else {
725  old = cellhds.north;
726  dsec = old * 10;
727  dsec2 = floor(dsec + 0.5);
728  new = dsec2 / 10;
729  diff = fabs(dsec2 - dsec) / (cellhds.ns_res * 10);
730  if (diff > 0 && diff < llepsilon) {
731  G_lat_format(old / 3600, buf);
732  G_lat_format(new / 3600, buf2);
733  if (strcmp(buf, buf2))
734  G_verbose_message(_("North rounded from %s to %s"), buf, buf2);
735  ll_adjust = 1;
736  cellhds.north = new;
737  }
738 
739  old = cellhds.south;
740  dsec = old * 10;
741  dsec2 = floor(dsec + 0.5);
742  new = dsec2 / 10;
743  diff = fabs(dsec2 - dsec) / (cellhds.ns_res * 10);
744  if (diff > 0 && diff < llepsilon) {
745  G_lat_format(old / 3600, buf);
746  G_lat_format(new / 3600, buf2);
747  if (strcmp(buf, buf2))
748  G_verbose_message(_("South rounded from %s to %s"), buf, buf2);
749  ll_adjust = 1;
750  cellhds.south = new;
751  }
752  }
753  cellhds.ns_res = (cellhds.north - cellhds.south) / cellhds.rows;
754 
755  /* E - W */
756  /* resolution */
757  res_adj = 0;
758  old = cellhds.ew_res;
759 
760  if (old > 0.4) {
761  /* round to nearest 0.1 sec */
762  dsec = old * 10;
763  dsec2 = floor(dsec + 0.5);
764  new = dsec2 / 10;
765  diff = fabs(dsec2 - dsec) / dsec;
766  if (diff > 0 && diff < llepsilon) {
767  G_llres_format(old / 3600, buf);
768  G_llres_format(new / 3600, buf2);
769  if (strcmp(buf, buf2))
770  G_verbose_message(_("EW resolution rounded from %s to %s"), buf,
771  buf2);
772  ll_adjust = 1;
773  res_adj = 1;
774  cellhds.ew_res = new;
775  }
776  }
777 
778  if (res_adj) {
779  double w_off, e_off;
780 
781  old = cellhds.west;
782  dsec = old * 10;
783  dsec2 = floor(dsec + 0.5);
784  diff = fabs(dsec2 - dsec) / (cellhds.ew_res * 10);
785  w_off = diff;
786 
787  old = cellhds.east;
788  dsec = old * 10;
789  dsec2 = floor(dsec + 0.5);
790  diff = fabs(dsec2 - dsec) / (cellhds.ew_res * 10);
791  e_off = diff;
792 
793  if (w_off < llepsilon || w_off <= e_off) {
794  old = cellhds.west;
795  dsec = old * 10;
796  dsec2 = floor(dsec + 0.5);
797  new = dsec2 / 10;
798  diff = w_off;
799  if (diff > 0 && diff < llepsilon) {
800  G_lon_format(old / 3600, buf);
801  G_lon_format(new / 3600, buf2);
802  if (strcmp(buf, buf2))
803  G_verbose_message(_("West rounded from %s to %s"), buf,
804  buf2);
805  cellhds.west = new;
806  }
807 
808  old = cellhds.east;
809  new = cellhds.west + cellhds.ew_res *cellhds.cols;
810  diff = fabs(new - old) / cellhds.ew_res;
811  if (diff > 0) {
812  G_lon_format(old / 3600, buf);
813  G_lon_format(new / 3600, buf2);
814  if (strcmp(buf, buf2))
815  G_verbose_message(_("East adjusted from %s to %s"), buf,
816  buf2);
817  }
818  cellhds.east = new;
819  }
820  else {
821  old = cellhds.east;
822  dsec = old * 10;
823  dsec2 = floor(dsec + 0.5);
824  new = dsec2 / 10;
825  diff = e_off;
826  if (diff > 0 && diff < llepsilon) {
827  G_lon_format(old / 3600, buf);
828  G_lon_format(new / 3600, buf2);
829  if (strcmp(buf, buf2))
830  G_verbose_message(_("East rounded from %s to %s"), buf,
831  buf2);
832  cellhds.east = new;
833  }
834 
835  old = cellhds.west;
836  new = cellhds.east - cellhds.ew_res *cellhds.cols;
837  diff = fabs(new - cellhds.west) / cellhds.ew_res;
838  if (diff > 0) {
839  G_lon_format(old / 3600, buf);
840  G_lon_format(new / 3600, buf2);
841  if (strcmp(buf, buf2))
842  G_verbose_message(_("West adjusted from %s to %s"), buf,
843  buf2);
844  }
845  cellhds.west = new;
846  }
847  }
848  else {
849  old = cellhds.west;
850  dsec = old * 10;
851  dsec2 = floor(dsec + 0.5);
852  new = dsec2 / 10;
853  diff = fabs(dsec2 - dsec) / (cellhds.ew_res * 10);
854  if (diff > 0 && diff < llepsilon) {
855  G_lon_format(old / 3600, buf);
856  G_lon_format(new / 3600, buf2);
857  if (strcmp(buf, buf2))
858  G_verbose_message(_("West rounded from %s to %s"), buf, buf2);
859  ll_adjust = 1;
860  cellhds.west = new;
861  }
862 
863  old = cellhds.east;
864  dsec = old * 10;
865  dsec2 = floor(dsec + 0.5);
866  new = dsec2 / 10;
867  diff = fabs(dsec2 - dsec) / (cellhds.ew_res * 10);
868  if (diff > 0 && diff < llepsilon) {
869  G_lon_format(old / 3600, buf);
870  G_lon_format(new / 3600, buf2);
871  if (strcmp(buf, buf2))
872  G_verbose_message(_("East rounded from %s to %s"), buf, buf2);
873  ll_adjust = 1;
874  cellhds.east = new;
875  }
876  }
877  cellhds.ew_res = (cellhds.east - cellhds.west) / cellhds.cols;
878 
879  cellhd->ns_res = cellhds.ns_res / 3600;
880  cellhd->ew_res = cellhds.ew_res / 3600;
881  cellhd->north = cellhds.north / 3600;
882  cellhd->south = cellhds.south / 3600;
883  cellhd->west = cellhds.west / 3600;
884  cellhd->east = cellhds.east / 3600;
885 
886  return ll_adjust;
887 }
void 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:163
#define LL_TOLERANCE
Definition: adj_cellhd.c:19
void G_adjust_Cell_head(struct Cell_head *cellhd, int row_flag, int col_flag)
Adjust cell header.
Definition: adj_cellhd.c:51
int G_adjust_window_ll(struct Cell_head *cellhd)
Adjust window for lat/lon.
Definition: adj_cellhd.c:552
void void void void G_fatal_error(const char *,...) __attribute__((format(printf
void G_warning(const char *,...) __attribute__((format(printf
void G_lat_format(double, char *)
Definition: ll_format.c:40
void G_lon_format(double, char *)
Definition: ll_format.c:55
int G_lon_scan(const char *, double *)
Definition: ll_scan.c:50
void void G_verbose_message(const char *,...) __attribute__((format(printf
void G_llres_format(double, char *)
Definition: ll_format.c:70
void void void G_important_message(const char *,...) __attribute__((format(printf
int G_lat_scan(const char *, double *)
Definition: ll_scan.c:45
int G_llres_scan(const char *, double *)
Definition: ll_scan.c:55
int G_debug(int, const char *,...) __attribute__((format(printf
#define PROJECTION_LL
Projection code - Latitude-Longitude.
Definition: gis.h:130
#define _(str)
Definition: glocale.h:10
if(!(yy_init))
Definition: sqlp.yy.c:775
2D/3D raster map header (used also for region)
Definition: gis.h:440
int cols3
Number of columns for 3D data.
Definition: gis.h:461
double ew_res
Resolution - east to west cell size for 2D data.
Definition: gis.h:476
double north
Extent coordinates (north)
Definition: gis.h:486
double bottom
Extent coordinates (bottom) - 3D data.
Definition: gis.h:496
int depths
number of depths for 3D data
Definition: gis.h:463
double east
Extent coordinates (east)
Definition: gis.h:490
double ew_res3
Resolution - east to west cell size for 3D data.
Definition: gis.h:478
double ns_res
Resolution - north to south cell size for 2D data.
Definition: gis.h:480
double ns_res3
Resolution - north to south cell size for 3D data.
Definition: gis.h:482
double top
Extent coordinates (top) - 3D data.
Definition: gis.h:494
int rows3
Number of rows for 3D data.
Definition: gis.h:457
int rows
Number of rows for 2D data.
Definition: gis.h:455
int cols
Number of columns for 2D data.
Definition: gis.h:459
int proj
Projection code.
Definition: gis.h:472
double south
Extent coordinates (south)
Definition: gis.h:488
double tb_res
Resolution - top to bottom cell size for 3D data.
Definition: gis.h:484
double west
Extent coordinates (west)
Definition: gis.h:492