[GRASS-SVN] r66593 - grass/trunk/raster/r.in.lidar

svn_grass at osgeo.org svn_grass at osgeo.org
Sat Oct 24 19:52:29 PDT 2015


Author: wenzeslaus
Date: 2015-10-24 19:52:29 -0700 (Sat, 24 Oct 2015)
New Revision: 66593

Added:
   grass/trunk/raster/r.in.lidar/rast_segment.c
   grass/trunk/raster/r.in.lidar/rast_segment.h
Modified:
   grass/trunk/raster/r.in.lidar/local_proto.h
   grass/trunk/raster/r.in.lidar/main.c
   grass/trunk/raster/r.in.lidar/support.c
Log:
r.in.lidar: separate code for base raster

But leaving the two options (row array and segment) in the main code.


Modified: grass/trunk/raster/r.in.lidar/local_proto.h
===================================================================
--- grass/trunk/raster/r.in.lidar/local_proto.h	2015-10-25 02:29:41 UTC (rev 66592)
+++ grass/trunk/raster/r.in.lidar/local_proto.h	2015-10-25 02:52:29 UTC (rev 66593)
@@ -52,5 +52,8 @@
 int update_sum(void *, int, int, int, RASTER_MAP_TYPE, double);
 int update_sumsq(void *, int, int, int, RASTER_MAP_TYPE, double);
 
+/* raster reading */
+int row_array_get_value_row_col(void *array, int arr_row, int arr_col,
+                                int cols, RASTER_MAP_TYPE rtype, double *value);
 
 #endif /* __LOCAL_PROTO_H__ */

Modified: grass/trunk/raster/r.in.lidar/main.c
===================================================================
--- grass/trunk/raster/r.in.lidar/main.c	2015-10-25 02:29:41 UTC (rev 66592)
+++ grass/trunk/raster/r.in.lidar/main.c	2015-10-25 02:52:29 UTC (rev 66593)
@@ -28,7 +28,9 @@
 #include <grass/gprojects.h>
 #include <grass/glocale.h>
 #include <liblas/capi/liblas.h>
+
 #include "local_proto.h"
+#include "rast_segment.h"
 
 struct node
 {
@@ -671,36 +673,21 @@
         use_base_raster_res = 1;
     if (base_raster_opt->answer && (res_opt->answer || use_base_raster_res))
         use_segment = 1;
-    if (base_raster_opt->answer) {
+    if (base_raster_opt->answer && !use_segment) {
+        /* TODO: do we need to test existence first? mapset? */
+        base_raster = Rast_open_old(base_raster_opt->answer, "");
+        base_raster_data_type = Rast_get_map_type(base_raster);
+        base_array = G_calloc((size_t)rows * (cols + 1), Rast_cell_size(base_raster_data_type));
+    }
+    if (base_raster_opt->answer && use_segment) {
         if (use_base_raster_res) {
             /* TODO: how to get cellhd already stored in the open map? */
             Rast_get_cellhd(base_raster_opt->answer, "", &input_region);
             /* TODO: make it only as small as the output is or points are */
             Rast_set_input_window(&input_region);  /* we have split window */
         }
-        /* TODO: do we need to test existence first? mapset? */
-        base_raster = Rast_open_old(base_raster_opt->answer, "");
-        base_raster_data_type = Rast_get_map_type(base_raster);
+        rast_segment_open(&base_segment, base_raster_opt->answer, &base_raster_data_type);
     }
-    if (base_raster_opt->answer && use_segment) {
-        if (!use_base_raster_res)
-            Rast_get_input_window(&input_region);  /* we have split window */
-        /* TODO: these numbers does not fit with what we promise about percentage */
-        int segment_rows = 64;
-        /* writing goes row by row, so we use long segments as well */
-        int segment_cols = cols;
-        int segments_in_memory = 4;
-        if (Segment_open(&base_segment, G_tempfile(), input_region.rows, input_region.cols,
-                segment_rows, segment_cols,
-                Rast_cell_size(base_raster_data_type), segments_in_memory) != 1)
-            G_fatal_error(_("Cannot create temporary file with segments of a raster map"));
-        for (row = 0; row < input_region.rows; row++) {
-            raster_row = Rast_allocate_input_buf(base_raster_data_type);
-            Rast_get_row(base_raster, raster_row, row, base_raster_data_type);
-            Segment_put_row(&base_segment, raster_row, row);
-        }
-        Rast_close(base_raster);  /* we won't need the raster again */
-    }
 
     if (!scan_flag->answer) {
 	/* check if rows * (cols + 1) go into a size_t */
@@ -726,9 +713,6 @@
 	if (bin_index)
 	    index_array =
 		G_calloc((size_t)rows * (cols + 1), Rast_cell_size(CELL_TYPE));
-	if (base_raster_opt->answer && !use_segment)
-	    base_array = G_calloc((size_t)rows * (cols + 1), Rast_cell_size(base_raster_data_type));
-	/* we don't free the raster, we just use it again */
 	/* TODO: perhaps none of them needs to be freed */
 
 	/* and then free it again */
@@ -909,61 +893,22 @@
             arr_col = (int)((x - region.west) / region.ew_res);
 
             if (base_array) {
-                ptr = base_array;
-                ptr =
-                    G_incr_void_ptr(ptr,
-                                    ((arr_row * cols) +
-                                     arr_col) * Rast_cell_size(base_raster_data_type));
-
                 double base_z;
-                if (Rast_is_null_value(ptr, base_raster_data_type)) {
+                if (row_array_get_value_row_col(base_array, arr_row, arr_col,
+                                                cols, base_raster_data_type,
+                                                &base_z))
+                    z -= base_z;
+                else
                     continue;
-                }
-                else {
-                    if (base_raster_data_type == DCELL_TYPE)
-                        base_z = *(DCELL *) ptr;
-                    else if (base_raster_data_type == FCELL_TYPE)
-                        base_z = (double) *(FCELL *) ptr;
-                    else
-                        base_z = (double) *(CELL *) ptr;
-                }
-                z -= base_z;
             }
             else if (use_segment) {
                 double base_z;
-                /* Rast gives double, Segment needs off_t */
-                off_t base_row = Rast_northing_to_row(y, &input_region);
-                off_t base_col = Rast_easting_to_col(x, &input_region);
-                /* skip points which are outside the base raster
-                 * (null propagation) */
-                if (base_row < 0 || base_col < 0 ||
-                        base_row >= input_region.rows ||
-                        base_col >= input_region.cols)
+                if (rast_segment_get_value_xy(&base_segment, &input_region,
+                                              base_raster_data_type, x, y,
+                                              &base_z))
+                    z -= base_z;
+                else
                     continue;
-                if (base_raster_data_type == DCELL_TYPE) {
-                    DCELL value;
-                    Segment_get(&base_segment, &value, base_row, base_col);
-                    if (Rast_is_null_value(&value, base_raster_data_type))
-                        continue;
-                    base_z = value;
-                }
-                else if (base_raster_data_type == FCELL_TYPE) {
-                    
-                    
-                    FCELL value;
-                    Segment_get(&base_segment, &value, base_row, base_col);
-                    if (Rast_is_null_value(&value, base_raster_data_type))
-                        continue;
-                    base_z = value;
-                }
-                else {
-                    CELL value;
-                    Segment_get(&base_segment, &value, base_row, base_col);
-                    if (Rast_is_null_value(&value, base_raster_data_type))
-                        continue;
-                    base_z = value;
-                }
-                z -= base_z;
             }
 
             if (zrange_opt->answer) {

Added: grass/trunk/raster/r.in.lidar/rast_segment.c
===================================================================
--- grass/trunk/raster/r.in.lidar/rast_segment.c	                        (rev 0)
+++ grass/trunk/raster/r.in.lidar/rast_segment.c	2015-10-25 02:52:29 UTC (rev 66593)
@@ -0,0 +1,86 @@
+
+#include <grass/gis.h>
+#include <grass/glocale.h>
+#include <grass/raster.h>
+#include <grass/segment.h>
+
+#include "rast_segment.h"
+
+static void rast_segment_load(SEGMENT * segment, int rowio,
+                              RASTER_MAP_TYPE map_type)
+{
+    void *raster_row = Rast_allocate_input_buf(map_type);
+    int row;
+
+    for (row = 0; row < Rast_input_window_rows(); row++) {
+        /* TODO: free mem */
+        Rast_get_row(rowio, raster_row, row, map_type);
+        Segment_put_row(segment, raster_row, row);
+    }
+}
+
+/* TODO: close function */
+
+void rast_segment_open(SEGMENT * segment, const char *name,
+                       RASTER_MAP_TYPE * map_type)
+{
+    /* TODO: check if not passing the mapset is OK */
+    int rowio = Rast_open_old(name, "");
+
+    *map_type = Rast_get_map_type(rowio);
+    int segment_rows = 64;
+
+    /* we use long segments because this is how the values a binned */
+    int segment_cols = Rast_input_window_cols();
+    int segments_in_memory = 4;
+
+    if (Segment_open(segment, G_tempfile(), Rast_input_window_rows(),
+                     Rast_input_window_cols(), segment_rows, segment_cols,
+                     Rast_cell_size(*map_type), segments_in_memory) != 1)
+        G_fatal_error(_("Cannot create temporary file with segments of a raster map"));
+    rast_segment_load(segment, rowio, *map_type);
+    Rast_close(rowio);          /* we won't need the raster again */
+}
+
+
+/* 0 on out of region or NULL, 1 on success */
+int rast_segment_get_value_xy(SEGMENT * base_segment,
+                              struct Cell_head *input_region,
+                              RASTER_MAP_TYPE rtype, double x, double y,
+                              double *value)
+{
+    /* Rast gives double, Segment needs off_t */
+    off_t base_row = Rast_northing_to_row(y, input_region);
+    off_t base_col = Rast_easting_to_col(x, input_region);
+
+    /* skip points which are outside the base raster
+     * (null propagation) */
+    if (base_row < 0 || base_col < 0 ||
+        base_row >= input_region->rows || base_col >= input_region->cols)
+        return 0;
+    if (rtype == DCELL_TYPE) {
+        DCELL tmp;
+
+        Segment_get(base_segment, &tmp, base_row, base_col);
+        if (Rast_is_d_null_value(&tmp))
+            return 0;
+        *value = (double)tmp;
+    }
+    else if (rtype == FCELL_TYPE) {
+        FCELL tmp;
+
+        Segment_get(base_segment, &tmp, base_row, base_col);
+        if (Rast_is_f_null_value(&tmp))
+            return 0;
+        *value = (double)tmp;
+    }
+    else {
+        CELL tmp;
+
+        Segment_get(base_segment, &tmp, base_row, base_col);
+        if (Rast_is_c_null_value(&tmp))
+            return 0;
+        *value = (double)tmp;
+    }
+    return 1;
+}


Property changes on: grass/trunk/raster/r.in.lidar/rast_segment.c
___________________________________________________________________
Added: svn:mime-type
   + text/x-csrc
Added: svn:eol-style
   + native

Added: grass/trunk/raster/r.in.lidar/rast_segment.h
===================================================================
--- grass/trunk/raster/r.in.lidar/rast_segment.h	                        (rev 0)
+++ grass/trunk/raster/r.in.lidar/rast_segment.h	2015-10-25 02:52:29 UTC (rev 66593)
@@ -0,0 +1,14 @@
+#ifndef GRASS_RAST_SEGMENT_H
+#define GRASS_RAST_SEGMENT_H
+
+#include <grass/raster.h>
+#include <grass/segment.h>
+
+void rast_segment_open(SEGMENT * segment, const char *name,
+                       RASTER_MAP_TYPE * map_type);
+int rast_segment_get_value_xy(SEGMENT * base_segment,
+                              struct Cell_head *input_region,
+                              RASTER_MAP_TYPE rtype, double x, double y,
+                              double *value);
+
+#endif /* GRASS_RAST_SEGMENT_H */


Property changes on: grass/trunk/raster/r.in.lidar/rast_segment.h
___________________________________________________________________
Added: svn:mime-type
   + text/x-chdr
Added: svn:eol-style
   + native

Modified: grass/trunk/raster/r.in.lidar/support.c
===================================================================
--- grass/trunk/raster/r.in.lidar/support.c	2015-10-25 02:29:41 UTC (rev 66592)
+++ grass/trunk/raster/r.in.lidar/support.c	2015-10-25 02:52:29 UTC (rev 66593)
@@ -128,3 +128,23 @@
 
     return 0;
 }
+
+/* 0 on NULL, 1 on success */
+int row_array_get_value_row_col(void *array, int arr_row, int arr_col,
+                                int cols, RASTER_MAP_TYPE rtype, double *value)
+{
+    void *ptr = array;
+    ptr =
+        G_incr_void_ptr(ptr,
+                        ((arr_row * cols) +
+                         arr_col) * Rast_cell_size(rtype));
+    if (Rast_is_null_value(ptr, rtype))
+        return 0;
+    if (rtype == DCELL_TYPE)
+        *value = (double) *(DCELL *) ptr;
+    else if (rtype == FCELL_TYPE)
+        *value = (double) *(FCELL *) ptr;
+    else
+        *value = (double) *(CELL *) ptr;
+    return 1;
+}



More information about the grass-commit mailing list