histb2b2 = function(x1, x2, 
        names=NULL, breaks=NULL, 
        xlab="",ylab="",
        col=NULL, border="#444444", clinecol="#444444",
        xlim=NULL, ylim=NULL,
        ...) {
### Comments:
#     x1 and x2: vectors of values to plot
#     names: lables for the two variables
#     breaks: usual breakpoints for the histograms
#             may be the number of breaks, or the breaks themselves
#     col: two colors for the histograms
#     border: the border collor for the histogram rectangles
#     clinecol: color for the center line
#
#
#


# Check for defaults
if( length(col)<1 )   col=c(2,4)
if( length(col)>2 )   col=col[c(1,2)]

if( length(names)<1 ) names=c("x1","x2")
if( length(names)>2 ) names=names[c(1,2)]


if(length(breaks)==0) {
  breaks=seq( min(x1,x2), max(x1,x2), length=31 )
}

if(length(breaks)==1) {
  breaks=seq( min(x1,x2), max(x1,x2), length=breaks )
} 



# Create the two histograms
h1 = hist(x1, breaks=breaks, plot=FALSE)
h2 = hist(x2, breaks=breaks, plot=FALSE)


widest=max( c(h1$counts,h2$counts)*1.01 )



  if( length(xlim)!=2 ) xlim=c( min(breaks),max(breaks) )

  # Start plot
  plot.new()
  plot.window( xlim=xlim, ylim=c(-1.1,1.1) )

cAdj = 0.07

  # Plot the histogram boxes
  for(i in 1:(length(breaks)-1)) {
    rect( h1$breaks[i],+cAdj, h1$breaks[i+1],+h1$counts[i]/widest+cAdj, col=col[1], border=border)
    rect( h2$breaks[i],-cAdj, h2$breaks[i+1],-h2$counts[i]/widest-cAdj, col=col[2], border=border)
  }

# Print x-axis values

# ticks
  for(i in seq(1,length(breaks),4) ) {
    segments( h1$breaks[i],-cAdj, h1$breaks[i],-0.6*cAdj )
    segments( h1$breaks[i],+cAdj, h1$breaks[i],+0.6*cAdj )

    text(h1$breaks[i],0,label=h1$breaks[i])
  }

# END function

}