{"word_count": 264, "character_count": 3744}
Pour aller un peu plus loin par rapport au précédent billet, on peut faire une carte dynamique, pour visualiser une tempête, sur une journée. Je vais partir du fait qu'on a encore les objets du précédant billet en mémoire.On commence par garder uniquement les données relatives a la journée qui nous intéresse,date_sel = "2016-10-13"date_sel_OK = paste0(substr(date_sel,9,10),"/",substr(date_sel,6,7),"/",substr(date_sel,1,4))donnees_carte_0 = subset(data, data$date==as.POSIXct(date_sel))horaires = sort(unique(donnees_carte_0$heure))(mais on pourrait prendre une semaine). On recupere ensuite le contour de la France,download.file("http://biogeo.ucdavis.edu/data/gadm2.8/rds/FRA_adm0.rds","FRA_adm0.rds")FR0=readRDS("FRA_adm0.rds")P1=FR0@polygons[[1]]@Polygons[[355]]@coordsP2=FR0@polygons[[1]]@Polygons[[27]]@coordscar seule la France métropolitaine nous intéressera (incluant la Corse, soit 2 gros polygones) dans l'esprit du précédant billet sur le zonier. On va ainsi definir un maillage sur lesquel on va lisser la pluviometrie,grille = expand.grid(seq(min(donnees_carte_0$longitude),max(donnees_carte_0$longitude),length=101),seq(min(donnees_carte_0$latitude),max(donnees_carte_0$latitude),length=101))paslong=(max(donnees_carte_0$longitude)-min(donnees_carte_0$longitude))/100paslat=(max(donnees_carte_0$latitude)-min(donnees_carte_0$latitude))/100f=function(i){ (point.in.polygon (grille[i, 1]+paslong/2 , grille[i, 2]+paslat/2 , P1[,1],P1[,2])&#62;0)+(point.in.polygon (grille[i, 1]+paslong/2 , grille[i, 2]+paslat/2 , P2[,1],P2[,2])&#62;0) }indic=unlist(lapply(1:nrow(grille),f))grille=grille[which(indic==1),]Pour lisser, on utilise les [latex]k[/latex]-plus proches voisinsknn=function(i,k=10){  d=distHaversine(grille[i,1:2],donnees_carte[,    c("longitude","latitude")], r=6378.137)   r=rank(d)  ind=which(r&#60;=k)  weighted.mean(donnees_carte[ind,"rr3"],(1/d[ind])/sum(1/d[ind]))}On a ensuite la fonction suivante, pour faire une cartecarto_prec&#60;-function(){  grille2&#60;-grille  grille2$rr=Vectorize(knn)(i=1:nrow(grille2))  bk=seq(0,50,length=21)  grille2$cuty=cut(grille2$rr,breaks=bk,labels=1:20)  cols = rev(carto.pal(pal1 = "blue.pal", n1=20, pal2 = "white.pal", n2=1))  plot(FR0,border=NA)  polygon(P1)  polygon(P2)  points(grille2[,1]+paslong/2,grille2[,2]+paslat/2,col=cols[grille2$cuty],pch=19)  points(donnees_carte$longitude,donnees_carte$latitude, col="black",pch=19,cex=.5)  title(main = paste0(date_sel_OK," à ",heure_sel,"H"),        sub = "Précipitations des 3 dernières heures",        cex.main = 1.5,   font.main= 4,        cex.sub = 1, font.sub = 3)  legend(8.2, 50, legend=seq(0,50,length=5), title='en mm',         fill=cols[seq(1,20,length=5)], cex=0.8)}On va ensuite utiliser 8 fichiers - car on regarde toutes les 3 heuresfor (hh in horaires) {  heure_sel&#60;-hh  donnees_carte &#60;- subset(donnees_carte_0, donnees_carte_0$heure==heure_sel)  png(paste0("CartePrec",date_sel,"_",hh,".png"))  carto_prec()  dev.off()}(mais ça peut se raffiner, avec un lissage spatio-temporel). Et on va ensuite créer notre animation en concaténant les imagesframes &#60;- image_morph(  c(image_scale(image_read("CartePrec2016-10-13_1.png")),image_scale(image_read("CartePrec2016-10-13_3.png")),image_scale(image_read("CartePrec2016-10-13_6.png")),image_scale(image_read("CartePrec2016-10-13_9.png")),image_scale(image_read("CartePrec2016-10-13_12.png")),image_scale(image_read("CartePrec2016-10-13_15.png")),image_scale(image_read("CartePrec2016-10-13_18.png")),image_scale(image_read("CartePrec2016-10-13_21.png"))), frames = 16)frames.anim &#60;- image_animate(frames)image_write(frames.anim, paste0("CartePrecDay",date_sel,".gif"))Il faut avouer que ça a de gueule... On voit clairement la tempête arriver...