#SOME UTILITIES 
##JONES GRIDCELL LOOKUP
 jones<-function(lat,long) {
	i<-18-floor(lat/5);
	j<- 36+ceiling(long/5);
	jones<-72*(i-1) + j
	jones
	}

##jones.inv
   jones.inv<-function(k) {
	j<- k%%72
	temp<-k-j
	i<-(temp/72)+1
	long<-(j-36)*5-2.5
	lat<-(18-i)*5+2.5
	jones.inv<-c(lat,long)
	jones.inv
	}

  ts.annavg=function(X,M=6) { 
        if(class(X)=="ts") {start0<-floor(tsp(X)[1]);X<-ts.array(X)} else 
        if( ((class(X)=="matrix")|(class(X)=="data.frame") ) &(ncol(X)==12)) start0<-as.numeric(row.names(X)[1])  else {X
                start0<-X[1,1]; X<-X[,2:13]}
        count<-apply(!is.na(X),1,sum);temp<-(count<M)
        ts.annavg<-apply(X,1,mean,na.rm=TRUE)
        ts.annavg[temp]<-NA
        ts.annavg<-ts(as.matrix(ts.annavg),start=start0)
        ts.annavg
	}
   
	
     fill.array<-function(X,method="column") {
	if (method=="row.name") {xyear<-as.numeric(row.names(X));X<-X[,c(13,1:12)] }
	if(class(X)=="matrix") {xyear<-X[,1];X<-X[,2:13]} 
	if (class(X)=="data.frame") {xyear<-X$year;X<-as.matrix(X[,2:13])}
	years<-min(xyear,na.rm=TRUE):max(xyear,na.rm=TRUE)
	fill.array<-array(NA,dim=c(length(years),12))
	temp<-!is.na( match(years,xyear) )
	fill.array[temp,]<-X
	row.names(fill.array)<-years
	dimnames(fill.array)[[2]]<-month0
	fill.array
	}

	ts.array<-function(x) {
	year1<- floor(tsp(x)[1]):floor(tsp(x)[2])
	month1<- round(12*(tsp(x)[1:2]%%1),0)
	x<-c(rep(NA,month1[1]),x ,rep(NA, 11-month1[2]))
	ts.array<-t(array(x,dim=c(12,length(x)/12)) )
	row.names(ts.array)<-year1
	ts.array
	}


##annavg
	#CALCULATE ANNUAL AVERAGE FROM MONTHLY TIME SERIES
	#even number of months
   annavg<-function(x) {
	n<-length(x)
	m<-n-12*floor(n/12)
	if(m>0) x<-c(x, rep(NA,12-m))
	years<-length(x)/12
	x<-array(x,dim=c(12,years))
	annavg<-apply(x,2,mean,na.rm=TRUE)
	return(annavg)
	}

##strip
	#FUNCTION TO LISTFILES AND REMOVE SUFFIXES
   strip<-function(h,n) {
	h<-unlist(strsplit(h, "\\."))
	h<-array(h,dim=c(n,length(h)/n))
	h<-t(h)
	strip<-h[,1]
	strip
	}

##split
   split<-function(h,char) {
	y<-unlist(strsplit(h,paste("\\",char,sep="")))
	y<-array(y,dim=c(2,length(y)/2))
	split<-t(y)
	split<-split[,1]
	split
	}




filter.combine.pad=function(x,a,M=NA,method="norm") {
	temp<-!is.na(x)  
   	w<-x[temp]
	N<-length(w)
	if(is.na(M)) M<-trunc (length(a)/2)+1
	if(method=="norm")  w<-c (rep(mean(w[1:M]),M),w,rep(mean(w[(N-M+1):N]),M)  )
	#if(method=="reflect") w<-c (  w[M:1],w,w[M:(N-M+1):N])
	if(method=="mannomatic") {w<-c( c (w[1]-(w[(M+1):2]-w[1])),w, c(w[N]- (w[(N-1):(N-M)]-w[N]) ) ) }
	  if(method=="mann2")  w<-c (rep(w[1],M),w,rep(w[N],M)  )
	   y<-filter(w,a,method="convolution",sides=2)
	y<-y[(M+1):(length(y)-M)]
	z<- rep(NA,length(x))
	z[temp]<-y
	if(class(x)=="ts") z<-ts(z,start=tsp(x)[1],freq=tsp(x)[3])
	   filter.combine.pad<-cbind(x,z)
	   #filter.combine<-y[,2]
	filter.combine.pad
	   }

##gaussian.filter.weights
   gaussian.filter.weights<-function(N,year) {
		#N number of points; year - is number of years
	sigma<-year/6
	i<-( (-(N-1)/2):((N-1)/2)   )  /sigma
	w<- ((2*pi)^-0.5) * exp(-0.5* i^2)
	gaussian.filter.weights<-w/sum(w)
	gaussian.filter.weights
	}

	#GAUSSIAN FILTER
	#see http://www.ltrr.arizona.edu/~dmeko/notes_8.pdf
	#The appropriate Gaussian distribution therefore has standard deviation
	# s = ? / 6,
	# where ? is the desired wavelength at which the amplitude of frequency response is 0.5. The
	# filter weights are obtained by sampling the pdf of the standard normal distribution at
	# values  i/s  where i={-(N-1)/s, ..., 0, .. (N-1)/s }
	# MEKo has values to verify for case of N=21, year=10 #this is verified

##truncated.gauss.weights
   truncated.gauss.weights<-function(year) {
	a<-gaussian.filter.weights(2*year+1,year)
	temp<-(a>0.05*max(a))
	a<-a[temp]
	a<-a/sum(a)
	truncated.gauss.weights<-a
	truncated.gauss.weights
	}

 maketable=function(X,align=-999,widths=-999){
	N=nrow(X);M=ncol(X); out= "d:/temp/temp.txt"
	if (align[1]== -999) right=rep(TRUE,M) else right=align
	code=rep("",M);code[right]=" <p align=&quotright&quot>"
	if(widths[1]== -999) pct= floor(100/M) else pct=widths

	write.table("<table border=&quot1&quot style=&quotborder-collapse: collapse&quot bordercolor=&quot#111111&quot width=&quot80%&quot id=&quotAutoNumber1&quot>",file=out,col.names=FALSE,row.names=FALSE,quote=FALSE)
   	for (j in 1:N){
	 write.table("<tr>",append=TRUE,file=out,row.names=FALSE,col.names=FALSE,quote=FALSE)
	  for(i in 1:M){
	    write.table( paste("<td width=&quot",pct[i],"%&quot>",code[i],X[j,i],"</td>",sep=""),append=TRUE,file=out,row.names=FALSE,col.names=FALSE,quote=FALSE)
   	    }
	 write.table("</tr>",append=TRUE,file=out,row.names=FALSE,col.names=FALSE,quote=FALSE)
	 }
	write.table("</table>", append=TRUE,file=out,row.names=FALSE,col.names=FALSE,quote=FALSE)
     }
	
#TO EXTEND SERIES BY PERSISTENCE
extend.persist<-function(tree) {
	extend.persist<-tree
	for (j in 1:ncol(tree) )  {
		test<-is.na(tree[,j])
		end1<-max ( c(1:nrow(tree)) [!test])
		test2<-( c(1:nrow(tree))>end1) & test
		extend.persist[test2,j]<-tree[end1,j]
		}
	extend.persist
	}