1#include <grass/config.h>
6#include <grass/lidar.h>
53int P_set_regions(
struct Cell_head *Elaboration,
struct bound_box *General,
54 struct bound_box *Overlap,
struct Reg_dimens dim,
int type)
58 struct Cell_head orig;
66 Elaboration->south = Elaboration->north - dim.
sn_size;
67 General->N = Elaboration->north - dim.
edge_h;
68 General->S = Elaboration->south + dim.
edge_h;
69 Overlap->N = General->N - dim.
overlap;
70 Overlap->S = General->S + dim.
overlap;
74 Elaboration->west = Elaboration->east - dim.
overlap - (2 * dim.
edge_v);
75 Elaboration->east = Elaboration->west + dim.
ew_size;
76 General->W = Elaboration->west + dim.
edge_v;
77 General->E = Elaboration->east - dim.
edge_v;
78 Overlap->W = General->W + dim.
overlap;
79 Overlap->E = General->E - dim.
overlap;
83 Elaboration->north = orig.north + 2 * dim.
edge_h;
84 Elaboration->south = Elaboration->north - dim.
sn_size;
85 General->N = orig.north;
86 General->S = Elaboration->south + dim.
edge_h;
87 Overlap->N = General->N;
88 Overlap->S = General->S + dim.
overlap;
92 Elaboration->south = orig.south - 2 * dim.
edge_h;
93 General->S = orig.south;
94 Overlap->S = General->S;
98 Elaboration->west = orig.west - 2 * dim.
edge_v;
99 Elaboration->east = Elaboration->west + dim.
ew_size;
100 General->W = orig.west;
101 General->E = Elaboration->east - dim.
edge_v;
102 Overlap->W = General->W;
103 Overlap->E = General->E - dim.
overlap;
107 Elaboration->east = orig.east + 2 * dim.
edge_v;
108 General->E = orig.east;
109 Overlap->E = General->E;
120 int total_splines, edge_splines, n_windows;
121 int lastsplines, lastsplines_min, lastsplines_max;
122 double E_extension, N_extension, edgeE, edgeN;
123 struct Cell_head orig;
128 E_extension = orig.east - orig.west;
129 N_extension = orig.north - orig.south;
140 total_splines = ceil(E_extension / pe);
141 edge_splines = edgeE / pe;
142 n_windows = E_extension / edgeE;
152 lastsplines = total_splines - edge_splines * n_windows;
153 while (lastsplines > lastsplines_max || lastsplines < lastsplines_min) {
158 edge_splines = edgeE / pe;
159 n_windows = E_extension / edgeE;
160 lastsplines = total_splines - edge_splines * n_windows;
166 total_splines = ceil(N_extension / pn);
167 edge_splines = edgeN / pn;
168 n_windows = N_extension / edgeN;
178 lastsplines = total_splines - edge_splines * n_windows;
179 while (lastsplines > lastsplines_max || lastsplines < lastsplines_min) {
184 edge_splines = edgeN / pn;
185 n_windows = N_extension / edgeN;
186 lastsplines = total_splines - edge_splines * n_windows;
223 return (2 * nsplines + 1);
226 return (4 * nsplines + 3);
234 int i, mean_count = 0;
236 struct bound_box mean_box;
238 Vect_region_box(Elaboration, &mean_box);
244 for (i = 0; i < npoints; i++) {
245 if (Vect_point_in_box(obs[i].coordX, obs[i].coordY, obs[i].coordZ,
254 mean /= (double)mean_count;
261 int type, npoints = 0;
262 double xmin = 0, xmax = 0, ymin = 0, ymax = 0;
264 struct line_pnts *points;
265 struct line_cats *categories;
266 struct bound_box region_box;
267 struct Cell_head orig;
270 Vect_region_box(&orig, ®ion_box);
272 points = Vect_new_line_struct();
273 categories = Vect_new_cats_struct();
276 while ((type = Vect_read_next_line(Map, points, categories)) > 0) {
277 if (!(type & GV_POINT))
282 if (points->z !=
NULL)
288 if (Vect_point_in_box(
x, y, z, ®ion_box)) {
307 Vect_destroy_cats_struct(categories);
308 Vect_destroy_line_struct(points);
312 *dist = sqrt(((xmax - xmin) * (ymax - ymin)) / npoints);
314 *dens = npoints / ((xmax - xmin) * (ymax - ymin));
323 struct Cell_head *Elaboration,
324 int *num_points,
int dim_vect,
int layer)
326 int line_num, pippo, npoints,
cat, type;
329 struct line_pnts *points;
330 struct line_cats *categories;
331 struct bound_box elaboration_box;
334 obs = (
struct Point *)G_calloc(pippo,
sizeof(
struct Point));
336 points = Vect_new_line_struct();
337 categories = Vect_new_cats_struct();
340 Vect_region_box(Elaboration, &elaboration_box);
346 while ((type = Vect_read_next_line(Map, points, categories)) > 0) {
348 if (!(type & GV_POINT))
355 if (points->z !=
NULL)
361 if (Vect_point_in_box(
x, y, z, &elaboration_box)) {
363 if (npoints >= pippo) {
365 obs = (
struct Point *)G_realloc(
366 (
void *)obs, (
signed int)pippo *
sizeof(
struct Point));
371 obs[npoints - 1].
coordY = y;
372 obs[npoints - 1].
coordZ = z;
376 Vect_cat_get(categories, layer, &
cat);
377 obs[npoints - 1].
cat =
cat;
380 Vect_destroy_line_struct(points);
381 Vect_destroy_cats_struct(categories);
383 *num_points = npoints;
388 struct Cell_head *Elaboration,
389 struct Cell_head *Original,
390 int *num_points,
int dim_vect)
392 int col, row, startcol, endcol, startrow, endrow, nrows, ncols;
396 struct bound_box elaboration_box;
399 obs = (
struct Point *)G_calloc(pippo,
sizeof(
struct Point));
402 Vect_region_box(Elaboration, &elaboration_box);
405 nrows = Original->rows;
406 ncols = Original->cols;
408 if (Original->north > Elaboration->north)
410 (Original->north - Elaboration->north) / Original->ns_res - 1;
413 if (Original->north > Elaboration->south) {
414 endrow = (Original->north - Elaboration->south) / Original->ns_res + 1;
420 if (Elaboration->west > Original->west)
421 startcol = (Elaboration->west - Original->west) / Original->ew_res - 1;
424 if (Elaboration->east > Original->west) {
425 endcol = (Elaboration->east - Original->west) / Original->ew_res + 1;
432 for (row = startrow; row < endrow; row++) {
433 for (col = startcol; col < endcol; col++) {
437 if (!Rast_is_d_null_value(&z)) {
439 x = Rast_col_to_easting((
double)(col) + 0.5, Original);
440 y = Rast_row_to_northing((
double)(row) + 0.5, Original);
442 if (Vect_point_in_box(
x, y, 0, &elaboration_box)) {
444 if (npoints >= pippo) {
446 obs = (
struct Point *)G_realloc(
448 (
signed int)pippo *
sizeof(
struct Point));
453 obs[npoints - 1].
coordY = y;
454 obs[npoints - 1].
coordZ = z;
460 *num_points = npoints;
467 dbTable *auxiliar_tab;
470 auxiliar_tab = db_alloc_table(2);
471 db_set_table_name(auxiliar_tab, tab_name);
472 db_set_table_description(auxiliar_tab,
"Intermediate interpolated values");
474 column = db_get_table_column(auxiliar_tab, 0);
475 db_set_column_name(column,
"ID");
476 db_set_column_sqltype(column, DB_SQL_TYPE_INTEGER);
478 column = db_get_table_column(auxiliar_tab, 1);
479 db_set_column_name(column,
"Interp");
480 db_set_column_sqltype(column, DB_SQL_TYPE_REAL);
482 if (db_create_table(
driver, auxiliar_tab) == DB_OK) {
483 G_debug(1, _(
"<%s> created in database."), tab_name);
487 G_warning(_(
"<%s> has not been created in database."), tab_name);
495 dbTable *auxiliar_tab;
498 auxiliar_tab = db_alloc_table(4);
499 db_set_table_name(auxiliar_tab, tab_name);
500 db_set_table_description(auxiliar_tab,
"Intermediate interpolated values");
502 column = db_get_table_column(auxiliar_tab, 0);
503 db_set_column_name(column,
"ID");
504 db_set_column_sqltype(column, DB_SQL_TYPE_INTEGER);
506 column = db_get_table_column(auxiliar_tab, 1);
507 db_set_column_name(column,
"Interp");
508 db_set_column_sqltype(column, DB_SQL_TYPE_REAL);
510 column = db_get_table_column(auxiliar_tab, 2);
511 db_set_column_name(column,
"X");
512 db_set_column_sqltype(column, DB_SQL_TYPE_DOUBLE_PRECISION);
514 column = db_get_table_column(auxiliar_tab, 3);
515 db_set_column_name(column,
"Y");
516 db_set_column_sqltype(column, DB_SQL_TYPE_DOUBLE_PRECISION);
518 if (db_create_table(
driver, auxiliar_tab) == DB_OK) {
519 G_debug(1, _(
"<%s> created in database."), tab_name);
523 G_warning(_(
"<%s> has not been created in database."), tab_name);
533 db_init_string(&drop);
534 db_append_string(&drop,
"drop table ");
535 db_append_string(&drop, tab_name);
536 return db_execute_immediate(
driver, &drop);
542 int ncols, col, nrows, row;
545 nrows = Rast_window_rows();
546 ncols = Rast_window_cols();
548 raster = Rast_allocate_buf(DCELL_TYPE);
550 for (row = 0; row < nrows; row++) {
553 Rast_set_d_null_value(raster, ncols);
555 for (col = 0, ptr = raster; col < ncols;
557 Rast_set_d_value(ptr, (DCELL)(matrix[row][col]), DCELL_TYPE);
559 Rast_put_d_row(fd, raster);
566 dbDriver *
driver,
char *tab_name)
572 struct line_pnts *point;
573 struct line_cats *cat;
582 point = Vect_new_line_struct();
583 cat = Vect_new_cats_struct();
585 db_init_string(&sql);
586 db_zero_string(&sql);
588 sprintf(buf,
"select ID, X, Y, sum(Interp) from %s group by ID, X, Y",
591 db_append_string(&sql, buf);
592 db_open_select_cursor(
driver, &sql, &cursor, DB_SEQUENTIAL);
594 while (db_fetch(&cursor, DB_NEXT, &more) == DB_OK && more) {
595 table = db_get_cursor_table(&cursor);
597 column = db_get_table_column(table, 0);
598 type = db_sqltype_to_Ctype(db_get_column_sqltype(column));
599 if (type == DB_C_TYPE_INT)
600 value = db_get_column_value(column);
604 column = db_get_table_column(table, 1);
605 type = db_sqltype_to_Ctype(db_get_column_sqltype(column));
606 if (type == DB_C_TYPE_DOUBLE)
607 value = db_get_column_value(column);
610 coordZ = db_get_value_double(value);
612 column = db_get_table_column(table, 2);
613 type = db_sqltype_to_Ctype(db_get_column_sqltype(column));
614 if (type == DB_C_TYPE_DOUBLE)
615 value = db_get_column_value(column);
618 coordX = db_get_value_double(value);
620 column = db_get_table_column(table, 3);
621 type = db_sqltype_to_Ctype(db_get_column_sqltype(column));
622 if (type == DB_C_TYPE_DOUBLE)
623 value = db_get_column_value(column);
626 coordY = db_get_value_double(value);
628 Vect_copy_xyz_to_pnts(point, &coordX, &coordY, &coordZ, 1);
629 Vect_reset_cats(cat);
630 Vect_cat_set(cat, 1, 1);
631 Vect_write_line(Out, GV_POINT, point, cat);
void * G_incr_void_ptr(const void *ptr, size_t size)
Advance void pointer.
int G_debug(int level, const char *msg,...)
Print debugging message.
void G_get_window(struct Cell_head *window)
Get the current region.
void G_warning(const char *msg,...)
Print a warning message to stderr.
void G_get_set_window(struct Cell_head *window)
Get the current working window (region)
void G_percent(long n, long d, int s)
Print percent complete messages.
int Segment_get(SEGMENT *SEG, void *buf, off_t row, off_t col)
Get value from segment file.
void P_Aux_to_Raster(double **matrix, int fd)
double P_Mean_Calc(struct Cell_head *Elaboration, struct Point *obs, int npoints)
double P_estimate_splinestep(struct Map_info *Map, double *dens, double *dist)
void P_Aux_to_Vector(struct Map_info *Map UNUSED, struct Map_info *Out, dbDriver *driver, char *tab_name)
struct Point * P_Read_Vector_Region_Map(struct Map_info *Map, struct Cell_head *Elaboration, int *num_points, int dim_vect, int layer)
int P_get_edge(int interpolator, struct Reg_dimens *dim, double pe, double pn)
void P_zero_dim(struct Reg_dimens *dim)
int P_get_BandWidth(int interpolator, int nsplines)
struct Point * P_Read_Raster_Region_Map(SEGMENT *in_seg, struct Cell_head *Elaboration, struct Cell_head *Original, int *num_points, int dim_vect)
int P_Create_Aux2_Table(dbDriver *driver, char *tab_name)
int P_set_dim(struct Reg_dimens *dim, double pe, double pn, int *nsplx, int *nsply)
int P_Drop_Aux_Table(dbDriver *driver, char *tab_name)
int P_Create_Aux4_Table(dbDriver *driver, char *tab_name)
int P_set_regions(struct Cell_head *Elaboration, struct bound_box *General, struct bound_box *Overlap, struct Reg_dimens dim, int type)