[GRASS-user] v.in.csv available?
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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'
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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]
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]
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 SherpaDate: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]
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 ShepardDate: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
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
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 AdamsDate: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
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
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
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
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
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,
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
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
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
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?
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?
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
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 DELAYDate: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
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
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
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
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
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 Sinhawrote: > > > 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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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
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