##MANN's CO2DETREND.f EMULATION
	
#Mann's program: see http://www/climateaudit.org/data/mbh99/COMPARE/co2detrend.f
#MAnn comments:

# [8] "c      regress out co2-correlated trend (r=0.9 w/ co2)"        
# [9] "c      after 1800 from pc1 of ITRDB data"              

# [58] "c      remove co2-correlated portion (r=0.9) of 1800-1980"     
# [59] "c "                                                            
# [60] "       corr= 0.9"                                              

	# open (unit=3,file='pc01.out',status='old')
       
	# open (unit=4,file='ghg-lowf.dat',status='old')
       
	# open (unit=7,file='pc01-fixed.dat',status='unknown')
       do i=1,nmax
          read (3,*,end=99) iyear,pc1(iyear)
       end do
99     nx = i-1
       do i=1,nmax
          read (4,*,end=199) adum,val2
    
	
	ser26= read.table("http://www.climateaudit.org/data/mbh99/COMPARE/pc01.out") # mannomatic pc1 1000-1980
	ser26=ts(ser26[,2:ncol(ser26)],start=ser26[1,1])
	ser12= read.table("http://www.climateaudit.org/data/mbh99/COMPARE/ghg-lowf.dat") # CO2 values less 1610-1995 average
	ser12=ts(ser12[,2],start=1610) #this is the starting date based on other info

	pc1= -ser26; #invert orienation
	m1=mean(pc1[(1800:1980)-999]);sd1=sd(pc1[(1800:1980)-999]); c(m1,sd1) #[1] -0.01463  0.01858
	standard1= (pc1-m1)/sd1 #1000 1980 

	m2=mean(co2[(1800:1980)-1609]);sd2=sd(co2[(1800:1980)-1609]); c(m2,sd2) #[1]   7.219 14.059
	standard2= (co2-m2)/sd2 #1610 1980

	newpc=pc1; corr=.9
	newpc[(1610:1980)-999]=sd1*( standard1[(1610:1980)-999]- corr*standard2) + m1
	m2= mean(standard2[1:200]); m2 # -1.31
	newpc[(1000:1609)-999]=sd1*( standard1[(1000:1609)-999]- corr*m2) + m1
	

	fm=lm(ser34~pc1);coef(fm) #  -0.03043     0.06151  # correlation - 1
	test=coef(fm)[1]+coef(fm)[2]*newpc   #MAnnian rescaling emulation

	X=ts.union(test,pc1f[,1]); dimnames(X)[[2]]=c("CO2-detrended","Mann-adjusted")
	
	nf=layout(array(1:2,dim=c(2,1)),heights=c(1.1,1.3))
	par(mar=c(0,3,2,1))
	plot(1000:1980,X[,1],type="l",xlab="",ylab="",axes=FALSE)
	axis(side=1,labels=FALSE);axis(side=2,las=1);box();abline(h=0,lty=2)
	title("Rescaled Mannomatic PC1 with Adjustments")	
	text(1000,.8,font=2,pos=4,"co2detrend.f")
	par(mar=c(3,3,0,1))
	plot(1000:1980,X[,2],type="l",xlab="",ylab="",axes=FALSE)
	axis(side=1);axis(side=2,las=1);box();abline(h=0,lty=2)
	text(1000,.65,font=2,pos=4,"Jacoby coercion")