1 #include <string.h>
2 #include <grass/gis.h>
3 #include <grass/raster.h>
4 #include <grass/glocale.h>
5 #include "global.h"
6 
read_cells(void)7 void read_cells(void)
8 {
9     int fd, i, j;
10     RASTER_MAP_TYPE data_type;
11     CELL *ccell = NULL;
12     FCELL *fcell = NULL;
13     struct Cell_head inhead;
14     char buf_wrns[32], buf_wrew[32], buf_mrns[32], buf_mrew[32];
15 
16     fd = Rast_open_old(input, "");
17 
18     data_type = Rast_get_map_type(fd);
19     Rast_get_cellhd(input, "", &inhead);
20 
21     if (data_type == CELL_TYPE)
22 	ccell = (CELL *) G_malloc(sizeof(CELL) * window.cols);
23     else if (data_type == FCELL_TYPE)
24 	fcell = (FCELL *) G_malloc(sizeof(FCELL) * window.cols);
25 
26     cell = (DCELL **) G_malloc(sizeof(DCELL *) * window.rows);
27     atb = (DCELL **) G_malloc(sizeof(DCELL *) * window.rows);
28     a = (DCELL **) G_malloc(sizeof(DCELL *) * window.rows);
29 
30     if (window.ew_res < inhead.ew_res || window.ns_res < inhead.ns_res) {
31 	G_format_resolution(window.ew_res, buf_wrew, G_projection());
32 	G_format_resolution(window.ns_res, buf_wrns, G_projection());
33 	G_format_resolution(inhead.ew_res, buf_mrew, G_projection());
34 	G_format_resolution(inhead.ns_res, buf_mrns, G_projection());
35 	G_fatal_error(_("The current region resolution [%s x %s] is finer "
36 			"than the input map's resolution [%s x %s]. "
37 			"The current region resolution must be identical "
38 			"to, or coarser than, the input map's resolution."),
39 		      buf_wrew, buf_wrns, buf_mrew, buf_mrns);
40     }
41 
42     G_message(_("Reading elevation map..."));
43 
44     for (i = 0; i < window.rows; i++) {
45 	G_percent(i, window.rows, 2);
46 
47 	cell[i] = (DCELL *) G_malloc(sizeof(DCELL) * window.cols);
48 	atb[i] = (DCELL *) G_malloc(sizeof(DCELL) * window.cols);
49 	a[i] = (DCELL *) G_malloc(sizeof(DCELL) * window.cols);
50 
51 	if (data_type == CELL_TYPE) {
52 	    Rast_get_c_row(fd, ccell, i);
53 	    for (j = 0; j < window.cols; j++) {
54 		if (Rast_is_c_null_value(&ccell[j]))
55 		    Rast_set_d_null_value(&cell[i][j], 1);
56 		else
57 		    cell[i][j] = (DCELL) ccell[j];
58 	    }
59 	}
60 	else if (data_type == FCELL_TYPE) {
61 	    Rast_get_f_row(fd, fcell, i);
62 	    for (j = 0; j < window.cols; j++) {
63 		if (Rast_is_f_null_value(&fcell[j]))
64 		    Rast_set_d_null_value(&cell[i][j], 1);
65 		else
66 		    cell[i][j] = (DCELL) fcell[j];
67 	    }
68 	}
69 	else
70 	    Rast_get_d_row(fd, cell[i], i);
71     }
72     if (data_type == CELL_TYPE)
73 	G_free(ccell);
74     else if (data_type == FCELL_TYPE)
75 	G_free(fcell);
76     G_percent(i, window.rows, 2);
77     Rast_close(fd);
78 }
79 
write_cells(void)80 void write_cells(void)
81 {
82     int fd, i;
83     struct History history;
84 
85     fd = Rast_open_new(output, DCELL_TYPE);
86 
87     G_message(_("Writing topographic index map..."));
88 
89     for (i = 0; i < window.rows; i++) {
90 	G_percent(i, window.rows, 2);
91 	Rast_put_d_row(fd, atb[i]);
92     }
93     G_percent(i, window.rows, 2);
94     Rast_close(fd);
95 
96     Rast_short_history(output, "raster", &history);
97     Rast_set_history(&history, HIST_DATSRC_1, input);
98     Rast_write_history(output, &history);
99 }
100