Automatically downloading and processing TRMM rainfall data

TRMM rainfall data maybe is the most accurate rainfall data derived from satellite measurements and a valuable source in regions with scarse weather-stations. It has a good spatial (0.25°) and temporal (daily) resolution and is available since 1998. However, downloading and processing may be a lot of work, if not scripted. The following script may be badly coded, but it works. All you need is the open source software R. This script is written for a Linux environment.

First, load some packages:

library(raster)
 library(rgdal)
 library(RCurl)

Then, download the data from the FTP server: ftp://disc2.nascom.nasa.gov/data/TRMM/Gridded/. This example downloads 3B43 monthly data version 7 for 2012. The first line has to be changed for the desired data:

url = "ftp://disc2.nascom.nasa.gov/data/TRMM/Gridded/3B43_V7/2012/"
filenames = getURL(url, ftp.use.epsv = FALSE, ftplistonly=TRUE, crlf=TRUE)
 filenames = paste(url, strsplit(filenames, "\r*\n") [[1]], sep="")
filenames = filenames[grep(filenames, pattern="*precipitation.accum")] # this line allows to download only files with a certain pattern, e.g. only certain months or days. "*precipitation.accum" means monthly accumulated rainfall here.
filenames # list all files
mapply(download.file, filenames, basename(filenames)) # download files

Now I create a virtual file as the downloaded TRMM data come as binarys. The VRT-file contains all the downloaded binary files with the appropriate geo-informations. To automate the process, the following script generates the VRT-file (TRMM.vrt) for all 2012 data. Change “3B43.12” accordingly to the downloaded data. Save this example as trmm2012.sh.

echo '<VRTDataset rasterXSize="1440" rasterYSize="400"> 
 <Geotransform>-180,0.25,0,50, 0, -0.25</Geotransform> 
 <SRS>WGS84</SRS>' >TRMM.vrt;
for i in 3B43.12*
do echo '<VRTRasterBand dataType="Float32" band="1" subClass="VRTRawRasterBand"> 
 <SourceFilename relativeToVRT="1">'$i'</SourceFilename> 
 <ByteOrder>MSB</ByteOrder>
 <ImageOffset>0</ImageOffset>
 <PixelOffset>4</PixelOffset>
 <LineOffset>5760</LineOffset>
 </VRTRasterBand>' >>TRMM.vrt
done;
echo '</VRTDataset>' >>TRMM.vrt

Within R, the script (trmm2012.sh) is executed and the virtual-file (TRMM.vrt) loaded as a rasterbrick. This is flipped in y direction and the files written as multilayer Geotiff. This Geotiff contains all the layers for 2012 and can be opened in every GIS software.

system("./trmm2012.sh")
 b=brick("TRMM.vrt")
 b
 trmm = flip(b, direction='y')
 writeRaster(trmm, filename="trmm_acc_2012.tif", format="GTiff", overwrite=TRUE)