inst/main.R

In martindurocher/emdCwa: Hybrid method for analyzing teleconnection

#!/usr/bin/Rscript

source('../R/libEmdCwa.R')
library(EMD)
library(biwavelet)

library(SpatialExtremes)
library(sqldf)

arg <- list()
arg$naoFileName <- '../data/nao.Rdata'
arg$pImf <- 4
arg$kNoise <- 2
arg$cacheFile <- 'cache.Rdata'
arg$useCache <- TRUE

main <- function(arg)
{	
	load(arg$naoFileName)
	data(rainfall)
	
	cmd <- "select year, avg(value) as value from nao
		where (year between 1962 and 2008)
			and (month in (6,7,8))
		group by year"
			
	naoSummer <- sqldf(cmd)
	yvec <- naoSummer[,1]
	
	mydata <- newEmdCwa(naoSummer[,1],
		index = scale(naoSummer[,2]),
		site = scale(log(rain)))
	
	if(!arg$useCache)
	{
		mydata <- eemd(mydata, p = arg$pImf, k = arg$kNoise, max.sift = 5)
		save(mydata, file = arg$cacheFile)
	}
	else
		load(arg$cacheFile)
	
	overlayPlot(mydata, m = 1, k = 1, type = 'l', col = 2)

	plot(test0 <- cwt(mydata,m=1,1,1),type = 'power.corr.norm',
		plot.cb = TRUE)

	#test <- multiCwt(mydata, 1)
	#print(test$sigma)
	 print(names(test0))
	 print(test0$d1.sigma)
	
	plotWt(log2(c(test0$power.corr)),test0$t,test0$period,
		test0$d1.sigma,test0$d2.sigma)

	
}

pdf()
par(mfrow = c(2,2))
main(arg)
dev.off()