##Retreive AMSU data and make monthlies


#http://discover.itsc.uah.edu/amsutemps/data/amsu_daily_85N85S_ch04.r000.txt
#http://discover.itsc.uah.edu/amsutemps/data/amsu_daily_85N85S_ch07.r004.txt
#near surface - channel 4 

#850 mb values end in March 2008 

r=data.frame(channel=c(0,4:14),r=c(0,0,2:11),mb=c(1000,850,600,400,250,150,90,50,25,10,5,2.5),
	km=c(0,2,4.3,7.5,11,14,17,21,25,31,36,41) ) 

numbertochar=function(x) if (nchar(x)==1) paste("0",x,sep="") else x

f=function(day) { f=sapply(amsu, function(A) A[day,"2010"])
	dimnames(f)[[2]]=r$mb;f}

get.amsu=function(channel){
	loc=paste("http://discover.itsc.uah.edu/amsutemps/data/amsu_daily_85N85S_ch",numbertochar(channel),
	  ".r00",r$r[match(channel,r$channel)],".txt",sep="")
	if(channel>=13) loc=gsub("r00","r0",loc) 
	if(channel=="LT") loc="http://discover.itsc.uah.edu/amsutemps/data/amsu_daily_85N85S_chLT.r001.txt"
	test=readLines(loc)
	fred=test[(length(test)-365):length(test)]
	fred=gsub("\\t"," ",fred)
	fred=gsub("/"," ",fred)
	writeLines(fred,"temp")
	x=scan("temp")
	x=t(array(x,dim=c(length(x)/366,366) ) )
	K=ncol(x)
	dimnames(x)[[2]]=c("month","day",paste( (2011-(K-3)) :2011) )
	# x=x[,1:(ncol(x)-4)] #rem max, min, mean
	x[x==-999]=NA
	x[x==0]=NA
	dimnames(x)[[1]]=1:366
	#x=x[,3:ncol(x)]
	return(x)
} 

makef=function() {
 amsu=list()
 amsu[["LT"]]=get.amsu("LT")
 for (channel in 4:14) amsu[[channel-2]]=get.amsu(channel)
 names(amsu)=r$mb
 return(amsu)
}

make_anom=function(A) {
	K=ncol(A)  
        X= data.frame(year= gl(K-2,366,labels=dimnames(A)[[2]][3:K]),month=rep(A[,1],K-2),day=rep(A[,2],K-2),
		x=c(A[,3:K]) )
	X$mm=paste(X$year,sapply(X$month,numbertochar),sep="")
	monthly=unlist( tapply(X$x,X$mm,mean,na.rm=T) )
	monthly=ts( c(monthly),start=min(as.numeric(as.character(X$year)) ),freq=12) 
	return(monthly)
}

make_monthly=function(amsu) {
	K=length(amsu)
	work=NULL	
	for (i in 1:K ) work=ts.union(work, make_anom(amsu[[i]]) )
 	dimnames(work)[[2]]=names(amsu)
 	window(work[,2],2008.25,2008.3)=NA # patch to eliminate tail value of dying satellite
	return(work)
}


# amsu = makef()
#  amsum=make_monthly(amsu)