Hola, mi trozo de c�digo es:
# Primero creamos el mapa de toda la instrumentaci�n
rm(list=ls()) # remove all the variables from the workspace
setwd("C:/Users/Marcos/Documents/TFM/Situacion_pluviometros/")
library(ggmap)
library(rgdal)
library(ggplot2)
library(scales)
library(maptools)
# Convertimos coordenadas de UTM a geogr�ficas:
# Leemos el archivo de texto donde est�n las coordenadas UTM de los Pluvi�metros
coord <-
read.delim("C:/Users/Marcos/Documents/TFM/Situacion_pluviometros/pluviometrosVenero.txt",
header=TRUE)[,4:5]
# A�adimos coordenadas del TDR, los dos limn�metros y el radar meteorol�gico
(figura 3.1 TFM_Carlos)
coord.todo <- rbind(coord,c(356489,4474915), c(360744,4470212),
c(360040,4471782),
c(359088,4474084))
rownames(coord.todo)=make.names(c(rep("Pluvi�metros",7), "Radar", "TDR",
rep("Limn�metros",2)), unique=TRUE)
coord.todo.spdf <- coord.todo # ".spdf": Class for spatial attributes that
have spatial point locations
colnames(coord.todo.spdf) <- c("x","y")
coordinates(coord.todo.spdf) <- ~x+y
# To retrieve the CRS for a spatial object: proj4string(x)
proj4string(coord.todo.spdf) <- CRS("+proj=utm +zone=30 ellps=WGS84")
# Ahora pasamos de un CRS a otro (pasamos de coordenadas UTM a geogr�ficas)
coord.geog <- as.data.frame(spTransform(coord.todo.spdf, CRS("+proj=longlat
+datum=WGS84")))
colnames(coord.geog) <- c("lon", "lat")
# Creamos un data.frame llamado "datos" para poder pintar bien la leyenda del
mapa
datos <- coord.todo
datos$lon <- coord.geog$lon
datos$lat <- coord.geog$lat
datos$INSTRUMENTACI�N <-
c("Pluvi�metro","Pluvi�metro","Pluvi�metro","Pluvi�metro","Pluvi�metro","Pluvi�metro","Pluvi�metro",
"Radar
meteorol�gico","TDR","Limn�metro","Limn�metro")
# Ahora con la funci�n "get_map" consultamos en Google Maps el mapa
# que queremos crear de nuestra zona de estudio
mapa <- get_map(location=c(lon=-4.66, lat=40.39), zoom=13, maptype="terrain",
source="google")
nombres.pluvio <- c("Arromoro", "Atalaya", "Dehesa", "Gasolinera", "Collado
Morales",
"Pe�a Parda", "Trampalones")
# Cercamos el �rea de la cuenca
area <- readShapePoly("venero_Project.shp") # Lee el shape con el area de
la cuenca
proj4string(area) # Describes data's current
coordinate reference system
# To change to correct projection:
area <- spTransform(area, CRS("+proj=longlat +datum=WGS84"))
proj4string(area) <- "+proj=utm +zone=30 +north +ellps=WGS84 +datum=WGS84
+units=m +no_defs"
png("mapa_zona_estudio.png", units="in", width=9, height=8.5, res=300)
# Vamos a pintar los distintos instrumentos de medida sobre el mapa (junto con
�ste) y a�adimos leyenda
ggmap(mapa, extend='device', legend="left", base_layer=ggplot(datos, aes(x=lon,
y=lat))) +
geom_point(data=datos, aes(shape = INSTRUMENTACI�N, color=INSTRUMENTACI�N,
fill=INSTRUMENTACI�N), size=7.5) +
geom_text(data=coord.geog[1:7,], aes(x=coord.geog$lon[1:7],
y=coord.geog$lat[1:7], label=nombres.pluvio, hjust=0.2, vjust=1), size=5.2) +
scale_shape_manual(values=c(21,24,22,23))
geom_polygon(aes(x=long, y=lat, group=id), data=area, color="white",
alpha=.4, size=.2)
dev.off()
El error que me da la consola de R es:
Error in spTransform(xSP, CRSobj, ...) :
No transformation possible from NA reference system
�Alguien me puede decir lo que hay que hacer para que ya no me de este error?
Muchas gracias de antemano.
Un saludo,
Marcos
[[alternative HTML version deleted]]
_______________________________________________
R-help-es mailing list
R-help-es@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-help-es
