[GRASS-stats] Re: Problem with v.distance in spgrass6 package in R

Markus Metz markus.metz.giswork at googlemail.com
Thu Nov 3 03:19:39 EDT 2011


2011/11/3 Roger Bivand <Roger.Bivand at nhh.no>:
> 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?

The spearfish example works for me too, exactly the same result like yours.

Markus M

>
> 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 at PERMANENT',fields1")
>
> Copy vector <fields at PERMANENT> to current mapset as <fields1>
>>
>> execGRASS("g.copy", vect="'archsites at PERMANENT',archsites1")
>
> Copy vector <archsites at 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, at alaska, at ntv2_0.gsb, at 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, at alaska, at ntv2_0.gsb, at 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 at PERMANENT',fields1 # adding  polygon to mapset
>>
>> g.copy vect='archsites at 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
>> grass-stats at lists.osgeo.org
>> 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: Roger.Bivand at nhh.no
>
> _______________________________________________
> grass-stats mailing list
> grass-stats at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-stats
>
>


More information about the grass-stats mailing list