On Wed, 2 Nov 2011, toke wrote:
Hi Roger and everybody else
I finally found time to translate the problem into spearfish data.
To recapitulate I would like to use the v.distance command to transfer
information from polygons to point data.
The problem is that the v.distance command do not update the column
specified to receive the calculation done by the v.distance command.
I do not see the same problem, and the command works for me - can others
please try to reproduce it?
Roger
My output in R:
library(spgrass6)
Loading required package: sp
Loading required package: rgdal
Geospatial Data Abstraction Library extensions to R successfully loaded
Loaded GDAL runtime: GDAL 1.8.1, released 2011/07/09
Path to GDAL shared files: /usr/local/share/gdal
Loaded PROJ.4 runtime: Rel. 4.7.1, 23 September 2009, [PJ_VERSION: 470]
Path to PROJ.4 shared files: (autodetected)
Loading required package: XML
GRASS GIS interface loaded with GRASS version: GRASS 6.4.2svn47586 (2011)
and location: spearfish60
execGRASS("g.copy", vect="'fields@PERMANENT',fields1")
Copy vector <fields@PERMANENT> to current mapset as <fields1>
execGRASS("g.copy", vect="'archsites@PERMANENT',archsites1")
Copy vector <archsites@PERMANENT> to current mapset as <archsites1>
execGRASS("v.db.addcol", map="archsites1",
+ columns="\"fields double precision\"")
fields1<-readVECT6("fields1")
WARNING: The map contains islands. To preserve them in the output map,
use
the -c flag
Exporting 65 areas (may take some time)...
100%
v.out.ogr complete. 67 features written to <fields1> (ESRI_Shapefile).
OGR data source with driver: ESRI Shapefile
Source: "/home/rsb/topics/grassdata/spearfish60/rsb/.tmp/reclus2", layer:
"fields1"
with 67 features and 2 fields
Feature type: wkbPolygon with 2 dimensions
archsites1<-readVECT6("archsites1")
Exporting 25 geometries...
100%
v.out.ogr complete. 25 features written to <archsite> (ESRI_Shapefile).
OGR data source with driver: ESRI Shapefile
Source: "/home/rsb/topics/grassdata/spearfish60/rsb/.tmp/reclus2", layer:
"archsite"
with 25 features and 3 fields
Feature type: wkbPoint with 2 dimensions
summary(archsites1)
Object of class SpatialPointsDataFrame
Coordinates:
min max
coords.x1 589860 608355
coords.x2 4914479 4926490
Is projected: TRUE
proj4string :
[+proj=utm +zone=13 +datum=NAD27 +units=m +no_defs +ellps=clrk66
+nadgrids=@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat]
Number of points: 25
Data attributes:
cat str1 fields
Min. : 1 No Name :12 Min. : NA
1st Qu.: 7 Bob Miller : 1 1st Qu.: NA
Median :13 Boulder Creek Cabin: 1 Median : NA
Mean :13 Canyon Station : 1 Mean :NaN
3rd Qu.:19 Cole Creek Mine : 1 3rd Qu.: NA
Max. :25 Elkhorn Peak : 1 Max. : NA
(Other) : 8 NA's : 25
set.echoCmdOption(TRUE)
[1] FALSE
execGRASS("v.distance", flags="p", from="archsites1", to="fields1",
from_type="point", to_type="point,line,area", dmax=as.integer(1),
upload="to_attr", column="fields", to_column="cat")
GRASS command: v.distance -p from=archsites1 to=fields1 from_type=point
to_type=point,line,area dmax=1 upload=to_attr column=fields to_column=cat
100%
100%
100%
from_cat|fields
1|63
2|63
3|63
4|63
5|null
6|null
7|null
8|null
9|25
10|63
11|null
12|null
13|63
14|63
15|null
16|63
17|null
18|63
19|63
20|63
21|null
22|63
23|63
24|null
25|63
v.distance complete.
set.echoCmdOption(FALSE)
[1] TRUE
archsites1<-readVECT6("archsites1")
Exporting 25 geometries...
100%
v.out.ogr complete. 25 features written to <archsite> (ESRI_Shapefile).
OGR data source with driver: ESRI Shapefile
Source: "/home/rsb/topics/grassdata/spearfish60/rsb/.tmp/reclus2", layer:
"archsite"
with 25 features and 3 fields
Feature type: wkbPoint with 2 dimensions
summary(archsites1)
Object of class SpatialPointsDataFrame
Coordinates:
min max
coords.x1 589860 608355
coords.x2 4914479 4926490
Is projected: TRUE
proj4string :
[+proj=utm +zone=13 +datum=NAD27 +units=m +no_defs +ellps=clrk66
+nadgrids=@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat]
Number of points: 25
Data attributes:
cat str1 fields
Min. : 1 No Name :12 Min. :25.00
1st Qu.: 7 Bob Miller : 1 1st Qu.:63.00
Median :13 Boulder Creek Cabin: 1 Median :63.00
Mean :13 Canyon Station : 1 Mean :60.47
3rd Qu.:19 Cole Creek Mine : 1 3rd Qu.:63.00
Max. :25 Elkhorn Peak : 1 Max. :63.00
(Other) : 8 NA's :10.00
table(archsites1$fields)
25 63
1 14
## alternatively:
execGRASS("v.distance", from="archsites1", to="fields1",
from_type="point", to_type="point,line,area", dmax=as.integer(1),
upload="to_attr", column="fields", to_column="cat")
GRASS command: v.distance from=archsites1 to=fields1 from_type=point
to_type=point,line,area dmax=1 upload=to_attr column=fields to_column=cat
100%
100%
100%
25 categories read from the map
25 categories exist in the table
25 categories read from the map exist in the table
25 records updated
v.distance complete.
## also OK (overwriting field values)
sessionInfo()
R version 2.14.0 (2011-10-31)
Platform: x86_64-unknown-linux-gnu (64-bit)
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=C LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] spgrass6_0.7-4 XML_3.4-3 rgdal_0.7-1 sp_0.9-91
loaded via a namespace (and not attached):
[1] grid_2.14.0 lattice_0.20-0
Here is my script
grass
g.gisenv
g.mapset mapset=user1
g.copy vect='fields@PERMANENT',fields1 # adding polygon to mapset
g.copy vect='archsites@PERMANENT',archsites1 # adding points to mapset
v.db.addcol map=archsites1 columns="fields double precision" # adding
column to points data
v.distance from=archsites1 to=fields1 from_type=point
to_type=point,line,area dmax=1 upload=to_attr column=fields to_column=cat #
adding cat information from polygon to point data
### v.distance works in GRASS ###
### trying to do the same thing in R ###
v.db.dropcol archsites1 column=fields ## removing column
v.db.addcol map=archsites1 columns="fields double precision" ## adding a
clean column
R # starting R within grass
library(spgrass6) # calling library
gmeta6()
fields1<-readVECT6("fields1") # reading polygon into R
archsites1<-readVECT6("archsites1") # reading points into R
archsites1 # looking at data
# Trying to do the same v.distance calculation in R.
### The specified column is not updated ###
execGRASS("v.distance", from="archsites1",
to="fields1", from_type="point", to_type="point,line,area",
dmax=as.integer(1),
upload="to_attr", column="fields", to_column="cat")
### However if the flags “p” is add to the v.distance calculation it is
possible to see that v.distance
### actually do the calculations.
execGRASS("v.distance", flags="p", from="archsites1",
to="fields1", from_type="point", to_type="point,line,area",
dmax=as.integer(1),
upload="to_attr", column="fields", to_column="cat")
### script finished ###
There seems to be a bug in spgrass6.
Cheers Toke
--
View this message in context:
http://osgeo-org.1803224.n2.nabble.com/Problem-with-v-distance-in-spgrass6-package-in-R-tp6939286p6955833.html
Sent from the Grass - Stats mailing list archive at Nabble.com.
_______________________________________________
grass-stats mailing list
[email protected]
http://lists.osgeo.org/mailman/listinfo/grass-stats
--
Roger Bivand
Department of Economics, NHH Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: [email protected]
_______________________________________________
grass-stats mailing list
[email protected]
http://lists.osgeo.org/mailman/listinfo/grass-stats
