P1: Basic visualisation of radar scans in R

1. Preparation

Radar data needed for this practical can be downloaded from Dropbox. Don’t download everything right now, but when you are asked to in the exercises (so we spread our Wifi use). https://www.dropbox.com/sh/8dl30htd9zxdbx3/AAAMOHk4_UlCsXOYb6LDFONOa?dl=0

Throughout the practical we will use the bioRad R package. It’s currently not (yet) in the default CRAN repository, but can be downloaded from Github. See https://github.com/adokter/bioRad for install instructions. Go through steps 1-3 of these install instructions if you haven’t done so already.

Execute each of the code examples provided below in Rstudio, and try to complete the exercises.

# make sure you start with a fresh R session
# load the bioRad package
library(bioRad)
# check the package version
packageVersion("bioRad")
# make sure you have the latest version (0.2.6). If you have an older version, download and reinstall it as follows:
library(devtools)
install_github("adokter/bioRad")

All the functions of the bioRad package have fairly elaborated manual pages, that you can refer to at any moment:

# bring up the package general help page:
?bioRad

In the help window, you can click on each of the links to pull up help pages of specific functions. Click on the Index link at the bottom of each help page to get an overview of all functions and available manual pages.

Start by making a new directory on your local machine that you will use for this practical

# make a new local directory on your machine where to download data for this practical
# replace the string below with the path of that directory:
HOME="your/personal/working/directory/"
# check that the directory exists. If the next statement evaluates to FALSE, something went wrong: the directory does not exist or you didn't specify its path correctly
file.exists(HOME)
# we will make HOME our work directory, the default folder in which R will look
# for files, and where it will output newly generated files.
setwd(HOME)
# Finally, we set the local time zone to UTC, so all plotted time axes will be in UTC
Sys.setenv(TZ="UTC")

Your R session is now properly set up

To work with US NEXRAD data, bioRad currently requires a working installation of Docker (linux / mac) or Docker for Windows (not Docker Toolbox, this is an older implementation of Docker for Windows operating systems that is not supported). If you managed to successfully install Docker, test whether it works

# start your local Docker installation
# we first test whether R can communicate with Docker:
checkDocker()

If you get a “Hello from Docker!” welcome message, everything is working.

If you did not manage to install Docker, you will not be able to load NEXRAD data into R at this time, but you will be able to continue the exercises with European radar data that can be loaded without Docker.

2. Downloading and loading polar volume data

US NEXRAD data can be accessed in the Amazon cloud

A useful site for inspecting premade movies of the US composite is http://www.pauljhurtado.com/US_Composite_Radar/

The names of the radars can be found on this figure: NEXRAD radar network

  1. Inspect the radar movie for 11 Sep 2015 at http://www.pauljhurtado.com/US_Composite_Radar/2015-9-11/. Pick a radar where you see clear signatures of bird migration, and download one polar volume file for this radar at a timestamp where you expect a lot of migration.

  2. Navigate to https://www.dropbox.com/sh/8dl30htd9zxdbx3/AAAMOHk4_UlCsXOYb6LDFONOa?dl=0 and also download a European radar data file, volume.h5, and put it in your working directory.

3. The structure of polar volumes

# the bioRad package also includes an example radar volume file, that we will inspect first
# first locate this example file:
pvol.path <- system.file("extdata", "volume.h5", package="bioRad")
# print the local path of the volume file:
pvol.path
# load the file into R:
pvol=read.pvol(pvol.path)
## print some information about the polar volume
pvol
# print information about the polar scans contained in this polar volume:
pvol$scans
  1. Load the polar volumes you downloaded in exercise 1 and 2 into R. Which is the minimum and maximum scan elevation contained in these files? And which radar quantities are available? Check the help file of the read.pvol function for the nomenclature of various available quantities.

4. Plotting radar scans

# (if you haven't done so already) load the polar volume data from the volume.h5 file you just downloaded
pvol=read.pvol("volume.h5")
# let's extract the third scan, which was collected at 1.5 degree elevation:
pscan = pvol$scans[[3]]
# print some information about this scan:
pscan
# before we can plot the scan, we need to project it on a Cartesian grid,
# i.e. we need to make a Plan Position Indicator (PPI)
my_ppi = ppi(pscan)
# print some information about this ppi:
my_ppi
# you can see we projected it on a 500 metre grid.
# Check the manual of the ppi function to see how you can change the projection
# Now we are ready to plot the ppi
# plot the reflectivity factor image:
plot(my_ppi, param="DBZH")
  1. The European radar file contains convective precipitation, characterised by localized but intense thunderstorms, as well as biological scattering. Make also a ppi plot of the radial velocity and correlation coefficient. Verify how the texture of the radial velocity, and the values of correlation coefficient of the precipitive areas differs from the areas with biological scattering.

  2. Based on the radial velocity image, are the biological scatterers birds or insects? Why?

Navigate to https://www.dropbox.com/sh/8dl30htd9zxdbx3/AAAMOHk4_UlCsXOYb6LDFONOa?dl=0 and download a set of European radar volume files, volumes.zip. Put the file in your working directory and unzip it.

# first we load the file names of the files you just downloaded in a list
# make sure the h5 files you just downloaded are in your working directory.
files=list.files(paste(HOME, "/20150927", sep=""), pattern="seang_pvol_20150927",full.names=TRUE)
# print the file names
files
# we will plot the file at 4:30 UTC first, which is the fourth one:
pvol=read.pvol(files[4])
# let's extract the 2 degree elevation scan, in this case the fifth scan
pscan = pvol$scans[[5]]
# project it as a Plan Position Indicator (PPI)
my_ppi = ppi(pscan,cellsize=200,range.max=25000)
# plot the reflectivity factor (we adjust the plotting scale to -20 to 15 dBZ):
plot(my_ppi, param="DBZH", zlim=c(-20,15))
  1. Check the manual of the ppi function to see how you can change the projection (argument cellsize) and the maximum distance from the radar up to where to plot data (argument range.max). Make a ppi for the lowest elevation scan at 4:30 UTC, up to 20 km distance from the radar, on a 200x200 m resolution grid. Make a plot for the radial velocity (VRAD).

  2. (for those with working Docker installation) Make a PPI plot for reflectivity one of your selected NEXRAD polar volume files. Verify that reflectivity factors during bird migration at S-band (NEXRAD radars US) are much higher than at C-band (most European radars).

5. Overlaying radar scans on maps

# It is often informative to plot radar data on a base layer, such as google earth maps.
# first download the background image:
satelliteImage=basemap(my_ppi,maptype="satellite")
# then overlay the PPI on the satellite image:
map(my_ppi,param="DBZH",satelliteImage, zlim=c(-20,15))
  1. Plot the same ppi but now using a cartographic map instead of satellite maps as baselayer. See help of bioRad’s map function for available types of basemaps.
# Note that in R, spatial data is usually contained in class objects of packag 'sp'
# bioRad also uses these objects in the background, and they can be extracted if you want to.
# The spatial data is stored in the data slot, as in:
my_spatialgrid=my_ppi$data
# you can use the sp package to save the spatial data to all kinds GIS formats, for example ArcGis:
library(sp)
write.asciigrid(my_spatialgrid,"PPI_in_arcgis_format.asc")

6. Plotting a series of scans for multiple volumes

In this section we will make images for all the polar volumes you downloaded.

The Swedish case includes the transition from nocturnal migration into the start of daytime migration, when several cohorts of daytime migrants are departing. There is a lot to be seen in this fairly complex case. Sunrise was at 04:10 UTC, which is when the nocturnal migrants have stopped, and daytime migrants start taking off.

  1. Make ppi maps for the 2.5 degree polar scan for all the volume files you downloaded for this migration event, and save them to disk. If you feel confident with R - use below example code, which make use of a for-loop by which multiple images can be generated and saved automatically.
# let's store all the polar volumes we want to plot in a vector:
# first, get all the files with 'seang_pvol_20150927'
files=list.files(paste(HOME, "/20150927", sep=""), pattern="seang_pvol_20150927",full.names=TRUE)
# second, make sure we only use the files that end with 'h5', so hdf5 files only
files=files[grep("h5$",files)]
# load the ggplot2 library, because we will use its ggsave function later on to save images
library(ggplot2)
# we will use a for-loop to process the images one by one
# this simple loop prints the contents of variable 'i' iterating over a numeric vector
for (i in c(1,2,3)){
  print(paste("in this iteration variable i equals",i))
}
# let's regenerate the satellite image basemap
vol=read.pvol(files[1])
satelliteImage=basemap(ppi(vol$scans[[3]]),maptype="satellite")
# now let's make a more advanced for-loop that makes and save the images one by one:
for (file in files){
  # 1) read the polar volume 
  vol=read.pvol(file)
  # 2) extract the fifth scan (which is at 2 degree in this case)
  scan=vol$scans[[5]]
  # 3) make a ppi
  ppi = ppi(scan)
  # 4) plot the ppi, and store it in the object 'myPlot'.
  myPlot=map(ppi,param="DBZH",satelliteImage,zlim=c(-20,15))
  # 5) make a string with the filename in which to store the image
  myPlotFilename=paste(file,".jpeg",sep="")
  # 6) print the filename to the console
  print(paste("saving file",myPlotFilename,"..."))
  # 7) save the image
  ggsave(myPlotFilename,myPlot)
  # after ending up here, the code will be executed again from the start for the next image
}
  1. Cycle through the images and inspect where you see birds emerging. Not all daytime migrants depart synchronously. Explain what you see in the series of images in terms of the timing of departure and flight altitudes of the daytime migrants.

We will continue with this case on Thursday, and work on the quantification of bird densities, speeds and altitudes from these images.