peaks <- function(series, span = 3, do.pad = TRUE) {
    if((span <- as.integer(span)) %% 2 != 1) stop("'span' must be odd")
    s1 <- 1:1 + (s <- span %/% 2)
    if(span == 1) return(rep.int(TRUE, length(series)))
    z <- embed(series, span)
    v <- apply(z[,s1] > z[, -s1, drop=FALSE], 1, all)
    if(do.pad) {
        pad <- rep.int(FALSE, s)
        c(pad, v, pad)
    } else v
}

peaksign <- function(series, span = 3, do.pad = TRUE)
{
    ## Purpose: return (-1 / 0 / 1) if series[i] is ( trough / "normal" / peak )
    ## ----------------------------------------------------------------------
    ## Author: Martin Maechler, Date: 25 Nov 2005

    if((span <- as.integer(span)) %% 2 != 1 || span == 1)
        stop("'span' must be odd and >= 3")
    s1 <- 1:1 + (s <- span %/% 2)
    z <- embed(series, span)
    d <- z[,s1] - z[, -s1, drop=FALSE]
    ans <- rep.int(0:0, nrow(d))
    ans[apply(d > 0, 1, all)] <- as.integer(1)
    ans[apply(d < 0, 1, all)] <- as.integer(-1)
    if(do.pad) {
        pad <- rep.int(0:0, s)
        c(pad, ans, pad)
    } else ans
}


check.pks <- function(y, span = 3)
    stopifnot(identical(peaks( y, span), peaksign(y, span) ==  1),
              identical(peaks(-y, span), peaksign(y, span) == -1))

for(y in list(1:10, rep(1,10), c(11,2,2,3,4,4,6,6,6))) {
    for(sp in c(3,5,7))
        check.pks(y, span = sp)
    stopifnot(peaksign(y) == 0)
}

## load some sample data
d <- read.csv("tioga1_35-65.csv", header=FALSE)

## name the columns
names(d) <- c('x', 'MG', 'MG+GLY', 'K', 'K 350', 'K 550')

## locate peaks in the 'K' signal
## the second argument is the "sensitivity" of the peak finding algorithm
d.peaks <- peaks(d$K, 35)

## save a vector of the positions in the K signal where the peaks were identified
peak_idx <- which(d.peaks)

## simple plot of raw K signal
plot(K ~ x, data=d, type="l", cex=.25, main="Tioga1 35-65cm\nK Treatment", xlab="Deg. 2 Theta", ylab="Intensity", ylim=c(0,max(d$K) + 50))

## add peaks: note that we are sub0setting the original data by the peak location index
points(K ~ x, data=d[peak_idx, ], col = 2, cex = 1.5)

## compute peak d-spacings
peak_d_spacings <- signif( (0.154/ (sin(d$x[peak_idx] * pi/180))), 3)

## annotate d-spacings / or peak index
## text(d$x[peak_idx], d$K[peak_idx] + 20, peak_d_spacings, col='blue', cex=0.75 )
text(d$x[peak_idx], d$K[peak_idx] + 20, 1:length(peak_d_spacings), col='blue', cex=0.75 )

## print a simple table of d-spacings by index
 data.frame(peak_d_spacings)
   peak_d_spacings
1            1.370
2            1.040
3            0.751
4            0.638
5            0.563
6            0.488
7            0.452
8            0.435
9            0.367
10           0.347