# TODO: Add comment
# 
# Author: E.Korsching 15.5.2009/21.11.2021, 2024
###############################################################################


bar.point.plot.err.ind <- function(x=NULL, y, z=NULL, labx=NULL, plot.type="b", whisker="u", bar.0=0, fit=F,
	xlim=NULL, ylim=NULL, numticks=NULL, axis.onoff=c(1,1), x.axis.l=NULL, width=1,
	label.y="", label.x="", off.x.lab=0.7, cex=1, cex.xy=1, srt=45, mar.i=c(5.1, 4.1, 4.1, 1.1),
	lwd=1, font.type=1, family="sans", col.p="black", col.b=NA, col.f="blue")
{
	# plot bars - with or without error/var indicator on top - with or without trend regression line
	# can plot bars  or  points  plus/minus regression line
	# x: position vector for y values or NULL (an index will be created)
	# y: vec with bar height values
	# z: vec with error indicator values or NULL,
	# labx: x label vector or NULL or "",   fit: Bezier line through y values
	# plot.type: b: bar, p: point,   bar.0: start of bar plot,
	#  xlim,ylim: manually adjust the proportion of the graph
	#  whisker: u: upper, l: lower, b: both whisker to bar, point
	# col.p: symbol color,   col.b: border col or NA,   col.f: fill color (one or as many as bars)
	# axis.onoff: c(x,y): 1 on, 0 off,   x.axis.l: NULL or one or more values for horizontal axis lines
	#  width: width of bar,   numticks: NULL or give a vector of numbers,   off.x.lab: offset of x scale labels
	#  cex: for x-y main labels,   cex.xy: for x-y scale labels,   mar.i: tune margins
	#  lwd: line width of axis,   font.type: 1: normal default, 2:bold, 3: italic, 4:bold+italic, 5:symbol
	#  family: "sans" default,"serif","mono","symbol"
	
	# functions
	source("../0functions/prox01/fn_test_norm.R")	# Bezier fit curve
	
	is.even <- function(y){ y %% 2 == 0 }
	
	sd.ind <- function(xvals, xmean, xsd, whisker=whisker, width=1, errcol="black"){
		# whisker: u: plot upper, l: plot lower, b: plot both error indicator(s)
		# width: from bar width to create a fraction for the horizontal err line
		#ini
		width <- width - (width/1.8)
		if(length(xmean)!=length(xsd)){ stop("mean and sd vector not of same size") }
		#
		err.U <- xmean+xsd		# upper err bar
		err.L <- xmean-xsd		# lower err bar
		# Draw error bars
		if(whisker=="b" | whisker=="u"){
			segments(xvals, xmean, xvals, err.U, col=errcol)   #line
			segments(xvals-width, err.U, xvals+width, err.U, col=errcol)   #top
		}
		if(whisker=="b" | whisker=="l"){
			segments(xvals, xmean, xvals, err.L, col=errcol)
			segments(xvals-width, err.L, xvals+width, err.L, col=errcol)
		}
		#
		return()
	}
	
	# ini
	mar <- par()$mar	# 5.1 4.1 4.1 2.1 lines
	par(mar=mar.i)		# mai <- par()$mai	# 1.02 0.82 0.82 0.42 inch
	
	len <- length(y)
	max.y <- max(y)
	min.y <- min(y)
	if(is.null(z)){	# err indicator
		max.dev <- 0
	}else{
		lenz <- length(z)
		if(lenz != len){ stop(paste("\n length of y ",len," != length of z ",lenz,sep="")) }
		max.dev <- max(z)
	}
	if(is.null(x)){
		vec.x <- c(1:len)
	}else{
		vec.x <- x
	}
	if(is.null(labx)){
		labx <- paste(vec.x,sep="")	# numbers as labels
	}else{
		vec.lab.l <- length(labx)
		if(vec.lab.l>1 & vec.lab.l!=len){
			stop(paste("\n number of x positions ",len," != length of labx ",vec.lab.l))
		}
	}
	
	# linear regression line
#	if(len>1 & fit){
#		# regression based on least square estimation - ssq in y
#		xbar <- mean(vec.x)
#		ybar <- mean(y)
#		b <- sum((vec.x-xbar)*(y-ybar))/sum((vec.x-xbar)^2)		# estimator for b
#		a <- ybar-b*xbar		# estimator for a
#		ei <- y-a-b*vec.x	# estimator for the theoretical residuals
#		# ei : distance to line - vertical up or respectively down
#		ei <- ei * (-1)		# customize for calculation
#		fit.y <- y + ei
#	}
	
	#plot
	if(is.null(xlim)){
		xlim <- c(min(vec.x)-0.5, max(vec.x)+0.5)
	}
	if(is.null(ylim)){
		ylim <- c(min.y-max.dev, max.y+max.dev)
	}else{
		ylim <- c(ylim[1]-max.dev, ylim[2]+max.dev)
	}
	plot(0,0, type="n", xlim=xlim, ylim=ylim, axes=F, xlab="", ylab="", cex=cex)		#empty plot window
	
	if(plot.type=="b"){
		for(i in 1:len){
			if(length(col.f)==1){
				selcol <- rep(col.f,len)
			}else{
				selcol <- col.f
			}
			polygon(x=c(vec.x[i]-width/2, vec.x[i]+width/2, vec.x[i]+width/2, vec.x[i]-width/2),
				y=c(y[i], y[i], bar.0, bar.0),
				density=-1, angle=45, border=col.b, col=selcol[i])		# -1:fill
		}
	}
	
	if(plot.type=="p"){
		points(x=vec.x, y=y, col=col.p)
	}
	
	if(len>1 & fit){
		#segments(x0=vec.x[1], y0=fit.y[1], x1=vec.x[len], y1=fit.y[len], col=col.p)
		lines(bezierCurve(vec.x, y, 40), lwd=2, col=col.f)
	}
	
	if(!is.null(z)){
		sd.ind(x=vec.x, xmean=y, xsd=z, whisker=whisker, width=width, errcol="black")
	} 
	
	if(axis.onoff[1]==1){	# x
		if(!is.null(x.axis.l)){		# axis line(s)
			abline(h=x.axis.l, lwd=lwd)
		}
		if(min(y)>=0){	# labels
			text(x=vec.x, y=bar.0, labels=labx, pos=1, offset=off.x.lab, srt=srt, cex=cex.xy, font=font.type, family=family)
		}else{
			text(x=vec.x, y=min(y), labels=labx, pos=1, offset=off.x.lab, srt=srt, cex=cex.xy, font=font.type, family=family)
		}
	}
	if(axis.onoff[2]==1){	# y
		axis(side=2, at=numticks, cex=cex.xy, lwd=lwd, font.axis=font.type, family=family)
	}
	title(xlab=label.x, ylab=label.y, cex.lab=cex, font.lab=font.type, family=family)
	
	par(mar=mar)
	return()
}