[GRASS-SVN] r66094 - in grass/trunk/raster/r.in.lidar: . test
svn_grass at osgeo.org
svn_grass at osgeo.org
Thu Sep 3 19:40:46 PDT 2015
Author: wenzeslaus
Date: 2015-09-03 19:40:46 -0700 (Thu, 03 Sep 2015)
New Revision: 66094
Added:
grass/trunk/raster/r.in.lidar/test/
grass/trunk/raster/r.in.lidar/test/sample_test.sh
Modified:
grass/trunk/raster/r.in.lidar/main.c
grass/trunk/raster/r.in.lidar/r.in.lidar.html
Log:
r.in.lidar: add base raster to get height above ground
Modified: grass/trunk/raster/r.in.lidar/main.c
===================================================================
--- grass/trunk/raster/r.in.lidar/main.c 2015-09-03 22:06:20 UTC (rev 66093)
+++ grass/trunk/raster/r.in.lidar/main.c 2015-09-04 02:40:46 UTC (rev 66094)
@@ -99,7 +99,7 @@
int main(int argc, char *argv[])
{
- int out_fd;
+ int out_fd, base_raster;
char *infile, *outmap;
int percent;
int method = -1;
@@ -107,11 +107,11 @@
double zrange_min, zrange_max, d_tmp;
unsigned long estimated_lines;
- RASTER_MAP_TYPE rtype;
+ RASTER_MAP_TYPE rtype, base_raster_data_type;
struct History history;
char title[64];
void *n_array, *min_array, *max_array, *sum_array, *sumsq_array,
- *index_array;
+ *index_array, *base_array;
void *raster_row, *ptr;
struct Cell_head region;
int rows, cols; /* scan box size */
@@ -146,7 +146,7 @@
struct GModule *module;
struct Option *input_opt, *output_opt, *percent_opt, *type_opt, *filter_opt, *class_opt;
- struct Option *method_opt, *zrange_opt, *zscale_opt;
+ struct Option *method_opt, *base_raster_opt, *zrange_opt, *zscale_opt;
struct Option *trim_opt, *pth_opt, *res_opt;
struct Flag *print_flag, *scan_flag, *shell_style, *over_flag, *extents_flag, *intens_flag;
@@ -197,6 +197,12 @@
type_opt->answer = "FCELL";
type_opt->description = _("Storage type for resultant raster map");
+ base_raster_opt = G_define_standard_option(G_OPT_R_INPUT);
+ base_raster_opt->key = "base_raster";
+ base_raster_opt->required = NO;
+ base_raster_opt->label = _("Subtract raster values from the z coordinates");
+ base_raster_opt->description = _("The scale for z is applied beforehand, the filter afterwards");
+
zrange_opt = G_define_option();
zrange_opt->key = "zrange";
zrange_opt->type = TYPE_DOUBLE;
@@ -507,6 +513,7 @@
sum_array = NULL;
sumsq_array = NULL;
index_array = NULL;
+ base_array = NULL;
if (strcmp(method_opt->answer, "n") == 0) {
method = METHOD_N;
@@ -641,6 +648,12 @@
npasses = (int)ceil(1.0 * region.rows / rows);
+ if (base_raster_opt->answer) {
+ /* 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);
+ }
+
if (!scan_flag->answer) {
/* check if rows * (cols + 1) go into a size_t */
if (sizeof(size_t) < 8) {
@@ -665,6 +678,10 @@
if (bin_index)
index_array =
G_calloc((size_t)rows * (cols + 1), Rast_cell_size(CELL_TYPE));
+ if (base_raster_opt->answer)
+ 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 */
if (bin_n)
@@ -719,6 +736,13 @@
pass_south = region.south; /* exact copy to avoid fp errors */
}
+ if (base_array) {
+ G_debug(2, "filling base raster array");
+ for (row = 0; row < rows; row++) {
+ Rast_get_row(base_raster, base_array + (row * cols * Rast_cell_size(base_raster_data_type)), row, base_raster_data_type);
+ }
+ }
+
G_debug(2, "pass=%d/%d pass_n=%f pass_s=%f rows=%d",
pass, npasses, pass_north, pass_south, rows);
@@ -832,19 +856,41 @@
z = z * zscale;
- if (zrange_opt->answer) {
- if (z < zrange_min || z > zrange_max) {
- continue;
- }
- }
+ /* find the bin in the current array box */
+ arr_row = (int)((pass_north - y) / region.ns_res);
+ 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)) {
+ 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;
+ }
+
+ if (zrange_opt->answer) {
+ if (z < zrange_min || z > zrange_max) {
+ continue;
+ }
+ }
+
count++;
/* G_debug(5, "x: %f, y: %f, z: %f", x, y, z); */
- /* find the bin in the current array box */
- arr_row = (int)((pass_north - y) / region.ns_res);
- arr_col = (int)((x - region.west) / region.ew_res);
-
if (bin_n)
update_n(n_array, cols, arr_row, arr_col);
if (bin_min)
@@ -1207,7 +1253,8 @@
max_nodes = 0;
nodes = NULL;
}
-
+ if (base_array)
+ Rast_close(base_raster);
} /* passes loop */
G_percent(1, 1, 1); /* flush */
Modified: grass/trunk/raster/r.in.lidar/r.in.lidar.html
===================================================================
--- grass/trunk/raster/r.in.lidar/r.in.lidar.html 2015-09-03 22:06:20 UTC (rev 66093)
+++ grass/trunk/raster/r.in.lidar/r.in.lidar.html 2015-09-04 02:40:46 UTC (rev 66094)
@@ -200,7 +200,7 @@
<em>r.neighbors</em> to smooth the stddev map before further use.]
-<h2>EXAMPLE</h2>
+<h2>EXAMPLES</h2>
Import of a LAS file into an existing location/mapset (metric):
@@ -211,8 +211,9 @@
r.univar lidar_dem_mean
</pre></div>
+<h3>Serpent Mound dataset</h3>
+
<p>
-Serpent Mound dataset:
This example is analogous to the example used in the GRASS wiki page for
<a href="http://grasswiki.osgeo.org/wiki/LIDAR#Import_LAS_as_raster_DEM">importing LAS as raster DEM</a>.
<p>The sample LAS data are in the file "Serpent Mound Model LAS Data.las",
@@ -240,6 +241,21 @@
output=Serpent_Mound_Model_LAS_Data method=mean
</pre></div>
+<h3>Height above ground</h3>
+
+<p>
+Compute mean height above ground of the points for each raster cell
+(ground elevation is given by the raster map <tt>elevation</tt>):
+
+<div class="code"><pre>
+g.region raster=elevation
+r.in.lidar input=points.las output=height_above_ground base_raster=elevation
+</pre></div>
+
+In this type of computation, it might be advantageous to change the resolution
+to match the precision of the points rather then derive it from the base raster.
+
+
<h2>NOTES</h2>
The typical file extensions for the LAS format are .las and .laz (compressed).
The compressed LAS (.laz) format can be imported only if libLAS has been compiled
Added: grass/trunk/raster/r.in.lidar/test/sample_test.sh
===================================================================
--- grass/trunk/raster/r.in.lidar/test/sample_test.sh (rev 0)
+++ grass/trunk/raster/r.in.lidar/test/sample_test.sh 2015-09-04 02:40:46 UTC (rev 66094)
@@ -0,0 +1,11 @@
+# GNU GPL
+# Vaclav Petras
+# this is an idea for a test
+export GRASS_OVERWRITE=1
+r.in.lidar in=points.las out=test_pure
+r.in.lidar in=points.las out=test_hag base=elevation
+r.mapcalc "test_difference = (elevation + test_hag) - test_pure"
+r.univar test_difference
+
+
+g.region -p -a raster=test_difference at manual res=0.1 zoom=test_difference at manual
Property changes on: grass/trunk/raster/r.in.lidar/test/sample_test.sh
___________________________________________________________________
Added: svn:executable
+ *
Added: svn:mime-type
+ text/x-sh
Added: svn:eol-style
+ native
More information about the grass-commit
mailing list