pathsim<-function(PathLen=200,SST=SST, filepath)
{
    #  PathLen = number of location observations in pathway, SST = ocean surface temperature field (matrix) [included in supplementary
    #           material as a separate file (SST), needs to be placed in working S+ dir, e.g., source("c:/dirpath/SST")]
    #           filepath = directory path, including filename for dataset, in quotes e.g., "c:/dirpath/rw200"
    #
    #  Model assumes that variation in animal movement (distances and direction) is a negative exponential function of sea surface temperature with no
    #  directional bias
    #
    #  sigeta = standard dev. of process errors (ie the variance in movement not explained by our biological model) 
    sigeta <- 1
    # obs.err = standard dev. of measurement errors (ie our instruments are imprecise)
    sigeps <- 2
    # beta = parameter to describe relationship between movement variance and water temperature
    beta <- 0.9
    # keeps track of sea surface temperatures
    tempC <- numeric(PathLen)
    # Initializations
    x.true <- numeric(PathLen)
    y.true <- numeric(PathLen)
    obs.error <- numeric(PathLen)
    #  Start animal out in cold water
    #  Will need to change these values for other random fields
    x.true[1] <- 20
    y.true[1] <- 35
    tempC[1] <- SST[x.true[1], y.true[1]]
    #  Iterate temperature-dependent movements
    i <- 2
    repeat {
        if(i > PathLen)
            break
            x.true[i] <- x.true[i - 1] + rnorm(1, 0, sigeta * exp( - beta * SST[x.true[i - 1], y.true[i - 1]]))
            y.true[i] <- y.true[i - 1] + rnorm(1, 0, sigeta * exp( - beta * SST[x.true[i - 1], y.true[i - 1]]))
            if(x.true[i] > dim(SST)[1] || x.true[i] < 1 || y.true[i] > dim(SST)[1] || y.true[i] < 1) {
                break
            }
            tempC[i] <- SST[x.true[i], y.true[i]]
            i <- i + 1
    }
    #  Assume observation errors constant over range of temperatures
    obs.error <- rnorm(PathLen, 0, sigeps)
    # No observation error in first observation since this is when animal released
    obs.error[1] <- 0
    #  Add observation error onto true locations to generate 'observed' locations
    x.obs <- x.true + obs.error
    y.obs <- y.true + obs.error
    #  Truncate pathway at SST boundary
    x.true[x.true < 0.51 | x.true > dim(SST)[1] + 0.5 | y.true < 0.51 | y.true > dim(SST)[1] + 0.5] <- NA
    y.true[y.true < 0.51 | y.true > dim(SST)[1] + 0.5 | is.na(x.true)] <- NA
    x.obs[is.na(x.true) | is.na(y.true)] <- NA
    y.obs[is.na(x.true) | is.na(y.true)] <- NA
    tempC[is.na(x.true) | is.na(y.true)] <- NA
    print(length(x.true[!is.na(x.true)]))
    #  Plot pathways on top of SST
    par(pty = "s")
    image(SST, xlim = c(0, dim(SST)[1]), ylim = c(0, dim(SST)[1]))
    lines(x.obs, y.obs, lty = 1, col = 5)
    lines(x.true, y.true, lty = 1, col = 1)
    points(x.true, y.true, lty = 1, col = 1, cex = 0.2, pch = 16)
    write.table(data.frame(x.obs, y.obs), file = filepath, sep = ",", append = F, quote.strings = F, dimnames.write = F, na = NA, end.of.row = ",", 
      justify.format = "decimal")
    write.table(data.frame(tempC), file = paste(filepath,".tmpC", sep = ""), sep = ",", append = F, quote.strings = F, dimnames.write = F, na = NA,
      end.of.row = ",", justify.format = "decimal")
    return(list(beta, x.true, y.true, x.obs, y.obs))
}