Y<- scan("C:/CLIM/HoPrWi.txt")
n<- length(Y)

#----------Morlet Wavelet Spectrum----------------------------
Wspectr<- function(Y,n,s,J,t0,om){
spec<- 1:J                                #spec vector of dim J
for (j in 1:J){ sumc<- 0; sums<- 0
for (i in 1:n){ tij<- (t0 - i)/s[j]
sumc<- sumc + Y[i]*cos(om*tij)*exp(-0.5*(tij^2)) 
sums<- sums + Y[i]*sin(om*tij)*exp(-0.5*(tij^2)) }
spec[j]<- (sumc^2 + sums^2)/s[j]}
return(spec)
}
#-------------------------------------------------------------
om<- 2*pi; dj<- 0.50; s0<- 2               #Wavelet parameters
J<- trunc((1/dj)*log2(n/s0)) 
jot<- 0:(J-1); s<-s0*2^(jot*dj)        #jot,s vectors of dim J
c("om"=om,"dj"=dj,"s0"=s0,"J"=J)
"Vector s of length J"; s 
       
sink("C:/CLIM/WaveOut.txt")    
options(digits=4)                                          ###
        
w<- 1:((n-1)*J); dim(w)<- c((n-1),J)  #w matrix of dim (n-1)xJ
for(k in 1:(n-1)){ spec<- Wspectr(Y,n,s,J,k,om)
w[k,]<- spec
write(w[k,],ncolumns=6,file="")
}  

rm(list=objects())                                         ###