Hace unos días AEMET cerró el acceso público a sus servidores FTP. Esta decisión contrasta con las estrategias de open data defendidas desde la Comisión Europea, o con las posiciones adoptadas por otros organismos (por ejemplo, CM-SAF o NREL-MIDC).
ACTUALIZACIÓN: La web tiempodiario ha cambiado su estructura. El código que explico ya no es de utilidad para esta web. Para descargar los datos hay que seguir el formato que indica el webmaster de la página en uno de los comentarios de este artículo.
Afortunadamente, los datos publicados por AEMET hasta el 28 de Octubre de este año han sido recopilados en la página web tiempodiario.com. Con mis limitados conocimientos de HTML y la inestimable ayuda del paquete XML, he aprendido a acceder a la información de radiación solar de manera automática con R. Expongo a continuación lo que he conseguido. Debo avisar de que no siempre funciona bien el acceso a las páginas, y es posible que se produzca un error en el bucle. Una posible solución es reiniciar el bucle a partir del enlace problemático, o añadir controles de error con try:
library(XML)
library(zoo)
mainURL <- 'http://www.tiempodiario.com'
radURL <- 'http://www.tiempodiario.com/radiacion/'
## Días disponibles en tiempodiario.com
days <- seq(as.POSIXct('2010-12-05'), as.POSIXct('2012-10-28'), by='day')
for (i in seq_along(days)){
day <- days[i]
dayURL <- paste(radURL, format(day, '%Y/%m/%d'), '/', sep='')
## Obtengo el contenido HTML de la página del día en cuestión
content <- htmlParse(dayURL, encoding='utf8')
## y extraigo los enlaces a cada estación para ese día
links <- xpathSApply(content, "//a/@href")
linksData <- links[grep('radiacion-estacion', links)]
## Los nombres de las estaciones están escritos en los encabezados H2
h2 <- xpathSApply(content, "//h2")
namesStations <- sapply(h2[-1], function(s){
val <- xmlValue(s)
substr(val, 10, nchar(val))
})
## pero sus códigos de identificación están en los enlaces correspondientes
idStation <- strsplit(linksData, '/')
idStations <- sapply(idStation, function(x)x[4])
dayStationURL <- paste(mainURL, linksData, sep='')
## Accedo a la página de cada estación para ese día
zList <- lapply(dayStationURL, function(ds){
print(ds)
## Obtengo el contenido HTML de la página
content <- htmlParse(ds, encoding='utf8')
## y extraigo las tablas
tabs <- readHTMLTable(content, header=TRUE)
## Cada tabla codifica el índice temporal en la primera
## columna. A pesar de usar header=TRUE, la primera fila
## contiene el nombre las columnas. De ahí que tenga que usar
## -1 en el indexado.
time <- tabs[[1]][-1,1]
dayTime <- format(as.POSIXct(paste(day, time)), '%Y-%m-%d %H:%M')
## El nombre de la variable de cada tabla está codificado en
## los encabezados H2
h2 <- getNodeSet(content, '//h2')
tabLabels <- sapply(h2, xmlValue)
## Sólo me quedaré con las radiación difusa, directa y global
radVars <- c('Radiación Difusa (DF)',
'Radiación Directa (DT)',
'Radiación Global (GL)')
whichTabs <- which(tabLabels %in% radVars)
## Extraigo la información de las tablas seleccionadas. Dado
## que la primera fila contiene los nombres de las columnas,
## todo el contenido es leído como un "factor". De ahí que
## tenga que usar la conversión as.character y después
## as.numeric.
dat <- sapply(tabs[whichTabs],
function(tab)as.numeric(as.character(tab[-1, 2])))
## Finalmente, si todo ha ido bien, construyo un objeto zoo
## para almacenar la serie temporal de esta estación para este
## día. Si algo no ha funcionado obtendré NULL.
if (length(dat) > 0 & !is.list(dat)) {
z <- zoo(dat, as.POSIXct(dayTime))
names(z) <- tabLabels[whichTabs]
} else {
z <- NULL
}
z
})
## El resultado de lapply será una lista de objetos zoo, uno
## para cada estación. Pongo los nombres de las estaciones a
## cada componente de esta lista y grabo el resultado como un
## fichero RData.
names(zList) <- make.names(namesStations)
save(zList, file=paste(day, 'RData', sep='.'))
## Este código se repetirá para cada día, de forma que obtendré
## un conjunto de ficheros RData, uno por día. Cada uno de ellos
## contiene una lista de zoo's, correspondiendo cada componente
## de la lista a una estación de medida.
}
Una vez ejecutado el código anterior obtendremos una colección de ficheros diarios. Si nos interesa tener el registro completo de una estación determinada, es posible obtenerlo recorriendo todos los ficheros:
## Ficheros que he obtenido del bucle anterior.
files <- dir(pattern='RData')
## Por ejemplo, para obtener el registro de la estación de
## Valladolid, recorro todos los ficheros con lapply, y extraigo
## de cada uno el componente que me interesa.
valladolid <- lapply(files, FUN=function(f){
load(f)
dat <- zList$Valladolid
})
## El resultado de lapply será una lista de objetos zoo, que puedo
## transformar en un único objeto zoo con do.call
valladolid <- do.call(rbind, valladolid)
Con otro bucle adicional se puede automatizar el paso anterior y obtener toda la colección de estaciones:
stations <- names(zList)
lapply(stations, FUN=function(name){
datStation <- lapply(files, FUN=function(f){
load(f)
dat <- zList[[name]]
})
datStation <- do.call(rbind, datStation)
write.zoo(datStation, sep=';', file=paste(name, '.txt', sep=''))
})
