##MAKE WEIGHTS MAP

library(fields)
library(signal)
library("R.matlab");library("Rcompression")
library(xlsReadWrite)
 	source("d:/climate/scripts/mann.2008/utilities.v6.txt")
  	period_MBH=seq(1800,200,-100);M=length(period_MBH);M
     		#source("http://www.climateaudit.info/scripts/mann.2008/collation.cps.txt") 	
 	use0="pairwise.complete.obs"	
	 #Basics$details[899,c("lat","long")]=Basics$details[899,c("long","lat")] 
	notilj =  is.na(match( 1:1209,1061:1064) )

	K=17
	name0=c("latem","earlym","whole")
	emulation=rep(list(NA),K); names(emulation)=period_MBH
	criterion=notilj&noluter##;sum(criterion) #1134
	## criterion=notilj&noluter&nobristle
	   #used to get NOtilj_nobristle

#emulation
	for (k in 1:17) {
	emulation[[k]]=rep(list(NA),3); names(emulation[[k]])=name0
	for (j in 1:3) {
	    emulation[[k]][[j]]=manniancps(k,case=name0[j],criterion,verbose="verbose")  #verbose gives interim steps for benchmark
	   }#j
	}
	length(emulation) #17

#cps recons and stats for 4 cases 
	name0=c("latem","earlym","whole")
	target_name=c("CRU","iCRU","HAD","iHAD")
	hemi_name=c("NH","SH")
	name1=c(outer(target_name,hemi_name,function(x,y) paste(x,y,sep="_") ) )
	Target= data.frame(id= name1,target=rep(target_name,2), hemi=as.character(gl(2,4,labels=hemi_name)) )
	cps=rep(list(NA),4); names(cps)=target_name
	for (n in 1:4) {
	  cps[[n]]=rep(list(NA),17)
	  for(k in 1:17) {
	    cps[[n]][[k]]= rep(list(NA),3); names(cps[[n]][[k]])=name0
	      for (j in 1:3) {
	        x=Basics$cal[[ name0[j] ]]
	        cps[[n]][[k]][[ name0[j] ]]= cpsf( emulation[[k]][[ name0[j] ]]$recon,target=Target$target[n], index= x[1]:x[2])
	  }
	} 
	}

	Stat=rep(list(NA),8); names(Stat)=Target$id
	name1=c("century", c(t( outer(c("earlym","latem","avg","adj"),c("RE","CE","r2"), function(x,y) paste(x,y,sep="_") ) )) )[1:12]
	for (j in 1:8) {
	 Stat[[j]]=statf(emulationx=cps[[Target$target[j] ]], target_id=Target$target[j], hemi=Target$hemi[j])
  	}


Notilj=list(emulation=emulation,cps=cps,stat=Stat)
save(Notilj,file="d:/climate/scripts/mann.2008/notes/Notilj.tab")


#Notilj_nobristle=list(emulation=emulation,cps=cps,stat=Stat)
#save(Notilj_nobristle,file="d:/climate/scripts/mann.2008/notes/NOtilj_nobristle.tab")