Saturday, 26 September 2015

ProTrackR: chiptunes in R (part one)

ProTracker is a popular music tracker to sequence music on a Commodore Amiga machine. As I still own a Commodore Amiga 500, I thought it would be a nice challenge to develop an R package that could handle ProTracker module files, and so I did. The first version of this package was just released on CRAN.

In this post I will briefly explain what ProTracker modules are and how they work. Then I will give a short demonstrations of the current sweetness the ProTrackR package has to offer. On this blog I will try to keep you updated on new versions of this package and its new features.

The module (MOD) is a computer file format used primarily to represent music. A MOD file contains a set of instruments in the form of samples, a number of patterns indicating how (which notes and effects) and when the samples are to be played, and a list of what patterns to play in what order.

The first release of the ProTrackR package contains the bare bones required to import and export such module files, and modify them any way you'd like (and the file format restrictions allow you to). This is a good time to start with an example. Let's start by loading the module 'Elekfunk' from The Mod Archive.

# required packages:
require(plot3D)
require(ProTrackR)
# Read the module 'elekfunk' from the Mod Archive:
con <- url("http://api.modarchive.org/downloads.php?moduleid=41528", "rb")
mod <- read.module(con)
close(con)
rm(con)
# see how much memory the module object is eating:
format(object.size(mod), "Kb")
# calculate the ratio between the size of the object
# in memory and the size of the mod.file:
print(as.numeric(object.size(mod)/moduleSize(mod)))
# note that the object is roughly four times the size as the original file.
# This is probably because the Wave object (which the ProTrackerSample objec
# enherits), stores 8 bit wave data as integer values. In R, integers are
# actually 32-bit long, which explains the factor (32 / 8 =) 4 difference.

This module was used in a demo by Sanity to show off the Amiga's capabilities, which were impressive ate the time.

Playing routines in the ProTrackR package are currently limited to just playing the samples at specific notes (more fancy routines are in development). The examples below show how to play all samples in a module, a specific sample at different notes and applying different fine tunes. But at least you can export them as wave files.

# play all samples in the module at note "C-3"
playSample(mod)
# play chromatic scales of octaves 2 and 3 using the 23rd sample
data(period_table)
for (i_note in as.vector(outer(names(period_table)[-1:-2],
c(2:3), function(x, y) paste(x, y, sep = ""))))
{
print(i_note)
playSample(PTSample(mod, 23), 0.1, T, i_note, 0.5)
}
# play all finetunes of C-3
for (ft in -8:7)
{
print(ft)
playSample(PTSample(mod, 23), 0.1, T, "C-3", 0.5, finetune = ft)
}
# export one of the samples as *.wav file
writeWave(PTSample(mod, 25), "guitar slide.wav")

Combined with the power of R, all kinds of fancy analyses become possible. The module for instance contains a sample of a guitar slide. In R we can now generate a power spectrum (which would have been tricky to do real time on an original chipset Amiga).

# select a sample:
plot_samp <- PTSample(mod, 25)
# generate power spectrum
spec <- powspec(waveform(plot_samp),
wintime = 0.1,
steptime = 0.005)
spec_cutoff <- quantile(log10(spec), 0.01)
spec[log10(spec) < spec_cutoff] <- 10^spec_cutoff
# fancy 3D plot of waveform and powerspectrum
png("powerspec.png", 600, 300, type = "cairo",
bg = "transparent", pointsize = 14)
par(oma = c(0,0,0,0), mar = c(0.5,0.5,0.5,0.5))
lines3D(x = 60*(0:(length(plot_samp) - 1))/(length(plot_samp) - 1),
y = rep(15, length(plot_samp)),
z = 3 + 11*waveform(plot_samp)/255,
phi = 40, theta = 21, lwd = 0.15,
xlim = c(0,60), ylim = c(0,15), zlim = c(-3, 14.5),
xlab = "time", ylab = "Frequency", zlab = "Amplitude",
col = "black", lwd = 0.5, scale = F, bty = "b2")
persp3D(x = 60*((1:ncol(spec)) - 1)/(ncol(spec) - 1),
y = 15*((1:(nrow(spec) - 3)) - 1)/(nrow(spec) - 4),
z = t(log10(spec))[,-(1:3)],
xlab = "time", ylab = "Frequency", zlab = "Amplitude",
add = T, colkey = F, shade = 0.1)
dev.off()

The power spectrum shows the shift in frequency which you would expect for a guitar slide. There are plenty other examples I could show you here but leave it at this for the moment. Hopefully you enjoyed these examples of the ProTrackR package, please stay tuned for updates on this package.

Saturday, 5 September 2015

Geotags from JPEG images

Most modern camera's and telephones can geotag their pictures. Kind of neat that you can use your camera as GPS logger. Even better if we could do something useful with it in R. In this post I will show how GPS coordinates can be extracted from jpeg images and utilized in R. Although it sounds pretty simple, and it is, it is not very straight forward to achieve.

In this example I will use pictures I made on a day out. We made a small round trip by steamtram from Hoorn to Medemblik where we visited castle Radboud. I will show how to plot the locations where I took the pictures can be plotted onto a map. For privacy reasons I will not share the pictures. However, the example should also work with your own geotagged pictures.

# required packages
require(rgdal)
require(OpenStreetMap)
# select files for which geotags need to
# be extracted and store filenames in data.frame
file_info <- data.frame(file.name = choose.files(filters = Filters["jpeg",]))
view raw geotag01.r hosted with ❤ by GitHub

Geotags and other meta-information is stored in jpeg files in the so-called exchangeable image file format (or EXIF). Luckily it can be easily obtained using the function GDALinfo. This is where things get messy. The function returns the EXIF codes as character strings in the attributes of the object that is returned. So, we need to fetch the relevant attributes; find the relevant EXIF codes; and extract the information from the code and convert it into something workable.

The EXIF codes for GPS information will look something like this ‘EXIF_GPSLatitude=(52) (40) (57.45)’. I use regular expressions to extract the coordinate information from between the brackets. I'm no expert in the field of regular expressions, so there may be more efficient regular expressions to get the right information. If you have improvements, please feel free to leave them in the comments.

I've wrapped the routines to isolate the respective EXIF code and to extract the numerical coordinates from the literal EXIF character string in slick functions (‘get_EXIF_value’ for getting the literal string and ‘get_EXIF_coordinate’ to get the numerical coordinate from that string).

# function to get a parameter value from EXIF code,
# x = a vector of character strings representing the EXIF codes of a file
# parameter = character string representing the parameter
# for which you want to extract the value
# the literal character string is returned for the requested parameter
get_EXIF_value <- function(x, parameter)
{
line <- x[grepl(paste("EXIF_", parameter, "=", sep = ""), x)]
if (!is.null(line) && length(line) > 0)
{
pos <- gregexpr("=", line, fixed = T)[[1]]
pos <- pos[length(pos)]
line <- substr(line, pos + 1, nchar(line))
line <- gsub("[(]|[)]", "", unlist(strsplit(line, "([)]\\s+[(])|[)][()]")))
}
return(line)
}
# function to get the longitude and latitude from the
# character string obtained with the function above
# x = a character string with the coordinates
get_EXIF_coordinate <- function(x)
{
lat <- get_EXIF_value(x, "GPSLatitude")
lat <- as.numeric(lat)
lat <- lat[1] + lat[2]/60 + lat[3]/(60*60)
latref <- get_EXIF_value(x, "GPSLatitudeRef")
lat <- ifelse(toupper(latref) == "S", -lat, lat)
if (length(lat) == 0) lat <- NA
lon <- get_EXIF_value(x, "GPSLongitude")
lon <- as.numeric(lon)
lon <- lon[1] + lon[2]/60 + lon[3]/(60*60)
lonref <- get_EXIF_value(x, "GPSLongitudeRef")
lon <- ifelse(toupper(lonref) == "W", -lon, lon)
if (length(lon) == 0) lon <- NA
return (c(lon, lat))
}
# for each file obtain the EXIF codes and extract the coordinates
GPStags <- apply(array(file_info$file.name), 1, function(x) {
EXIF_code <- NULL
try (EXIF_code <- attr(GDALinfo(x), "mdata"), silent = T)
if (is.null(EXIF_code)) return(c(NA, NA)) else return(get_EXIF_coordinate(EXIF_code))
})
view raw geotag02.r hosted with ❤ by GitHub

After some rearranging of the data, we have a data.frame holding the filename and the longitude and latitude in decimal WGS84 degrees. I noticed that not all my pictures had a geotag (even though I had geotagging turned on, on my device). Apparently taking pictures is not a full proof method for logging locations (at least not on my device).

# rearrang the data
GPStags <- t(matrix(unlist(GPStags), 2))
# store in the data.frame
file_info$lon <- GPStags[,1]
file_info$lat <- GPStags[,2]
# check if any files lack GPS longitudes
table(is.na(file_info$lon))
# clean up after us
rm(GPStags)
view raw geotag03.r hosted with ❤ by GitHub

Let's plot the coordinates on a map. I download some OpenStreetMap-based map tiles using the function ‘openmap’. This function conveniently allows you to specify the upperleft and lowerright coordinates of your map. That way, you can be certain that your locations will fit on the map.

As we cannot change the projection of the downloaded maptiles. We therefore have to transform the coordinates of our locations to the same projection as the maptiles. We do this by first creating a SpatialPoints object from our coordinates and then transform those coordinates to the proper projection by applying spTransform. Plotting is now simply a matter of calling some plot routines.

# Download nice-looking map-tiles for the relevant location
photo_map <- openmap(upperLeft = c(max(file_info$lat, na.rm = T) + 0.2, min(file_info$lon, na.rm = T) - 0.2),
lowerRight = c(min(file_info$lat, na.rm = T) - 0.2, max(file_info$lon, na.rm = T) + 0.2),
type = "stamen-watercolor")
# we need to turn our coordinates into an sp object in order to
# properly project the points onto the map.
photo_points <- SpatialPoints(na.omit(file_info[,c("lon", "lat")]),
CRS("+proj=longlat +ellps=WGS84 +datum=WGS84"))
# transform the lon lat coordinates to the same projection as the map tiles.
photo_points <- spTransform(photo_points, photo_map$tiles[[1]]$projection)
# plot the map and locations of the pictures
png("jpg_geotags01.png", photo_map$tiles[[1]]$yres/3, photo_map$tiles[[1]]$xres/3, type = "cairo")
plot(photo_map)
plot(photo_points, add = T, pch = 19, cex = 0.8)
plot(photo_points, add = T, pch = 19, cex = 0.5, col = rainbow(length(photo_points)))
dev.off()
view raw geotag04.r hosted with ❤ by GitHub

To make things a bit more fancy, I used the rainbow palette to color the locations in the order in which the pictures were taken. Note that the use of geotagged images is limited. Proper GPS loggers for instance also log the precision of the position and your speed. But I will save that for a future post.