The density is represented as a polyline, which is a pair of parallel arrays, one for $x$, one for $y$, forming vertices along the graph of the density (with equal spacings in the $x$ direction). As such it is a discrete approximation to the idealized continuous density and we can use discrete versions of the relevant integrals to compute statistics. Because the spacing is typically so close, there's probably little need to interpolate between successive points: we can use simple algorithms.
x <- seq(-2.5, 10, length=1000000)
hx5 <- rnorm(x,0,1) + rexp(x,1/5) # tau=5 (rate = 1/tau)
# Compute the density.
dens <- density(hx5)
# Compute some measures of location.
n <- length(dens$y) #$
dx <- mean(diff(dens$x)) # Typical spacing in x $
y.unit <- sum(dens$y) * dx # Check: this should integrate to 1 $
dx <- dx / y.unit # Make a minor adjustment
x.mean <- sum(dens$y * dens$x) * dx
y.mean <- dens$y[length(dens$x[dens$x < x.mean])] #$
x.mode <- dens$x[i.mode <- which.max(dens$y)]
y.mode <- dens$y[i.mode] #$
y.cs <- cumsum(dens$y) #$
x.med <- dens$x[i.med <- length(y.cs[2*y.cs <= y.cs[n]])] #$
y.med <- dens$y[i.med] #$
# Plot the density and the statistics.
plot(dens, xlim=c(-2.5,10), type="l", col="green",
xlab="x", main="ExGaussian curve",lwd=2)
temp <- mapply(function(x,y,c) lines(c(x,x), c(0,y), lwd=2, col=c),
c(x.mean, x.med, x.mode),
c(y.mean, y.med, y.mode),
c("Blue", "Gray", "Red"))