covfield <- function(array3d,ts, ...) {

# Compute covariance between each point of a three-dimensional array, with 
# first two space dimensions and third time dimension, and a given vector of
# the same length as the time dimension of the three-dimentional array.
#
# Description:
#
#      Returns a space matrix with the covariances 
#
# Usage:
#
#      covfield(array3d,ts, ...)
#
# Input:
#
#      array3d: three-dimensional array with first two space dimensions
#               (e.g. longitude and latitude) and third time dimension 
#
#      ts: a vector containing the time series to be covariate with the
#          three-dimensional array
#
#     ...: Additional argument passed to the covariance function 
#          (e.g. method = "pearson", method = "kendall", method = "spearman",
#           use = "all.obs", use = "complete.obs", use = "pairwise.complete.obs")
#          Default computes pearson covariance (method = "pearson")
#          using all data available (use = "all.obs")
#
# Output:
#
#      Matrix of covariances
#
# Authors:
#
#      Dag Johan Steinskog <dag.johan.steinskog@nersc.no> 25 October 2005
#      Caio Coelho <c.a.d.s.coelho@reading.ac.uk> 
#
# Examples:
#
#      x <- seq(-20, 20, 5)
#      y <- seq(30, 60, 5)
#      dim <- c(length(x), length(y), 100)
#      data <- array(rnorm(prod(dim)), dim)
#      timeseries<-rnorm(100)
#      covfield(data,timeseries)
#      covfield(data,timeseries,method = "pearson")
#      covfield(data,timeseries,method = "kendall")
#      covfield(data,timeseries,method = "spearman")


covar <- function(z) {cov(z,ts,use="complete.obs", ...)}
apply(array3d, c(1, 2), covar)


}