When I started using R, one of the first things I wanted to do is make reliable colourplots. Colourplots are a way of representing 3D data as a 2D image, with the x and y coordinates representing variables and the z values represented as a colour.
For simple colorplots, the image() function from the base package, and the image.plot() function from the fields package are excellent choices. However, in my work, there is often a need to plot non-rectangular datasets, and this can become problematic.
The most common way around this, and the method used in graphing packages like Origin, is to interpolate a non-rectangular dataset to a rectangular grid. The result of this method can be seen in R too, using the disp.plot() function from my photonMonkey package.
For this example, I’ll be using the example dataset for a surface plasmon dispersion diagram available in the photonMonkey package.
library(photonMonkey)data(SPPdispersion) head(SPPdispersion)wavelength angle reflection 1 4.00e-07 7 0.4687 2 4.02e-07 7 0.4734 3 4.04e-07 7 0.4413 4 4.06e-07 7 0.4499 5 4.08e-07 7 0.4469 6 4.10e-07 7 0.4604
This is a standard rectangular dataset with reflection from a sample (in this case a zigzag grating) measured as a function of equally spaced steps of wavelength of the illuminating light and polar angle.
If we were to plot this as a non-rectangular dataset (such as a dispersion plot) using the method of interpolation to a rectangular grid, it would look like this:
angle <- SPPdispersion$angle
wavelength <- SPPdispersion$wavelength
reflection <- SPPdispersion$reflection
kx <- (2 * pi/wavelength) * sin(angle * pi/180)
omega <- 2 * pi * 3e+08/wavelength
disp.plot(kx, omega, reflection)
Ghastly.
- There is a horrible jagged edge to the dataset.
- There are missing values from poor interpolation.
- Overall resolution of the data has been lost due to interpolation.
The solution to all these problems is to plot the data as a non-rectangular dataset, with each element in the plot not being required to be the same shape, and certainly not requiring them to be square. We are going to fill the space using arbitrary trapezia, coloured based on the central value for that (x,y) coordinate.
The function we use for this is the disp.plot3() function from photonMonkey.
disp.plot3(x = wavelength, y = angle, z = reflection, fx = kx, fy = omega)
This is already far better. Absolutley no resolution from the dataset has been lost, and no values are missing. There is also no horrible jagged edge.
To use disp.plot3 the original and transformed datasets have to be specified. See ?disp.plot3 for more details.
Finally, we can take the best of both worlds and use a function that interpolates the rectangular dataset, and then transforms and plots using non-rectangular elements. The function for this is, again, in photonMonkey and is called disp.plot4. This function works a little differently, as now the transformed coordinates must be supplied as functions.
wavelength<-wavelength*1e9
kx<-function(x,y) (2*pi/(x*1e-9))*sin(y*pi/180)
omega<-function(x,y) 2*pi*3e8/(x*1e-9)
disp.plot4(x=wavelength,y=angle,z=reflection,
fx=kx,fy=omega,nx=300,ny=300)
Note, for the interpolation to work we needed to make the x and y scales similar. See ?disp.plot4 for details.
We can now use the utility functions in photonMonkey to tart this plot up nicely,
library(scales)
disp.plot4(x=wavelength,y=angle,z=reflection,
fx=kx,fy=omega,nx=300,ny=300,
xlab=label_in_plane_wavevector(),
ylab=label_angular_frequency(),
xaxt='n',yaxt='n')
lines(c(0,0),c(0,0)) # This is a silly work-a-round for a bug in R.
axis(1,at=axTicks(1),labels=label_scientific10(axTicks(1)))
axis(2,at=axTicks(2),labels=label_scientific10(axTicks(2)))
add.diffraction(kgx=2*pi/600e-9,kgv=2*pi/150e-9)
This isn’t just for dispersion plots, any functions can be used. disp.plot3 or 4 are also a great options for polar colourplots:
px<-function(x,y) x*sin(y*pi/180)
py<-function(x,y) x*cos(y*pi/180)
disp.plot4(x=wavelength,y=angle,z=reflection,
fx=px,fy=py,nx=300,ny=300,asp=1)
