#realizat : Alexandru Dumitrescu #data : 2012-01-03 #descriere : prelucrarea produsului SC MODIS #email : alexandru.dumitrescu [at] gmail.com library(RCurl) library(raster) library(rgdal) #setwd("~/DR/doc/text/R/Modis") # variabilele de selectie a produselor MODIS url <- "ftp://n4ftl01u.ecs.nasa.gov/MOSA/MYD10A1.005/" an<-"2011" luna<-"12" zi<-"26" data<-paste(an,luna,zi,sep=".") # creeaza adresa ftp unde se gasesc imaginile modis ftp<-paste(url,data,"/",sep="") fileNames <- getURL(ftp) #variabila cu indecșii imaginilor pentru Romania p<-c("h19v04","h20v04") # descarcarea imaginilor MODIS pentru Romania for (i in p) {print(i) files = strsplit(fileNames,"\r*\n")[[1]] # selectie imagini pentru Romania ff<-files[grep(files,pattern=paste(i,"*.hdf",sep="."))] ff<-ff[nchar(ff)==83] ff<-unlist(strsplit(ff," ")) ff<-ff[length(ff)] download.file(paste(url,data, "/", ff,sep=""),destfile=paste(getwd(), ff, sep="/"),mode='wb') } #lista cu fisierele *.hdf descarcate files<-list.files(pattern=".hdf$",full.names=T) #GDALinfo(files[1]) #citirea datelor MODIS in R ls1<-readGDAL(paste('HDF4_EOS:EOS_GRID:\"',files[1],'\":MOD_Grid_Snow_500m:Snow_Cover_Daily_Tile',sep="")) ls2<-readGDAL(paste('HDF4_EOS:EOS_GRID:\"',files[2],'\":MOD_Grid_Snow_500m:Snow_Cover_Daily_Tile',sep="")) # mozaicarea datelor MODIS rm<-merge(raster(ls1),raster(ls2)) rm<-crop(rm,extent(1550000,2400000,4818268.68831,55000000)) #proiecteaza din proiectie Sinusoidala in WGS84 srm<-projectRaster(rm, crs="+init=epsg:4326", method="ngb") #comaseaza Lake or inland water cu Open water (ocean) srm[srm==39]<-37 #comaseaza Snow-covered lake ice cu Snow srm[srm==200]<-100 #vizualizarea rezultatelor cu frontiere plot(srm,xlim=c(20,30),ylim=c(43.5,48.3),breaks=c(0,25,37,50,100),col=terrain.colors(4), axis.args=list(at=c(10,30,42,75), labels=c('fara zapada','ape','nori','strat zapada'))) download.file('http://earth.unibuc.ro/file_download/1251','fro.zip',mode = "wb") unzip('fro.zip') ro<-readOGR(dsn=getwd(), layer="frontiera_ro") plot(ro,add=T) # scrie datele MODIS in format GEOTIFF, cu numele corespunzator zilei srm<-crop(srm,extent(20,30,43.5,48.3)) writeRaster(srm, paste("SC_",gsub("\\.","_",data),".tif",sep=""),overwrite=T) ##stergerea fisierelor HDF4 descarcate unlink(files)