Re: [R-sig-Geo] sp::spDists*

2017-01-13 Thread Edzer Pebesma 
A new version of gstat, fixing the same bug, is now also on CRAN.

On 09/01/17 11:39, Roger Bivand wrote:
> New versions of spdep and spgwr fixing the same bug are now on CRAN.
> 
> Roger
> 
> On Thu, 8 Dec 2016, Edzer Pebesma wrote:
> 
>> While building support for geosphere distance functions in sf and
>> comparing results from geosphere::distMeeuws with those of sp::spDists,
>> I found a typo in the C code underlying the great circle distance
>> computation functions in sp (sp::spDists, sp::spDistsN1) [*]:
>>
>> https://github.com/edzer/sp/commit/d8374ff7efc6735cba9a054748c602bed0672f23
>>
>>
>> Before:
>>> library(sf)
>> Linking to GEOS 3.5.0, GDAL 2.1.0, proj.4 4.9.2
>>> library(sp)
>>> library(units)
>>>
>>> x = st_sfc(
>> + st_point(c(0,0)),
>> + st_point(c(1,0)),
>> + st_point(c(2,0)),
>> + st_point(c(3,0)),
>> + crs = 4326
>> + )
>>>
>>> y = st_sfc(
>> + st_point(c(0,10)),
>> + st_point(c(1,0)),
>> + st_point(c(2,0)),
>> + st_point(c(3,0)),
>> + st_point(c(4,0)),
>> + crs = 4326
>> + )
>>>
>>> st_distance(x, y)
>> Units: m
>>[,1] [,2] [,3] [,4] [,5]
>> [1,] 1105855 111319.5 222639.0 333958.5 445278.0
>> [2,] 387  0.0 111319.5 222639.0 333958.5
>> [3,] 1127822 111319.5  0.0 111319.5 222639.0
>> [4,] 1154693 222639.0 111319.5  0.0 111319.5
>>>
>>> d.sf = st_distance(x, y, geosphere::distMeeus)
>>> d.sp = spDists(as(x, "Spatial"), as(y, "Spatial"))
>>> units(d.sp) = make_unit("km")
>>> d.sf - d.sp
>> Units: m
>> [,1] [,2] [,3] [,4] [,5]
>> [1,] 1851.990 0.00e+00 0.00e+00 0.00e+000
>> [2,] 1842.938 0.00e+00 0.00e+00 2.910383e-110
>> [3,] 1816.564 0.00e+00 0.00e+00 1.455192e-110
>> [4,] 1775.032 2.910383e-11 1.455192e-11 0.00e+000
>>>
>>
>>
>>
>> After:
>>
>>> library(sf)
>> Linking to GEOS 3.5.0, GDAL 2.1.0, proj.4 4.9.2
>>> library(sp)
>>> library(units)
>>>
>>> x = st_sfc(
>> + st_point(c(0,0)),
>> + st_point(c(1,0)),
>> + st_point(c(2,0)),
>> + st_point(c(3,0)),
>> + crs = 4326
>> + )
>>>
>>> y = st_sfc(
>> + st_point(c(0,10)),
>> + st_point(c(1,0)),
>> + st_point(c(2,0)),
>> + st_point(c(3,0)),
>> + st_point(c(4,0)),
>> + crs = 4326
>> + )
>>>
>>> st_distance(x, y)
>> Units: m
>>[,1] [,2] [,3] [,4] [,5]
>> [1,] 1105855 111319.5 222639.0 333958.5 445278.0
>> [2,] 387  0.0 111319.5 222639.0 333958.5
>> [3,] 1127822 111319.5  0.0 111319.5 222639.0
>> [4,] 1154693 222639.0 111319.5  0.0 111319.5
>>>
>>> d.sf = st_distance(x, y, geosphere::distMeeus)
>>> d.sp = spDists(as(x, "Spatial"), as(y, "Spatial"))
>>> units(d.sp) = make_unit("km")
>>> d.sf - d.sp
>> Units: m
>>  [,1] [,2] [,3] [,4] [,5]
>> [1,] -2.328306e-10 0.00e+00 0.00e+00 0.00e+000
>> [2,]  2.328306e-10 0.00e+00 0.00e+00 2.910383e-110
>> [3,]  0.00e+00 0.00e+00 0.00e+00 1.455192e-110
>> [4,]  0.00e+00 2.910383e-11 1.455192e-11 0.00e+000
>>>
>>
>> [*] assuming the source code of
>> http://www.abecedarical.com/javascript/script_greatcircle.html is
>> correct; I don't have access to the Meeuws book mentioned there under
>> (1).
>>
>> The same C function is used by gstat; I corrected that as well.
>>
> 

-- 
Edzer Pebesma
Institute for Geoinformatics  (ifgi),  University of Münster
Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
Journal of Statistical Software:   http://www.jstatsoft.org/
Computers & Geosciences:   http://elsevier.com/locate/cageo/



signature.asc
Description: OpenPGP digital signature
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: [R-sig-Geo] sp::spDists*

2017-01-09 Thread Roger Bivand 

New versions of spdep and spgwr fixing the same bug are now on CRAN.

Roger

On Thu, 8 Dec 2016, Edzer Pebesma wrote:


While building support for geosphere distance functions in sf and
comparing results from geosphere::distMeeuws with those of sp::spDists,
I found a typo in the C code underlying the great circle distance
computation functions in sp (sp::spDists, sp::spDistsN1) [*]:

https://github.com/edzer/sp/commit/d8374ff7efc6735cba9a054748c602bed0672f23

Before:

library(sf)

Linking to GEOS 3.5.0, GDAL 2.1.0, proj.4 4.9.2

library(sp)
library(units)

x = st_sfc(

+ st_point(c(0,0)),
+ st_point(c(1,0)),
+ st_point(c(2,0)),
+ st_point(c(3,0)),
+ crs = 4326
+ )


y = st_sfc(

+ st_point(c(0,10)),
+ st_point(c(1,0)),
+ st_point(c(2,0)),
+ st_point(c(3,0)),
+ st_point(c(4,0)),
+ crs = 4326
+ )


st_distance(x, y)

Units: m
   [,1] [,2] [,3] [,4] [,5]
[1,] 1105855 111319.5 222639.0 333958.5 445278.0
[2,] 387  0.0 111319.5 222639.0 333958.5
[3,] 1127822 111319.5  0.0 111319.5 222639.0
[4,] 1154693 222639.0 111319.5  0.0 111319.5


d.sf = st_distance(x, y, geosphere::distMeeus)
d.sp = spDists(as(x, "Spatial"), as(y, "Spatial"))
units(d.sp) = make_unit("km")
d.sf - d.sp

Units: m
[,1] [,2] [,3] [,4] [,5]
[1,] 1851.990 0.00e+00 0.00e+00 0.00e+000
[2,] 1842.938 0.00e+00 0.00e+00 2.910383e-110
[3,] 1816.564 0.00e+00 0.00e+00 1.455192e-110
[4,] 1775.032 2.910383e-11 1.455192e-11 0.00e+000






After:


library(sf)

Linking to GEOS 3.5.0, GDAL 2.1.0, proj.4 4.9.2

library(sp)
library(units)

x = st_sfc(

+ st_point(c(0,0)),
+ st_point(c(1,0)),
+ st_point(c(2,0)),
+ st_point(c(3,0)),
+ crs = 4326
+ )


y = st_sfc(

+ st_point(c(0,10)),
+ st_point(c(1,0)),
+ st_point(c(2,0)),
+ st_point(c(3,0)),
+ st_point(c(4,0)),
+ crs = 4326
+ )


st_distance(x, y)

Units: m
   [,1] [,2] [,3] [,4] [,5]
[1,] 1105855 111319.5 222639.0 333958.5 445278.0
[2,] 387  0.0 111319.5 222639.0 333958.5
[3,] 1127822 111319.5  0.0 111319.5 222639.0
[4,] 1154693 222639.0 111319.5  0.0 111319.5


d.sf = st_distance(x, y, geosphere::distMeeus)
d.sp = spDists(as(x, "Spatial"), as(y, "Spatial"))
units(d.sp) = make_unit("km")
d.sf - d.sp

Units: m
 [,1] [,2] [,3] [,4] [,5]
[1,] -2.328306e-10 0.00e+00 0.00e+00 0.00e+000
[2,]  2.328306e-10 0.00e+00 0.00e+00 2.910383e-110
[3,]  0.00e+00 0.00e+00 0.00e+00 1.455192e-110
[4,]  0.00e+00 2.910383e-11 1.455192e-11 0.00e+000




[*] assuming the source code of
http://www.abecedarical.com/javascript/script_greatcircle.html is
correct; I don't have access to the Meeuws book mentioned there under (1).

The same C function is used by gstat; I corrected that as well.



--
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; e-mail: roger.biv...@nhh.no
http://orcid.org/-0003-2392-6140
https://scholar.google.no/citations?user=AWeghB0J&hl=en
http://depsy.org/person/434412

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: [R-sig-Geo] sp::spDists*

2016-12-28 Thread Edzer Pebesma 
The problem described below has been corrected ins sp 1.2-4 which is now
on CRAN; new releases of gstat, spgwr and spdep will follow to correct
for the same bug.

On 08/12/16 15:48, Edzer Pebesma wrote:
> While building support for geosphere distance functions in sf and
> comparing results from geosphere::distMeeuws with those of sp::spDists,
> I found a typo in the C code underlying the great circle distance
> computation functions in sp (sp::spDists, sp::spDistsN1) [*]:
> 
> https://github.com/edzer/sp/commit/d8374ff7efc6735cba9a054748c602bed0672f23
> 
> Before:
>> library(sf)
> Linking to GEOS 3.5.0, GDAL 2.1.0, proj.4 4.9.2
>> library(sp)
>> library(units)
>>
>> x = st_sfc(
> + st_point(c(0,0)),
> + st_point(c(1,0)),
> + st_point(c(2,0)),
> + st_point(c(3,0)),
> + crs = 4326
> + )
>>
>> y = st_sfc(
> + st_point(c(0,10)),
> + st_point(c(1,0)),
> + st_point(c(2,0)),
> + st_point(c(3,0)),
> + st_point(c(4,0)),
> + crs = 4326
> + )
>>
>> st_distance(x, y)
> Units: m
> [,1] [,2] [,3] [,4] [,5]
> [1,] 1105855 111319.5 222639.0 333958.5 445278.0
> [2,] 387  0.0 111319.5 222639.0 333958.5
> [3,] 1127822 111319.5  0.0 111319.5 222639.0
> [4,] 1154693 222639.0 111319.5  0.0 111319.5
>>
>> d.sf = st_distance(x, y, geosphere::distMeeus)
>> d.sp = spDists(as(x, "Spatial"), as(y, "Spatial"))
>> units(d.sp) = make_unit("km")
>> d.sf - d.sp
> Units: m
>  [,1] [,2] [,3] [,4] [,5]
> [1,] 1851.990 0.00e+00 0.00e+00 0.00e+000
> [2,] 1842.938 0.00e+00 0.00e+00 2.910383e-110
> [3,] 1816.564 0.00e+00 0.00e+00 1.455192e-110
> [4,] 1775.032 2.910383e-11 1.455192e-11 0.00e+000
>>
> 
> 
> 
> After:
> 
>> library(sf)
> Linking to GEOS 3.5.0, GDAL 2.1.0, proj.4 4.9.2
>> library(sp)
>> library(units)
>>
>> x = st_sfc(
> + st_point(c(0,0)),
> + st_point(c(1,0)),
> + st_point(c(2,0)),
> + st_point(c(3,0)),
> + crs = 4326
> + )
>>
>> y = st_sfc(
> + st_point(c(0,10)),
> + st_point(c(1,0)),
> + st_point(c(2,0)),
> + st_point(c(3,0)),
> + st_point(c(4,0)),
> + crs = 4326
> + )
>>
>> st_distance(x, y)
> Units: m
> [,1] [,2] [,3] [,4] [,5]
> [1,] 1105855 111319.5 222639.0 333958.5 445278.0
> [2,] 387  0.0 111319.5 222639.0 333958.5
> [3,] 1127822 111319.5  0.0 111319.5 222639.0
> [4,] 1154693 222639.0 111319.5  0.0 111319.5
>>
>> d.sf = st_distance(x, y, geosphere::distMeeus)
>> d.sp = spDists(as(x, "Spatial"), as(y, "Spatial"))
>> units(d.sp) = make_unit("km")
>> d.sf - d.sp
> Units: m
>   [,1] [,2] [,3] [,4] [,5]
> [1,] -2.328306e-10 0.00e+00 0.00e+00 0.00e+000
> [2,]  2.328306e-10 0.00e+00 0.00e+00 2.910383e-110
> [3,]  0.00e+00 0.00e+00 0.00e+00 1.455192e-110
> [4,]  0.00e+00 2.910383e-11 1.455192e-11 0.00e+000
>>
> 
> [*] assuming the source code of
> http://www.abecedarical.com/javascript/script_greatcircle.html is
> correct; I don't have access to the Meeuws book mentioned there under (1).
> 
> The same C function is used by gstat; I corrected that as well.
> 
> 
> 
> ___
> R-sig-Geo mailing list
> R-sig-Geo@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 

-- 
Edzer Pebesma
Institute for Geoinformatics  (ifgi),  University of Münster
Heisenbergstraße 2, 48149 Münster, Germany; +49 251 83 33081
Journal of Statistical Software:   http://www.jstatsoft.org/
Computers & Geosciences:   http://elsevier.com/locate/cageo/



signature.asc
Description: OpenPGP digital signature
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo