[Live-demo] [OSGeo] #998: R demo portion fails
OSGeo
trac_osgeo at osgeo.org
Sun Aug 19 00:13:36 PDT 2012
#998: R demo portion fails
-----------------------+----------------------------------------------------
Reporter: darkblueb | Owner: live-demo@…
Type: defect | Status: new
Priority: normal | Milestone: OSGeoLive6.5
Component: LiveDVD | Keywords: R
-----------------------+----------------------------------------------------
due to a missing dependency ''maps''
{{{
> library(gstat)
> demo(wind)
demo(wind)
---- ~~~~
Type <Return> to start :
> #pdf("wind.pdf")
> library(gstat)
> library(rgdal)
Geospatial Data Abstraction Library extensions to R successfully loaded
Loaded GDAL runtime: GDAL 1.9.1, released 2012/05/15
Path to GDAL shared files: /usr/share/gdal/1.9
Loaded PROJ.4 runtime: Rel. 4.7.1, 23 September 2009, [PJ_VERSION: 470]
Path to PROJ.4 shared files: (autodetected)
> library(maptools)
Loading required package: foreign
Checking rgeos availability: FALSE
Note: when rgeos is not available, polygon geometry
computations in maptools depend on gpclib,
which has a restricted licence. It is disabled by default;
to enable gpclib, type gpclibPermit()
> # load wind data, run test:
> example(wind)
wind> data(wind)
wind> summary(wind)
year month day RPT
Min. :61.0 Min. : 1.000 Min. : 1.00 Min. : 0.67
1st Qu.:65.0 1st Qu.: 4.000 1st Qu.: 8.00 1st Qu.: 8.12
Median :69.5 Median : 7.000 Median :16.00 Median :11.71
Mean :69.5 Mean : 6.523 Mean :15.73 Mean :12.36
3rd Qu.:74.0 3rd Qu.:10.000 3rd Qu.:23.00 3rd Qu.:15.92
Max. :78.0 Max. :12.000 Max. :31.00 Max. :35.80
VAL ROS KIL SHA
Min. : 0.21 Min. : 1.50 Min. : 0.000 Min. : 0.13
1st Qu.: 6.67 1st Qu.: 8.00 1st Qu.: 3.580 1st Qu.: 6.75
Median :10.17 Median :10.92 Median : 5.750 Median : 9.96
Mean :10.65 Mean :11.66 Mean : 6.306 Mean :10.46
3rd Qu.:14.04 3rd Qu.:14.67 3rd Qu.: 8.420 3rd Qu.:13.54
Max. :33.37 Max. :33.84 Max. :28.460 Max. :37.54
BIR DUB CLA MUL
Min. : 0.000 Min. : 0.000 Min. : 0.000 Min. : 0.000
1st Qu.: 4.000 1st Qu.: 6.000 1st Qu.: 5.090 1st Qu.: 5.370
Median : 6.830 Median : 9.210 Median : 8.080 Median : 8.170
Mean : 7.092 Mean : 9.797 Mean : 8.494 Mean : 8.496
3rd Qu.: 9.670 3rd Qu.:12.960 3rd Qu.:11.420 3rd Qu.:11.210
Max. :26.160 Max. :30.370 Max. :31.080 Max. :25.880
CLO BEL MAL
Min. : 0.040 Min. : 0.13 Min. : 0.67
1st Qu.: 5.330 1st Qu.: 8.71 1st Qu.:10.71
Median : 8.290 Median :12.50 Median :15.00
Mean : 8.707 Mean :13.12 Mean :15.60
3rd Qu.:11.630 3rd Qu.:16.88 3rd Qu.:19.83
Max. :28.210 Max. :42.38 Max. :42.54
wind> wind.loc
Station Code Latitude Longitude MeanWind
1 Valentia VAL 51d56'N 10d15'W 5.48
2 Belmullet BEL 54d14'N 10d00'W 6.75
3 Claremorris CLA 53d43'N 8d59'W 4.32
4 Shannon SHA 52d42'N 8d55'W 5.38
5 Roche's Point RPT 51d48'N 8d15'W 6.36
6 Birr BIR 53d05'N 7d53'W 3.65
7 Mullingar MUL 53d32'N 7d22'W 4.38
8 Malin Head MAL 55d22'N 7d20'W 8.03
9 Kilkenny KIL 52d40'N 7d16'W 3.25
10 Clones CLO 54d11'N 7d14'W 4.48
11 Dublin DUB 53d26'N 6d15'W 5.05
12 Roslare ROS 52d16'56.791"N 6d21'25.056"W 6.00
wind> wind.loc$y =
as.numeric(char2dms(as.character(wind.loc[["Latitude"]])))
wind> wind.loc$x =
as.numeric(char2dms(as.character(wind.loc[["Longitude"]])))
wind> coordinates(wind.loc) = ~x+y
wind> # fig 1:
wind> ## No test:
wind> if (require(mapdata)) {
wind+ map("worldHires", xlim = c(-11,-5.4), ylim = c(51,55.5))
wind+ plot(wind.loc, add=TRUE, pch=16)
wind+ text(coordinates(wind.loc), pos=1, label=wind.loc$Station)
wind+ }
Loading required package: mapdata
wind> ## End(No test)
wind> wind$time = ISOdate(wind$year+1900, wind$month, wind$day)
wind> # time series of e.g. Dublin data:
wind> plot(DUB~time, wind, type= 'l', ylab = "windspeed (knots)", main =
"Dublin")
Hit <Return> to see next plot:
wind> # fig 2:
wind> #wind = wind[!(wind$month == 2 & wind$day == 29),]
wind> wind$jday = as.numeric(format(wind$time, '%j'))
wind> windsqrt = sqrt(0.5148 * as.matrix(wind[4:15]))
wind> Jday = 1:366
wind> windsqrt = windsqrt - mean(windsqrt)
wind> daymeans = sapply(split(windsqrt, wind$jday), mean)
wind> plot(daymeans ~ Jday)
Hit <Return> to see next plot:
wind> lines(lowess(daymeans ~ Jday, f = 0.1))
wind> # subtract the trend:
wind> meanwind = lowess(daymeans ~ Jday, f = 0.1)$y[wind$jday]
wind> velocity = apply(windsqrt, 2, function(x) { x - meanwind })
wind> # match order of columns in wind to Code in wind.loc:
wind> pts = coordinates(wind.loc[match(names(wind[4:15]),
wind.loc$Code),])
wind> # fig 3, but not really yet...
wind> dists = spDists(pts, longlat=TRUE)
wind> corv = cor(velocity)
wind> sel = !(as.vector(dists) == 0)
wind> plot(as.vector(corv[sel]) ~ as.vector(dists[sel]),
wind+ xlim = c(0,500), ylim = c(.4, 1), xlab = "distance (km.)",
wind+ ylab = "correlation")
Hit <Return> to see next plot:
wind> # plots all points twice, ignores zero distance
wind>
wind> # now really get fig 3:
wind> ros = rownames(corv) == "ROS"
wind> dists.nr = dists[!ros,!ros]
wind> corv.nr = corv[!ros,!ros]
wind> sel = !(as.vector(dists.nr) == 0)
wind> plot(as.vector(corv.nr[sel]) ~ as.vector(dists.nr[sel]), pch = 3,
wind+ xlim = c(0,500), ylim = c(.4, 1), xlab = "distance (km.)",
wind+ ylab = "correlation")
Hit <Return> to see next plot:
wind> # add outlier:
wind> points(corv[ros,!ros] ~ dists[ros,!ros], pch=16, cex=.5)
wind> xdiscr = 1:500
wind> # add correlation model:
wind> lines(xdiscr, .968 * exp(- .00134 * xdiscr))
> m = map2SpatialLines(
+ map("worldHires", xlim = c(-11,-5.4), ylim = c(51,55.5), plot=F))
Loading required package: maps
Error in cbind(map$x, map$y) : could not find function "map"
In addition: Warning messages:
1: In library(package, lib.loc = lib.loc, character.only = TRUE,
logical.return = TRUE, :
there is no package called ‘mapdata’
2: In library(package, lib.loc = lib.loc, character.only = TRUE,
logical.return = TRUE, :
there is no package called ‘maps’
}}}
--
Ticket URL: <http://trac.osgeo.org/osgeo/ticket/998>
OSGeo <http://www.osgeo.org/>
OSGeo committee and general foundation issue tracker.
More information about the Osgeolive
mailing list