# TODO: Add comment
# 
# Author: E.Korsching 10.9.2009
###############################################################################



plot.differential <- function(
	x,		# spec(significant): xbar,ybar,Gene.Symbol,diff.means,fold.change,rawtp,rawUp,ci,t.p,U.p,SIG,sw.n1,sw.n2,sw.wratio,rawFp
	t.or.U,				# chosen value
	alpha,				# " 
	fold.change.k,		# "
	diff.k,				# "
	gin,				# chosen gin
	test.corr,			# chosen value
	label.col,			# "
	main = "",
	rows.on.page = 30,
	lwd = 2,
	col = c(1, 10, 16),	# 1. confidence plot, 2. group1, 3. group2
	cex = 0.7,
	data.name="test",
	relative.path="/results/"
)
{
	# plot analysis results from the 't.diffential' function using the 'Significant' object
	#  always save in file(s)
	
	# check if something to plot
	if(nrow(x)==0){ cat("\n no rows no differential plot \n"); return() }
    # check in col 4:10 on NA
	tmp <- apply(x[ ,c(4:10),drop=F], 1, function(x){any(is.na(x))} )
	tmp2 <- sum(tmp)
    if(tmp2>0){
		cat("\n NA values ",tmp2,"/",nrow(x)," - data set will be returned (aa.tmp...) \n")
		assign(x=paste("aa.tmp",data.name,sep="."), value=x, envir=.GlobalEnv)
	}
    
    # ini
	rnames <- row.names(x)
	nr <- nrow(x)		#number of positions on the y axis in the graph
	cext <- cex+0.3		#adjust text size in graph
	
	# extract interesting x cols in separate variables
	xbar <- x[, "xbar"]
	ybar <- x[, "ybar"]
	diff.means <- x[, "diff.means"]
	fold.change <- x[, "fold.change"]
	if(t.or.U=="t"){
		if(test.corr=="none"){ p.value <- x[, "rawtp"] }else{ p.value <- x[, "t.p"] }
 	}
	if(t.or.U=="U"){
		if(test.corr=="none"){ p.value <- x[, "rawUp"] }else{ p.value <- x[, "U.p"] }
 	}
	ci <- x[, "ci"]
	SIG <- x[, "SIG"]
	
	# customize
    xy <- cbind(xbar, ybar)
	ci1 <- diff.means - ci		# ci line left start point
	ci2 <- diff.means + ci		# ci line right end point
	
#	fc1 <- ifelse(fold.change > 1,fold.change,-1/fold.change)		# symmetric appearance -- already in t.differential
	
	# customize page-line counter - sp indicators   111......222......333... for rows.on.page 1,2,3...
	sp <- rep(1:ceiling(nr/rows.on.page), each = rows.on.page)
	# line counter on page i : 1,2,3...
	spi <- 1:nr

	npage <- length(unique(sp))		#number of pages
	n1 <- 0			# counter for processed rows - see figure legend
	
	pdf(file=paste(getwd(),relative.path,data.name,".",format(Sys.time(), "%Y%m%d%H%M"),".pdf", sep=""),
			width=11, height=7, onefile=T, title=paste(data.name,".diff plot", sep=""), pointsize=12)

    for(ipage in sort(unique(sp))) {
		
#		png(filename = paste(getwd(),relative.path,data.name,".diff.",format(Sys.time(), "%Y%m%d%H%M."),ipage,".png", sep=""),
#				width = 270, height = 170, units = "mm", res = 200, pointsize = 12, bg = "white")
		
		par(cex=cex)
		layout(matrix(c(1,2,3,4,5,6,7,8,9,rep(10,9)), 2, 9, byrow = TRUE), widths=c(2,0.8,2,0.8,0.8,0.8,0.8,0.9,5), heights=c(4,1))

		# revert plot order 1,2,3 to 3,2,1 , first on top
		ip <- rev(spi[sp == ipage])
		# remove NA values from incomplete filled pages
		ipNA <- !is.na(ip)
		ip <- ip[ipNA]

		nr1 <- length(ip)		#adjusted number of lines on present page
		
		# defining max x in confidence plot - per page
		tmp1 <- ci1[ip]
		tmp2 <- ci2[ip]
#		cat("\n tmp1 ",tmp1)
#		cat("\n tmp2 ",tmp2)
		m <- max(abs(range( tmp1[tmp1!=Inf], tmp2[tmp2!=Inf], na.rm=T) ) )		# filter on Inf and NA !
#		cat("\n m ",m)
		if(m==Inf){ m <- 1 }	# all values are NA --> Inf , one value in range() will return that value
		
		# 1/1 - t test
		par(mar=c(2, 6, 1, 1))
		plot(range(ci1[ip], ci2[ip]), c(0, rows.on.page), xlim=c(-m, m), type="n", xlab="p", ylab="", cex=cex,
            	adj=0, axes=F)
		
		axis(1, cex=cex, line=-1)
		axis(2, at=1:nr1, labels=signif(p.value[ip], 3), adj=1, cex=cex, las=2)			# p values t-test
		lines(c(0, 0), c(0, nr1), lty=2)

		y <- cbind(ci1[ip], ci2[ip], 1:nr1)
		diff1 <- diff.means[ip]
		for(i in 1:nr1) {
			points(diff1[i], i, cex=cex, col=palette()[col[1]], pch=19)		# sig points  + confidence interval
			lines(c(y[i, 1], y[i, 2]), c(y[i, 3], y[i, 3]), col=palette()[col[1]])
		}
		
		# 1/2 - significance labels : n, *, **, ***, ****
		par(mar=c(2, 2, 1, 1))
		plot(c(0, 2), c(0, rows.on.page), xlim=c(-m, m), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)
		text(x=1.8, y=1:nr1, labels=SIG[ip], cex=cex, adj=1, srt=0)	#cext
		text(x=0, y=0, labels="significance", cex=cex, font=2, adj=1, srt=90, xpd=TRUE)
		
		# 1/3 - plot the xbar ybar pairs
		xy1 <- xy[ip,,drop=F]
		maxy <- max(xy1) > 1		#flag: controls number of digits after the point
		if(maxy){ digi.num <- 1 } else { digi.num <- 2 }
		
		plot(c(0, max(xy1)), c(0, rows.on.page), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)

		for(i in 1:nr1) {			# bar plot of pairs of xbar ybar
		    my.barplot(x=0.25 + i, y=xy1[i, 1],		# xbar
				    col=c("black",palette()[col[2]]), width=0.4, xlab="", horiz=T, tick=F)
		    my.barplot(x=-0.25 + i, y=xy1[i, 2],		# ybar
				    col=c("black",palette()[col[3]]), width=0.4, xlab="", horiz=T, tick=F)
		}
		axis(1, cex=cex, line=-1)
		
		# 1/4 - xbar
		plot(c(0, 2), c(0, rows.on.page), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)
		text(x=1, y=1:nr1, labels=round(xbar[ip], digi.num), cex=cext, adj=1, srt=0, xpd=TRUE)
		text(x=1, y=0, labels="xbar", cex=cext, font=2, adj=1, srt=0, xpd=TRUE)
		
		# 1/5 - ybar
		plot(c(0, 2), c(0, rows.on.page), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)
		text(x=1, y=1:nr1, labels=round(ybar[ip], digi.num), cex=cext, adj=1, srt=0, xpd=TRUE)
		text(x=1, y=0, labels="ybar", cex=cext, font=2, adj=1, srt=0, xpd=TRUE)
		
		# 1/6 - diff.means
		plot(c(0, 2), c(0, rows.on.page), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)
		if(maxy){
			text(x=1, y=1:nr1, labels=round(diff.means[ip], 1), cex=cext, adj=1, srt=0, xpd=TRUE)
		}else{
			text(x=1, y=1:nr1, labels=round(diff.means[ip], 2), cex=cext, adj=1, srt=0, xpd=TRUE)
        }
		text(x=1, y=0, labels="delta", cex=cext, font=2, adj=1, srt=0, xpd=TRUE)
		
		# 1/7 - fc
		plot(c(0, 2), c(0, rows.on.page), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)
		text(x=1, y=1:nr1, labels=round(fold.change[ip],2), cex=cext, adj=1, srt=0, xpd=TRUE)
		text(x=1, y=0, labels="fc", cex=cext, font=2, adj=1, srt=0, xpd=TRUE)
		
		# 1/8 - row names
		plot(c(0, 1), c(0, rows.on.page), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)
		text(x=1, y=1:nr1, labels=rnames[ip], cex=cext, adj=1, srt=0, xpd=TRUE)
		text(x=1, y=0, labels="row ID", cex=cext, font=2, adj=1, srt=0, xpd=TRUE)
		
		# 1/9 - gin
		if(!missing(gin)){ label.g = gin[rnames[ip], label.col] }else{ label.g = rnames[ip] }	# genome information
		plot(c(0, 10), c(0, rows.on.page), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)
		text(x=0, y=1:nr1, labels=label.g, cex=cext, adj=0, srt=0, xpd=TRUE)
		text(x=0, y=0, labels="description", cex=cext, font=2, adj=0, srt=0, xpd=TRUE)
		
		# 2/10 - gin
		#cext <- cext + 0.2		# adjust text size for layout box at the bottom
		plot(c(0, 20), c(-4, 0), type="n", xlab="", ylab="", cex=cex,
				adj=0, axes=F)
		txt1 <- paste("\n treshold values: av.diff >= ", diff.k,", fold change >= ", fold.change.k,", ",t.or.U,"-test: p <= ", alpha,
					",  significance: <0.0001: ****, <0.001: ***, <0.01: **, <0.05: *, >0.05: n ",
					"   --   page ", ipage, " of ", npage, ",  results ", n1+1, ":", n1+nr1, " of ", nr, sep = "")
		
	    legend(x=1, y=0, legend=c("average value in group 1", "average value in group 2"), density=-1, fill=palette()[c(col[2], col[3])], yjust=1, border="black", bty="n", cex=cext)
		text(x=1, y=-2.5, labels=main, cex=cext, adj=0, srt=0, xpd=TRUE)
		text(x=1, y=-3.5, labels=txt1, cex=cext, adj=0, srt=0, xpd=TRUE)
		
		# row counter
		n1 <- n1 + nr1
		
		# png off
#		dev.off()
	}
	# pdf off
	dev.off()
}