# MBH98 DATE CONFIGURATION

# Directories where region data is stored
	mbhdatatree<-"c:/climate/data/new/TREE"   #tree data storage
	mbhdata<-"c:/climate/data/mann" #directory where Mann data is stored
	pcstore<-"c:/climate/data/pcs"# Directory were PC will be stored during this collation

##TEMPERATURE PRINCIPAL COMPONENT INFORMATION USED IN PROXY CALIBRATION AND RECONSTRUCTION
  # temperature principal components
	load(  file.path(mbhdata,"pc.tab"))

  # temperature eofs used in reconstruction
	load(  file.path(mbhdata,"eof.ordermann.tab"))

  # temperature eigen-values and choose first 16 as only 16 PCs are used
	load(file.path(mbhdata,"tpca-eigenvals.tab"))
	lambda<-lambda[1:16,2]

# STEPWISE RECONSTRUCTION PARAMETERS
  #there are 11 steps and different temperature PCs used in each step
	period.m<-c(1981,1820,1800,1780,1760,1750,1730,1700,1600,1500,1450,1400) #sub-interval boundary 
	period<-period.m-1399
	#eigenvectors selected in each sub-interval- one less than "period"
	select<-list(c(1:5,7,9,11,14:16),c(1:5,7,9,11,14:16),c(1:5,7,9,11,14:16),
		c(1:5,7,9,11,15),c(1:3,5,7,9,11,15),c(1,2,5,11,15),c(1,2,5,11,15),
		c(1,2,11,15),c(1,2),c(1,2),c(1))

##LOAD PROXY
   #this is a collation of the non-tree ring PC proxies 
   #columns are preserved for 31 non-tree ring proxies from columns 68-99	
   #for each reconstruction step, the specified tree ring PCs will be pasted in 
	load("c:/climate/data/mann/proxy.tab")
	proxy.MBH<-ts( proxy[1:581,],start=1400,end=1980)

   # weights used in proxy regression
	load( file.path(mbhdata,"weights1.tab"))

##TREE RING PC INFORMATION
   # there are 6 different tree ring regions with varying number of PCs used and varying calls to subdirectories

   # Names of tree ring regions
	region<-c("STAHLE/OK","STAHLE/SWM","ITRDB/NOAMER","ITRDB/SOAMER","ITRDB/AUSTRAL","VAGANOV")
	regionid<-c("STAHLE.OK","STAHLE.SWM","NOAMER","SOAMER","AUSTRAL","VAGANOV") # for storage

   # Parameters for tree ring PC collection
	mbh98.backto<-list(c(1700),c(1400,1450,1500,1600,1700),c(1400,1450,1600,1750),c(1600),c(1600,1750),c(1450,1600,1750)) #Labels on PC files
	mannpc<-c(3,9,9,3,4,3) # Number of PC pc*.out files to read
	mannid<-c(69,72,84,93,96,81) # Location of first column of regional PC in proxy collation

   # Schedules for the number of PCs by calculation step and region; and call subdirectories by calculation step and region
	series.MBH<-read.table ( file.path(mbhdata,"treepc.byregion.count.SI.2004.txt"),sep="\t",header=TRUE) #table of retained PCs by period and region
	directory.MBH<-read.table( file.path(mbhdata,"treepc.byregion.directory.SI.2004.txt"),sep="\t",header=TRUE) #table of look-up directories for above
	series.MBH<-series.MBH[2:ncol(series.MBH)]
	directory.MBH<-directory.MBH[2:ncol(directory.MBH)]

##GRIDCELL INFORMATION
    # MBH uses 1082 gridcells in their "dense" reconstruction and 219 gridcells in their "sparse" reconstruction specified here
    # load list of MBH "dense" grid-cells
	load(file.path(mbhdata,"gridpoints.tab")) #this table includes the "Jones number" of the gridcell
	tempmann<-sort(grid2[,5])
	manngrid<-!is.na(match(c(1:2592),tempmann)) 
	ordermann<-order(grid2[,5])
	grid2<-grid2[ordermann,]
	nhtemp<-(grid2[,4]>0)   #logical vector for NH gridcells
	mu<-cos(grid2[,4]*pi/180)  #used for areal weighting

    # allow for gridcells with 0 values #this is no longer used after July 2004 Corrigendum information 
	#load(file.path(mbhdata,"nil.tab"))#load list of 4 grid-cells identified as having no values. See read.jones for discussion.
	nil<-rep(FALSE,1082) #this uses new SI information (May 2004)
	nhtemp.nil<-nhtemp&!nil

    # load list of grid-cell standard deviations # truncate to exclude nil cells (allowed for but not applicable here) 
	load(file.path(mbhdata,"gridcell.stdev.ordermann.tab")) #R version of archived information
	mann.scale<-gridcell.stdev.ordermann #no nil values
	nhscale<-1/(mann.scale[nhtemp.nil]*sum(nhtemp.nil)) #
	nhscale<-array(  nhscale,  dim= c(length(nhscale),1) ) #dim 

    # mask for sparse subset of 219 cells #mask nowhere specified but determiend by careful inspection of diagram
	load("c:/climate/data/mann/mask.edited.tab") #defined in toolbox/notes.2 #list of 219 jones cells
	mask.logical<-!is.na(match(grid2[,5],mask.edited)) #sum(mask.logical) #219 #length(mask.logical)
	mask<-list(mask.edited,mask.logical)

## TEMPERATURE AND RECONSTRUCTION DATA
	load(file.path(mbhdata,"nhmann.tab")) #this is R version of archived information
	NH.scale.factor<-sd(nhmann[503:581,3])
	NH.MBH<-ts(nhmann[,2],start=1400)

     # archived "sparse" data  	#no information at NOAA on zero but dense is zero on 1902-1980 by experiment
	url<-"ftp://ftp.ngdc.noaa.gov/paleo/paleocean/by_contributor/mann1998/nhem-sparse.dat"
	g<-read.table(url,skip=1)
	sparse<-ts(g[,2:3],start=g[1,1],end=g[nrow(g),1])
	dimnames(sparse)[[2]]<-c("instr","recon")
	sparse<-sparse[,1]-mean(sparse[(1902-1853):(1980-1853),1])


##VERIFICATION TESTS
   	#these are done to compare against diagrams to make sure spacings are correct etc.
	#plot.ts(mann.scale) #this is OK

	#library(fields)

	##PLOTTING FUNCTION
	#	plotworld <-function (h,title,col2)   {
	#	x<-5*(-35:36);
	#	y1<-5*(-17:18);
	#	h<-array(h,dim=c(72,36));
	#	h<-h[,36:1];
	#	image.plot(x,y1,h,main=paste(title),col=col2,xlab="LONGITUDE", ylab="LATITUDE");
	#	world(xlim=c(-87.5,87.5),ylim=c(-180,180),add=TRUE)
	#	}

	#	h<-rep(NA,2592)
	#	h[manngrid]<-eof[,1]
	#	col9<-c("grey30",rep("grey80",2),rep("yellow",5),rep("red",8))
	#	plotworld(h,"",col9) #this is correct

	#	h<-rep(NA,2592)
	#	h[manngrid][nhtemp]<-rnorm(707)
	#	col9<-c("grey30",rep("grey80",2),rep("yellow",5),rep("red",8))
	#	plotworld(h,"",col9) #this is screwed up

	#	h<-rep(NA,2592)
	#	h[manngrid]<-gridcell.stdev.ordermann
	#	col9<-c("grey30",rep("grey80",2),rep("yellow",5),rep("red",8))
	#	plotworld(h,"",col9) #this is correct