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