[GRASS-user] v.in.csv available?

2021-10-28 Thread Johannes Radinger

Hi all,

I was wondering if the addon v.in.csv is available or not? I want to 
import a simple point file (*.csv) which is unfortunately in another 
projection then my location. So the tool v.in.csv looks handy.


This tool is listed here https://grass.osgeo.org/grass78/manuals/addons/

but not here https://grasswiki.osgeo.org/wiki/AddOns/GRASS7/vector and 
seems not be available via g.extension(?)


Any suggestion? What would be the most convenient way to do such an 
import with coordinate transformation to match my current location?


Cheers,

Johannes

___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user


Re: [GRASS-user] Merge spatially connected features

2020-03-16 Thread Johannes Radinger
On Thu, Mar 12, 2020 at 10:48 PM Markus Metz 
wrote:

>
>
> On Thu, Mar 12, 2020 at 12:34 PM Johannes Radinger <
> johannesradin...@gmail.com> wrote:
> >
> > Thank you Markus,
> > indeed your approach looks like what I need..The hint with
> v.net.components was the part that I was missing;
>
> Note that v.net.components does not need a network prepared with v.net,
> you can use the extract of all lines with the same stream order as it is.
>
> I just tested your suggested approach using v.net.components
(without running v.net before)...works really smoothly and very fast and I
get to the desired results of extracting subnets. Thanks a lot Markus for
this hint!

/Johannes



> Markus M
>
> > I'll try as soon as possible and will report back on how this works.
> > cheers,
> > Johannes
> >
> > On Wed, Mar 11, 2020 at 10:16 PM Markus Metz <
> markus.metz.gisw...@gmail.com> wrote:
> >>
> >> Hi Johannes,
> >>
> >> IIUC, what you want to do is an operation that involves topological
> relations of vector geometries (connected lines) and a common attribute.
> There is no easy common recipe for this.
> >>
> >> Just a suggestion:
> >> for each stream order:
> >>   extract all lines with this stream order (v.extract)
> >>   identify connected lines (v.net + v.net.components)
> >>   update a new attribute of the original lines with the comp attribute
> of the output of v.net.components plus some offset to separate different
> stream orders
> >>
> >> HTH,
> >>
> >> Markus M
> >>
> >>
> >> On Tue, Mar 10, 2020 at 5:20 PM Johannes Radinger <
> johannesradin...@gmail.com> wrote:
> >> >
> >> > So...no also with GRASS-user as recipient...
> >> >
> >> > On 05.03.20 16:21, Micha Silver wrote:
> >> > >
> >> > > On 3/5/20 10:47 AM, Johannes Radinger wrote:
> >> > >>
> >> > >> Hi Micha, hi all,
> >> > >>
> >> > >> sorry for my late response...however, just today I managed to try
> >> > >> your approach of building polylines to connect "touching stream
> >> > >> lines"...but...
> >> > >>
> >> > >> On 24.02.20 16:48, Micha Silver wrote:
> >> > >>>
> >> > >>> On 24/02/2020 10:45, Johannes Radinger wrote:
> >> > >>>> Hi all,
> >> > >>>> I have a large river network dataset (lines). Now I'd to assign
> >> > >>>> unique categories to each group of connected lines that have an
> >> > >>>> attribute in common.
> >> > >>>>
> >> > >>>> For example, my rivers are categorized based on some kind of
> stream
> >> > >>>> order. I want to group all rivers that belong to stream order 2
> and
> >> > >>>> are spatially connected; each group should get a unique category
> >> > >>>> value. I thought that I could first extract all rivers with a
> >> > >>>> particular attribute (e.g. stream order = 2) which will provide
> me
> >> > >>>> some scattered pattern of lines. Then I need a spatial join tool
> to
> >> > >>>> make subgroups of lines that are connected. How can I achieve the
> >> > >>>> latter? Any idea?
> >> > >>>
> >> > >
> >> > >
> >> > >>>
> >> > >>> Here's a procedure that might work for you. Somewhat clunky, but I
> >> > >>> think it gets what you want.
> >> > >>>
> >> > >>> It's based on the v.build.polylines module to connect all touching
> >> > >>> stream reaches. First extract each order from the stream vector
> into
> >> > >>> a new vector. Then build polylines. Patch them all together. Now
> you
> >> > >>> have a polyline vector with a single cat value for each set of
> >> > >>> original stream reaches that had the same order and that were
> touching.
> >> > >>
> >> > >> Unfortunately, the v.build.polylines tool does not work as it only
> >> > >> does not connect multiple (intersecting) lines like in a river
> >> > >> network. As an example I tried to build polylines from the stream
> >> > >> network of the NC dataset. Yous suggested approach 

Re: [GRASS-user] v.edit break lines at coordinates

2020-03-14 Thread Johannes Radinger
Hi all,
FYI: I just issued a new ticket in github (
https://github.com/OSGeo/grass/issues/415) reporting that issue and
accompanied with some more information.
/Johannes

On Fri, Mar 13, 2020 at 11:53 AM Johannes Radinger <
johannesradin...@gmail.com> wrote:

> Hi GRASS users,
>
> maybe someone of you can explain me following problem (of understanding) I
> am facing:
> I want to break a single line into several. More specifically, I am using
> v.edit with the option break and two coordinate pairs to break a single
> line. I'd expect that breaking a single line at two points would result in
> three lines; However this results in 4 categories?!
>
> Here a small example using the NC dataset:
> #
> # Extract single line for testing
> v.extract --o input=streams@PERMANENT cats=92551 output=selected_stream
>
> # Break line at two point locations and clean (recommended after v.edit
> break)
> v.edit map=selected_stream tool=break threshold=30
> coords='635861,224469,635793,224553'
> v.clean --o input=selected_stream output=selected_stream2 tool=rmdupl
>
> # Add categories to layer 2 to check
> v.category --o input=selected_stream2 layer=2 output=selected_stream3
> option=add
> v.db.addtable map=selected_stream3 layer=2
> d.vect -c map=selected_stream3@test4 layer=2 label_layer=2
> attribute_column=cat
> 
>
> So the display shows three colors which makes sense as these would be the
> three lines that I'd expect from breaking a single line at two points. But
> there are four categories (that are also listed in the attribute table).
> Can anyone explain that to me, what I am doing wrong; or even better, what
> would you recommend to do if I want to break a line at two coordinates into
> three lines?
>
> Thank you for your help!
> Johannes
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] v.edit break lines at coordinates

2020-03-13 Thread Johannes Radinger
Hi GRASS users,

maybe someone of you can explain me following problem (of understanding) I
am facing:
I want to break a single line into several. More specifically, I am using
v.edit with the option break and two coordinate pairs to break a single
line. I'd expect that breaking a single line at two points would result in
three lines; However this results in 4 categories?!

Here a small example using the NC dataset:
#
# Extract single line for testing
v.extract --o input=streams@PERMANENT cats=92551 output=selected_stream

# Break line at two point locations and clean (recommended after v.edit
break)
v.edit map=selected_stream tool=break threshold=30
coords='635861,224469,635793,224553'
v.clean --o input=selected_stream output=selected_stream2 tool=rmdupl

# Add categories to layer 2 to check
v.category --o input=selected_stream2 layer=2 output=selected_stream3
option=add
v.db.addtable map=selected_stream3 layer=2
d.vect -c map=selected_stream3@test4 layer=2 label_layer=2
attribute_column=cat


So the display shows three colors which makes sense as these would be the
three lines that I'd expect from breaking a single line at two points. But
there are four categories (that are also listed in the attribute table).
Can anyone explain that to me, what I am doing wrong; or even better, what
would you recommend to do if I want to break a line at two coordinates into
three lines?

Thank you for your help!
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] v.in.ogr layer name

2020-03-13 Thread Johannes Radinger
Hi Markus,

On Thu, Mar 12, 2020 at 10:19 PM Markus Metz 
wrote:

>
>
> On Thu, Mar 12, 2020 at 12:43 PM Johannes Radinger <
> johannesradin...@gmail.com> wrote:
> >
> > Hi all,
> > I am using v.in.ogr to import a shape file (line vector) that has only
> one layer. Although I am specifying the output parameter in v.in.ogr the
> layer name of the imported map still remains the original; it is only the
> output map and the associated table that gets a new name as defined by
> 'output'.
> >
> > When checking the connections of that map using v.db.connect:
> > v.db.connect -p map=my_new_map@analysis2
>
> > Vector map  is connected by:
> > layer <1/Watercourse> table  in database
>  through driver  with key 
> >
> > Is that the intended behavior of v.in.ogr that the layer name does not
> get renamed based on the output parameter?
>
> Not sure how much intention is in this behaviour, but the GRASS-internal
> layer name is mostly informative. GRASS vector layers are usually addressed
> by number, not by name. It is probably safe to ignore the layer name of
> native GRASS vectors.
>

I agree that the layer name is probably mostly informative (at least I
didn't realize so far how/where it is involved in daily work with GRASS).
However, during some recent analysis I discovered that the layer name
affected how tables get renamed when new additional tables in new layers
are created... Then it is weird when the old initial layer name pops up in
a table name.
To better illustrate the issue I generated a working example to reproduce
what I mean. This is based on an INSPIRE GML file where the layer name did
not get renamed during import and where it affected the (re)naming of a
table:


# Create location with EPSG:4258
# Download GML file from
https://geoportal.bafg.de/inspire/download/reporting_units/WFDDesignatedWatersFish/WFDDesignatedWatersFish_DE.gml

# Import vector data with user-defined output name
v.in.ogr input=/home/jradinger/Downloads/WFDDesignatedWatersFish_DE.gml
output=Fishdata
g.region vector=Fishdata

# Check connection parameters
v.db.connect -p map=Fishdata@PERMANENT

#Vector map  is connected by:
#layer <1/WFDDesignatedWatersFish_LIN> table  in database

through driver  with key 
### Layer name got not renamed!

# To reproduce: Add new cat and table on layer 2 and layer 3
v.category input=Fishdata layer=2 output=Fishdata2 option=add
v.db.addtable map=Fishdata2 layer=2
v.category input=Fishdata2 layer=3 output=newFishdata option=add
v.db.addtable map=newFishdata layer=3

db.tables -p
# A original table in layer 1 was renamed to
newFishdata_WFDDesignatedWatersFish_LIN using the initial layer name
#

Don't know the underlying mechanisms how tables get renamed

/Johannes


>
> Markus M
>
> >
> > cheers,
> > Johannes
> >
> >
> > ___
> > grass-user mailing list
> > grass-user@lists.osgeo.org
> > https://lists.osgeo.org/mailman/listinfo/grass-user
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Meaning of flag 'i' in d.vect

2020-03-12 Thread Johannes Radinger
Hi all,

I am just wondering what exactly the i-flag in d.vector means? It is
described as "Use values from 'cats' option as feature id". What is the
effect of this in the display?

Furthermore, I was wondering how the query tool in the map display works.
For example, I display a vector map with 3 layers (I chose to display only
layer 3 with d.vect layer=3). However, the query tool selects more features
than I clicked on; it seems it rather selects the features based on the
category value of the first layer? So the behaviour I'd expect is that when
only one layer is selected in the display (and for example colorized based
on cat, c-flag) than the query click would select only the feature of this
specific category?!

Is there a mistake I am making or some missunderstanding of the concept of
the query tool? Is is this something that could be improved?

cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] v.in.ogr layer name

2020-03-12 Thread Johannes Radinger
Hi all,
I am using v.in.ogr to import a shape file (line vector) that has only one
layer. Although I am specifying the output parameter in v.in.ogr the layer
name of the imported map still remains the original; it is only the output
map and the associated table that gets a new name as defined by 'output'.

When checking the connections of that map using v.db.connect:
v.db.connect -p map=my_new_map@analysis2

Vector map  is connected by:
layer <1/Watercourse> table  in database
 through driver  with key 

Is that the intended behavior of v.in.ogr that the layer name does not get
renamed based on the output parameter?

cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Merge spatially connected features

2020-03-12 Thread Johannes Radinger
Thank you Markus,
indeed your approach looks like what I need..The hint with v.net.components
was the part that I was missing; I'll try as soon as possible and will
report back on how this works.
cheers,
Johannes

On Wed, Mar 11, 2020 at 10:16 PM Markus Metz 
wrote:

> Hi Johannes,
>
> IIUC, what you want to do is an operation that involves topological
> relations of vector geometries (connected lines) and a common attribute.
> There is no easy common recipe for this.
>
> Just a suggestion:
> for each stream order:
>   extract all lines with this stream order (v.extract)
>   identify connected lines (v.net + v.net.components)
>   update a new attribute of the original lines with the comp attribute of
> the output of v.net.components plus some offset to separate different
> stream orders
>
> HTH,
>
> Markus M
>
>
> On Tue, Mar 10, 2020 at 5:20 PM Johannes Radinger <
> johannesradin...@gmail.com> wrote:
> >
> > So...no also with GRASS-user as recipient...
> >
> > On 05.03.20 16:21, Micha Silver wrote:
> > >
> > > On 3/5/20 10:47 AM, Johannes Radinger wrote:
> > >>
> > >> Hi Micha, hi all,
> > >>
> > >> sorry for my late response...however, just today I managed to try
> > >> your approach of building polylines to connect "touching stream
> > >> lines"...but...
> > >>
> > >> On 24.02.20 16:48, Micha Silver wrote:
> > >>>
> > >>> On 24/02/2020 10:45, Johannes Radinger wrote:
> > >>>> Hi all,
> > >>>> I have a large river network dataset (lines). Now I'd to assign
> > >>>> unique categories to each group of connected lines that have an
> > >>>> attribute in common.
> > >>>>
> > >>>> For example, my rivers are categorized based on some kind of stream
> > >>>> order. I want to group all rivers that belong to stream order 2 and
> > >>>> are spatially connected; each group should get a unique category
> > >>>> value. I thought that I could first extract all rivers with a
> > >>>> particular attribute (e.g. stream order = 2) which will provide me
> > >>>> some scattered pattern of lines. Then I need a spatial join tool to
> > >>>> make subgroups of lines that are connected. How can I achieve the
> > >>>> latter? Any idea?
> > >>>
> > >
> > >
> > >>>
> > >>> Here's a procedure that might work for you. Somewhat clunky, but I
> > >>> think it gets what you want.
> > >>>
> > >>> It's based on the v.build.polylines module to connect all touching
> > >>> stream reaches. First extract each order from the stream vector into
> > >>> a new vector. Then build polylines. Patch them all together. Now you
> > >>> have a polyline vector with a single cat value for each set of
> > >>> original stream reaches that had the same order and that were
> touching.
> > >>
> > >> Unfortunately, the v.build.polylines tool does not work as it only
> > >> does not connect multiple (intersecting) lines like in a river
> > >> network. As an example I tried to build polylines from the stream
> > >> network of the NC dataset. Yous suggested approach should result that
> > >> each sub-network (i.e. river network that is not connected to another
> > >> one) should get its own ID/cat...however, v.build.polylines results
> > >> in a connected stream network that consists of multiple cats:
> > >>
> > > Maybe I misunderstood your question. The steps I tried use a
> > > stream_order column to group stream segments, then apply a new
> > > attribute "merged_id" to those stream orders that touch. i.e. that
> > > connect to the same confluence point.
> > >
> > >
> > > Here's what I get using the nc_basic_spm mapset:
> > >
> > >
> > > r.watershed elev=elevation accum=nc_facc drain=nc_fdir bas=nc_bas
> > > stream=nc_str thresh=1000
> > > r.stream.order stream_rast=nc_str direct=nc_fdir elev=elevation
> > > accum=nc_facc stream_vect=nc_streams
> > > ORDERS=`v.db.select -c nc_streams group=strahler column=strahler`
> > > echo $ORDERS
> > >
> > > # Create a new stream vector for each stream order
> > >
> > > for o in $ORDERS; do
> > >
> > > v.extract input=nc_streams output=streams_${o}
> where="strahler=$

Re: [GRASS-user] Merge spatially connected features

2020-03-10 Thread Johannes Radinger

So...no also with GRASS-user as recipient...

On 05.03.20 16:21, Micha Silver wrote:


On 3/5/20 10:47 AM, Johannes Radinger wrote:


Hi Micha, hi all,

sorry for my late response...however, just today I managed to try 
your approach of building polylines to connect "touching stream 
lines"...but...


On 24.02.20 16:48, Micha Silver wrote:


On 24/02/2020 10:45, Johannes Radinger wrote:

Hi all,
I have a large river network dataset (lines). Now I'd to assign 
unique categories to each group of connected lines that have an 
attribute in common.


For example, my rivers are categorized based on some kind of stream 
order. I want to group all rivers that belong to stream order 2 and 
are spatially connected; each group should get a unique category 
value. I thought that I could first extract all rivers with a 
particular attribute (e.g. stream order = 2) which will provide me 
some scattered pattern of lines. Then I need a spatial join tool to 
make subgroups of lines that are connected. How can I achieve the 
latter? Any idea?







Here's a procedure that might work for you. Somewhat clunky, but I 
think it gets what you want.


It's based on the v.build.polylines module to connect all touching 
stream reaches. First extract each order from the stream vector into 
a new vector. Then build polylines. Patch them all together. Now you 
have a polyline vector with a single cat value for each set of 
original stream reaches that had the same order and that were touching.


Unfortunately, the v.build.polylines tool does not work as it only 
does not connect multiple (intersecting) lines like in a river 
network. As an example I tried to build polylines from the stream 
network of the NC dataset. Yous suggested approach should result that 
each sub-network (i.e. river network that is not connected to another 
one) should get its own ID/cat...however, v.build.polylines results 
in a connected stream network that consists of multiple cats:


Maybe I misunderstood your question. The steps I tried use a 
stream_order column to group stream segments, then apply a new 
attribute "merged_id" to those stream orders that touch. i.e. that 
connect to the same confluence point.



Here's what I get using the nc_basic_spm mapset:


r.watershed elev=elevation accum=nc_facc drain=nc_fdir bas=nc_bas 
stream=nc_str thresh=1000
r.stream.order stream_rast=nc_str direct=nc_fdir elev=elevation 
accum=nc_facc stream_vect=nc_streams

ORDERS=`v.db.select -c nc_streams group=strahler column=strahler`
echo $ORDERS

# Create a new stream vector for each stream order

for o in $ORDERS; do

    v.extract input=nc_streams output=streams_${o} where="strahler=${o}"

    # Give each polyline it's own cat value

    v.build.polylines input=streams_${o} output=streams_${o}_polyline 
type=line cat=first


done


# patch the stream orders back together

POLYLINES=`g.list vect pattern="streams*polyline" separator=comma`

v.patch input=$POLYLINES output=streams_polylines

v.db.addcolumn map=streams column="merged_id INTEGER"


# And use v.distance to update that merged_id column from cat values 
in polylines vector

v.distance from=streams to=streams_polylines upload=cat column=merged_id
v.db.addcolumn map=nc_streams column="merged_id INTEGER"
v.distance from=nc_streams to=streams_polylines upload=cat 
column=merged_id


Now, all stream reaches that have the same order and are "touching" 
have the same merged_id. See the attached image.



If that's not your purpose, then just ignore...

Micha thank you for your help and of course, you're fully correct! 
Merging lines that belong to the same stream order works in this case 
well...but this is because of the definition of the Strahler ordering 
system, where there is only one "touching node" (i.e. river junction) of 
two rivers of the same stream order (i.e. when two 2nd order streams 
meet, the become a 3rd order stream). Thus your solution works because 
of this specifics and might not work if streams are grouped based on a 
different (ordering) system.


I was already thinking of the next step (beyond simple Strahler): As 
mentioned in my initial post I am dealing with "some kind" of stream 
order. It is similar to grouped stream orders (e.g. stream order 1-2 = 
"headwater streams"). I tried to somehow reproduce my situation based on 
your example of the NC dataset. What I basically did was to reassign a 
new stream order "99" to all former 1st and 2nd order streams. Then I 
did exactly what you did in your example, and of course I don't unique 
merged_ids for the subnetworks of touching lines (see attached Figs) 
that all belong the the same "order" 99 (the original strahler order 3 
works of course, see Fig.)...So is there a more general way (as said 
something like v.dissolve but for lines/networks?):


#
g.region raster=elevation

r.watershed --o e

Re: [GRASS-user] Merge spatially connected features

2020-03-05 Thread Johannes Radinger
Hi Micha, hi all,

sorry for my late response...however, just today I managed to try your
approach of building polylines to connect "touching stream lines"...but...
On 24.02.20 16:48, Micha Silver wrote:


On 24/02/2020 10:45, Johannes Radinger wrote:

Hi all,
I have a large river network dataset (lines). Now I'd to assign unique
categories to each group of connected lines that have an attribute in
common.

For example, my rivers are categorized based on some kind of stream order.
I want to group all rivers that belong to stream order 2 and are spatially
connected; each group should get a unique category value. I thought that I
could first extract all rivers with a particular attribute (e.g. stream
order = 2) which will provide me some scattered pattern of lines. Then I
need a spatial join tool to make subgroups of lines that are connected. How
can I achieve the latter? Any idea?



Here's a procedure that might work for you. Somewhat clunky, but I think it
gets what you want.

It's based on the v.build.polylines module to connect all touching stream
reaches. First extract each order from the stream vector into a new vector.
Then build polylines. Patch them all together. Now you have a polyline
vector with a single cat value for each set of original stream reaches that
had the same order and that were touching.

Unfortunately, the v.build.polylines tool does not work as it only does not
connect multiple (intersecting) lines like in a river network. As an
example I tried to build polylines from the stream network of the NC
dataset. Yous suggested approach should result that each sub-network (i.e.
river network that is not connected to another one) should get its own
ID/cat...however, v.build.polylines results in a connected stream network
that consists of multiple cats:

v.clean --overwrite input=streams@PERMANENT output=streams_break tool=break
v.build.polylines --overwrite input=streams_break@test output=streams_poly
cats=first type=line
d.vect -c map=streams_poly

So what would be needed here is some kind of tool that connects all
touching lines and assigns a common category value, similar to the
v.dissolve tool for polygon features. I can imagine that such a task might
be not that uncommon also in another context? Any suggestions how to
achieve this in GRASS?

A workaround that came into my mind was to create buffers around lines in
order to make areas out of lines. Subsequently these touching areas can be
merged using v.dissolve and the information about the common category can
be queried using v.distance. Nevertheless, a rather cumbersome way to just
assign a common category value to all lines that are touching...

Any further ideas?

cheers,

Johannes


Finally, with v.distance you can upload that cat value to the original
streams.


# Get a list of stream orders
ORDERS=`v.db.select -c streams group=strahler column=strahler`
echo $ORDERS
#1 2 3 4 5 6
# How many stream segments in original
v.info -t streams | grep lines
# lines=1420

# Now loop thru list of stream orders and extract stream segments for each
order
for o in $ORDERS; do
v.extract input=streams output=streams_${o} where="strahler=${o}"
# Create polyline for each stream order
# Line "connects" all touching stream segments
v.build.polylines input=streams_${o} output=streams_${o}_polyline
type=line cat=first
done

# Patch stream order polylines together
POLYLINES=`g.list vect pattern="streams*polyline" separator=comma`
echo $POLYLINES
v.patch input=$POLYLINES output=streams_polylines
v.info -t streams_polylines | grep lines
# lines=738

# Add a new column to the original streams for new ID value
v.db.addcolumn map=streams column="merged_id INTEGER"
# And use v.distance to update that column from cat values in polylines
vector
v.distance from=streams to=streams_polylines upload=cat column=merged_id

HTH


Cheers,
Johannes

___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Merge spatially connected features

2020-03-02 Thread Johannes Radinger

Hi Micha, hi all,

sorry for my late response...however, just today I managed to try your 
approach of building polylines to connect "touching stream lines"...but...


On 24.02.20 16:48, Micha Silver wrote:


On 24/02/2020 10:45, Johannes Radinger wrote:

Hi all,
I have a large river network dataset (lines). Now I'd to assign 
unique categories to each group of connected lines that have an 
attribute in common.


For example, my rivers are categorized based on some kind of stream 
order. I want to group all rivers that belong to stream order 2 and 
are spatially connected; each group should get a unique category 
value. I thought that I could first extract all rivers with a 
particular attribute (e.g. stream order = 2) which will provide me 
some scattered pattern of lines. Then I need a spatial join tool to 
make subgroups of lines that are connected. How can I achieve the 
latter? Any idea?



Here's a procedure that might work for you. Somewhat clunky, but I 
think it gets what you want.


It's based on the v.build.polylines module to connect all touching 
stream reaches. First extract each order from the stream vector into a 
new vector. Then build polylines. Patch them all together. Now you 
have a polyline vector with a single cat value for each set of 
original stream reaches that had the same order and that were touching.


Unfortunately, the v.build.polylines tool does not work here as it does 
not connect multiple (intersecting) lines like in a river network. For 
example, I tried to build polylines from the stream network of the NC 
dataset. The suggested approach should result that each sub-network 
(i.e. river network that is not connected to another one) should get its 
own ID/cat...however, v.build.polylines results in a  stream sub-network 
that consists of multiple cats:


v.clean --overwrite input=streams@PERMANENT output=streams_break tool=break
v.build.polylines --overwrite input=streams_break@test 
output=streams_poly cats=first type=line

d.vect -c map=streams_poly

So is there any other way to join all lines that are touching...This 
would be something similar to v.dissolve but based on lines rather than 
on polygons.


Any ideas for a simple work around?

cheers,

Johannes




Finally, with v.distance you can upload that cat value to the original 
streams.



# Get a list of stream orders
ORDERS=`v.db.select -c streams group=strahler column=strahler`
echo $ORDERS
#1 2 3 4 5 6
# How many stream segments in original
v.info -t streams | grep lines
# lines=1420

# Now loop thru list of stream orders and extract stream segments for 
each order

for o in $ORDERS; do
    v.extract input=streams output=streams_${o} where="strahler=${o}"
    # Create polyline for each stream order
    # Line "connects" all touching stream segments
    v.build.polylines input=streams_${o} 
output=streams_${o}_polyline type=line cat=first

done

# Patch stream order polylines together
POLYLINES=`g.list vect pattern="streams*polyline" separator=comma`
echo $POLYLINES
v.patch input=$POLYLINES output=streams_polylines
v.info -t streams_polylines | grep lines
# lines=738

# Add a new column to the original streams for new ID value
v.db.addcolumn map=streams column="merged_id INTEGER"
# And use v.distance to update that column from cat values in 
polylines vector

v.distance from=streams to=streams_polylines upload=cat column=merged_id

HTH


Cheers,
Johannes

___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user



___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Splitting lines by other lines that overlap

2020-02-28 Thread Johannes Radinger
Hi all,

I have two line vectors: LV1) a complete stream network LV2) sections of
the same stream network representing impoundments.

LV2 is a subset of LV1 and fully spatially overlapping. However, the full
stream network consists of polylines with start/end nodes that do not
correspond to the impoundments vector:

LV1 (complete network):
+---+--+
LV2 (subset with different start/end nodes)
 +---+   +-+

So how can I get the information of impoundments (LV2) into the attribute
table of the first line vector (LV1). A way that came into my mind but
which is somehow cumbersome:

1) Extract the start/end node coordinates each impoundment (LV2) using
v.to.db
2) Use these coordinates in v.edit to break the river network (LV1)
3) Add new categories in layer 2 of LV1 and add an attribute table
4) Use v.distance to query the information from LV2 and copy it to the
corresponding new reaches in LV1
..This is not yet tested and not the most straight forward way.

All this I guess would be easier if there would be something like v.overlay
that works with two line vectors...

Has anybody done similar analysis in GRASS so far?

cheers,
Johannes

tart/end nodes of the
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Merge spatially connected features

2020-02-24 Thread Johannes Radinger
Hi all,
I have a large river network dataset (lines). Now I'd to assign unique
categories to each group of connected lines that have an attribute in
common.

For example, my rivers are categorized based on some kind of stream order.
I want to group all rivers that belong to stream order 2 and are spatially
connected; each group should get a unique category value. I thought that I
could first extract all rivers with a particular attribute (e.g. stream
order = 2) which will provide me some scattered pattern of lines. Then I
need a spatial join tool to make subgroups of lines that are connected. How
can I achieve the latter? Any idea?

Cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Error when importing gml file

2020-02-21 Thread Johannes Radinger
Hi all,
I am trying to import a *.gml file which I received from a governmental
agency. However when trying v.import I get an error message (see below). I
can successfully import the file in QGIS (and display the attribute table).
Other gml files from the same region and from the same provider work also.

I can imagine it might be related to some attribute names?! Is the
attribute name 'order' not a valid column name as ORDER is also an sqlite
command? Anybody any ideas what's wrong with the file or the import?

cheers,
Johannes

#
Importing  ...
Check if OGR layer  contains polygons...
Creating attribute table for layer ...
Default driver / database set to:
driver: sqlite
database: $GISDBASE/$LOCATION_NAME/$MAPSET/sqlite/sqlite.db
Column name  renamed to

Column name  renamed to

DBMI-SQLite driver error:
Error in sqlite3_prepare():
near "order": syntax error
DBMI-SQLite driver error:
Error in sqlite3_prepare():
near "order": syntax error
ERROR: Unable to create table: 'create table Watercourse (cat integer,
gml_id text, language varchar ( 3 ), nativeness text, sourceOfName text,
pronunciation text, text varchar ( 90 ), script varchar ( 4 ),
transliterationScheme text, grammaticalGender text, grammaticalNumber text,
classificationScheme varchar ( 8 ), localId integer, namespace varchar ( 83
), beginLifespanVersion varchar ( 25 ), inspireId_Identifier_localId
integer, inspireId_Identifier_namespace varchar ( 83 ), versionId varchar (
5 ), origin text, persistence text, tidal text, drainsBasin text,
delineationKnown integer, length double precision, length_uom varchar ( 1
), level text, order varchar ( 36 ), orderScheme varchar ( 83 ), scope
varchar ( 8 ), width text)'
ERROR: Unable to import OGR datasource 


My system:
GRASS version: 7.9.dev

Code revision: ac8bd2777

Build date: 2020-01-21

Build platform: x86_64-pc-linux-gnu

GDAL: 2.2.3

PROJ: 4.9.3

GEOS: 3.6.2

SQLite: 3.22.0

Python: 3.6.9

wxPython: 4.0.1
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] r.sim.water simulation considering dams

2020-01-23 Thread Johannes Radinger
Dear GRASS users,

maybe someone of you has already worked with r.sim.water while considering
dams as impermeable structures that retain water, i.e. create
impoundments/reservoirs.

I am wondering how to best specify dams within r.sim.water: There is the
parameter flow_control which is a raster map that indicates the
permeability of the landscape. The manual mentions here a permeability
ratio (0-1) and I guess this should be 0 for all cells that represent a dam
and 1 for all other cells of the computational region?!
But what about the elevation data then; Do I need to provide an
elevation raster that includes the dam height or should this be the
elevation as without dams?
I'd like to know whether I need to provide a flow_control map or an
elevation map that includes dams or both? I couldn't find any example of
r.sim.water that uses the flow_control parameter so far... any
ideas/experiences from your side?

cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Create line perpendicular to other line at specific point

2019-11-08 Thread Johannes Radinger
Hi all,
I've vector lines representing rivers and a set of points. Now, I'd like to
create short lines (e.g. 150 m) that are perpendicular to the river lines
and cross at the coordinates of the points vector. This is basically to
create transects similar to the tool v.transect. However v.transect creates
the perpendicular lines at regular intervals and I'd like to create them
only for the locations of the points on the river line. Is there any
straight forward tool to accomplish this? I could retrieve the azimuth of
the river lines at each point (using v.to.db), but how can I extrude a
point into a line with a given length and azimuth? Any suggestions?

cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Import CORDEX data in GRASS

2019-10-23 Thread Johannes Radinger
Hi Markus, hi all,

thanks for your help!

I just uploaded two example data files (modelled precipitation at 0.11
degrees) of EURO CORDEX data with rotated-pole grids to play with:
https://we.tl/t-ZaZWaPyIaA

Attached also a description of the EURO CORDEX data and the grid domains.

Typically such data are processed and manipulated with the CDO toolset
(Climate Data Operators, https://code.mpimet.mpg.de/projects/cdo, see
attached in the shared folder). CDO can easily load the netcdf data and
even remap it to other grids (although I am not sure how that works
exactly); So maybe another way (if rotated pole is not possible with
GRASS/gdal) I could try to transform using CDO and then import into a GRASS
location

/Johannes

On Mon, Oct 21, 2019 at 10:39 PM Markus Neteler  wrote:

> Hi Johannes,
>
> On Fri, Oct 18, 2019 at 2:57 PM Johannes Radinger
>  wrote:
> >
> > Hi all,
> >
> > I was wondering if anybody is also working with CORDEX data (Regionally
> downscaled climate data - Europe) that are provided via the ESGF platforms
> as netcdf files?
> >
> > I want to import such files in GRASS but first want to create a location
> with the projection of these data. AFAIK the data are in a rotated-pole
> coordinate system (
> https://is-enes-data.github.io/cordex_archive_specifications.pdf) which
> should be supported by GDAL (https://trac.osgeo.org/gdal/ticket/4285)
> however I am unsure in which version of gdal (here I am using 2.2.2 (Ubuntu
> 16LTS, GRASS76)
> >
> > I tried to use the location wizard to read the projection information
> from the netcdf file. However, here only a XY location (unprojected) is
> suggested. So how (and what) do I need to manually specify to define a
> Location that is able to correctly display CORDEX data? Any suggestions?
>
> Is there anywhere a file to play with?
>
> best
> Markus
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Import CORDEX data in GRASS

2019-10-18 Thread Johannes Radinger
Hi all,

I was wondering if anybody is also working with CORDEX data (Regionally
downscaled climate data - Europe) that are provided via the ESGF platforms
as netcdf files?

I want to import such files in GRASS but first want to create a location
with the projection of these data. AFAIK the data are in a rotated-pole
coordinate system (
https://is-enes-data.github.io/cordex_archive_specifications.pdf) which
should be supported by GDAL (https://trac.osgeo.org/gdal/ticket/4285)
however I am unsure in which version of gdal (here I am using 2.2.2 (Ubuntu
16LTS, GRASS76)

I tried to use the location wizard to read the projection information from
the netcdf file. However, here only a XY location (unprojected) is
suggested. So how (and what) do I need to manually specify to define a
Location that is able to correctly display CORDEX data? Any suggestions?

cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] v.stream.* questions

2019-10-04 Thread Johannes Radinger

Hi Rich,

AFAIK, the tool v.stream.order was developed to work with grass vector 
files only (no raster input needed). So, the input needed for 
v.stream.order is the vector lines representing the streams (lines must 
be connected) and points that define the outlet of each river network 
that should be processed. I haven't tested v.stream.order for a while, 
so I am not sure if it still works with newer versions of GRASS?! 
However, tool comes with a testsuite (see source code) and a working 
example based on the NC-dataset to test operability of v.stream.order.


HTH,

Johannes

On 03.10.19 22:56, Rich Shepard wrote:

Reading the manual pages for v.stream.network, v.stream.order, and
v.stream.profiler in the extensions doc I see details for the input 
options

but all descriptions refer to raster files.

Are there more detailed docs that show examples using the vector modules
along with the raster modules?

I've a complex, extensive, river network as a vector file. I've not
extracted streams, basins, or watersheds from the DEM and I'd like to
run the v.stream.* modules on my vector stream network. Possible?

TIA,

Rich
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Transparent points in ps.map output

2019-05-09 Thread Johannes Radinger
Hi all,

I'm using ps.map to produce nice GRASS maps that combine some vector lines
(vlines) and points (vpoints).

For the points I want to distinguish between two types (presence and
absence of species at a location) by using different colouring/filling: One
set of vpoints should be printed as black-filled circles the other set
should be transparent circles (i.e. just the outer line, no filling). If
possible, I'd like to provide this information via the 'rgbcolumn' argument
of vpoints in ps.map.

While using a value of '0:0:0' in the rgbcolum for black-filled circles
works well, using "none" in the rgbcolumn does not work to produce
transparent circles, i.e. circles with no filling. This rather raises an
error "Invalid RGB colour definition in column". I also couldn't find any
answer here:
https://grass.osgeo.org/grass70/manuals/ps.map.html#NAMED_COLORS

How can I set a point transparent (i.e. no filling) in ps.map?

Thanks for yours answers!
/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] ps.map with interactive map instructions using python

2019-03-04 Thread Johannes Radinger
Hi all,

I want to use ps.map in python to create *.eps maps from a GRASS vector
line using some interactive mapping instructions, i.e. I don't want to save
my mapping instructions into a *.txt file but rather provide the
instructions as a python string.

For example my mapping instructions:
map_instructions = """
# my instructions
vlines my_lines
  type line
  layer 2
  color color
  rgbcolumn rgb_column
  width 2
  end
"""
# And now I'd like to use this string as input to ps.map within a python
command:
grass.run_command("ps.map",
  flags="e",
  input="-",
  stdin_=map_instructions,
  output="myoutput.eps")

However, stdin is not the correct way to specify this...what would be the
correct way to avoid saving the instructions as text but rather using the
python text-string as direct input to ps.map within python? At least the
GUI of ps.map allows to use interactive input for mapping instructions...
Any suggestions?

Cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] v.extract and layer/table settings with layer=2

2018-11-28 Thread Johannes Radinger
Thank you Moritz for clarification...
...I also issued an enhancement ticket for v.extract to include a new
output_layer parameter (https://trac.osgeo.org/grass/ticket/3700#ticket)

/J

On Wed, Nov 28, 2018 at 8:26 AM Moritz Lennert 
wrote:

>
>
> Am 27. November 2018 20:51:37 MEZ schrieb Johannes Radinger <
> johannesradin...@gmail.com>:
> >Hi all,
> >
> >this is maybe I quite easy routine for somebody that is using v.extract
> >regularly.
> >
> >Here my challenge: I have a vector map (initially created by v.net
> >tools)
> >that has several layers and associated tables to each layer. Now for
> >further processing I want to extract only some features from layer=2
> >and
> >the corresponding attribute table.
> >
> >I can do that using "v.extract input=mymap1 layer=2 where='mycol=x'
> >output=mymap2". The resulting 'mymap2' has an attribute table connected
> >to
> >layer 2 and there are no other layers (e.g. there is no layer 1, which
> >is
> >considered in most grass tools as default layer). So I wanted to change
> >that and to shift layer 2 (and the attribute table) to layer 1 using
> >"v.category option=chlayer layer='2,1' input=mymap2 output=mymap3"
> >
> >...however, this does not change the attribute table connection (as
> >checked
> >by "v.db.connect -p "). So I guess I also have to use v.db.connect in
> >some
> >way (??) to get the attribute table shifted from connection at former
> >layer
> >2 to layer 1 (unused).
>
> Yes:
>
> v.db.connect -d layer=2
> v.db.connect layer=1 table=
>
>
> >
> >In my opinion this is quite a bit complicated, considering that I just
> >wanted to extract some features from a 3-layered map?? Is there a
> >easier
> >way or a general description if one want to extract some feature that
> >are
> >not stored in layer 1?
>
> The best would probably be some form of output_layer parameter in
> v.extract, probably with a default of 1. Worth a trac ticket I think.
>
> Moritz
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] v.extract and layer/table settings with layer=2

2018-11-27 Thread Johannes Radinger
Hi all,

this is maybe I quite easy routine for somebody that is using v.extract
regularly.

Here my challenge: I have a vector map (initially created by v.net tools)
that has several layers and associated tables to each layer. Now for
further processing I want to extract only some features from layer=2 and
the corresponding attribute table.

I can do that using "v.extract input=mymap1 layer=2 where='mycol=x'
output=mymap2". The resulting 'mymap2' has an attribute table connected to
layer 2 and there are no other layers (e.g. there is no layer 1, which is
considered in most grass tools as default layer). So I wanted to change
that and to shift layer 2 (and the attribute table) to layer 1 using
"v.category option=chlayer layer='2,1' input=mymap2 output=mymap3"

...however, this does not change the attribute table connection (as checked
by "v.db.connect -p "). So I guess I also have to use v.db.connect in some
way (??) to get the attribute table shifted from connection at former layer
2 to layer 1 (unused).

In my opinion this is quite a bit complicated, considering that I just
wanted to extract some features from a 3-layered map?? Is there a easier
way or a general description if one want to extract some feature that are
not stored in layer 1?

All the best,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Cut a floodplain polygon with 1meter wide polygons perpendicular to its river channel

2018-10-15 Thread Johannes Radinger
It seems that would would like to have your floodplain cut based on the
columns of the underlying raster-resolution (and not perpendicular to the
river line as initially stated). You could use mapcalc to generate a column
raster map (r.mapcalc expression="colmap = col()") which you could then
transform to a vector map and use to intersect with your floodplain vector.

best
Johannes

On Tue, Oct 16, 2018 at 1:24 AM Shane Carey  wrote:

> Hi All,
>
> I don't think v.transects is going to work for me as lines appear to be
> cutting across each other. Is there another tool in Grass that can cut the
> floodplain polygon into 5m lengths (as green lines in image attached)
>
> Thanks
>
>
>
> On Luan 15 DFómh 2018 at 15:38, Shane Carey  wrote:
>
>> Great thanks, is there anything available that will snap these lines to
>> the polygon edges?
>>
>> Thanks
>> Le gach dea ghui,
>> *Shane Carey*
>> *GIS and Data Solutions Consultant*
>>
>>
>> On Mon, Oct 15, 2018 at 3:30 PM Markus Neteler  wrote:
>>
>>> On Sun, Oct 14, 2018 at 11:28 PM Shane Carey 
>>> wrote:
>>> >
>>> > Hi,
>>> >
>>> > I need to cut a floodplain polygon into one meter strips which are
>>> perpendicular to the river. This needs to be an automated process. Just
>>> wondering has anyone ever done this or similar in Grass GIS through python
>>> or a bash script?
>>> >
>>> > The link attached shows two slides. Stage one is the first slide and
>>> slide two is what I hope to end up with.
>>> >
>>> >
>>> https://docs.google.com/presentation/d/1x7IQ3BJlLLSb6R2atxIIiaeWzyCqCZQ3uWaGPZyo6ds/edit?usp=sharing
>>>
>>> This addon seems to do (similar) perpendicular line creation:
>>>   v.civil - Generates a alignment for designing roads, channels, and
>>> ports in civil engineering
>>>   https://grass.osgeo.org/grass7/manuals/addons/v.civil.html
>>>
>>> Best
>>> Markus
>>>
>>> --
>>> Markus Neteler, PhD
>>> http://www.mundialis.de - free data with free software
>>> http://grass.osgeo.org
>>> http://courses.neteler.org/blog
>>>
>> --
> Le gach dea ghui,
> *Shane Carey*
> *GIS and Data Solutions Consultant*
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/grass-user
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Cut a floodplain polygon with 1meter wide polygons perpendicular to its river channel

2018-10-14 Thread Johannes Radinger
Shane,

have a look at the GRASS v.transect add-on:
https://grass.osgeo.org/grass74/manuals/addons/v.transects.html
The manual says: "Creates transect lines or quadrilateral areas at regular
intervals perpendicular to a polyline."

HTH
Johannes

On Sun, Oct 14, 2018 at 11:28 PM Shane Carey  wrote:

> Hi,
>
> I need to cut a floodplain polygon into one meter strips which are
> perpendicular to the river. This needs to be an automated process. Just
> wondering has anyone ever done this or similar in Grass GIS through python
> or a bash script?
>
> The link attached shows two slides. Stage one is the first slide and slide
> two is what I hope to end up with.
>
>
> https://docs.google.com/presentation/d/1x7IQ3BJlLLSb6R2atxIIiaeWzyCqCZQ3uWaGPZyo6ds/edit?usp=sharing
>
> Any help would be appreciated. Thanks in advance.
>
> Le gach dea ghui,
> *Shane Carey*
> *GIS and Data Solutions Consultant*
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/grass-user
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Identify duplicate points

2018-09-21 Thread Johannes Radinger
Hi all,

this is maybe a very easy question: But how is it possible to identify
duplicate points (that have the same pair of coordinates) of a set of
points. I'd like to get the categories of all non-single
(i.e.duplicate/triplicate etc.) points. I know that I could use the
'rmdupl' tool from v.clean to remove duplicate geometries but I'd rather
like to manually select the duplicate points I'd like to remove. Therefore,
I'd need to identify the points with similar geometries first. I thought of
using v.distance to get the all points that have 0 distance from another
point but the result is 0 for all points (due to the self-distance which is
0 of course). I could update the X,Y coordinates of each point and the use
a SQL command to identify points with similar coordinate pairs; however, is
there a more direct way or tool in GRASS to perform such a task?

cheers,

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] v.net manual describing 'arc_layer'

2018-09-20 Thread Johannes Radinger
Hi all,

this is probably a simple question, but it is unclear to me when reading
the manual of the tool v.net:
The parameter 'arc_layer' (same holds true for 'arc_node') is
described as "Vector
features can have category values in different layers. This number
determines which layer to use. When used with direct OGR access this is the
layer name." Here it is unclear if this refers to the layer of the input or
the output? So if I set 'arc_layer=2', does this mean that my arcs are
saved in the output network in layer 2, OR does this mean that v.net uses
layer 2 from the input vector to create arcs from?
I have the feeling that the descriptions of 'arc_layer' and 'node_layer'
are default descriptions of the parameter 'layer' (e.g. see manual of
v.to.rast); IMHO this could be improved in the manual of v.net.

cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Elevation above a river

2018-09-16 Thread Johannes Radinger
That is were the third step should follow, i.e. using r.mapcalc to identify
all the cells that are 1m higher than the grown river (irrespective in
which distance from the river these cells are located). The rasterized area
could then be translated into a vector format using r.to.vect.

Here a small example of how the working flow could be using the North
Caroline example dataset:
##
# Set region
g.region raster=elevation@PERMANENT

# Extract elevation of the streams
r.mapcalc --o expression="streams_elevation = if( streams_derived@PERMANENT,
elevation@PERMANENT,null())"

# Grow stream_elevation map by a some hundred meters
r.grow -m --overwrite input=streams_elevation output=streams_elevation_grow
radius=500

# Calculate difference between original elevatoin and stream channel
elevation
r.mapcalc expression="stream_elevation_diff = elevation@PERMANENT -
streams_elevation_grow" --overwrite

# Extract areas that are than x meters higher than stream channel elevation
r.mapcalc expression="stream_elevation_diff_smaller1 = if(
stream_elevation_diff < 1,1,null())" --overwrite

# Maybe here it needs some cleaning of spurious cells etc. (using some kind
of neighbourhood filtering)

# Raster areas to a vector area (with smoothed corners; s-flag)
r.to.vect -s --overwrite input=stream_elevation_diff_smaller1
output=stream_elevation_diff_smaller1_area type=area
##

HTH
/Johannes

On Sun, Sep 16, 2018 at 8:53 PM Shane Carey  wrote:

> Hey Johannes,
>
> Thanks for your reply. How does r.grow work if let's say the height above
> the river reaches 1m at 3meters away from the river. And in an other area
> it reaches the 1meter height at 2meters away from the river. Is it able to
> follow that line?
>
> Thanks
>
> On Domh 16 MFómh 2018 at 17:21, Johannes Radinger <
> johannesradin...@gmail.com> wrote:
>
>> To me this looks like a flooding-related question, i.e. to extract the
>> shore lines of a river if it's water level is raised by 1m or 3m?
>> Maybe (1) extract the raster cells of the elevation map that represents
>> the river channel, (2) then apply r.grow and (3) then r.mapcalc to subtract
>> the grown river channel from the original elevation map.
>> /Johannes
>>
>> On Sun, Sep 16, 2018 at 5:16 PM Markus Neteler  wrote:
>>
>>> Hi,
>>>
>>> Shane Carey  schrieb am Fr., 14. Sep. 2018, 23:02:
>>>
>>>> Hi All,
>>>>
>>>> Does anyone know is it possible to calculate the elevation above a
>>>> river channel (actual river network that was digitised as opposed to being
>>>> extracted from a DTM) from a DTM and create a polygon from it.
>>>>
>>>> I need to calculate heights of 1m and 3m above a river channel on both
>>>> sides of the channel and create a polygon from it.
>>>>
>>>
>>> This isn't clear to me. Could you elaborate?
>>>
>>> Best
>>> Markus
>>>
>>>
>>> Thanks all.
>>>>
>>>> Le gach dea ghui,
>>>> *Shane Carey*
>>>> *GIS and Data Solutions Consultant*
>>>> ___
>>>> grass-user mailing list
>>>> grass-user@lists.osgeo.org
>>>> https://lists.osgeo.org/mailman/listinfo/grass-user
>>>
>>> ___
>>> grass-user mailing list
>>> grass-user@lists.osgeo.org
>>> https://lists.osgeo.org/mailman/listinfo/grass-user
>>
>> --
> Le gach dea ghui,
> *Shane Carey*
> *GIS and Data Solutions Consultant*
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Elevation above a river

2018-09-16 Thread Johannes Radinger
To me this looks like a flooding-related question, i.e. to extract the
shore lines of a river if it's water level is raised by 1m or 3m?
Maybe (1) extract the raster cells of the elevation map that represents the
river channel, (2) then apply r.grow and (3) then r.mapcalc to subtract the
grown river channel from the original elevation map.
/Johannes

On Sun, Sep 16, 2018 at 5:16 PM Markus Neteler  wrote:

> Hi,
>
> Shane Carey  schrieb am Fr., 14. Sep. 2018, 23:02:
>
>> Hi All,
>>
>> Does anyone know is it possible to calculate the elevation above a river
>> channel (actual river network that was digitised as opposed to being
>> extracted from a DTM) from a DTM and create a polygon from it.
>>
>> I need to calculate heights of 1m and 3m above a river channel on both
>> sides of the channel and create a polygon from it.
>>
>
> This isn't clear to me. Could you elaborate?
>
> Best
> Markus
>
>
> Thanks all.
>>
>> Le gach dea ghui,
>> *Shane Carey*
>> *GIS and Data Solutions Consultant*
>> ___
>> grass-user mailing list
>> grass-user@lists.osgeo.org
>> https://lists.osgeo.org/mailman/listinfo/grass-user
>
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/grass-user
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Not matching categories/line IDs

2018-09-14 Thread Johannes Radinger
Hi Markus,
Hi all,

I created a zipped version of a test location including one PERMANENT
mapset that includes the vector map that causes these inconsistencies
(Ebro_river_vector_poly_clean_tmp). The test location can be downloaded
from here: https://we.tl/t-eSUZqfP0TQ. I also included the initial vector
map before running the v.clean tool (Ebro_river_vector_poly_tmp). More
specifically, I am using following code:

grass.run_command("v.clean",
overwrite=True,
flags="c",
input="Ebro_river_vector_poly_tmp",
output="Ebro_river_vector_poly_clean_tmp",
tool="break,rmdupl,rmline,rmsa,rmdangle",
threshold="0,0,0,0,0.001")

Maybe this helps to find out what happens here.

/Johannes

On Thu, Sep 13, 2018 at 10:26 PM Markus Metz 
wrote:

>
>
> On Thu, Sep 13, 2018 at 10:50 AM Johannes Radinger <
> johannesradin...@gmail.com> wrote:
> >
> > Hi all,
> >
> > I am a little bit puzzled about the actual number of lines/categories of
> a specific vector map (river network). As discussed yesterday, the output
> of v.category with the report option provides a column 'count' which is a
> feature count. Running v.category option=report on my vector provides me
> following:
> >
> > v.category input=Ebro_river_vector_poly_clean_tmp@Ebro2_alien_spread
> option=report
> > Layer/table: 1/Ebro_river_vector_poly_clean_tmp
> > type   countminmax
> > point  0  0  0
> > line 804  1   1053
> > 
> >
> > This would indicate that I have 804 line features. However, using the
> tool v.info on the same vector I get following different result regarding
> the number of lines:
> > v.info map=Ebro_river_vector_poly_clean_tmp@Ebro2_alien_spread
>
> >
>  
> ++
> >  | Name:Ebro_river_vector_poly_clean_tmp
>  |
> >  ...
> >  |   Number of points:   0   Number of centroids:  0
>  |
> >  |   Number of lines:799 Number of boundaries: 0
>  |
> > ...
> >
> > And even more puzzling, I then added a table to this vector map and the
> output of v.db.addtable is:
> > v.db.addtable map=Ebro_river_vector_poly_clean_tmp@Ebro2_alien_spread
>
> > Reading features...
> > Updating database...
> > 803 categories read from vector map (layer 1)
> > 803 categories read from vector map don't exist in selection from table
> > 803 records updated/inserted (layer 1)
>
> If there are 804 features, 2 features must have the same category value,
> resulting in 803 updates to the database
> >
> > So, is it 804,799 or 803??
>
> The output of v.info and v.category option=report must be identical with
> regard to the count of features.
> >
> > The map I am using is a polyline river network (cat=first) that has gone
> trough the tool v.clean with the options:
> break,rmdupl,rmline,rmsa,rmdangle. The number of lines and cats before the
> cleaning tool is 805 (consistent output among all tools).
> >
> > So I am wondering what happens during cleaning and how to get a clean
> map with a corresponding table that has entries for each cat, and were each
> cat represents one line feature (i.e. that the number of cats/line features
> is consistent over the multiple tools).
> >
> > Please let me know if I should share the vector map (and how, which
> format).
>
> Ideal would be a stripped down GRASS mapset with only this vector
> Ebro_river_vector_poly_clean_tmp. Exporting the vector would change the
> geometry representation and it would no longer be possible to reproduce the
> problem.
>
> Markus M
>
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Defect v.edit tool=break with multiple coordinates

2018-09-13 Thread Johannes Radinger
Hi all,

I presumably detected an problem/bug with v.edit and its tool "break". I
wanted to break a line at multiple coordinates. But the result is an
overlay of multiple lines instead of the original line that is broken into
segments. Here a small reproducible example using the NC dataset:

from grass.script import core as grass

# Copy data from PERMANENT mapset
grass.run_command("g.copy",
overwrite=True,
vector="railroads@PERMANENT,railroads")
grass.run_command("v.extract",
overwrite=True,
input="firestations@PERMANENT",
cat="58,66,29,35,31,32",
output="firestations")

# Connectors to find exact break points coords
grass.run_command("v.distance",
overwrite=True,
to="railroads",
from_="firestations",
output="firestations_connectors")
firestations_connectors_coors = grass.read_command("v.to.db",
flags="p",
map="firestations_connectors",
option="end",
separator="comma").splitlines()
firestations_connectors_coors = [[",".join(x.split(",")[1:3])] for x
in firestations_connectors_coors[1:]]

# Break at point coords
grass.run_command("v.edit",
  tool="break",
  map="railroads",
  coords=firestations_connectors_coors,
  threshold=50,
  layer=1,
  cats="1-")


My system:
GRASS version: 7.4.1

GRASS SVN revision: r72807M

Build date: 2018-06-13

Build platform: x86_64-apple-darwin17.6.0

GDAL: 2.0.0

PROJ.4: 5.1.0

GEOS: 3.6.2

SQLite: 3.19.3

Python: 2.7.15

wxPython: 4.0.0

Platform: Darwin-17.7.0-x86_64-i386-64bit

I seems that this issue has already been detected long time ago (see
https://trac.osgeo.org/grass/ticket/2903). I extended the bug ticket with
my experience and details using the new GRASS versions. Interestingly, all
people that reported about this problem were working on MacOS?!

Can anybody else reproduce this problem? Or does someone know a solution?

Cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Not matching categories/line IDs

2018-09-13 Thread Johannes Radinger
Hi all,

I am a little bit puzzled about the actual number of lines/categories of a
specific vector map (river network). As discussed yesterday, the output of
v.category with the report option provides a column 'count' which is a
feature count. Running v.category option=report on my vector provides me
following:

v.category input=Ebro_river_vector_poly_clean_tmp@Ebro2_alien_spread
option=report
Layer/table: 1/Ebro_river_vector_poly_clean_tmp
type   countminmax
point  0  0  0
line 804  1   1053


This would indicate that I have 804 line features. However, using the tool
v.info on the same vector I get following different result regarding the
number of lines:
v.info map=Ebro_river_vector_poly_clean_tmp@Ebro2_alien_spread

 ++
 | Name:Ebro_river_vector_poly_clean_tmp
  |
 ...
 |   Number of points:   0   Number of centroids:  0
  |
 |   Number of lines:799 Number of boundaries: 0
  |
...

And even more puzzling, I then added a table to this vector map and the
output of v.db.addtable is:
v.db.addtable map=Ebro_river_vector_poly_clean_tmp@Ebro2_alien_spread

Reading features...
Updating database...
803 categories read from vector map (layer 1)
803 categories read from vector map don't exist in selection from table
803 records updated/inserted (layer 1)

So, is it 804,799 or 803??

The map I am using is a polyline river network (cat=first) that has gone
trough the tool v.clean with the
options: break,rmdupl,rmline,rmsa,rmdangle. The number of lines and cats
before the cleaning tool is 805 (consistent output among all tools).

So I am wondering what happens during cleaning and how to get a clean map
with a corresponding table that has entries for each cat, and were each cat
represents one line feature (i.e. that the number of cats/line features is
consistent over the multiple tools).

Please let me know if I should share the vector map (and how, which format).

Cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] v.category report cats

2018-09-12 Thread Johannes Radinger
Thanks for your quick answer, Moritz!
I filed a ticket for this (https://trac.osgeo.org/grass/ticket/3643)

/Johannes

On Wed, Sep 12, 2018 at 2:58 PM Moritz Lennert 
wrote:

> On 12/09/18 14:40, Johannes Radinger wrote:
> > Hi all, I was wondering about the report output of 'v.category
> > option=report'. This output contains a column for each layer that is
> > called 'count'. Does this 'count' refer to the number of vector
> > objects or the number of unique category numbers in that layer (I
> > guess it is the first)? This might differ if two objects share the
> > same category. I think both would be useful, however which one is
> > provided here? Maybe this could be added to the manual, as one could
> > assume that 'count' refers to categories numbers?!
>
> The first examples in the man page do implicitely give the answer, i.e.
> that it is feature count, not unique category numbers count. Actually,
> you cannot get this latter information from v.category op=report. AFAIK,
> you would have to do something like this (in *nix):
>
> v.category YourMap op=print | sort -n | uniq | wc -l
>
> But I agree that the count info is not clearly defined in the man page.
> Please file a ticket for this, so we don't forget.
>
> Moritz
>
>
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/grass-user
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] v.category report cats

2018-09-12 Thread Johannes Radinger
Hi all,
I was wondering about the report output of 'v.category option=report'. This
output contains a column for each layer that is called 'count'. Does this
'count' refer to the number of vector objects or the number of unique
category numbers in that layer (I guess it is the first)? This might differ
if two objects share the same category. I think both would be useful,
however which one is provided here? Maybe this could be added to the
manual, as one could assume that 'count' refers to categories numbers?!

Cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] GRASS-way to delete points from a map and attribute table

2017-11-09 Thread Johannes Radinger
Indeed, finally I also used a combination of v.edit to remove the point
from the map and db.execute ("DELETE FROM my_point_map WHERE cat IN
(1,2,3)") to delete the entries from the attribute table.

/Johannes

On Thu, Nov 9, 2017 at 10:27 AM, Stefan Blumentrath <
stefan.blumentr...@nina.no> wrote:

> So, removing records from the attribute table would be a job for
> db.execute of db.select I guess…
>
>
>
> *From:* Johannes Radinger [mailto:johannesradin...@gmail.com]
> *Sent:* torsdag 9. november 2017 10.17
> *To:* Stefan Blumentrath <stefan.blumentr...@nina.no>
> *Cc:* GRASS user list <grass-user@lists.osgeo.org>
> *Subject:* Re: [GRASS-user] GRASS-way to delete points from a map and
> attribute table
>
>
>
> Dear Stefan,
>
> thanks for the tip, but as said I tried already v.edit but somehow this
> tools does not modify the attribute table (or do I miss something here?)
>
> /Johannes
>
>
>
> On Thu, Nov 9, 2017 at 10:10 AM, Stefan Blumentrath <
> stefan.blumentr...@nina.no> wrote:
>
> Did you try v.edit:
>
> https://grass.osgeo.org/grass72/manuals/v.edit.html
>
> Cheers
>
> Stefan
>
>
>
> *From:* grass-user [mailto:grass-user-boun...@lists.osgeo.org] *On Behalf
> Of *Johannes Radinger
> *Sent:* torsdag 9. november 2017 10.04
> *To:* GRASS user list <grass-user@lists.osgeo.org>
> *Subject:* [GRASS-user] GRASS-way to delete points from a map and
> attribute table
>
>
>
> Hi,
>
>
>
> what would be the most GRASS-like way to delete several points from a
> points vector map.
>
>
>
> I know the specific cats of the points I would like to delete (10 of
> 20) and I want to delete the points from the map and, of course, also
> from the corresponding attribute table. I already tried v.edit but here the
> points are not deleted from the attribute table. I could also use v.extract
> with the reverse flag, however this requires to define a new output map
> (which I would like to avoid if possible). Any straight forward way in
> GRASS?
>
>
>
> /Johannes
>
>
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] GRASS-way to delete points from a map and attribute table

2017-11-09 Thread Johannes Radinger
Dear Stefan,
thanks for the tip, but as said I tried already v.edit but somehow this
tools does not modify the attribute table (or do I miss something here?)
/Johannes

On Thu, Nov 9, 2017 at 10:10 AM, Stefan Blumentrath <
stefan.blumentr...@nina.no> wrote:

> Did you try v.edit:
>
> https://grass.osgeo.org/grass72/manuals/v.edit.html
>
> Cheers
>
> Stefan
>
>
>
> *From:* grass-user [mailto:grass-user-boun...@lists.osgeo.org] *On Behalf
> Of *Johannes Radinger
> *Sent:* torsdag 9. november 2017 10.04
> *To:* GRASS user list <grass-user@lists.osgeo.org>
> *Subject:* [GRASS-user] GRASS-way to delete points from a map and
> attribute table
>
>
>
> Hi,
>
>
>
> what would be the most GRASS-like way to delete several points from a
> points vector map.
>
>
>
> I know the specific cats of the points I would like to delete (10 of
> 20) and I want to delete the points from the map and, of course, also
> from the corresponding attribute table. I already tried v.edit but here the
> points are not deleted from the attribute table. I could also use v.extract
> with the reverse flag, however this requires to define a new output map
> (which I would like to avoid if possible). Any straight forward way in
> GRASS?
>
>
>
> /Johannes
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] t.register: register 4 consecutive rasters maps per month for long time series

2017-11-09 Thread Johannes Radinger
Hi Veronica,

thank you for your input regarding the registering maps with 0.25*months
time steps.
Here the solution I came up with was to import the maps in groups of 4
bands from the netCDF file. Then I calculate average over the 4 monthly map
and register the average map as a single map for the respective months:

# Register raster maps to assign correct date
scpdsi_maps_no = range(1,2688+1) # Total number of bands/maps
scpdsi_maps_groups_month = [scpdsi_maps_no[i:i + 4] for i in xrange(0,
len(scpdsi_maps_no), 4)]
scpdsi_maps_groups_monthyear = [scpdsi_maps_groups_month[i:i + 12] for i in
xrange(0, len(scpdsi_maps_groups_month), 12)]

year=1961 # First year
for idx, i in enumerate(scpdsi_maps_groups_monthyear):
month=1 # First month
for j in i:
print(str(year) + "-" + str(month).zfill(2) + "-01")
grass.run_command("r.import",
flags="o",
overwrite=True,
input="/path/to/my/scpdsi.nc",
band=j,
extent="region",
output="scPDSI")
scpdsi_maps = ["scPDSI."+str(x) for x in j]
grass.run_command("r.series",
overwrite=True,
input=scpdsi_maps,
output="scPDSI"+"-"+str(year)+"-"+str(month).zfill(2),
method="average")
grass.run_command("g.remove",
flags="f",
type="raster",
name=scpdsi_maps)
grass.run_command("t.register",
overwrite=True,
#flags="i",
type="raster",
input="scPDSI_monthly",
maps="scPDSI"+"-"+str(year)+"-"+str(month).zfill(2),
start=str(year) + "-" + str(month).zfill(2) + "-01",
increment="1 months")
month=month+1
year=year+1

This solution works fine for me.

/Johannes

On Tue, Oct 24, 2017 at 7:04 PM, Veronica Andreo <veroand...@gmail.com>
wrote:

> Hello Johanes,
>
> there's and example in t.register manual page with netCDF files [0] that
> might help. You would have 48 maps per year (IMHO, 0.25*month is a bit
> strange as time step, but well). I would say that you need to write a
> script to get dates for each 1/4 of a month (maybe someone else know of a
> function to do a similar task, dunno) in your time series and then make a
> file with mapname|start_time|end_time to pass to t.register.
>
> However, if you only need to estimate yearly means and long term means,
> you can use r.series [1] with every 48 maps for yearly means and, the whole
> list of maps for the long term mean.
>
> hth,
> Vero
>
> [0] https://grass.osgeo.org/grass72/manuals/t.register.html
> [1] https://grass.osgeo.org/grass72/manuals/r.series.html
>
> 2017-10-24 11:56 GMT+02:00 Johannes Radinger <johannesradin...@gmail.com>:
>
>> Hi,
>> I am very new to the temporal functionalities of GRASS. Specifically I
>> have a netCDF file that contains many layers where each layer represents a
>> raster map of a drought index (http://monitordesequia.csic.es/map/) for
>> a timepoint of a timeseries. The time series is from 1/1961 to 12/2015 with
>> 4 maps for each month. So these are not really weekly maps but maps with an
>> interval of 0.25*months. The layers of the netCDF are only numbered
>> consecutively so there is no indication of a specific date etc. (I just now
>> that they start with 1/1961 and then the interval of 0.25*months). What
>> would be the procedure to register these raster maps in GRASS so that I can
>> calculate later only e.g. yearly means and a long-time year average?
>>
>> Best regards,
>> Johannes
>>
>> ___
>> grass-user mailing list
>> grass-user@lists.osgeo.org
>> https://lists.osgeo.org/mailman/listinfo/grass-user
>>
>
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] GRASS-way to delete points from a map and attribute table

2017-11-09 Thread Johannes Radinger
Hi,

what would be the most GRASS-like way to delete several points from a
points vector map.

I know the specific cats of the points I would like to delete (10 of
20) and I want to delete the points from the map and, of course, also
from the corresponding attribute table. I already tried v.edit but here the
points are not deleted from the attribute table. I could also use v.extract
with the reverse flag, however this requires to define a new output map
(which I would like to avoid if possible). Any straight forward way in
GRASS?

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] t.register: register 4 consecutive rasters maps per month for long time series

2017-10-24 Thread Johannes Radinger
Hi,
I am very new to the temporal functionalities of GRASS. Specifically I have
a netCDF file that contains many layers where each layer represents a
raster map of a drought index (http://monitordesequia.csic.es/map/) for a
timepoint of a timeseries. The time series is from 1/1961 to 12/2015 with 4
maps for each month. So these are not really weekly maps but maps with an
interval of 0.25*months. The layers of the netCDF are only numbered
consecutively so there is no indication of a specific date etc. (I just now
that they start with 1/1961 and then the interval of 0.25*months). What
would be the procedure to register these raster maps in GRASS so that I can
calculate later only e.g. yearly means and a long-time year average?

Best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Create a point along a line in a specified distance from another point

2017-10-17 Thread Johannes Radinger
Hi,

is there an easy straight forward command to create a new point on a line
that is in a specified distance (e.g. 1000 m) from another point (on the
same line)?

Just an example (of course from river science): I have a sampling site on a
river line. I'd like to use these new point(s) to break the original river
lines in order to extract a line segment that is exactly e.g. 2000 m long
(1000 up- and downstream the sampling site).

I had a look in v.to.points but this tool only creates points in
equidistant steps irrespective of any other points (e.g. sampling points).
Another option would be using v.buffer around the sampling points and then
v.overlay. However, this approach would "cut" the river lines in euclidean
distance from the sampling site rather then in line/river distance.

Any suggestions?

/J
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Extract subnet from grass vector network based on selected nodes

2017-10-17 Thread Johannes Radinger
Dear Stefan,
Dear all,

thank you for your tipps, and especially the recommendations regarding the
use of the igraph library. In fact, I decided to use do most of the network
analysis now in igraph/R. I use GRASS and the report option of v.net to get
a table of all vertices and edges my graph/network consists of.
Subsequently, I load these tables in igraph and conduct specific analysis.
That seems a very convenient way to obtain most of the information I am
looking for.

/J

On Fri, Oct 13, 2017 at 3:17 PM, Stefan Blumentrath <
stefan.blumentr...@nina.no> wrote:

> Hi Johannes,
>
>
>
> Sorry for the delayed response.
>
>
>
> You could try v.net and then use igraph.
>
>
>
> See: https://github.com/NINAnor/gudbrand_hydro/blob/master/v.
> igraph.order.py
>
> as an example for network analysis using a GRASS python script in
> combination with igraph.
>
>
>
> What you would need is to compute the paths between all possible
> combinations of sampling sites (using e.g. http://igraph.org/python/doc/
> igraph.GraphBase-class.html#get_shortest_paths).
>
> Once you got a unique list of edges you can use v.extract I guess to get
> your subnetwork…
>
>
>
> Or if you continue with network analysis you could just use:
>
> http://igraph.org/python/doc/igraph.GraphBase-class.html#subgraph_edges
>
>
>
> Hope that is somehow useful…
>
>
>
> Cheers
>
> Stefan
>
>
>
> *From:* grass-user [mailto:grass-user-boun...@lists.osgeo.org] *On Behalf
> Of *Johannes Radinger
> *Sent:* torsdag 28. september 2017 13.04
> *To:* Markus Metz <markus.metz.gisw...@gmail.com>
> *Cc:* GRASS user list <grass-user@lists.osgeo.org>
> *Subject:* Re: [GRASS-user] Extract subnet from grass vector network
> based on selected nodes
>
>
>
> Hi all,
>
>
>
> I just tried two different tools that both work to extract a subnetwork
> that connects a set of selected nodes:
>
>
>
> 1) v.net.allpairs works fine and a subnetwork is extracted. However, this
> definitely needs some clean up and layer/attribute operations to get the
> connected tables/cats from the original network
>
>
>
> 2) v.net.steiner also works fine. However, it seems that this approach
> takes slightly longer than v.net.allpairs.
>
>
>
> In general both approaches are rather slow if one wants to extract a
> subnetwork based on e.g. >1000 selected nodes from a very large network.
> For example, my initial (large) network consists of >10 lines which
> makes any further analysis rather slow. Thus, I wanted to minimize the
> network to one that still connects my target nodes but skips parts that are
> not needed.
>
> Thank you for you suggestions, anyway.
>
> /Johannes
>
>
>
>
>
> On Thu, Sep 28, 2017 at 10:06 AM, Markus Metz <
> markus.metz.gisw...@gmail.com> wrote:
>
>
>
> On Thu, Sep 28, 2017 at 9:43 AM, Moritz Lennert <
> mlenn...@club.worldonline.be> wrote:
> >
> > On 28/09/17 08:51, Markus Metz wrote:
> >>
> >>
> >>
> >> On Wed, Sep 27, 2017 at 11:55 PM, Moritz Lennert <
> mlenn...@club.worldonline.be <mailto:mlenn...@club.worldonline.be>> wrote:
> >>  >
> >>  > On 27/09/17 21:03, Markus Metz wrote:
> >>  >>
> >>  >>
> >>  >>
> >>  >> On Wed, Sep 27, 2017 at 4:07 PM, Moritz Lennert <
> mlenn...@club.worldonline.be <mailto:mlenn...@club.worldonline.be>
> <mailto:mlenn...@club.worldonline.be <mailto:mlenn...@club.worldonline.be>>>
> wrote:
> >>  >>  >
> >>  >>  >
> >>  >>  >
> >>  >>  > Le 27 septembre 2017 13:11:54 GMT+02:00, Johannes Radinger <
> johannesradin...@gmail.com <mailto:johannesradin...@gmail.com>  johannesradin...@gmail.com <mailto:johannesradin...@gmail.com>>> a écrit :
> >>  >>  > >Hi,
> >>  >>  > >
> >>  >>  > >I have a GRASS vector network that represents a river network
> (with
> >>  >>  > >many
> >>  >>  > >first order tributaries) and that has additional connected
> nodes that
> >>  >>  > >represent sampling sites.
> >>  >>  > >
> >>  >>  > >I'd like to extract a minimum subnetwork of the full network
> that still
> >>  >>  > >connects a set of selected nodes (e.g. identified by their cat).
> >>  >>  > >However,
> >>  >>  > >network edges (i.e. river segments) that are not necessary to
> connect
> >>  >>  > >the

[GRASS-user] Version GRASS for Mac OSX

2017-10-01 Thread Johannes Radinger
Hi all,

the binaries of GRASS for Mac OSX available via
http://grassmac.wikidot.com/downloads are already older than a year (from
14.09.2016). Are there any plans to upload and provide a newer version of
GRASS as binary download in near future?

Best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Extract subnet from grass vector network based on selected nodes

2017-09-28 Thread Johannes Radinger
Hi all,

I just tried two different tools that both work to extract a subnetwork
that connects a set of selected nodes:

1) v.net.allpairs works fine and a subnetwork is extracted. However, this
definitely needs some clean up and layer/attribute operations to get the
connected tables/cats from the original network

2) v.net.steiner also works fine. However, it seems that this approach
takes slightly longer than v.net.allpairs.

In general both approaches are rather slow if one wants to extract a
subnetwork based on e.g. >1000 selected nodes from a very large network.
For example, my initial (large) network consists of >10 lines which
makes any further analysis rather slow. Thus, I wanted to minimize the
network to one that still connects my target nodes but skips parts that are
not needed.

Thank you for you suggestions, anyway.

/Johannes



On Thu, Sep 28, 2017 at 10:06 AM, Markus Metz <markus.metz.gisw...@gmail.com
> wrote:

>
>
> On Thu, Sep 28, 2017 at 9:43 AM, Moritz Lennert <
> mlenn...@club.worldonline.be> wrote:
> >
> > On 28/09/17 08:51, Markus Metz wrote:
> >>
> >>
> >>
> >> On Wed, Sep 27, 2017 at 11:55 PM, Moritz Lennert <
> mlenn...@club.worldonline.be <mailto:mlenn...@club.worldonline.be>> wrote:
> >>  >
> >>  > On 27/09/17 21:03, Markus Metz wrote:
> >>  >>
> >>  >>
> >>  >>
> >>  >> On Wed, Sep 27, 2017 at 4:07 PM, Moritz Lennert <
> mlenn...@club.worldonline.be <mailto:mlenn...@club.worldonline.be>
> <mailto:mlenn...@club.worldonline.be <mailto:mlenn...@club.worldonline.be>>>
> wrote:
> >>  >>  >
> >>  >>  >
> >>  >>  >
> >>  >>  > Le 27 septembre 2017 13:11:54 GMT+02:00, Johannes Radinger <
> johannesradin...@gmail.com <mailto:johannesradin...@gmail.com>  johannesradin...@gmail.com <mailto:johannesradin...@gmail.com>>> a écrit :
> >>  >>  > >Hi,
> >>  >>  > >
> >>  >>  > >I have a GRASS vector network that represents a river network
> (with
> >>  >>  > >many
> >>  >>  > >first order tributaries) and that has additional connected
> nodes that
> >>  >>  > >represent sampling sites.
> >>  >>  > >
> >>  >>  > >I'd like to extract a minimum subnetwork of the full network
> that still
> >>  >>  > >connects a set of selected nodes (e.g. identified by their cat).
> >>  >>  > >However,
> >>  >>  > >network edges (i.e. river segments) that are not necessary to
> connect
> >>  >>  > >the
> >>  >>  > >sampling points should be excluded in the new subnetwork. Is
> there a
> >>  >>  > >function or process in GRASS GIS to extract such a subnetwork
> that
> >>  >>  > >fully
> >>  >>  > >connects a set of selected nodes?
> >>  >>  >
> >>  >>  > not sure but maybe v.net.spanningtree ?
> >>  >>
> >>  >> v.net.spanningtree calculates a tree covering all nodes in the
> network, not only selected nodes, therefore v.net.spanningtree does not
> apply here.
> >>  >
> >>  >
> >>  > If you connect only the selected nodes to the network, wouldn't that
> work ? Or does v.net.spanningtree consider all connections between lines as
> nodes ?
> >>
> >> v.net.spanningtree considers all internal nodes of the network. See also
> >> https://en.wikipedia.org/wiki/Spanning_tree
> >
> >
> > Ok, thanks. So, one would need to "disconnect" lines at non-selected
> nodes for this to work.
>
> or use v.net.steiner (see my previous reply)
> >
> > And maybe some clarification on what is meant by "nodes" in the sentence
> "A spanning tree is a minimum cost subnetwork connecting all nodes in an
> undirected network" in the man page might help future users.
>
> Yes, that would help. I needed to look at the library fn
> NetA_spanning_tree() to be sure.
>
> Markus M
> >
> > Moritz
>
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Extract subnet from grass vector network based on selected nodes

2017-09-27 Thread Johannes Radinger
Hi,

I have a GRASS vector network that represents a river network (with many
first order tributaries) and that has additional connected nodes that
represent sampling sites.

I'd like to extract a minimum subnetwork of the full network that still
connects a set of selected nodes (e.g. identified by their cat). However,
network edges (i.e. river segments) that are not necessary to connect the
sampling points should be excluded in the new subnetwork. Is there a
function or process in GRASS GIS to extract such a subnetwork that fully
connects a set of selected nodes?

Cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Merge lines into common categories based an attribute information

2017-09-17 Thread Johannes Radinger
Thank you Moritz!

Of course, v.reclass does the job. I just didn't see the column option...

Regarding the second (additional) part of this question: It seems that
v.to.db works well to report line lengths on a category-base. However,
sinuosity seems to be calculated based on features (i.e. feature ids).
Thus, I only get reasonable results if a cat only consists of one feature.
For connected lines that consist of more than one feature (lines) the
sinuosity index is not correct. So, I probably need to loop over the cats,
extract the lines that belong to one cat and build a polyline before
calculating sinuosity.

/J

On Sun, Sep 17, 2017 at 3:53 PM, Moritz Lennert <
mlenn...@club.worldonline.be> wrote:

>
>
> Le 17 septembre 2017 15:16:02 GMT+02:00, Johannes Radinger <
> johannesradin...@gmail.com> a écrit :
> >Hi,
> >
> >how can I merge multiple lines into one category based on some
> >attribute
> >information.
> >For example I have a line vector consisting of 10 lines and several
> >lines
> >belong to three different groups (lines belonging to the same group
> >share
> >the same value in a specific attribute column). Now I'd like to merge
> >them
> >(assign a common cats) and create an new attribute with the "group
> >value"
> >as new category (the table should then have 3 rows)? I think this
> >should be
> >easily possible, but I could not yet find the correct tool.
>
> Try v.reclass.
>
>
> Moritz
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Merge lines into common categories based an attribute information

2017-09-17 Thread Johannes Radinger
Hi,

how can I merge multiple lines into one category based on some attribute
information.
For example I have a line vector consisting of 10 lines and several lines
belong to three different groups (lines belonging to the same group share
the same value in a specific attribute column). Now I'd like to merge them
(assign a common cats) and create an new attribute with the "group value"
as new category (the table should then have 3 rows)? I think this should be
easily possible, but I could not yet find the correct tool.

Afterwards I'd like to use v.to.db to update the attribute table with some
info on the vector features. BTW, how is sinuosity (v.to.db) calculated if
several lines are grouped into one category?

Any suggestion is highly welcome! Thanks!


Best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Running GRASS within R (Rstudio) on MacOS

2017-09-07 Thread Johannes Radinger
Hi,

I try to run GRASS73 within R (actually R-studio) using the rgrass7 library
on MacOS 10.12.6. I followed the instructions stated in the wiki (
https://grasswiki.osgeo.org/wiki/R_statistics/rgrass7), however some of the
commands were not working for me. For instance, "grass70 --config path"
(also with grass73 etc) does not work here, however I found the path to the
application which is: /Applications/GRASS-7.3.app; i.e. using
"open /Applications/GRASS-7.3.app" opens my grass installation.

So I tried to initialize GRASS (North Carolina location - PERMANENT) from
within R using the initGRASS() command:

initGRASS(gisBase = "/Applications/GRASS-7.3.app", home = tempdir(),
  gisDbase = "/Users/Johannes_Radinger/Documents/GRASS-GIS/",
  location = "nc_spm_08_grass7", mapset = "PERMANENT",
SG="elevation",
  override=TRUE)

However I get following error(s):
sh: g.gisenv: command not found
sh: g.gisenv: command not found
sh: g.gisenv: command not found
sh: g.gisenv: command not found
sh: g.gisenv: command not found
sh: g.version: command not found
Error in system(paste("g.version", get("addEXE", envir = .GRASS_CACHE),  :
  error in running command

So what I am doing wrong? Do I miss some arguments in my command or is this
related to my Mac environment? Any suggestion I should try to solve that
problem?

Thanks and best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Image segmentation to separate object from background

2017-06-26 Thread Johannes Radinger
Hi GRASS users,

has anyone of you tried to use GRASS image tools (e.g. segmentation etc) to
identify an object in a picture. For example I have multiple photos of fish
and would like to separate the fish from its background in an automatized
way. The images look like:

http://fishbase.org/photos/PicturesSummary.php?ID=4730=species
http://fishbase.org/photos/PicturesSummary.php?StartRow=
0=4662=species=5

So what I already know apriori: I would like two final classes (fish and
background), the fish is more or less centred within each picture, and
there is only one fish in each picture. I played already around with
i.segment but the results are not yet satisfying.

Any ideas?

/johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Single centroid/point for multiple polygons sharing the same category

2017-06-19 Thread Johannes Radinger
Hi GRASS users,

I am working on a polygon map of France that shows the single areas
belonging to a specific postal code. I want to extract for each area the
centroid (using v.extract), transform it to a point map (using v.type) and
finally add the X and Y coordinates to the attribute table for each postal
code region.

However, now I have the problem that some "postal code areas" have more
than one centroid, i.e. two or more closeby but not connected polygons
belong the same postal region and share a common category and entry in the
attribute table. So for these cases I get two or more points (sharing the
same cat) and of course I cannot calculate a single XY pair for that postal
code region.

Consequently, I need to find a method to get centroids of all the areas
sharing a common cat or I need to calculate the mean position of points
afterwards (using v.centerpoint). As far as I understood v.centerpoint all
points are used to calculate a mean point but not stratified for points
sharing the same cat!? FYI: I have around 2000 postal code areas of which
200 consists of more than a single polygon.

Any suggestions?

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Obtain land use within 50 km buffers around vector points

2017-06-16 Thread Johannes Radinger
Hi Moritz,
Hi All...

here my python code snippets that worked for me. I think there might be
many ways to simplify the code...

import grass.script as grass
import sqlite3
import pandas

# Def variable names
respondents_coord="respondents_coord"
respondents_points="respondents_points"
respondents_buffer="respondents_buffer"
respondents_buffer_cat_2_layer="respondents_buffer_cat_2_layer"
respondents_buffer_cat_2_raster="respondents_buffer_cat_2_raster"
correspondance_table="correspondance_table"
clc_summary_table="clc_summary_table"
clc_summary_pivot="clc_summary_pivot"

# Import points from X-Y table
grass.run_command("v.in.db",
overwrite=True,
table=respondents_coord,
output=respondents_points)

 Calculate buffers around points
grass.run_command("v.buffer",
overwrite=True,
flags="t",
input=respondents_points,
type="point",
output=respondents_buffer,
distance=5)
grass.run_command("v.category",
overwrite=True,
input=respondents_buffer,
option="add",
layer=2,
out=respondents_buffer_cat_2_layer)

# Get CORINE land use special classes per buffer
correspondance_table_values = []
for line in grass.read_command('v.category',
input_=respondents_buffer_cat_2_layer,
layer='1,2',
option='print').splitlines():
layers=line.split('|')
l1 = layers[0].split('/')
l2 = layers[1]
for cat in l1:
correspondance_table_values.append((cat, l2))

grass.run_command("db.execute",
sql="CREATE TABLE correspondance_table (buffer_cat_1 INTEGER, buffer_cat_2
INTEGER)")

conn = sqlite3.connect('/path/to/sqlite/sqlite.db')
cur = conn.cursor()
cur.executemany('INSERT INTO correspondance_table (buffer_cat_1,
buffer_cat_2) values (?,?)',correspondance_table_values)
conn.commit()
conn.close()

grass.run_command("v.to.rast",
overwrite=True,
input=respondents_buffer_cat_2_layer,
layer=2,
output=respondents_buffer_cat_2_raster,
use="cat")

if clc_summary_pivot in
grass.read_command("db.tables",flags="p").split("\n"):
grass.run_command("db.droptable",
flags="f",
table=clc_summary_pivot)

grass.run_command("db.execute",
sql="CREATE TABLE {} (buffer_cat_2 INTEGER, CLC_built_up INTEGER,
CLC_arable INTEGER, CLC_permanent_crops INTEGER, CLC_grassland INTEGER,
CLC_forest INTEGER, CLC_others INTEGER, CLC_intertidal_coastal INTEGER,
CLC_water_bodies INTEGER, CLC_sea INTEGER)".format(clc_summary_pivot))

clc_values = []
for line in grass.read_command("r.stats",
flags="cn",
input="{},{}".format(respondents_buffer_cat_2_raster,"CLC_reclass@Corine_LandCover
")).splitlines():
clc_values.append((line.split(' ')[0], line.split(' ')[1],line.split('
')[2]))

df = pandas.DataFrame(clc_values, columns=['buffer_cat_2', 'CLC_class',
'count'])
CLC_class_cols = ["1","2","3","4","5","6","7","8","9"]
clc_summary_pivot_values = [tuple(x) for x in
df.pivot(index='buffer_cat_2', columns='CLC_class',
values='count').reindex(columns=CLC_class_cols).fillna(0).astype(int).to_records(index=True)]

conn = sqlite3.connect('/path/to/sqlite/sqlite.db')
cur = conn.cursor()
cur.executemany('INSERT INTO {} (buffer_cat_2, CLC_built_up, CLC_arable,
CLC_permanent_crops, CLC_grassland, CLC_forest, CLC_others,
CLC_intertidal_coastal, CLC_water_bodies, CLC_sea) values
(?,?,?,?,?,?,?,?,?,?)'.format(clc_summary_pivot),clc_summary_pivot_values)
conn.commit()
conn.close()

grass.run_command("db.execute",
sql="CREATE TABLE {} AS SELECT buffer_cat_1,SUM(CLC_built_up) AS
CLC_built_up, SUM(CLC_arable) AS CLC_arable, SUM(CLC_permanent_crops) AS
CLC_permanent_crops, SUM(CLC_grassland) AS CLC_grassland, SUM(CLC_forest)
AS CLC_forest, SUM(CLC_others) AS CLC_others, SUM(CLC_intertidal_coastal)
AS CLC_intertidal_coastal, SUM(CLC_water_bodies) AS CLC_water_bodies,
SUM(CLC_sea) AS CLC_sea  FROM {} AS A LEFT JOIN {} AS B ON
A.buffer_cat_2=B.buffer_cat_2 GROUP BY
buffer_cat_1".format(clc_summary_table,correspondance_table,clc_summary_pivot))

# Join information on land use back to original center points of buffers
grass.run_command("v.db.join",
map=respondents_points,
column="id_INT",
other_table=clc_summary_table,
other_column="resp_id")




cheers, Johannes

On Wed, Jun 14, 2017 at 11:10 AM, Moritz Lennert <
mlenn...@club.worldonline.be> wrote:

> Hi Johannes,
>
> On 07/06/17 10:20, Johannes Radinger wrote:
>
>> Thank you Moritz,
>>
>> your suggestion using r.stats and the rasterized layer 2 of buffers
>> works really nice. It took me just a while to summarize and join all
>> data and get them back into the db in the right format. However, I think
>> I managed this task now. Thank you for your help!
>>
>
> Would you be willing to share the final version of your approach ? This
> might be a nice seed for a module that would offer this function.
>
> Moritz
>
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Obtain land use within 50 km buffers around vector points

2017-06-01 Thread Johannes Radinger
Hi,

I have a vector map of 50 km buffers surrounding 4000 points in entire
Europe. The single buffers where calculated using v.buffer with the -t
flag. Thus, many of my buffers overlap (i.e. geometries of buffers are not
merged but split up). Now, I want to obtain for each 50 km buffer the
underlying land use (from a CORINE land use map). For example, I'd like to
know the relative share of "forest" within each buffer and I want to update
this information in the attribute table of the vector buffer map.

So what would be the GRASS way for such type of analysis involving
overlapping buffers? Do I need to loop over the single buffers, extract
each buffer, and run something like v.rast.stats using maps for the single
land use classes as raster input?

Any suggestions for a computationally efficient way for such type of
analyis?

Best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Indentify cell with maximum value per category of a cover map

2017-05-17 Thread Johannes Radinger
Hi,

I'd like to identify the cell of maximum flow accumulation for each
subbasin.

For example, I have a cover raster map representing all my subbasins (each
subbasin has its own cat value). Additionally, I have a flow accumulation
map. Now, I'd like to get a map with cells that are actually those with the
largest flow accumulation per subbasin (all other cell should be NULL). I
guess, using r.stats.zonal/r.statistics I can identify the maximum flow
accumulation value for each subbasin, but this doesn't provide me a map
with the corresponding cells!? Is there a straight forward way to do that
in GRASS with e.g. the mapcalc? Or do I need to loop over each subbasin?

Any suggestions?

cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Column names of sqlite attribute tables while using v.import

2017-05-17 Thread Johannes Radinger
Thank you Moritz,

that's what I thought. So I'll try the 'columns' option in v.in.ogr to
rename the columns during import.

/johannes

On Tue, May 16, 2017 at 11:50 PM, Moritz Lennert <
mlenn...@club.worldonline.be> wrote:

>
>
> Le 16 mai 2017 22:16:03 GMT+02:00, Johannes Radinger <
> johannesradin...@gmail.com> a écrit :
> >Hi all,
> >
> >I tried to import a polygon shape, specifically the HydroBASINS from
> >hydrosheds.org however I got following error:
> >
> >v.import
> >input=/.../Hydrosheds/hybas_na_lev00_v1c/hybas_na_lev00_v1c.shp
> >layer=hybas_na_lev00_v1c output=hybas_na_lev00_v1c
> >WARNING: All available OGR layers will be imported into vector map
> >
> >Check if OGR layer  contains polygons...
> >DBMI-SQLite driver error:
> >Error in sqlite3_prepare():
> >near "ORDER": syntax error
> >DBMI-SQLite driver error:
> >Error in sqlite3_prepare():
> >near "ORDER": syntax error
> >ERROR: Unable to create table: 'create table hybas_na_lev00_v1c (cat
> >integer, HYBAS_ID double precision, NEXT_DOWN double precision,
> >NEXT_SINK
> >double precision, MAIN_BAS double precision, DIST_SINK double
> >precision,
> >DIST_MAIN double precision, SUB_AREA double precision, UP_AREA double
> >precision, ENDO integer, COAST integer, ORDER integer, SORT double
> >precision, PFAF_1 integer, PFAF_2 integer, PFAF_3 integer, PFAF_4
> >integer,
> >PFAF_5 integer, PFAF_6 integer, PFAF_7 integer, PFAF_8 integer, PFAF_9
> >integer, PFAF_10 double precision, PFAF_11 double precision, PFAF_12
> >double
> >precision)'
> >ERROR: Unable to import
> >
> >
> >So I guess ORDER is a sqlite specific word that is reserved for SQL
> >operations and is not allowed as column name? Is this a general thing
> >of
> >sqlite or is this just related to v.import/v.in.ogr and the creation of
> >the
> >sqlite table?
>
> AFAIK, ORDER is a general SQL reserved keyword, not only in Sqlite. See
> [1] for example.
>
> Moritz
>
> [1] http://www.sql.org/sql-database/postgresql/manual/
> sql-keywords-appendix.html
>
>
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Column names of sqlite attribute tables while using v.import

2017-05-16 Thread Johannes Radinger
Hi all,

I tried to import a polygon shape, specifically the HydroBASINS from
hydrosheds.org however I got following error:

v.import input=/.../Hydrosheds/hybas_na_lev00_v1c/hybas_na_lev00_v1c.shp
layer=hybas_na_lev00_v1c output=hybas_na_lev00_v1c
WARNING: All available OGR layers will be imported into vector map

Check if OGR layer  contains polygons...
DBMI-SQLite driver error:
Error in sqlite3_prepare():
near "ORDER": syntax error
DBMI-SQLite driver error:
Error in sqlite3_prepare():
near "ORDER": syntax error
ERROR: Unable to create table: 'create table hybas_na_lev00_v1c (cat
integer, HYBAS_ID double precision, NEXT_DOWN double precision, NEXT_SINK
double precision, MAIN_BAS double precision, DIST_SINK double precision,
DIST_MAIN double precision, SUB_AREA double precision, UP_AREA double
precision, ENDO integer, COAST integer, ORDER integer, SORT double
precision, PFAF_1 integer, PFAF_2 integer, PFAF_3 integer, PFAF_4 integer,
PFAF_5 integer, PFAF_6 integer, PFAF_7 integer, PFAF_8 integer, PFAF_9
integer, PFAF_10 double precision, PFAF_11 double precision, PFAF_12 double
precision)'
ERROR: Unable to import


So I guess ORDER is a sqlite specific word that is reserved for SQL
operations and is not allowed as column name? Is this a general thing of
sqlite or is this just related to v.import/v.in.ogr and the creation of the
sqlite table?

cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] [snap pour point for r.water.outlet]

2017-03-30 Thread Johannes Radinger
Sorry I ment r.stream.snap...

 Original message From: Ang Sherpa 
<angsherpa...@gmail.com> Date:30/03/2017  07:32  (GMT+01:00) 
To: Johannes Radinger <johannesradin...@gmail.com> Cc: 
Helmut Kudrnovsky <hel...@web.de>, GRASS user list <grass-user@lists.osgeo.org> 
Subject: Re: [GRASS-user] [snap pour point for r.water.outlet] 

Thanks Johannes,

Installed. No I have not installed r.stream.order, is it mandatory or ease the 
process?


Best,
Ang Dawa Sherpa


On Thu, Mar 30, 2017 at 10:23 AM, Johannes Radinger 
<johannesradin...@gmail.com> wrote:
Dear Ang Sherpa,

r.stream.snap is an add-on. So you first need to install it using e.g 
g.extension. Have you installed r.stream.order before calling it?

/j


 Original message 
From: Ang Sherpa
Date:30/03/2017 05:34 (GMT+01:00)
To: Helmut Kudrnovsky , GRASS user list
Subject: Re: [GRASS-user] [snap pour point for r.water.outlet]

Thanks for the link Helmut,

However, it throws an error stating:
(Thu Mar 30 09:14:55 2017)  
r.stream.snap  
'r.stream.snap' is not recognized as an internal or external
command,
operable program or batch file.

I am using  standalone version Grass gis 7.2.0

Regards,
Ang Dawa Sherpa
GIS technician - Irrigation Master Plan
WRPPF - DOI, Nepal Government
Lalitpur
contact: 984 007 3861

On Thu, Mar 30, 2017 at 1:03 AM, Helmut Kudrnovsky <hel...@web.de> wrote:
Ang Sherpa wrote
> Hi users,
>
> While using "r.water.outlet" to delineate watershed basin, although the
> coordinates of the stream was noted from google earth and fed into
> "r.water.outlet" module, it produces plain raster.
>
> Is there any solution to make sure that the coordinates of outlet point
> automatically snaps to the nearest line of stream in drainage direction
> map?

what about?:

https://grass.osgeo.org/grass72/manuals/addons/r.stream.snap.html

r.stream.snap - Snap point to modelled stream network.
Input can be stream network, point vector map with outlets or outlet
coordinates.




-
best regards
Helmut
--
View this message in context: 
http://osgeo-org.1560.x6.nabble.com/snap-pour-point-for-r-water-outlet-tp5314835p5314856.html
Sent from the Grass - Users mailing list archive at Nabble.com.
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user


___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] [snap pour point for r.water.outlet]

2017-03-29 Thread Johannes Radinger
Dear Ang Sherpa,

r.stream.snap is an add-on. So you first need to install it using e.g 
g.extension. Have you installed r.stream.order before calling it?

/j

 Original message From: Ang Sherpa 
 Date:30/03/2017  05:34  (GMT+01:00) 
To: Helmut Kudrnovsky , GRASS user list 
 Subject: Re: [GRASS-user] [snap pour 
point for r.water.outlet] 
Thanks for the link Helmut,

However, it throws an error stating:
(Thu Mar 30 09:14:55 2017)  
r.stream.snap  
'r.stream.snap' is not recognized as an internal or external
command,
operable program or batch file.

I am using  standalone version Grass gis 7.2.0

Regards,
Ang Dawa Sherpa
GIS technician - Irrigation Master Plan
WRPPF - DOI, Nepal Government
Lalitpur
contact: 984 007 3861

On Thu, Mar 30, 2017 at 1:03 AM, Helmut Kudrnovsky  wrote:
Ang Sherpa wrote
> Hi users,
>
> While using "r.water.outlet" to delineate watershed basin, although the
> coordinates of the stream was noted from google earth and fed into
> "r.water.outlet" module, it produces plain raster.
>
> Is there any solution to make sure that the coordinates of outlet point
> automatically snaps to the nearest line of stream in drainage direction
> map?

what about?:

https://grass.osgeo.org/grass72/manuals/addons/r.stream.snap.html

r.stream.snap - Snap point to modelled stream network.
Input can be stream network, point vector map with outlets or outlet
coordinates.




-
best regards
Helmut
--
View this message in context: 
http://osgeo-org.1560.x6.nabble.com/snap-pour-point-for-r-water-outlet-tp5314835p5314856.html
Sent from the Grass - Users mailing list archive at Nabble.com.
___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Change elevation value for entire DEM raster cell [UPDATE 2]

2017-03-23 Thread Johannes Radinger
Rich,

Maybe also this conversation on setting a specific raster value using r.mapcalc 
might help: 
http://osgeo-org.1560.x6.nabble.com/Set-raster-value-for-specific-cell-td5052853.html

/j

 Original message From: Rich Shepard 
 Date:23/03/2017  22:07  (GMT+01:00) 
To: grass-user@lists.osgeo.org Subject: Re: [GRASS-user] 
Change elevation value for entire DEM raster cell
  [UPDATE 2] 
On Thu, 23 Mar 2017, Rich Shepard wrote:

> Thanks. I overlooked entering the column formats while entering the x, y,
> and z column positions.

   Turns out the problem is the pwd: it was ~/data/grassdata and I needed to
cd to the topography subdirectory where the text file is located. Sigh.

Mea culpa!

Rich

___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user___
grass-user mailing list
grass-user@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Auto-detection of GCP for rectifying image

2016-10-06 Thread Johannes Radinger
Hi GRASS users!

I have an image of a small object (e.g. photograph of blue square 15x5cm)
which I want to automatically georeference. In addition, the image contains
6 points where I know the exact coordinates, i.e. these points can be use
as ground control points (GCP) for rectifying the image. This is basically
a experimental setup, and I could use different colors or shapes for the
GCP; the image itself will then be rather a photograph of a fish on a white
background than a blue square; I just want to test with a simpler setup
before.

However, as I want to automatize the step of rectifying/georeferencing I am
looking for a way to autodetect these six points in the image. I am
thinking of tools like pattern/face recognition that are able to autodetect
objects (e.g. eyes, points etc.) and extract their position (coordinates)
within that image. I assume these coordinates together with their "true"
position could then be used for rectifying the picture using e.g. i.rectify.

Has anyone done this or a similar exercise before and can recommend tools
and approaches to auto-detect GCP from an image?

Any hints are welcome!

Best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] r.flow: define contributing area

2016-10-03 Thread Johannes Radinger
Rich,

Maybe have a look at r.lake: 
https://grass.osgeo.org/grass70/manuals/r.lake.html. This module fills the area 
upstream a dam or a blocked stream.

Cheers,
Johannes



 Original message From: Thomas Adams 
 Date:04/10/2016  04:32  (GMT+01:00) 
To: Rich Shepard  Cc: 
grass-user@lists.osgeo.org Subject: Re: [GRASS-user] r.flow: define 
contributing area 
Rich,

The only way to -correctly- do this is with hydrodynamic modeling, such as with 
HEC-RAS. It can be very crudely approximated with r.damflood 
(https://grass.osgeo.org/grass70/manuals/addons/r.damflood.html).

Tom

On Mon, Oct 3, 2016 at 7:39 PM, Rich Shepard  wrote:
  Reading the r.flow manual page suggests that the use of the '-u'
(upstream) option allows determination of flowlines and lengths that can be
used to delineate the area that would drain to a specific point. Is this
correct?

  To determine the area flooded if an outlet is plugged would require the
inverse of r.drain. That module calculates the area flooded with surface
water originating from a specific point; e.g., a hole in a dike or dam. The
closest I've seen to do the opposite seems to be r.flow. The quantity of
surface water and the rate of accumulation and infiltration into the vadoz
zone are not of interest, only the area flooded is needed.

  The modules listed under the Natural Hazards application page all appear
to delineate the area flooded from a point source. Can any of these be used
backwards to delineate the area flooded if an outlet is blocked?

  I'd appreciate pointers and suggestions based on cumulative experience
here for appropriate module(s).

TIA,

Rich
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user




___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Update categories/attribute table after v.clean with option break

2016-09-07 Thread Johannes Radinger
Hi,

I have a vector line map (river network) which I want to clean and break
the lines at intersections. Therefore, I use the tool v.clean with the
option break. This of course increases the number of vector lines as
several lines got broken into two or more separate lines. Now there is a
attribute table linked to the original map. I want to update this attribute
table by inserting a new column that contains a new unique ID (cat) for
each line but still keeps the old cat in another column (E.g. old_cat).
Thus the attribute table should also increase in number of rows (=number of
lines).

Now, most likely the tool to use now is v.category to get a column of new
unique categories (in another table/layer?) and the use probably a join
statement to get back the information from the original table (with doubled
entries for those lines that got broken during v.clean).

What I tried so far:

#
# Cleaning vector and breaking lines
v.clean -c --overwrite input=river_network_modified
output=river_network_modified_clean tool=break,snap,rmline threshold=0,50,0

# Check cats
v.category input=river_network_modified_clean option=report layer=-1
# type   countminmax
# line   16592  1  16448

# Update cats
v.category -t --overwrite input=river_network_modified_clean type=line
output=river_network_modified_clean2 option=del
v.category -t --overwrite input=river_network_modified_clean2 type=line
output=river_network_modified_clean3 option=add

# Check cats
v.category input=river_network_modified_clean2 option=report layer=-1
# Layer: 1
# type   countminmax
# line   16591  2  16448

# Check cats
v.category input=river_network_modified_clean3 option=report layer=-1
# Layer: 1
# type   countminmax
# line   16592  1  16448


However, I expected the max values of the latest (updated) vector to be
16592 to have really unique values for all 16592 lines (from 1 to 16592)?
But the output shows the same information as the initial input map?
Probably I am still mixing up things with categories/layers... and maybe
someone has a quick hint how to get an updated attribute table for a map
that has been processed using v.clean and where lines were broken.

/johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Select specific layer when copying vector map

2016-08-09 Thread Johannes Radinger
Indeed, v.extract using the layer option did the job!
Thank you Moritz!

Probably, I can then use v.category option="chlayer" to change the layer
Nr. 3 to Nr. 1!?


On Tue, Aug 9, 2016 at 2:23 PM, Moritz Lennert <mlenn...@club.worldonline.be
> wrote:

> On 09/08/16 13:48, Johannes Radinger wrote:
>
>> Hi,
>>
>> I'd like to create from an existing GRASS vector map that has 3 layers
>> (with three attribute tables) a new vector map that has only one
>> layer and its associated attribute table (here the information related
>> to layer 3).
>>
>> What I am doing so far (in python):
>>
>> # First copy the map to a new one and copy table (that of layer 3)
>> # that I want to re-connect to new map
>> grass.run_command("g.copy",
>> overwrite=True,
>> vector="oldmap,newmap")
>> grass.run_command("db.copy",
>> overwrite=True,
>> from_table="table_layer3",
>> to_table="table_layer3_copy")
>>
>> # Second delete all tables of the newmap
>> for i in [1,2,3]:
>> grass.run_command("v.db.droptable",
>> flags="f",
>> map="newmap",
>> layer=i)
>>
>> # Third reconnect a table to the newmap
>> grass.run_command("v.db.connect",
>> overwrite=True,
>> map="newmap",
>> table="table_layer3_copy")
>>
>> When I check in the attribute table manager, there is only one layer
>> listed (layer 1).
>> Similarly, v.db.connect reports:
>> v.db.connect -p map=newmap
>> Vector map  is connected by layer <1/table_layer3_copy> table
>>  in database [...]
>>
>> However, in d.vect I can still select from all three layers that were in
>> the old map (-1,1,2,3),
>> and all of them can be displayed (with -1 specified)?
>> So does may newmap have three layers or only one?
>>
>
> 3. Layers are defined by category values, not by attribute tables.
>
>
>> How can I make the third layer the default layer 1 and remove all other
>> layers?
>> Or is there a way to copy a vector map with just a single layer and its
>> associated table specified as default?
>>
>
> Try v.extract instead of g.copy.
>
> Moritz
>
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Select specific layer when copying vector map

2016-08-09 Thread Johannes Radinger
Hi,

I'd like to create from an existing GRASS vector map that has 3 layers
(with three attribute tables) a new vector map that has only one layer and
its associated attribute table (here the information related to layer 3).

What I am doing so far (in python):

# First copy the map to a new one and copy table (that of layer 3)
# that I want to re-connect to new map
grass.run_command("g.copy",
overwrite=True,
vector="oldmap,newmap")
grass.run_command("db.copy",
overwrite=True,
from_table="table_layer3",
to_table="table_layer3_copy")

# Second delete all tables of the newmap
for i in [1,2,3]:
grass.run_command("v.db.droptable",
flags="f",
map="newmap",
layer=i)

# Third reconnect a table to the newmap
grass.run_command("v.db.connect",
overwrite=True,
map="newmap",
table="table_layer3_copy")

When I check in the attribute table manager, there is only one layer listed
(layer 1).
Similarly, v.db.connect reports:
v.db.connect -p map=newmap
Vector map  is connected by layer <1/table_layer3_copy> table
 in database [...]

However, in d.vect I can still select from all three layers that were in
the old map (-1,1,2,3),
and all of them can be displayed (with -1 specified)?
So does may newmap have three layers or only one?

How can I make the third layer the default layer 1 and remove all other
layers?
Or is there a way to copy a vector map with just a single layer and its
associated table specified as default?

Thank you very much!

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Updating v.net-report results directly to attribute table

2016-06-29 Thread Johannes Radinger
Hi all,

is there an easy way to transfer the information reported by v.net
(operation=report)
directly to the attribute table of the vector network (arc-table)?

The approach I am using so far is to create an temporary database table
(using python)
where I store the output from v.net operation=report, and then use a SQL
statement to update the original arc-table from the temporary table:

###
fidimo_db.execute('''CREATE TEMP TABLE arcs_temp
(cat INTEGER, from INTEGER, to INTEGER)''')

# Get for each arc the orig cat for the start (from) and end point (to)
e = [(int(x.split()[0]),int(x.split()[1]),int(x.split()[2])) for x in
grass.read_command("v.net",
quiet=True,
input="my_net",
operation="report",
arc_layer=3).splitlines()]
my_db.executemany("INSERT INTO arcs_tmp (cat, from, to) VALUES (?,?,?)", e)

my_db.execute('''UPDATE arcs SET
from = (SELECT from FROM arcs_tmp WHERE cat=arcs.cat),
to = (SELECT to FROM arcs_tmp WHERE cat=arcs.cat)
WHERE EXISTS (SELECT cat FROM arcs_tmp WHERE cat=arcs.cat)''')
##

That approach works for me, but I was wondering if there is something
easier/more direct? However it seems a direct update of the attribute table
is not included in v.net. The module v.db.update has an option for adding
start/end points of lines, however this refers to coordinate pairs rather
than to category values of nodes in a network.

Any other ideas?

Best,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] New GRASS GIS addon: v.stream.order

2016-04-22 Thread Johannes Radinger
Dear GRASS Users,

As some of you might already have noticed, there is a new GRASS GIS add-on
available, called v.stream.order:
https://grass.osgeo.org/grass70/manuals/addons/v.stream.order.html

This module computes various types of stream order (Strahler, Shreve,
Scheidegger Drwal) of stream network vector maps. Specifically, stream
order is calculated for each sub-network as identified by their respective
outlet vector points. At its current state, v.stream.order works on GRASS
GIS 7.1svn.

The module v.stream.order includes a comprehensive manual and a testsuite
and is available via from the GRASS add-on SVN repository
https://svn.osgeo.org/grass/grass-addons/grass7/vector/v.stream.order/

Maybe a future extension of v.stream.order will also allow to consider more
complex networks that include e.g. loops and interconnections.

The module v.stream.order has been developed in the scope of the BiodiERsA
project FISHCON (Project partner: IGB Berlin) and was programmed and
implemented by Geoinformatikbüro Dassau GmbH, Soeren Gebbert.

I hope the tool is also of use for some of your future analysis.

Cheers,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Average raster value of the upslope contributing area,

2015-12-17 Thread Johannes Radinger
Hi Pierluigi,

I think here you could use the tool r.watershed (
https://grass.osgeo.org/grass70/manuals/r.watershed.html)
with your underlying DEM as input. This tool will provide you an
accumulation map with increasing values in
downstream direction.

I am not exactly sure if this is your task, but if you provide an input map
for 'flow' with values for each cell that represents
the actual cell size than the output accumulation map output is your
upslope contributing area. If you
use another input map for 'flow' (e.g. sediment input per cell) than the
output will be the accumulated sediment
in downstream direction. Subsequently you can use r.mapcalc to divide both
maps to get the ratio
of sediment delivery per unit upslope contributing area.

I hope that helps?!

cheers,
Johannes


On Thu, Dec 17, 2015 at 11:29 AM, Ing. Pierluigi De Rosa <
pierluigi.der...@gfosservices.it> wrote:

> Dear users,
>
> my deal is to calculate a raster of the sediment delivery ratio like
> explained here:
>
> http://data.naturalcapitalproject.org/nightly-build/invest-users-guide/html/sdr.html#sediment-delivery-ratio
>
> moreover I need to calculate, for each cell, the average of a raster (C
> factor) of the upslope contributing area.
> take a look to these image to better explication:
>
> http://data.naturalcapitalproject.org/nightly-build/invest-users-guide/html/_images/connectivity_diagram.png
>
> How can I do that?
> Thanks
> Pierluigi
>
>
> --
>
>
> 
> Ing. Pierluigi De Rosa (PhD)
> Studio Associato GFOSSERVICES
>
>
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Unique IDs for network segments

2015-12-04 Thread Johannes Radinger
Moritz,
thank you! Your proposed approach solved my problem!


On Thu, Dec 3, 2015 at 3:04 PM, Moritz Lennert <mlenn...@club.worldonline.be
> wrote:

> On 03/12/15 14:19, Johannes Radinger wrote:
>
>> Hi Moritz,
>>
>> The two-step v.category approach (del + add) and an additional v.to.db()
>>   was the
>> right way. Now I have the problem that there are 3 categories reported for
>> one single line that got split into to (1 old cat + 2 new cats). The old
>> category is still
>> associated with attributes but the two new rows are empty. How can
>> I get an attribute table that contains only rows for the new cats with
>> the attribute information
>> of the old category?
>>
>> [...]
>
>> # Remove existing cats
>> v.category --overwrite input=test_net type=line cat=-1
>> output=test_net_nocat option=del
>>
> [...]
>
>> # Add new cats
>> v.category --overwrite input=test_net_nocat type=line
>> output=test_net_newcat option=add
>> v.to.db map=test_net_newcat type=line option=cat columns=cat
>> #
>>
>
> Instead of removing the existing cats, add new cats in a new layer and
> transfer the info into that layer's attribute table:
>
>
> v.category --overwrite input=test_net type=line output=test_net_newcat
> option=add layer=3
> v.db.addtable test_new_newcat layer=3
> v.db.addcolumn test_net_newcat lay=3 col="col1 int, old_cat int"
> v.to.db test_net_newcat layer=3 op=query col=col1 query_layer=1
> query_colum=col1
> v.to.db test_net_newcat layer=3 op=query col=old_cat query_layer=1
> query_colum=cat
>
> Just make sure to set arc_layer=3 in the network analysis modules.
>
> If you really want to have arcs in layer 1, you can use v.category
> op=transfer, but you will also have to change table connections of the
> layers with v.db.connect.
>
> Moritz
>
>
>
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Unique IDs for network segments

2015-12-03 Thread Johannes Radinger
Hi Moritz,

The two-step v.category approach (del + add) and an additional v.to.db()
 was the
right way. Now I have the problem that there are 3 categories reported for
one single line that got split into to (1 old cat + 2 new cats). The old
category is still
associated with attributes but the two new rows are empty. How can
I get an attribute table that contains only rows for the new cats with the
attribute information
of the old category?

Here a working example (in the NC location) for illustrating my problem:
### Import test line to GRASS GIS
echo "ORGANIZATION: GRASS Development Team
DIGIT DATE:   03/12/2015
DIGIT NAME:   -
MAP NAME: test_lines
MAP DATE: 2015
MAP SCALE:1
OTHER INFO:   test lines
ZONE:  0
MAP THRESH:   0.50
VERTI:
L  3 1
 641045.70521951 218398.66077254
 641049.6221644 218261.56770124
 641026.12049504 218155.8101891
 1 5" | v.in.ascii --o format=standard input=- output=test_lines

# Add attribute table to imported line
v.db.addtable map=test_lines

# Add some attributes to table including the original cat as attribute
v.db.addcolumn map=test_lines@fidimo_test columns="col1 INT, old_cat INT"
v.db.update map=test_lines@fidimo_test column=col1 value=42
v.db.update map=test_lines@fidimo_test column=old_cat query_column=cat

# Import Point
echo "100|641045.940579|218320.669742" | v.in.ascii --o input=-
output=mypoints x=2 y=3 cat=1

# Create network
v.net -s --overwrite input=test_lines \
points=mypoints output=test_net operation=connect threshold=25

# Remove existing cats
v.category --overwrite input=test_net type=line cat=-1
output=test_net_nocat option=del

# Add new cats
v.category --overwrite input=test_net_nocat type=line
output=test_net_newcat option=add
v.to.db map=test_net_newcat type=line option=cat columns=cat
#

The attribute table of the final vector (test_net_newcat) should contain
only two rows for
cat 1 and 2 with the original attribute information for both rows. Is that
somehow possible
in GRASS?

/Johannes


On Wed, Dec 2, 2015 at 4:55 PM, Moritz Lennert <mlenn...@club.worldonline.be
> wrote:

> On 02/12/15 15:40, Johannes Radinger wrote:
>
>> Hi all,
>>
>> I built a GRASS7 network using v.net <http://v.net>:
>>
>> First, I used the 'connect' option to snap some points to my network
>> (using the -s flag)
>> Second, I used that first network and also added all nodes with the
>> 'nodes' option to get a final network.
>>
>> With the 'report' option of v.net <http://v.net> I can report all edges
>> with their start and end end points. However, these edges of the final
>> network are based on the initial category values of the input vector.
>> How can I get the network split into single lines with unique IDs
>> between each node of the final network. In other words I want a report
>> of v.net <http://v.net> that contains only unique edge categories? I'd
>> like to write the output of v.net <http://v.net> into a new table and
>> assign lengths to each edge.
>>
>> For example, at the moment, a line with cat=300 is confined by a start
>> node cat=1 and an end node cat=2. For that case v.net <http://v.net>
>> reports:
>> 300 1 2
>>
>> If I connect a new node with cat=3 along that line and calculate v.net
>> <http://v.net>, it reports:
>> 300 1 3
>> 300 3 2
>>
>> Any ideas/suggestions? Maybe I need to 'break' the network lines at each
>> node before running v.net <http://v.net>?
>>
>
> How about v.category op=delete type=line followed by v.category op=add
> type=line ?
>
> Moritz
>
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Unique IDs for network segments

2015-12-02 Thread Johannes Radinger
Hi all,

I built a GRASS7 network using v.net:

First, I used the 'connect' option to snap some points to my network (using
the -s flag)
Second, I used that first network and also added all nodes with the 'nodes'
option to get a final network.

With the 'report' option of v.net I can report all edges with their start
and end end points. However, these edges of the final network are based on
the initial category values of the input vector. How can I get the network
split into single lines with unique IDs between each node of the final
network. In other words I want a report of v.net that contains only unique
edge categories? I'd like to write the output of v.net into a new table and
assign lengths to each edge.

For example, at the moment, a line with cat=300 is confined by a start node
cat=1 and an end node cat=2. For that case v.net reports:
300 1 2

If I connect a new node with cat=3 along that line and calculate v.net, it
reports:
300 1 3
300 3 2

Any ideas/suggestions? Maybe I need to 'break' the network lines at each
node before running v.net?

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Error starting GRASS71 trunk: module catalog missing?

2015-11-24 Thread Johannes Radinger
Thanks Martin,

It seems that deleting the installation directory helped!

/Johannes

On Tue, Nov 24, 2015 at 2:58 PM, Martin Landa <landa.mar...@gmail.com>
wrote:

> Hi,
>
> 2015-11-24 14:39 GMT+01:00 Johannes Radinger <johannesradin...@gmail.com>:
> > from datacatalog.catalog   import DataCatalog
> > ImportError: No module named catalog
>
> several times discussed, please delete installation directory
> (/usr/local/grass-7.1.svn) and compile from scratch (distclean,
> configure, make, make install). Ma
>
> --
> Martin Landa
> http://geo.fsv.cvut.cz/gwiki/Landa
> http://gismentors.cz/mentors/landa
>
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Error starting GRASS71 trunk: module catalog missing?

2015-11-24 Thread Johannes Radinger
Hi,

I tried to start the most recent (SVN up) self-compiled version of GRASS71
on a Ubuntu-machine. GRASS starts and I can selected a location/mapset from
the GUI but after the selection I get following error and the GUI does not
start:

Launching  GUI in the background, please wait...
GRASS 7.1.svn (nc_spm_08):~ > NumPy not found.

This module requires the NumPy module, which could not be
imported.  It probably is not installed (it's not part of the
standard Python distribution). See the Numeric Python site
(http://numpy.scipy.org) for information on downloading source or
binaries.

Traceback (most recent call last):
  File "/usr/local/grass-7.1.svn/gui/wxpython/wxgui.py", line 140, in

sys.exit(main())
  File "/usr/local/grass-7.1.svn/gui/wxpython/wxgui.py", line 132, in main
app = GMApp(workspaceFile)
  File "/usr/local/grass-7.1.svn/gui/wxpython/wxgui.py", line 46, in
__init__
wx.App.__init__(self, False)
  File "/usr/lib/python2.7/dist-packages/wx-3.0-gtk2/wx/_core.py", line
8628, in __init__
self._BootstrapApp()
  File "/usr/lib/python2.7/dist-packages/wx-3.0-gtk2/wx/_core.py", line
8196, in _BootstrapApp
return _core_.PyApp__BootstrapApp(*args, **kwargs)
  File "/usr/local/grass-7.1.svn/gui/wxpython/wxgui.py", line 79, in OnInit
from lmgr.frame import GMFrame
  File "/usr/local/grass-7.1.svn/gui/wxpython/lmgr/frame.py", line 73, in

from datacatalog.catalog   import DataCatalog
ImportError: No module named catalog

I am not sure why I get the errors concerning numpy as this is installed
and can be loaded:
radinger@grassgis:~$ python
Python 2.7.3 (default, Jun 22 2015, 19:43:34)
[GCC 4.6.3] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> import numpy
>>>

However, I am more worried about the error (or the failure of the GUI which
does not start) related to the  messages on the wx-module "No module named
catalog". There were no errors during configuration and compilation of the
source code of GRASS.

I have also compiled and installed GRASS 7.0.2svn, here the GUI starts
without any problems (however, I get the messages concerning numpy, but it
seems they don't have any effect on running GRASS)

So what is the reason for GRASS71 to fail here? I have installed wxpython
3.0.1:
>>> import wx
>>> wx.VERSION_STRING
'3.0.1.0'
>>> wx.version()
'3.0.1.0 gtk2 (classic)'
>>>

Any suggestions?
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] river and dams in netwok analysis

2015-10-26 Thread Johannes Radinger
Hi Etienne,

Maybe you could use v.net.allpairs with negative/positve costs for up- and 
downstream direction and then select the smallest positive and negative 
distance for each dam to get the two neighbours.

Here what I found in the manual of v.net.distance which might also be useful 
for your purpose: 

"In order to find nearest neighbors within a group of nodes, you can either 
loop through each node as to and all other nodes as from or create a complete 
distance matrix with v.net.allpairs and select the lowest non-zero distance for 
each node."

https://grass.osgeo.org/grass70/manuals/v.net.distance.html

Maybe a start...

/Johannes



 Original message From: Etienne DELAY 
 Date:26/10/2015  15:15  (GMT+01:00) 
To: Markus Neteler  Cc: GRASS user 
list  Subject: Re: [GRASS-user] river 
and dams in netwok analysis 
Dear all dear Markus,
Above all great job GRASS dev team ! Grass 7 is so nice and powerful. 
I'm always in the dark with my dams problem.

I would like to export dams points as node and for each node in the 
attribute table 2 new column : dam before | dam after. I'm looking for 
network analysis but it may be a wrong way ?

There is a way to do that in Grass?


Cordialement

Etienne DELAY (Skype : etienne.delay.tic)
Chaire: Capital environnemental et gestion durable des cours d'eau
laboratoire GEOLAB UMR 6042 CNRS
Université de Limoges, FLSH
39E rue Camille Guérin 87036 Limoges
blog : http://elcep.legtux.org

Le 18/10/2015 17:30, Markus Neteler a écrit :
> On Sat, Oct 17, 2015 at 11:17 AM, Etienne DELAY
> > wrote:
>  > Re,
>  > I have read carefully the v.out.ogr as will said for tring to export the
>  > network layer but I've got an error :
>  >
>  > v.out.ogr --overwrite input=network@PERMANENT type=point dsn=~/ layer=2
>  > format=ESRI_Shapefile
>  >
>  > WARNING: 111738 line(s) found, but not requested to be exported. Verify
>  >  'type' parameter.
> ...
>
> (I guess that you need to set "type=line" and not "type=point")
>
> ... this is why you may consider to upgrade to the current 7.0.x stable
> version. Note the differences:
>
> https://grass.osgeo.org/grass64/manuals/v.out.ogr.html
> type=string[,string,...] Feature type(s).
>Combinations not supported by all output formats. Default: first type
> found in input.
>Options: point,line,boundary,centroid,area,face,kernel,auto
>Default: line,boundary
>
> https://grass.osgeo.org/grass70/manuals/v.out.ogr.html
> type=string[,string,...] Feature type(s)
>Combination of types is not supported by all output formats.
>Default is to use first type found in input vector map.
>Options: point, line, boundary, centroid, area, face, kernel, auto
>Default: auto  <<==
>
> Less to think about, much faster vector data processing, improved
> network analysis along with new graphical frontend:
> https://grasswiki.osgeo.org/wiki/Vector_network_analysis#Screenshots
>
> Best
> Markus
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Problems running g.gui.gcp

2015-10-23 Thread Johannes Radinger
On Fri, Oct 23, 2015 at 10:17 AM, Moritz Lennert <
mlenn...@club.worldonline.be> wrote:

> On 22/10/15 16:25, Johannes Radinger wrote:
>
>> Dear GRASS users,
>>
>> Today I wanted to try the georectification tool g.gui.gcp in GRASS 7 for
>> the first time. I was basically following the description on the GRASS
>> wiki (https://grasswiki.osgeo.org/wiki/Georeferencing). In particular, I
>> have a photograph from my camera with known reference points which I
>> load into an unprojected XY location. My target location should also be
>> a XY location as the picture is taken from a "virtual" coordinate
>> system. Inside my target XY-location I try to open g.gui.gcp and select
>> raster, and the location and mapset where my photograph is stored. When
>> I want to click next I get following error/message: "You must select a
>> valid location and mapset in order to continue" and I can't proceed with
>> the rectification. Just this GUI message appears, no error/warning in
>> the grass console etc. I am puzzled how to interpret that? What is wrong
>> with the location/mapset? The mapset I've selected is existing (I can
>> select it from the dropdown menu). I tried to select any other
>> location/mapset but receive always the same message.
>>
>
>
> Georectification is for georeferencing to a projected system. It does not
> make sense to georeference to an XY-location.
>
> I see two options for you:
>
> - Import you image directely into the target region and use r.region to
> adjust its boundary coordinates
> - Define your own projection system to fit your needs.


So far as I understand georectification and orthorectification in
particular is to rectify a picture/photo/image based on spatially defined
ground control points which might not be done by just adjusting the
boundary coordinates. BTW I also tried the European LAEA coordinate system
(EPSG 3035) as a target location which is not working either here in my
machine/setup (same message as already mentioned). Here (
https://www.dropbox.com/s/snjzd8eu7b8115k/2015-10-22%2015.26.30.jpg?dl=0)
you can find a a photo (*.jpg) which is an "areal" photo with known 3d
coordinates for the single points (ground control points). I think this
example clarifies and illustrates what I intend to do and why the target
system can also be a XY location in my opinion. A rectification of this
image should make the lines of the grid more or less parallel. Or do I get
something completely wrong?

/Johannes



>
>
> Moritz
>
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Problems running g.gui.gcp

2015-10-23 Thread Johannes Radinger
On Fri, Oct 23, 2015 at 11:19 AM, Moritz Lennert <
mlenn...@club.worldonline.be> wrote:

> On 23/10/15 10:56, Johannes Radinger wrote:
>
>>
>>
>> On Fri, Oct 23, 2015 at 10:17 AM, Moritz Lennert
>> <mlenn...@club.worldonline.be <mailto:mlenn...@club.worldonline.be>>
>> wrote:
>>
>> On 22/10/15 16:25, Johannes Radinger wrote:
>>
>> Dear GRASS users,
>>
>> Today I wanted to try the georectification tool g.gui.gcp in
>> GRASS 7 for
>> the first time. I was basically following the description on the
>> GRASS
>> wiki (https://grasswiki.osgeo.org/wiki/Georeferencing). In
>> particular, I
>> have a photograph from my camera with known reference points
>> which I
>> load into an unprojected XY location. My target location should
>> also be
>> a XY location as the picture is taken from a "virtual" coordinate
>> system. Inside my target XY-location I try to open g.gui.gcp and
>> select
>> raster, and the location and mapset where my photograph is
>> stored. When
>> I want to click next I get following error/message: "You must
>> select a
>> valid location and mapset in order to continue" and I can't
>> proceed with
>> the rectification. Just this GUI message appears, no
>> error/warning in
>> the grass console etc. I am puzzled how to interpret that? What
>> is wrong
>> with the location/mapset? The mapset I've selected is existing
>> (I can
>> select it from the dropdown menu). I tried to select any other
>> location/mapset but receive always the same message.
>>
>>
>>
>> Georectification is for georeferencing to a projected system. It
>> does not make sense to georeference to an XY-location.
>>
>> I see two options for you:
>>
>> - Import you image directely into the target region and use r.region
>> to adjust its boundary coordinates
>> - Define your own projection system to fit your needs.
>>
>>
>> So far as I understand georectification and orthorectification in
>> particular is to rectify a picture/photo/image based on spatially
>> defined ground control points which might not be done by just adjusting
>> the boundary coordinates. BTW I also tried the European LAEA coordinate
>> system (EPSG 3035) as a target location which is not working either here
>> in my machine/setup (same message as already mentioned). Here
>> (https://www.dropbox.com/s/snjzd8eu7b8115k/2015-10-22%2015.26.30.jpg?dl=0
>> )
>> you can find a a photo (*.jpg) which is an "areal" photo with known 3d
>> coordinates for the single points (ground control points). I think this
>> example clarifies and illustrates what I intend to do and why the target
>> system can also be a XY location in my opinion. A rectification of this
>> image should make the lines of the grid more or less parallel. Or do I
>> get something completely wrong?
>>
>
> Actually, I was wrong: I don't get any complaint trying to use g.gui.gcp
> to georectify from XY to XY. So there seems to be an issue with your XY
> locations.
>
> Have you tried recreating the locations from scratch ?


Your results look really reasonable that's what I was looking for...

Actually, both XY locations are created from scratch (via the location
wizzard)
I get following output for both XY locations:

g.proj -p

XY location (unprojected)

I tried to start the GCP manager from within different other locations
(e.g. laea, datum: etrs89) and selecting different source location (e.g.
WGS84)however I always get the same message. So I guess this is not
specifically related to my location but rather to the GCP manager/WX GUI?
As I don't get an error message it is difficult to judge what is going
wrong in my case...

I just looked up the source code and the error message seems to be related
to follwoing line:
https://trac.osgeo.org/grass/browser/grass/branches/releasebranch_7_0/gui/wxpython/gcp/manager.py#L399
As far as I understand, this error message is only given if no
 location/mapset is provided; however that is not the case here. Maybe my
selection of the location/mapset via the GUI does not get parsed
further?

Strange that orthorectification from one XY to antother XY location works
on other machines but not on mine.

PS: As an important side note: I am working on Ubuntu which is installed on
a virtual machine (VMWare) on our server. I remotely connect to this
virtual machine via RDP.


>
> Moritz
>
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Problems running g.gui.gcp

2015-10-23 Thread Johannes Radinger
On Fri, Oct 23, 2015 at 1:37 PM, Anna Petrášová <kratocha...@gmail.com>
wrote:

>
> On Oct 23, 2015 6:04 AM, "Johannes Radinger" <johannesradin...@gmail.com>
> wrote:
> >
> >
> >
> > On Fri, Oct 23, 2015 at 11:19 AM, Moritz Lennert <
> mlenn...@club.worldonline.be> wrote:
> >>
> >> On 23/10/15 10:56, Johannes Radinger wrote:
> >>>
> >>>
> >>>
> >>> On Fri, Oct 23, 2015 at 10:17 AM, Moritz Lennert
> >>> <mlenn...@club.worldonline.be <mailto:mlenn...@club.worldonline.be>>
> wrote:
> >>>
> >>> On 22/10/15 16:25, Johannes Radinger wrote:
> >>>
> >>> Dear GRASS users,
> >>>
> >>> Today I wanted to try the georectification tool g.gui.gcp in
> >>> GRASS 7 for
> >>> the first time. I was basically following the description on
> the
> >>> GRASS
> >>> wiki (https://grasswiki.osgeo.org/wiki/Georeferencing). In
> >>> particular, I
> >>> have a photograph from my camera with known reference points
> which I
> >>> load into an unprojected XY location. My target location should
> >>> also be
> >>> a XY location as the picture is taken from a "virtual"
> coordinate
> >>> system. Inside my target XY-location I try to open g.gui.gcp
> and
> >>> select
> >>> raster, and the location and mapset where my photograph is
> >>> stored. When
> >>> I want to click next I get following error/message: "You must
> >>> select a
> >>> valid location and mapset in order to continue" and I can't
> >>> proceed with
> >>> the rectification. Just this GUI message appears, no
> >>> error/warning in
> >>> the grass console etc. I am puzzled how to interpret that? What
> >>> is wrong
> >>> with the location/mapset? The mapset I've selected is existing
> >>> (I can
> >>> select it from the dropdown menu). I tried to select any other
> >>> location/mapset but receive always the same message.
> >>>
> >>>
> >>>
> >>> Georectification is for georeferencing to a projected system. It
> >>> does not make sense to georeference to an XY-location.
> >>>
> >>> I see two options for you:
> >>>
> >>> - Import you image directely into the target region and use
> r.region
> >>> to adjust its boundary coordinates
> >>> - Define your own projection system to fit your needs.
> >>>
> >>>
> >>> So far as I understand georectification and orthorectification in
> >>> particular is to rectify a picture/photo/image based on spatially
> >>> defined ground control points which might not be done by just adjusting
> >>> the boundary coordinates. BTW I also tried the European LAEA coordinate
> >>> system (EPSG 3035) as a target location which is not working either
> here
> >>> in my machine/setup (same message as already mentioned). Here
> >>> (
> https://www.dropbox.com/s/snjzd8eu7b8115k/2015-10-22%2015.26.30.jpg?dl=0)
> >>> you can find a a photo (*.jpg) which is an "areal" photo with known 3d
> >>> coordinates for the single points (ground control points). I think this
> >>> example clarifies and illustrates what I intend to do and why the
> target
> >>> system can also be a XY location in my opinion. A rectification of this
> >>> image should make the lines of the grid more or less parallel. Or do I
> >>> get something completely wrong?
> >>
> >>
> >> Actually, I was wrong: I don't get any complaint trying to use
> g.gui.gcp to georectify from XY to XY. So there seems to be an issue with
> your XY locations.
> >>
> >> Have you tried recreating the locations from scratch ?
> >
> >
> > Your results look really reasonable that's what I was looking for...
> >
> > Actually, both XY locations are created from scratch (via the location
> wizzard)
> > I get following output for both XY locations:
> >
> > g.proj -p
>
> > XY location (unprojected)
> >
> > I tried to start the GCP manager from within different other locations
> 

[GRASS-user] Problems running g.gui.gcp

2015-10-22 Thread Johannes Radinger
Dear GRASS users,

Today I wanted to try the georectification tool g.gui.gcp in GRASS 7 for
the first time. I was basically following the description on the GRASS wiki
(https://grasswiki.osgeo.org/wiki/Georeferencing). In particular, I have a
photograph from my camera with known reference points which I load into an
unprojected XY location. My target location should also be a XY location as
the picture is taken from a "virtual" coordinate system. Inside my target
XY-location I try to open g.gui.gcp and select raster, and the location and
mapset where my photograph is stored. When I want to click next I get
following error/message: "You must select a valid location and mapset in
order to continue" and I can't proceed with the rectification. Just this
GUI message appears, no error/warning in the grass console etc. I am
puzzled how to interpret that? What is wrong with the location/mapset? The
mapset I've selected is existing (I can select it from the dropdown menu).
I tried to select any other location/mapset but receive always the same
message.

I am running GRASS 7.0.2 on Ubuntu:
GRASS version: 7.0.2svn

GRASS SVN Revision: 66568

Build Date: 2015-10-22

Build Platform: i686-pc-linux-gnu

GDAL/OGR: 1.10.0

PROJ.4: 4.8.0

GEOS: 3.4.2

SQLite: 3.7.9

Python: 2.7.3

wxPython: 2.9.4.1

Platform: Linux-3.2.0-92-generic-pae-i686-with-Ubuntu-12.04-precise

Any ideas?
Thanks!

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Need help in reprojecting raster images

2015-09-23 Thread Johannes Radinger
Hi Uttam,

if I understand your problem correctly, you can use the -p flag of r.proj
to print input map's bounds in the current projection. With this
information you can set the region boundaries appropriately using g.region.

Afterwards you can 'import' your raster into your laea-projection location
using the r.proj as mentioned before.

best,
Johannes

On Wed, Sep 23, 2015 at 8:21 PM, Uttam Sinha  wrote:

>
>
> Hi,
>
> I tried the way you told. I created a new location in laea without setting
> the default boundary. Then from this mapset I used
>
> r.proj location*=input_location* mapset=*input_mapset* input=*input_raster
> name* output=
>
> *output_raster*
> I get an error:
>
>
> ERROR: Input raster map is outside current region
>
> Any suggestions.
>
> Uttam.
>
>
> On Tue, Sep 22, 2015 at 11:59 PM, sajid pareeth 
> wrote:
>
>> Hi Uttam
>>
>> If I understand you correctly you want to reproject a landsat image from
>> sinusoidal to laea projection!!
>> You dont have to manually set the bounding box while reprojecting a
>> raster, as GRASS will take care of it.
>>
>> First you create a new location and mapset with desired output projection
>> , which in this case is Lamberts Azimuthal Equal Area.
>> Then get into this mapset, use the following r.proj command.
>> r.proj location*=input_location* mapset=*input_mapset* input=*input_raster
>> name* output=
>> *output_raster*
>>
>> More details here: https://grass.osgeo.org/grass70/manuals/r.proj.html
>>
>> If you then need to set your working area to a particular bounding box,
>> use g.region
>>
>> g.region raster=*output_raster*
>>
>> regards
>>
>> Sajid
>>
>> On Tue, Sep 22, 2015 at 11:07 PM, Uttam Sinha 
>> wrote:
>>
>>>
>>>
>>> Hi All,
>>>
>>> I have a set of Landsat images in Sinusoidal Projection, WGS84 Datum. I
>>> know the north, south, east and west extent of the boundary of the image.
>>>
>>> I want to reproject the image to Lamberts Azimuthal Equal Area, WGS84
>>> Datum. I do not know what will be the coordinates of the north, south, east
>>> and west bounds in this projection.
>>>
>>> 1.) How do I convert upper left coordinates and bottom right coordinates
>>> of Sinusoidal Projection to  Lamberts Azimuthal Equal Area to create a new
>>> destination Location?
>>>
>>> 2.) How do I reproject the images from Sinusoidal to Lamberts without
>>> georeferencing the image at source location (in sinusoidal projection)
>>>
>>> Appreciate any reply.
>>>
>>> Uttam.
>>>
>>> ___
>>> grass-user mailing list
>>> grass-user@lists.osgeo.org
>>> http://lists.osgeo.org/mailman/listinfo/grass-user
>>>
>>
>>
>
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user
>
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Stream extraction from pseudo-elevation map using r.watershed

2015-07-13 Thread Johannes Radinger
Dear GRASS users,

I further elaborated on this issue of r.watershed and I found what is
actually causing this strange behaviour. It seems that r.watershed cannot
handle (pseudo)elevation values data are larger than 210. My procedure
is working value if the values are smaller, but I get an error message (see
previous email) if my elevation map contains e.g. 220.

Here small working example using the NC dataset to reproduce the behaviour
of getting stream segements (and stream order) from a raster river network
without having the original elevation map available.

##
# Extract river network to get a working example
g.region raster=elev_ned_30m@PERMANENT
r.watershed -m --overwrite elevation=elev_ned_30m@PERMANENT threshold=250
drainage=drainage stream=stream
r.stream.basins direction=drainage coordinates=638985, 218135
basins=selected_basin
r.mapcalc selected_streams = if( selected_basin   stream,1,null())
--overwrite


# For the actual purpose only a raster stream network is available (no DEM)
# So we need to calculate a pseudoelevation map using r.cost
r.cost --overwrite input=selected_streams output=stream_accum
start_coordinates=638985, 218135

# Create a buffer around stream (without outlet getting buffered). See
difference if the value is 220
r.grow.distance --overwrite input=stream_accum value=stream_raster_nearest
r.grow --overwrite input=stream_accum output=stream_accum_grow radius=2.01
old=1 new=210
r.mapcalc  stream_accum_grow_start = if( stream_raster_nearest  ==
0,null(), stream_accum_grow) --overwrite
r.grow --overwrite input=stream_accum_grow_start
output=stream_raster_buffer old=210 new=210
r.patch --overwrite input=stream_accum,stream_raster_buffer
output=stream_accum_buffered

# Extract stream with buffer does  work only if buffer value =210
r.watershed -m --overwrite elevation=stream_accum_buffered threshold=3
stream=stream_test2
#

Of course I agree a value of 220 is pretty much extreme for a raster
elevation map. However, from the error message it is difficult to actually
know what is causing the problem. Is there any special reason for that
threshold? Maybe r.watershed can have a test before to check if this
threshold is not exceeded. Most probably this issue can be solved just by
dividing the elevation by a certain value (e.g. 1000) before using it with
r.watershed.

Cheers,
Johannes


On Thu, Jul 2, 2015 at 2:18 PM, Johannes Radinger 
johannesradin...@gmail.com wrote:

 Hi,

 I tried something that was working at least several months ago, but
 somehow r.watershed has changed. What I am trying to do is to extract a
 stream network from a pseudo elevation map. This pseudo elevation map is
 created by getting the distance of each stream cell from to outlet in
 upstream direction (using r.cost). In addition, the streams are buffered (3
 cells wide) with cells of a value that is much greater than the actual
 stream value. This map has several features: 1) The cell values for the
 stream are certainly decreasing in downstream direction (=elevation is
 decreasing) 2) The stream is carved into a valley (important at river
 junctions where two streams meet). So this should actually meet the
 requirements for r.watershed (at least it was working some time ago) to
 extract the stream and get the drainage directions (e.g. for other tools
 like r.stream.order etc.).

 Now, when I try to run r.watershed (GRASS7.1) I get following warning
 several times:
 WARNING: MFD: cumulative proportion of flow distribution not 1.0 but
  -0.00

 Although I get a drainage map as output, no output is created for the
 extracted stream network. I am wondering what changed and why and if this
 is a bug in r.watershed or if I need to change settings differently.

 Attached you can find the pseudo elevation map as geotiff to try it
 yourself!

 The command I am using is:
 r.watershed -m --o  elevation=pseudo_elevation drainage=drainage_test
 stream=stream_test threshold=3

 My system:
 GRASS version: 7.1.svn

 GRASS SVN revision: 65534

 Build date: 2015-07-02

 Build platform: i686-pc-linux-gnu

 GDAL: 1.10.0

 PROJ.4: 4.8.0

 GEOS: 3.4.2

 SQLite: 3.7.9

 Python: 2.7.3

 wxPython: 2.9.4.1

 Any ideas or suggestions?

 Best regards,
 Johannes

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Stream extraction from pseudo-elevation map using r.watershed

2015-07-02 Thread Johannes Radinger
Hi,

I tried something that was working at least several months ago, but somehow
r.watershed has changed. What I am trying to do is to extract a stream
network from a pseudo elevation map. This pseudo elevation map is created
by getting the distance of each stream cell from to outlet in upstream
direction (using r.cost). In addition, the streams are buffered (3 cells
wide) with cells of a value that is much greater than the actual stream
value. This map has several features: 1) The cell values for the stream are
certainly decreasing in downstream direction (=elevation is decreasing) 2)
The stream is carved into a valley (important at river junctions where two
streams meet). So this should actually meet the requirements for
r.watershed (at least it was working some time ago) to extract the stream
and get the drainage directions (e.g. for other tools like r.stream.order
etc.).

Now, when I try to run r.watershed (GRASS7.1) I get following warning
several times:
WARNING: MFD: cumulative proportion of flow distribution not 1.0 but
 -0.00

Although I get a drainage map as output, no output is created for the
extracted stream network. I am wondering what changed and why and if this
is a bug in r.watershed or if I need to change settings differently.

Attached you can find the pseudo elevation map as geotiff to try it
yourself!

The command I am using is:
r.watershed -m --o  elevation=pseudo_elevation drainage=drainage_test
stream=stream_test threshold=3

My system:
GRASS version: 7.1.svn

GRASS SVN revision: 65534

Build date: 2015-07-02

Build platform: i686-pc-linux-gnu

GDAL: 1.10.0

PROJ.4: 4.8.0

GEOS: 3.4.2

SQLite: 3.7.9

Python: 2.7.3

wxPython: 2.9.4.1

Any ideas or suggestions?

Best regards,
Johannes


pseudo_elevation.gtiff
Description: Binary data
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Dataset file for add-on

2015-05-02 Thread Johannes Radinger
Hi,

I'd like to use a HDF5 datafile (multidimensional matrix with
parameter values) to be used and called from an python add-on.
Thus, I thought about providing the data file together with
the python add-on an load the data via a relative path (relative
to the path of the add-on). Is it generally possible to include
a datafile in the respecitive add-on folder so that it is also installed
together with the add-on. And if so, how should it be loaded (of course
with the hd5py functionalities) from within the python script (so that
it works on Windows and Linux platforms)

Best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Copy raster map from one unprojected location to another

2015-03-19 Thread Johannes Radinger
Hi,

I'd like to import a raster map from an unprojected location (XY-location)
to my actual location (also unprojected XY). Because these are two
different locations I can't use g.copy. So I tried r.proj. Here I get the
error that the projection is unknown, which makes sense. Of course, I could
export the map to geotiff and import them in the other location; however I
thought there might be a simpler GRASS way?!

Any suggestions? Or do I miss anything?

best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Stream order for vector networks

2015-03-19 Thread Johannes Radinger
Stefan,
thanks for your recommendation of the igraph-package. This would make
indeed sense for more complex analysis.
At the moment I am just interested in an GRASS (native) way to get a stream
order for a vector network. For this, purpose, and I guess also for using
the igraph package, one needs the information of the start (from_node) and
endnode (to_node) of each arc. So I am still looking how this information
can be extracted in GRASS and e.g. written into the attribute table?!

Best regards,
Johannes

On Wed, Mar 18, 2015 at 10:43 AM, Blumentrath, Stefan 
stefan.blumentr...@nina.no wrote:

  Hi Johannes,



 For more complex network analysis I can warmly recommend the igraph
 package (http://igraph.org/, which has both Python and R bindings).

 It comes very handy for customized network analysis, is not too
 complicated to learn, well documented and quite efficient.

 I could provide you with some sample code in R if that would be of
 interest.



 Cheers

 Stefan



 *From:* grass-user-boun...@lists.osgeo.org [mailto:
 grass-user-boun...@lists.osgeo.org] *On Behalf Of *Johannes Radinger
 *Sent:* 18. mars 2015 10:30
 *To:* GRASS user list
 *Subject:* [GRASS-user] Stream order for vector networks



 Hi,



 I am interested in calculating stream order (Strahler, Shreve) for stream
 networks based on GRASS vector networks (e.g. created by v.net) without
 using the raster approach of the GRASS add-on r.stream.order. Therefore, I
 came across following paper:



 Gleyzer, A., Denisyuk, M., Rimmer, A. and Salingar, Y. (2004), A FAST
 RECURSIVE GIS ALGORITHM FOR COMPUTING STRAHLER STREAM ORDER IN BRAIDED AND
 NONBRAIDED NETWORKS. JAWRA Journal of the American Water Resources
 Association, 40: 937–946. doi: 10./j.1752-1688.2004.tb01057.x



 which also provides pseudocode to calculate Strahler order for non-braided
 and braided vector networks defined by arcs and nodes. Thus, I am wondering
 if it is possible to get the information from which node a specific arc
 starts (from_node) and where it drains to (to_node). BTW, I seems some kind
 of this information is already provided with the NC-dataset (FNODE_ and
 TNODE column) in the streams attribute table, however I don't know where
 this information comes from. Of course to extract this information the
 network needs to be directed to define what is the upstream and what the
 downstream node of each arc. However, having this information available, it
 should not be to complicated to apply the algorithm proposed by Gleyzer et
 al 2004. Thus I am interested to get 'from_node' and 'to_node' information
 e.g. stored in the attribute table of v.net output?



 Here the pseudocode provided by Gleyzer et al 2004:



 

 MakeDictionaries(Network)

 1 for each Arc ∈ Network

 2 do FromNodesPerArc[Arc's ID] ← Arc’s FromNodeID

 3 InflowingArcsPerNode[Arc's ToNodeID] ← InflowingArcsPerNode[Arc’s
 ToNodeID] ∪ Arc’s ID

 4 OriginatingNode[Arc's ID] ← Arc’s FromNodeID



 StreamOrdering(ArcID)

 1 Visited[ArcID] ← true

 2 if | InflowingArcsPerNode[ FromNodesPerArc[ArcID] ] | = 0

 3 then StreamOrders[ArcID] ← 1

 4 else

 5 for each Arc ∈ InflowingArcsPerNode[ FromNodesPerArc[ArcID] ]

 6 do if not Visited[Arc]

 7 then UpstreamOrders[Arc] ← (StreamOrdering(Arc), OriginatingNode[Arc])

 8 else UpstreamOrders[Arc] ← (StreamOrders[Arc], OriginatingNode[Arc])

 9 MaxOrder ← 0

 10 MaxOrderCount ← 0

 11 for each (Order, Origin) ∈ UpstreamOrders

 12 do if Order  MaxOrder

 13 then MaxOrder ← Order

 14 MaxOrderCount ← 1

 15 MaxOrderOrigin ← Origin

 16 else if Order = MaxOrder

 17 then if Origin ≠ MaxOrderOrigin

 18 then MaxOrderCount ← MaxOrderCount + 1

 19 if MaxOrderCount  1

 20 then StreamOrders[ArcID] ← MaxOrder + 1

 21 OriginatingNode[ArcID] ← FromNodesPerArc[ArcID]

 22 else StreamOrders[ArcID] ← MaxOrder

 23 OriginatingNode[ArcID] ← MaxOrderOrigin

 24 if SegmentIDsPerOriginatingNode[ OriginatingNode[ArcID] ] = nil

 25 then SegmentIDs[ StreamOrders[ArcID] ] ← SegmentIDs[
 StreamOrders[ArcID] ] + 1

 26 SegmentIDsPerOriginatingNode[ OriginatingNode[ArcID] ] ← SegmentIDs[
 StreamOrders[ArcID] ]

 27 Segments[ArcID] ← SegmentIDsPerOriginatingNode[ OriginatingNode[ArcID] ]

 28 return StreamOrders[ArcID]

 



 Is there anyone having experience with that, or has already tried to
 implement a vector-stream order approach in GRASS GIS. Specifically I am
 interested

 in Shreve stream order for simple (non-braided) stream networks and
 Strahler stream order for braided river networks (Gleyzer et al 2004).



 Best regards,

 Johannes

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Copy raster map from one unprojected location to another

2015-03-19 Thread Johannes Radinger
Thanks for your ideas.
I tried the pack-unpack approach. However, it seems that I can unpack and
import the raster map,
as there is no PROJ-file associated with that map. Here the error when
unpacking the *.pack that
was created in an unprojected location.

r.unpack --overwrite
input=/media/grassgis_data/tmp/tmp_Cele_raster_pack/my_raster.pack
output=my_raster
Traceback (most recent call last):
  File /usr/local/grass-7.0.1svn/scripts/r.unpack, line
146, in module
sys.exit(main())
  File /usr/local/grass-7.0.1svn/scripts/r.unpack, line
104, in main
proj=True):
  File
/usr/local/grass-7.0.1svn/etc/python/grass/script/core.py,
line 863, in compare_key_value_text_files
checkunits=units)
  File
/usr/local/grass-7.0.1svn/etc/python/grass/script/core.py,
line 790, in _text_to_key_value_dict
text = open(filename, r).readlines()
IOError: [Errno 2] No such file or directory: 'PROJ_INFO'

So maybe I need to go with the geotiff-approach. Is there any reason why
there is no PROJ for an
unprojected location? Wouldn't it also make sense to define a EPSG-code for
such special type of projection?

/Johannes

On Thu, Mar 19, 2015 at 11:30 AM, Benjamin Ducke bendu...@fastmail.fm
wrote:

 On Thu, Mar 19, 2015, at 11:25, Roy wrote:
  Hi,
 
 
 
  Il 19/03/2015 11:04, Johannes Radinger ha scritto:
Of course, I could export the map to geotiff and import them in the
   other location;
 
  IMHO this is the simpler way,

 Bad idea. You will lose all associated metadata: colour table,
 categories,
 editing history, maybe even the correct cell stats and no data code.

 GRASS raster maps are sets of files that are unfortunately
 a little spread out over several subdirectories in the mapset directory.

 To copy the entire collection of files that constitutes a GRASS raster
 map from one place to another, use:

 http://grass.osgeo.org/grass64/manuals/r.pack.html

 and the associated r.unpack.

 Best,

 Ben


 
  Roy.
  ___
  grass-user mailing list
  grass-user@lists.osgeo.org
  http://lists.osgeo.org/mailman/listinfo/grass-user
 ___
 grass-user mailing list
 grass-user@lists.osgeo.org
 http://lists.osgeo.org/mailman/listinfo/grass-user

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Copy raster map from one unprojected location to another

2015-03-19 Thread Johannes Radinger
Doesn't work either with the -o flag:

r.unpack --o
input=/media/grassgis_data/tmp/tmp_Cele_raster_pack/my_raster.pack
output=my_raster
Traceback (most recent call last):
  File /usr/local/grass-7.0.1svn/scripts/r.unpack, line 146, in module
sys.exit(main())
  File /usr/local/grass-7.0.1svn/scripts/r.unpack, line 104, in main
proj=True):
  File /usr/local/grass-7.0.1svn/etc/python/grass/script/core.py, line
863, in compare_key_value_text_files
checkunits=units)
  File /usr/local/grass-7.0.1svn/etc/python/grass/script/core.py, line
790, in _text_to_key_value_dict
text = open(filename, r).readlines()
IOError: [Errno 2] No such file or directory: 'PROJ_INFO'


On Thu, Mar 19, 2015 at 12:26 PM, Nikos Alexandris n...@nikosalexandris.net
wrote:

 On 19.03.2015 13:15, Johannes Radinger wrote:

  Thanks for your ideas.
 I tried the pack-unpack approach. However, it seems that I can unpack and
 import the raster map,
 as there is no PROJ-file associated with that map. Here the error when
 unpacking the *.pack that
 was created in an unprojected location.


  r.unpack --overwrite
 input=/media/grassgis_data/tmp/tmp_Cele_raster_pack/my_raster.pack
 output=my_raster


 Try to use the -o   Override projection check (use current location's
 projection) flag.

 Nikos

 [rest deleted]

 ___
 grass-user mailing list
 grass-user@lists.osgeo.org
 http://lists.osgeo.org/mailman/listinfo/grass-user

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Stream order for vector networks

2015-03-18 Thread Johannes Radinger
Hi,

I am interested in calculating stream order (Strahler, Shreve) for stream
networks based on GRASS vector networks (e.g. created by v.net) without
using the raster approach of the GRASS add-on r.stream.order. Therefore, I
came across following paper:

Gleyzer, A., Denisyuk, M., Rimmer, A. and Salingar, Y. (2004), A FAST
RECURSIVE GIS ALGORITHM FOR COMPUTING STRAHLER STREAM ORDER IN BRAIDED AND
NONBRAIDED NETWORKS. JAWRA Journal of the American Water Resources
Association, 40: 937–946. doi: 10./j.1752-1688.2004.tb01057.x

which also provides pseudocode to calculate Strahler order for non-braided
and braided vector networks defined by arcs and nodes. Thus, I am wondering
if it is possible to get the information from which node a specific arc
starts (from_node) and where it drains to (to_node). BTW, I seems some kind
of this information is already provided with the NC-dataset (FNODE_ and
TNODE column) in the streams attribute table, however I don't know where
this information comes from. Of course to extract this information the
network needs to be directed to define what is the upstream and what the
downstream node of each arc. However, having this information available, it
should not be to complicated to apply the algorithm proposed by Gleyzer et
al 2004. Thus I am interested to get 'from_node' and 'to_node' information
e.g. stored in the attribute table of v.net output?

Here the pseudocode provided by Gleyzer et al 2004:


MakeDictionaries(Network)
1 for each Arc ∈ Network
2 do FromNodesPerArc[Arc's ID] ← Arc’s FromNodeID
3 InflowingArcsPerNode[Arc's ToNodeID] ← InflowingArcsPerNode[Arc’s
ToNodeID] ∪ Arc’s ID
4 OriginatingNode[Arc's ID] ← Arc’s FromNodeID

StreamOrdering(ArcID)
1 Visited[ArcID] ← true
2 if | InflowingArcsPerNode[ FromNodesPerArc[ArcID] ] | = 0
3 then StreamOrders[ArcID] ← 1
4 else
5 for each Arc ∈ InflowingArcsPerNode[ FromNodesPerArc[ArcID] ]
6 do if not Visited[Arc]
7 then UpstreamOrders[Arc] ← (StreamOrdering(Arc), OriginatingNode[Arc])
8 else UpstreamOrders[Arc] ← (StreamOrders[Arc], OriginatingNode[Arc])
9 MaxOrder ← 0
10 MaxOrderCount ← 0
11 for each (Order, Origin) ∈ UpstreamOrders
12 do if Order  MaxOrder
13 then MaxOrder ← Order
14 MaxOrderCount ← 1
15 MaxOrderOrigin ← Origin
16 else if Order = MaxOrder
17 then if Origin ≠ MaxOrderOrigin
18 then MaxOrderCount ← MaxOrderCount + 1
19 if MaxOrderCount  1
20 then StreamOrders[ArcID] ← MaxOrder + 1
21 OriginatingNode[ArcID] ← FromNodesPerArc[ArcID]
22 else StreamOrders[ArcID] ← MaxOrder
23 OriginatingNode[ArcID] ← MaxOrderOrigin
24 if SegmentIDsPerOriginatingNode[ OriginatingNode[ArcID] ] = nil
25 then SegmentIDs[ StreamOrders[ArcID] ] ← SegmentIDs[ StreamOrders[ArcID]
] + 1
26 SegmentIDsPerOriginatingNode[ OriginatingNode[ArcID] ] ← SegmentIDs[
StreamOrders[ArcID] ]
27 Segments[ArcID] ← SegmentIDsPerOriginatingNode[ OriginatingNode[ArcID] ]
28 return StreamOrders[ArcID]


Is there anyone having experience with that, or has already tried to
implement a vector-stream order approach in GRASS GIS. Specifically I am
interested
in Shreve stream order for simple (non-braided) stream networks and
Strahler stream order for braided river networks (Gleyzer et al 2004).

Best regards,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] sampling points in a grid

2015-02-18 Thread Johannes Radinger
Maybe you can use a database selection approach for example in sqlite.
First update your points with your grid cell ID and then randomly select
7 rows from the selection where the grid cell ID == 42.

In sqlite there exists the function random(). Here you can find
an approach to select random rows from a sqlite table:
http://stackoverflow.com/questions/2279706/select-random-row-from-an-sqlite-table

I have not tested that approach, but maybe worth a try...

/Johannes

On Wed, Feb 18, 2015 at 11:55 AM, Margherita Di Leo direg...@gmail.com
wrote:



 On Tue, Feb 17, 2015 at 5:30 PM, Markus Neteler nete...@osgeo.org wrote:

 On Tue, Feb 17, 2015 at 4:49 PM, Margherita Di Leo direg...@gmail.com
 wrote:
  Hi,
 
  I have a vector of points. According to a regular grid, I need to
 (randomly)
  sample a certain number n of these points in each cell of the grid. Such
  number n is given in a column of the attribute table of the grid. How
 would
  you do this?

 With a small script this should be doable:
 - generate the grid as polygons (v.mkgrid)
 - loop over each polygon
  - fetch category number or v.extract
  - fetch corresponding number of points to be created from DB
  - v.random, using the current grid polygon for restricted creation
 feature

 http://grass.osgeo.org/grass70/manuals/v.random.html#restriction-to-vector-areas


 Here lies the problem. I need to randomly pick points that are already
 existing. Example: for grid cell that has ID 42, I have 50 points, and have
 to randomly pick 7 out of these 50.


  - save point map
 - patch all point maps into single resulting map
 - remove single point maps

 Something like this...

 Markus




 --
 Best regards,

 Dr. Margherita DI LEO
 Scientific / technical project officer

 European Commission - DG JRC
 Institute for Environment and Sustainability (IES)
 Via Fermi, 2749
 I-21027 Ispra (VA) - Italy - TP 261

 Tel. +39 0332 78 3600
 margherita.di-...@jrc.ec.europa.eu

 Disclaimer: The views expressed are purely those of the writer and may not
 in any circumstance be regarded as stating an official position of the
 European Commission.

 ___
 grass-user mailing list
 grass-user@lists.osgeo.org
 http://lists.osgeo.org/mailman/listinfo/grass-user

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Choosing mapset from command line

2015-02-13 Thread Johannes Radinger
Hi,

I tried to access my machine with GRASS GIS via SSH and wanted to get the
text based welcome screen of GRASS 7 (RC in Ubuntu) to select a mapset.
Usually I am starting the GUI selection menu so I am not familiar with the
text based startup. I thought that the command grass70 -gtext should
launch the text-based menu however here I get:

radinger@grassgis:~$ grass70 -gtext

Cleaning up temporary files...

Starting GRASS GIS...

WARNING: It appears that the X Windows system is not active.

A graphical based user interface is not supported.

Switching to text based interface mode.


Hit RETURN to continue


With -gtext I already indicated that I don't want to start a graphical
based user interface so I am wondering about the warning. With hitting
return I automatically start with the last mapset I used (like the command
grass70 -text) and don't get to a text based menu to e.g. create and start
into a new mapset of a given location.

What I am doing wrong?


Of course I can create and switch to a new the mapset with g.mapset as an
alternative.


/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Choosing mapset from command line

2015-02-13 Thread Johannes Radinger
On Fri, Feb 13, 2015 at 4:58 PM, Johannes Radinger 
johannesradin...@gmail.com wrote:

 Maybe I expressed myself unclear or mixed things up:

 I want to get the text based menu to select e.g. a location and a mapset,
 which I think is called text-based location wizzard [1]. I mixed that up
 with the welcome screen. And of course I don't need the welcome screen.

 I was just confused because the wiki [2] says:
 *grass70 -text:* Start GRASS using the text-based user interface. The
 user will be *prompted* to choose the appropriate location and mapset.


sorry thats the info for -gtext: *grass70 -gtext:* Start GRASS using the
text-based user interface



 However I get automatically started into my last session which is actually
 the definition of:


 *grass70 -text *Start GRASS using the text-based user interface.
 Appropriate location and mapset must be set by environmental variables (see
 examples bellow) otherwise taken from the last GRASS session.

 So is the location wizzard not exisiting anymore in GRASS70? If not, is
 g.mapset from within an existing mapset the recommended way to go to create
 and start into a new mapset? Another option would be only to create the
 mapset with g.mapset from inside an existing mapset and than to restart
 grass with grass70 newmapset.

 [1]
 http://grasswiki.osgeo.org/wiki/GRASS_Location_Wizard#Text_based_location_wizard
 [2] http://grass.osgeo.org/grass70/manuals/grass7.html


 On Fri, Feb 13, 2015 at 3:43 PM, Moritz Lennert 
 mlenn...@club.worldonline.be wrote:

 On 13/02/15 15:33, Johannes Radinger wrote:

 Hi,

 I tried to access my machine with GRASS GIS via SSH and wanted to get
 the text based welcome screen of GRASS 7 (RC in Ubuntu) to select a
 mapset. Usually I am starting the GUI selection menu so I am not
 familiar with the text based startup. I thought that the command
 grass70 -gtext should launch the text-based menu however here I get:

 radinger@grassgis:~$ grass70 -gtext

 Cleaning up temporary files...

 Starting GRASS GIS...

 WARNING: It appears that the X Windows system is not active.

 A graphical based user interface is not supported.

 Switching to text based interface mode.


 Hit RETURN to continue


 With -gtext I already indicated that I don't want to start a graphical
 based user interface so I am wondering about the warning. With hitting
 return I automatically start with the last mapset I used (like the
 command grass70 -text) and don't get to a text based menu to e.g.
 create and start into a new mapset of a given location.

 What I am doing wrong?


 '-gtext' asks to show the GUI welcome screen
 '-text' skips even that so that's probably what you need.

 Moritz




___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] Choosing mapset from command line

2015-02-13 Thread Johannes Radinger
Maybe I expressed myself unclear or mixed things up:

I want to get the text based menu to select e.g. a location and a mapset,
which I think is called text-based location wizzard [1]. I mixed that up
with the welcome screen. And of course I don't need the welcome screen.

I was just confused because the wiki [2] says:
*grass70 -text:* Start GRASS using the text-based user interface. The user
will be *prompted* to choose the appropriate location and mapset.
However I get automatically started into my last session which is actually
the definition of:
*grass70 -text *Start GRASS using the text-based user interface.
Appropriate location and mapset must be set by environmental variables (see
examples bellow) otherwise taken from the last GRASS session.

So is the location wizzard not exisiting anymore in GRASS70? If not, is
g.mapset from within an existing mapset the recommended way to go to create
and start into a new mapset? Another option would be only to create the
mapset with g.mapset from inside an existing mapset and than to restart
grass with grass70 newmapset.

[1]
http://grasswiki.osgeo.org/wiki/GRASS_Location_Wizard#Text_based_location_wizard
[2] http://grass.osgeo.org/grass70/manuals/grass7.html


On Fri, Feb 13, 2015 at 3:43 PM, Moritz Lennert 
mlenn...@club.worldonline.be wrote:

 On 13/02/15 15:33, Johannes Radinger wrote:

 Hi,

 I tried to access my machine with GRASS GIS via SSH and wanted to get
 the text based welcome screen of GRASS 7 (RC in Ubuntu) to select a
 mapset. Usually I am starting the GUI selection menu so I am not
 familiar with the text based startup. I thought that the command
 grass70 -gtext should launch the text-based menu however here I get:

 radinger@grassgis:~$ grass70 -gtext

 Cleaning up temporary files...

 Starting GRASS GIS...

 WARNING: It appears that the X Windows system is not active.

 A graphical based user interface is not supported.

 Switching to text based interface mode.


 Hit RETURN to continue


 With -gtext I already indicated that I don't want to start a graphical
 based user interface so I am wondering about the warning. With hitting
 return I automatically start with the last mapset I used (like the
 command grass70 -text) and don't get to a text based menu to e.g.
 create and start into a new mapset of a given location.

 What I am doing wrong?


 '-gtext' asks to show the GUI welcome screen
 '-text' skips even that so that's probably what you need.

 Moritz



___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Transform single raster cells to 3d bars

2015-02-06 Thread Johannes Radinger
Hi,

I have a 2D raster map with scattered single cells (e.g. created by
r.random.cells)
with lots of empty space (NaNs) between the cells. Each cell has a
value that represents actually the height of this cell. Now I'd like to
transform
the 2D raster cells into a 3D raster where each cell should grow in height
so that the cells become independent 3D bars with a their specific height.
The bottom of the bars (i.e. the base elevation level where the grow out)
should
be defined by another 2D raster map with base elevation levels.

The bars should than be exported to vtk to represent them in paraview.

I guess this can be done with r.to.rast3elev, but what is the input and what
the elevation map to use? What should be use for the upper and lower value?
And what is then the appropriate command to export  the bars
with r3.out.vtk?


Here what I am doing so far to create the random cells

r.random.cells output=random_cells distance=40

r.null map=random_cells setnull=0


Any suggestions?


/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] temporary files of grass.script.array

2015-02-04 Thread Johannes Radinger
Hi,

I just tested in trunk, and it works: The temporary file that is create
when assigning a raster to an python array is deleted once e.g.
 the python session is closed again. So this would be a nice
thing also to have implemented in the next RC of GRASS7.

/johannes

On Tue, Feb 3, 2015 at 11:35 PM, Markus Neteler nete...@osgeo.org wrote:

 On Tue, Feb 3, 2015 at 4:23 PM, Glynn Clements gl...@gclements.plus.com
 wrote:
  The temporary files were supposed to be deleted when the array object
  was destroyed. However, this relied upon the ._close() method, which
  is no longer used. r64426 (trunk) should fix this.

 Please let me know if that should be backported then.

 Markus

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] r.in.onearth in GRASS 7

2015-02-04 Thread Johannes Radinger
Hi,

I just try to follow some examples r3.out.vtk
to get some nice 3D plots in paraview.

http://grass.osgeo.org/grass70/manuals/r3.out.vtk.html

In the manual there is an example that
uses the r.in.onearth add-on
(http://trac.osgeo.org/grass/browser/grass-addons/grass6/raster/r.in.onearth
)
that is only available in GRASS 6.
Is there a similar simple tool available in GRASS 7?

Or is there any other easy way to drap a satelitte image over
a e.g. DEM for viewing in paraview?

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] temporary files of grass.script.array

2015-02-03 Thread Johannes Radinger
Hi,

In GRASS 7 (RC1) python it is possible to export a raster map to an array
which can be used e.g. by Numpy.

For example:

import grass.script.array as garray
x1 = garray.array()
x1.read(my_raster)

It seems that this creates a temporary file in the .tmp folder of the
respective
location. I am not sure if this has been existing before, but these
temporary
files are not deleted after e.g. x1 is not used anymore (i.e. the
calculation of
the script finished). So do I need to manually (in the python script) close
that x1 so
that the associated temporary file gets deleted (and if yes how?)? Otherwise
when using this procedure to do numpy calculations in loops the tmp-files
pile up quite fast.

Of course I can empty the tmp folder after every loop, but I am not sure
if this is really the right way.

Any suggestions?

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] how to change all cats in a vector line

2014-12-01 Thread Johannes Radinger
Hi Michael,

I have not read the entire thread, but maybe follwowing might work
if you want to extract this upstream-downstream watercourse in raster
format:

1) Use r.cost with your outlet as starting point and your cell resolution as
cost surface. So the output is a pseudo-elevation
map with increasing values in upstream direction

2) use r.drain to extract the main course in downstream direction based on a
defined upstream end and the pseudo-elevation map created by r.cost.

3) Than you can use r.to.vect to convert your raster cells to vectors and
extract your distance per point (cummulative cost) and your elevation.

Maybe this might help in on or the other way.

best,
Johannes

On Mon, Dec 1, 2014 at 7:21 PM, Michael Barton michael.bar...@asu.edu
wrote:

 v.build.polylines is indeed part of what is needed. But it is not enough
 to be able to use v.to.points in the data I have—and this is not a weirdly
 meandering stream, just a Mediterranean barranco.

 A major part of the problem is in the raster to vector conversion. I’ve
 been using r.to.vect. But it breaks the resulting line into multiple
 segments and assigns different cat numbers to the segments. If the segments
 all connect cleanly, v.build.polylines may combine the segments. But I’ve
 tried using cat=first and it continues to assign multiple cats to the
 segments. To get a single line, what I have to do currently is:

 v.clean
 v.build.polylines with cats=no
 v.category with option=del ## may not be necessary but I’ve had enough
 inconsistent results to want to make sure on this
 v.category with option=add cat=1 step=0 ## this gets me a vector with
 cat=1 for all parts.

 repeat the above a 2nd time

 This gives a sufficiently clean single-object line with a single cat
 number that I can run v.to.points on.

 Maybe r.stream.order does a cleaner job in raster to vector conversion.
 It’s worth a try. Unfortunately, I can’t use the stream order info to
 extract a stream since I want the stream from its headwaters (stream order
 = 1) to its outlet (stream order  1).

 Ideally, there would be a module where I could enter the coordinates of
 the upstream end of a stream, it would extract the entire watercourse as a
 single object to its downstream end, edge of map, or a defined outlet stop
 point. Then I could run v.to.points and v.what.rast to get profile data. OR
 in a a stream profile module, it would run these for me. I could do it in
 Python if I could extract a clean vector (one object, one cat, no loops or
 self-crossing) of a single stream from a stream network created with
 r.watershed or r.stream.extract. AFAICT, there is nothing available in
 GRASS that can do this. Kind of surprising.

 Michael


 
 C. Michael Barton
 Director, Center for Social Dynamics  Complexity
 Professor of Anthropology, School of Human Evolution  Social Change
 Head, Graduate Faculty in Complex Adaptive Systems Science
 Arizona State University

 voice:  480-965-6262 (SHESC), 480-965-8130/727-9746 (CSDC)
 fax: 480-965-7671 (SHESC),  480-727-0709 (CSDC)
 www: http://www.public.asu.edu/~cmbarton, http://csdc.asu.edu



 On Dec 1, 2014, at 9:31 AM, Moritz Lennert mlenn...@club.worldonline.be
 wrote:

  On 30/11/14 06:28, Michael Barton wrote:
 
  
  C. Michael Barton
  Director, Center for Social Dynamics  Complexity
  Professor of Anthropology, School of Human Evolution  Social Change
  Head, Graduate Faculty in Complex Adaptive Systems Science
  Arizona State University
 
  voice: 480-965-6262 (SHESC), 480-965-8130/727-9746 (CSDC)
  fax: 480-965-7671 (SHESC),  480-727-0709 (CSDC)
  www: http://www.public.asu.edu/~cmbarton, http://csdc.asu.edu
 
 
 
 
  Thanks César,
 
  This is pretty much what I’m doing. However, it is not that simple, it
  turns out.
 
  First, for a stream profile, you want to isolate a single watercourse.
  r.stream.extract creates a stream network. Getting stream order does not
  help with profiling, however.
 
  I looked at r.stream.distance. But this calculates distance from each
  stream junction. I want it for the entire course of the selected
  stream—from headwaters to outlet. That is the only way to graph a stream
  profile for the entire stream.
 
  Once I have a set of points with the distance from the beginning or the
  end of the line stored as an attribute for each point, the rest is easy.
  It is getting to this step that is hard.
 
  If a line is composed of multiple segments—as is the case for any line
  created with r.stream.extract or r.watershed, and also with r.drain
  surprisingly—there is nothing that will ignore all of these segments and
  just give me the distance along the line from one end to the other. This
  is the case with v.to.points
 
  As MarkusN said: you have to use v.build.polylines with cats=first, but
 then it is easy to get total length with v.to.db or v.to.points if that is
 what you need (or many other options (get end node coordinates and use
 v.distance 

Re: [GRASS-user] A question before I embark on a programming exercise

2014-11-12 Thread Johannes Radinger
Hi Tom,

what comes into my mind is following approach (not tested):

1) calculate flow direction using r.watershed for each pixel of your map
2) transforming raster to vector points (r.to.vect), adding some columns
for your results (X, Y, ID). This step should provide you already with
unique IDs for your cells/points.
3) With v.to.db you can update the X  and Y coordinate of your cell center
and together with the raster resolution you can get the lower left corner
as desired.
3) By knowing the flow direction for each point you should be able to
calculate the X and Y coordinate of the downstream cell (adding or
substracting your cell dimensions from the cell center, correct for
diagonal flow). When you know X and Y of your downstream cell then you
should be able to also get the ID of the downstream cell (most probably a
SQL query I guess)

Maybe there are easier ways

hope it works,
Johannes

On Wed, Nov 12, 2014 at 5:10 AM, Thomas Adams tea...@gmail.com wrote:

 All:

 I need to generate an ascii text file from a flow direction grid that
 consists of (among a couple other things that don't really matter at this
 point) for each pixel:

 (1) a unique integer identifier (1 -- N) for the pixel
 (2) the integer identifier of the downstream pixel (assuming there is ONLY
 one)
 (3) the x,y location of the pixel (presumably, the lower left corner of
 the pixel)

 has anyone already done something like this?

 It's needed (along with some header information) as an input file for a
 gridded distributed hydrologic model, identifying the flow connectivity
 from pixel to pixel for streamflow routing purposes. I have already
 formulated conceptually how to do this, but if someone has already done
 such a thing, why reinvent the wheel??

 Any help is appreciated.

 Regards,
 Tom


 ___
 grass-user mailing list
 grass-user@lists.osgeo.org
 http://lists.osgeo.org/mailman/listinfo/grass-user

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Overlapping points in ps.map

2014-10-21 Thread Johannes Radinger
Hi,

I am trying to use ps.map to create a map of many vector points which should be 
colored based on the values (0-1) in a specific column in the attribute table. 
If I understood the manual correctly, I need to create an rgbcolumn indicating 
the color for each point. Is there an example explaining this on detail (eg how 
to create such a column based on a specific color scale)?

Moreover, as many of my points are overlapping I am interested how this is 
handled in ps.map? Is there any method to plot the darkest point on top? Like 
the darken feature in qgis?

Any suggestion is welcome.

Best,
Johannes___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Upper case letters for defining columns via v.db.addtable

2014-10-13 Thread Johannes Radinger
Hi,

when I want to specify a column in v.db.addtable, all upper case letters in
the column name get converted to lower case. I'm working on GRASS 71
(trunk). Is this an intended behaviour?

v.db.addtable map=myMap table=test layer=3 columns=My_COL INT


cheers,

/Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Categories in new layers of vector maps

2014-10-13 Thread Johannes Radinger
Hi,

I am slightly confused about creating new layers for vector maps and adding
a category value.
I want to add a layer and subsequently I want to use v.what.rast to query a
raster map and update
that layer. Hence, I need categories in the new layer. However, it seems
that adding a new table and connect it to a specific layer (v.db.addtable)
does not create/update the category column (in GRASS 71, trunk):

Here an example for the NC-dataset:

v.db.addtable map=firestations@PERMANENT layer=2

DB settings already defined, nothing to do
Reading features...
Updating database...
0 categories read from vector map (layer 2)
0 records updated/inserted (layer 2)


And when I try to manually update the category I am also not successful:
v.to.db map=firestations@PERMANENT layer=2 option=cat columns=cat

Reading features...
Updating database...
0 categories read from vector map (layer 2)
0 records updated/inserted (layer 2)

Maybe I just got something wrong?


The map firestations@PERMANENT has categories (stored in layer 1)
v.category input=firestations@PERMANENT option=report

Layer/table: 1/firestations
type   countminmax
point 71  1 71
line   0  0  0
boundary   0  0  0
centroid   0  0  0
area   0  0  0
face   0  0  0
kernel 0  0  0
all   71  1 71
v.category complete. 0 features modified.



Any suggestions?


best,
Johannes
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] my first script

2014-10-12 Thread Johannes Radinger
Hi, 

It seems to raise an error related to the flag(s) you are specifying in oneof 
your grass modules you are using in your script. Flags in grass python scripts 
should be used like: flags=a or for two flags: flags=ab. Note this does not 
need a - like in the bash command.

Best,
/johannes

 Original message 
From: Giuliano Urgeghe giulian...@gmail.com 
Date: 13/10/2014  02:44  (GMT+01:00) 
To: grass-user@lists.osgeo.org 
Subject: [GRASS-user] my first script 
 
Hi all,
I have my first script in python for grass 
but when I launch it on windows vista for grass 6.4.2 in quantum gis
I have this output:
...
File C:\PROGRA~1\QUANTU~3\apps\grass\grass-6.4.2\etc\pyth
on\grass\script\core.py, line 189, in run_command
    ps = start_command(*args, **kwargs)
  File C:\PROGRA~1\QUANTU~3\apps\grass\grass-6.4.2\etc\pyth
on\grass\script\core.py, line 167, in start_command
    args = make_command(prog, flags, overwrite, quiet,
verbose, **options)
  File C:\PROGRA~1\QUANTU~3\apps\grass\grass-6.4.2\etc\pyth
on\grass\script\core.py, line 124, in make_command
    raise ScriptError('-' is not a valid flag)
grass.script.core.ScriptError: '-' is not a valid flag

someone know this problem?

thanks in advance
Giuliano

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

[GRASS-user] Null-value in grass python array

2014-10-10 Thread Johannes Radinger
Hi,

I am using the GRASS-numpy functionality [1] to read a python numpy.array
from a GRASS raster map (GRASS7):

import grass.script.array as garray

a = garray.array()

a.read(map)



Here [2] it says that also the null value can be specified. So what is the
actual null value in GRASS so that the raster map is correctly saved as
numpy.array with numpy.nan where appropriate?

And a second related question: Does grass.script.array.read() respect an
existing
MASK so that areas that are masked will get numpy.nan? Of course this can
also be done within numpy like map[numpy.isnan(mask)] but this again needs
to properly define the nan value when reading the raster mask into an
numpy.array.


Best regards,
Johannes

[1]
http://grasswiki.osgeo.org/wiki/GRASS_Python_Scripting_Library#Interfacing_with_NumPy
[2]
http://grass.osgeo.org/programming6/classpython_1_1array_1_1array.html#a493c541122108eca3fc4b441fe1a5487
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] raster exchange between GRASS and R with nodata

2014-09-26 Thread Johannes Radinger
I think I discovered the problem. This was a mistake on my side:
The summary() function in R applied on a raster does subsample the
data in case of very large raster maps. However, this seems to
cause that summary information about NA gets lost and no NA are
mentioned in the summary output although there are NA included
in the raster map.

Here an example output of a map with lots of NA's which are ignored
by summary(rastermap):

 summary(streams[[1]])
Elbe_distance_from_mouth
Min.  0.
1st Qu. 488.7067
Median  789.8560
3rd Qu.1025.4501
Max.   1876.6429
NA's  0.
Warning message:
In .local(object, ...) :
  summary is an estimate based on a sample of 1e+05 cells (10.7% of all
cells)

 length(bioclim[[1]][is.na(bioclim[[1]])])
[1] 12688


So there was actually no problem, just some confusion concerning the
information about
NA in the summary.

Best regards,
Johannes


On Thu, Sep 4, 2014 at 6:26 PM, Robert J. Hijmans r.hijm...@gmail.com
wrote:

 Johannes,

 raster uses GDAL to access GTiff data. So it seems that there is no,
 or at least GDAL does not see, a no data value in the file.
 I would suggest running GDALinfo to troubleshoot. Can you email me a
 small example file?

 By the way, the spgrass6 route has a shortcoming: it may not work with
 very large files as it attempts to load all values into RAM.

 Robert

 On Thu, Sep 4, 2014 at 6:47 AM, Johannes Radinger
 johannesradin...@gmail.com wrote:
 
 
 
  On Thu, Sep 4, 2014 at 3:38 PM, Nuno Sá nunocesard...@gmail.com wrote:
 
  Hello!
 
  Did you try this one?
 
  r.out.gdal etc nodata='NA'
 
 
  As mentioned in the manual of r.out.gdal, the no data parameter takes
 only
  float values
  and no strings like 'NA'. Without stating as specific value in GRASS,
 this
  nodata-value is
  automatically set to e.g. 65535 for DCELL rasters if I remember correctly
  and to 255
  for BYTE rasters. However, this seems not to be recognized when imported
  into R with
  the package 'raster'.
 
  /Johannes
 
 
 
 
 
  On 4 September 2014 14:27, Johannes Radinger 
 johannesradin...@gmail.com
  wrote:
 
  Hi all,
 
  of course it is possible to load the raster maps directly via spgrass6.
  However, we use this work
  flow also to exchange some of the maps between different users (e.g.
 via
  email) and to permanently
  store single files (geotiffs that contain the proj information within
 the
  file). So, I agree that using spgrass6 would be much more efficient,
 but
  I'll stick to exporting to geotiffs and so I need to solve the issues
 with
  NA's.
 
  /Johannes
 
 
  On Thu, Sep 4, 2014 at 2:31 PM, Thomas Adams tea...@gmail.com wrote:
 
  Johannes,
 
  If you want to read your file into R, there is no need to export your
  map from GRASS to do this. Simply install and use the R contributed
 package
  'spgrass6' (spgrass6 has R dependencies that need to be installed
 first); it
  works wonderfully:
 
  Within GRASS, at the GRASS terminal prompt...
 
   library(spgrass6)
  Loading required package: sp
  Loading required package: XML
  GRASS GIS interface loaded with GRASS version: GRASS 7.0.0beta3 (2014)
  and location: ohrfc_mpe
   dat-readRAST6(xmrg0101200306z)
   image(dat)
 
  This is far more efficient.
 
  Tom
 
 
  On Thu, Sep 4, 2014 at 5:32 AM, Johannes Radinger
  johannesradin...@gmail.com wrote:
 
  Hi,
 
  I want to export a raster map (FCELL) from GRASS70 to the geotiff
  format using r.out.gdal and to import it later on in R. The map
 contains
  many no data values.
 
  Here some details about the raster:
  Type of Map:  raster   Number of Categories: 0
  Data Type:FCELL
  Rows: 750
  Columns:  750
  Total Cells:  562500
  total null and non-null cells: 15105636
  total null cells: 15105047
 
  So when I export the map, r.out.gdal reports: Input raster map
  contains cells with NULL-value (no-data). The value -nan will be
 used to
  represent no-data values in the input map. You can specify a nodata
 value
  with the nodata option.
 
  When I subsequently try to import the geotiff into R (using the
 package
  'Raster') the nodata values are not recognised as NA's:
 
  a - raster(*.tif)
  summary(a)
  Min. 0.5294496
  1st Qu.  0.7171210
  Median   0.7871540
  3rd Qu.  1.1581826
  Max. 1.5494517
  NA's 0.000
 
  So I am wondering if I need to set any specific parameter during the
  export (r.out.gdal) or import (raster()).
 
  As I am not only exporting FCELL (Float32) raster but also multiple
  (N=500) other rasters to R I would be interested in a solution also
 for
  DCELL (Float64). Of course I can export all of as Float64 as the
 file size
  should not be a problem.
 
  Any suggestions or experiences of handling NA's during raster
 exchange
  between GRASS and R?
 
  /Johannes
 
  ___
  grass-user mailing list
  grass-user@lists.osgeo.org
  http

  1   2   3   4   >