I want to compare two raster maps for any differences. The r.mapcalc
expression I used yielde no results ... unless getting the GRASS shell
prompt means there are no differences.
Command and result:
r.mapcalc "mapdiff = blocked_culvert - open_culvert"
No differences? Wrong syntax?
TIA,
Rich
On Tue, 28 Mar 2017, Rich Shepard wrote:
No differences? Wrong syntax?
GRASS hung due to an initial syntax error; corrected.
Now I have the difference map which is a solid color using the elevation
color table. Ergo, no differences.
Rich
On 28/03/17 22:30, Rich Shepard wrote:
On Tue, 28 Mar 2017, Rich Shepard wrote:
No differences? Wrong syntax?
GRASS hung due to an initial syntax error; corrected.
Now I have the difference map which is a solid color using the elevation
color table. Ergo, no differences.
For the record: there is a "differences" color table which is useful for this.
Moritz
The difference is computed for each pixel.
Just visualize the new map "mapdiff" of run r.univar mapdiff.
2017-03-28 22:23 GMT+02:00 Rich Shepard <rshepard@appl-ecosys.com>:
I want to compare two raster maps for any differences. The r.mapcalc
expression I used yielde no results ... unless getting the GRASS shell
prompt means there are no differences.Command and result:
r.mapcalc "mapdiff = blocked_culvert - open_culvert"
No differences? Wrong syntax?
TIA,
Rich
_______________________________________________
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user
On Tue, Mar 28, 2017 at 4:30 PM, Rich Shepard <rshepard@appl-ecosys.com>
wrote:
Now I have the difference map which is a solid color...
...Ergo, no differences.
As Soeren mentioned in the other part of the thread, use r.univar to get
the statistics. Don't rely on just looking at the resulting raster map.
There can be things like outliers which may throw the visualization off.
Vaclav
On Wed, 29 Mar 2017, Vaclav Petras wrote:
As Soeren mentioned in the other part of the thread, use r.univar to get
the statistics. Don't rely on just looking at the resulting raster map.
There can be things like outliers which may throw the visualization off.
Vaclav,
Here are the results of r.univar on the output of r.mapcalc "output_diff =
blocked_h_max - open_h_max":
n=140652
null_cells=71271
cells=211923
min=-1.89342188835144
max=0.00093058729544282
range=1.89435247564688
mean=-5.03205256742558e-05
mean_of_abs=5.24979771199771e-05
stddev=0.00847518705371129
variance=7.18287955953955e-05
coeff_var=-16842.4056389523
sum=-7.07768257713542
first_quartile=-7.4142e-09
median=0
third_quartile=2.60843e-09
percentile_90=7.45058e-08
With no experience interpreting r.univar output I'm not sure what the
above tells me about the two maps. Yes, I see a small range and many null
cells, but I expect that because 4 cells on one input DEM are set about 0.6m
higher elevation than those on the other input DEM.
Rich
Rich Shepard wrote:
No differences? Wrong syntax?
GRASS hung due to an initial syntax error; corrected.
Now I have the difference map which is a solid color using the elevation
color table. Ergo, no differences.
Moritz:
For the record: there is a "differences" color table which is useful for
this.
Don't forget to try r.colors with -e or -a -g if the first attempt does not
show much differences.
Nikos
On Wed, Mar 29, 2017 at 8:46 PM, Rich Shepard <rshepard@appl-ecosys.com>
wrote:
n=140652
null_cells=71271
cells=211923
Relatively hight number of null cells.
min=-1.89342188835144
max=0.00093058729544282
range=1.89435247564688
Probably small differences, but there is at least one outlier.
mean=-5.03205256742558e-05
mean_of_abs=5.24979771199771e-05
sum=-7.07768257713542
median=0
The overall differences are small or cancel out as they are distributed
around zero.
Looking at histogram may also help you interpret these statistics.
On Thu, Mar 30, 2017 at 2:46 AM, Rich Shepard <rshepard@appl-ecosys.com> wrote:
On Wed, 29 Mar 2017, Vaclav Petras wrote:
As Soeren mentioned in the other part of the thread, use r.univar to get
the statistics. Don’t rely on just looking at the resulting raster map.
There can be things like outliers which may throw the visualization off.Vaclav,
Here are the results of r.univar on the output of r.mapcalc “output_diff =
blocked_h_max - open_h_max”:[…]
min=-1.89342188835144
max=0.00093058729544282
[…][…] I see a small range and many null
cells, but I expect that because 4 cells on one input DEM are set about 0.6m
higher elevation than those on the other input DEM.
Apparently something went wrong because min should be 0 and max should be about 0.6.
The point locations used to edit the original DEM are probably wrong (not close to the center), and the computational region might have been changed somewhere in the process of editing the original DEM.
Markus M
On Thu, 30 Mar 2017, Nikos Alexandris wrote:
For the record: there is a "differences" color table which is useful for
this.
Don't forget to try r.colors with -e or -a -g if the first attempt does not
show much differences.
I _had_ not noticed the 'differences' color table. Applying that produces
interesting results (see attached .png).
The four blue cells are ~0.6m higher in the blocked map than in the open
culvert map. I don't know the reason for the lighter blue cell to the
southwest of the group of four. I suspect, but don't know, that the pinkish
cells represent slightly higher water accumulation in the blocked map than
in the open map.
Germane to what Markus M. wrote, "The point locations used to edit the
original DEM are probably wrong (not close to the center), ..." I now think
that's where I made one mistake. When I isolated those four cells and
queried them I truncated the easting and northing coordinates to the integer
value assuming that all elevations in a 1m x 1m cell were the same. If this
assumption is incorrect then I'll re-create the change and start the model
running again with a slightly different set of elevations for those cells.
If anyone can suggest other changes to the model input to better show
where water accummulates behind a blockage please let me know.
Rich
On Thu, 30 Mar 2017, Rich Shepard wrote:
Germane to what Markus M. wrote, "The point locations used to edit the
original DEM are probably wrong (not close to the center), ..." I now
think that's where I made one mistake.
Nope. I zoom in to isolate those 4 cells in in the blocked map. No matter
where within the cell I click the elevation is the same.
However, ... I'm now wondering if the issue is the elevation in feet
rather than meters.
The source DEMs were run through r.proj to convert from feet to meters.
This converted the coordinates but apparently not the elevations. If Itzi
wants all input in SI units then having elevations in feet with rainfall in
mm/hr would certainly get the model's briefs in a twist.
If having elevation data range from ~147 to ~168 feet rather than ~45m to
~51m would produce incorrect results then I suppose I need to multiply both
input maps (open and blocked) by 0.305 using r.mapcalc, then run both models
again.
Comments? Thoughts?
Rich
Rich,
If I understand well, you have a 'raw' DEM and you want to open the
culvert in that DEM.
If it is an under-road culvert, your 'raw' DEM will show the road
elevation, and the culvert will not be shown. i.e the culvert is
blocked.
If you want to 'open' the culvert, you'll have to cut through the
road. You can do this with r.carve, for example.
Regards,
Laurent
2017-03-30 10:04 GMT-06:00 Rich Shepard <rshepard@appl-ecosys.com>:
On Thu, 30 Mar 2017, Nikos Alexandris wrote:
For the record: there is a "differences" color table which is useful for
this.Don't forget to try r.colors with -e or -a -g if the first attempt does
not
show much differences.I _had_ not noticed the 'differences' color table. Applying that produces
interesting results (see attached .png).The four blue cells are ~0.6m higher in the blocked map than in the open
culvert map. I don't know the reason for the lighter blue cell to the
southwest of the group of four. I suspect, but don't know, that the pinkish
cells represent slightly higher water accumulation in the blocked map than
in the open map.Germane to what Markus M. wrote, "The point locations used to edit the
original DEM are probably wrong (not close to the center), ..." I now think
that's where I made one mistake. When I isolated those four cells and
queried them I truncated the easting and northing coordinates to the integer
value assuming that all elevations in a 1m x 1m cell were the same. If this
assumption is incorrect then I'll re-create the change and start the model
running again with a slightly different set of elevations for those cells.If anyone can suggest other changes to the model input to better show
where water accummulates behind a blockage please let me know.Rich
_______________________________________________
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user
On Thu, 30 Mar 2017, Laurent C. wrote:
If I understand well, you have a 'raw' DEM and you want to open the
culvert in that DEM. If it is an under-road culvert, your 'raw' DEM will
show the road elevation, and the culvert will not be shown. i.e the
culvert is blocked. If you want to 'open' the culvert, you'll have to cut
through the road. You can do this with r.carve, for example.
Laurent,
Interesting observation. On the source DEM the elevation of the culvert
inlet is about 2 feet lower than the road and the elevation of the culvert
outlet is about 3 feet lower than the road.
The attached .png shows the results of a 40-minute simulation of rainfall
accumulation. The two red dots are the culvert inlet and outlet. The
accumulation south of the road is because water pools in the swale below the
inlet until it reaches that elevation, then flows through and accumulates in
the deeper portion of the ditch on the north side.
I expected to see something like this using the actual rainfall data I
have.
Then, of course, if the culvert is blocked I expect water to accumulate
south of the road until it reaches that elevation, then overflow the road
and areas to the south.
Rich
On Thu, 30 Mar 2017, Laurent C. wrote:
If you want to 'open' the culvert, you'll have to cut through the road.
You can do this with r.carve, for example.
Certainly worth doing. The r.carve manual page requires as input a vector
map containing streams. Having only two points (culvert inlet and outlet), I
suppose the idea is to create a stream between those two points.
Would the process be to use v.edit on a copy (called 'stream') of the
culvert point map to generate a straight line between the two points? That
would provide the vector input r.carve wants.
Rich
On Thu, 30 Mar 2017, Rich Shepard wrote:
Would the process be to use v.edit on a copy (called 'stream') of the
culvert point map to generate a straight line between the two points? That
would provide the vector input r.carve wants.
How do I add a line between two specific points using the vector digitizer
in the GUI map display window?
Clicking on one point, then on the other point does not show a line
between the two points.
Clicking on one point and dragging the cursor to the other point before
releasing the button also does not display a line between the two points.
Clue stick needed.
Rich
2017-03-30 12:00 GMT-06:00 Rich Shepard <rshepard@appl-ecosys.com>:
On Thu, 30 Mar 2017, Laurent C. wrote:
If you want to 'open' the culvert, you'll have to cut through the road.
You can do this with r.carve, for example.Certainly worth doing. The r.carve manual page requires as input a vector
map containing streams. Having only two points (culvert inlet and outlet), I
suppose the idea is to create a stream between those two points.Would the process be to use v.edit on a copy (called 'stream') of the
culvert point map to generate a straight line between the two points? That
would provide the vector input r.carve wants.
Rich,
You can create a new vector map and create a line that link the two points.
However, r.carve will interpolate the elevation between the start
point and the end point.
If you culvert invert have a lower elevation, this will not show up.
But considering the amount of simplification you're already applying,
it could be acceptable.
Regards,
Laurent
2017-03-30 12:09 GMT-06:00 Rich Shepard <rshepard@appl-ecosys.com>:
On Thu, 30 Mar 2017, Rich Shepard wrote:
Would the process be to use v.edit on a copy (called 'stream') of the
culvert point map to generate a straight line between the two points? That
would provide the vector input r.carve wants.How do I add a line between two specific points using the vector digitizer
in the GUI map display window?Clicking on one point, then on the other point does not show a line
between the two points.Clicking on one point and dragging the cursor to the other point before
releasing the button also does not display a line between the two points.Clue stick needed.
You need first to create a new vector map (or edit an existing one)
and select "digitize new line".
Laurent
On Thu, 30 Mar 2017, Laurent C. wrote:
You need first to create a new vector map (or edit an existing one) and
select "digitize new line".
Laurent,
Oh. Didn't notice the combo box on the extreme left side. Sigh.
Thanks,
Rich
On Thu, 30 Mar 2017, Laurent C. wrote:
You can create a new vector map and create a line that link the two
points. However, r.carve will interpolate the elevation between the start
point and the end point.
Laurent,
Have this line in a vector map called 'stream.'
If you culvert invert have a lower elevation, this will not show up. But
considering the amount of simplification you're already applying, it could
be acceptable.
The invert has a slightly higher elevation than the outfall. Both are
lower than the elevation of the road over the culvert.
Will try this and see how it compares to the base map.
Thanks,
Rich
On Thu, 30 Mar 2017, Laurent C. wrote:
You can create a new vector map and create a line that link the two
points. However, r.carve will interpolate the elevation between the start
point and the end point.
Laurent,
If my visual interpretation of the output (detail of culvert attached) is
correct, this new DEM (with metric elevations) is suitable for input to
Itzi. The culvert depth is clearly less than the road surface.
Thanks,
Rich
