Bonjour
Je vous écris parce que malgré le temps que j'y ai passé et mes recherches sur internet, je n'ai pas trouvé de solution à mon problème.
Je dispose de données géolocalisées au format longitude-latitude (WGS84, EPSG:4326). Je veux faire une analyse géostatistique et du krigeage et pour ça, je dois les convertir en UTM, ce que j'ai fait sans problème (code ci dessous).
1 2 3
| coordinates(Mesures)= ~LON+LAT
proj4string(Mesures) <- CRS("+proj=longlat +datum=WGS84")
Mesures <- spTransform(Mesures, CRS("+proj=utm +zone=17 +south +datum=WGS84")) |
Le problème, c'est que j'aimerais reconvertir ensuite le vecteur de krigeage au format WGS84, EPSG:4326 et c'est là que j'ai un problème. Le code est :
1 2 3 4
| Prediction <- krige(NDVI~1, Mesures, Grille, VariNdviFit)
proj4string(Prediction) <- CRS("+proj=utm +zone=17 +south +datum=WGS84")
Prediction <- spTransform(Prediction, CRS("+proj=longlat +datum=WGS84")) |
Voici le message d'erreur qui m'est renvoyé*:
Grid warping not available, coercing to points
Les packages que j'utilise sont :
1 2 3 4 5 6 7
| library(gstat)
library(sp)
library(automap)
library(akima)
library(fields) # Lie automatiquement les packages spam & grid
library(rgdal) # For coordinates conversion
library(raster) |
Apparemment, d'après mes recherches, ça serait du à la grille qui serait déformée lors de la première projection. Sauf que ça ne m'explique pas comment résoudre le problème…
Voici le résultat d'un summary prediction à la fin de l'éxécution du programe*:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
| > summary(Prediction)
Object of class SpatialPointsDataFrame
Coordinates:
min max
X -83.76920 -83.75626
Y 10.01426 10.02539
Is projected: FALSE
proj4string :
[+proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0]
Number of points: 17466
Data attributes:
var1.pred var1.var
Min. :-0.02216 Min. :0.5869
1st Qu.: 0.54258 1st Qu.:0.6608
Median : 0.72105 Median :0.7509
Mean : 0.83489 Mean :0.8027
3rd Qu.: 0.96535 3rd Qu.:0.9684
Max. : 2.69623 Max. :1.0394
> |
Est ce que quelqu'un saurait m'aider s'il vous plait ? Merci d'avance !
Partager