;***********************************************************************;
; Function : hist_columns                                               ;
;                   wks: workstation object                             ;
;                    xy: graphic                                        ;
;             binvalues: numeric                                        ;
;               barlocs: numeric                                        ;
;              barwidth: numeric                                        ;
;                colors                                                 ;
;               compare: logical                                        ;
;                 gsres: logical                                        ;
;                                                                       ;
;  xy        - xy plot id to draw columns on                            ;
;  bins      - the center of each bin range.                            ;
;  binvalues - the number of values in the corresponding bin            ;
;  barlocs   - the start of the first bar in each bin                   ;
;  width     - the width of the bar                                     ;
;  colors    - array of colors to use                                   ;
;  gsres     - optional primitive resources                             ;
;                                                                       ;
; This function creates the columns for a histogram plot. The Y axis    ;
; will represent the number of values in a bin and a percentage.        ;
;                                                                       ;
;***********************************************************************;
undef("hist_columns")
function hist_columns(wks[1]:graphic,xy[1]:graphic,binvalues:numeric, \
                      barlocs[*][*]:numeric, barwidth[*]:numeric, \
                      colors, compare:logical, gsres:logical)
local i, nbins, dims, nbinvalues, gsres, xpoints, ypoints, multibars,\
col_dims, col_rank
begin
  dims       = dimsizes(barlocs)
  nbarsinbin = dims(0)
  nbars      = dims(1)
  delete(dims)

  dims = dimsizes(binvalues)
  if(dimsizes(dims).eq.1) then
    nbinvalues = dims(0)
  else
    nbinvalues = dims(1)
  end if
  delete(dims)
  if(nbars.ne.nbinvalues) then
    print("Error: hist_columns: Dimension sizes of bins (" + nbars+ ") and binvalues (" + nbinvalues + ") must be the same")
    return
  end if

  if(nbarsinbin.ge.2)
    multibars   = True
  else
    multibars   = False
  end if

;
; Set up arrays to hold polygon points.
;
  if(multibars) then
    xpoints   = new((/nbarsinbin,nbars,5/),float)
    ypoints   = new((/nbarsinbin,nbars,5/),float)
    polygons  = new((/nbarsinbin,nbars/),graphic)
  else
    xpoints  = new((/nbars,5/),float)
    ypoints  = new((/nbars,5/),float)
    polygons = new(nbars,graphic)
  end if
;
; Set up variable to hold resources.
;
  gsres           = True
  set_attr(gsres,"gsEdgesOn",True)
;
; Begin assigning polygon points.
;
  if(multibars)
    do i=0,nbarsinbin-1
      ypoints(i,:,0) = (/1/)               ;;; Changed 0 to 1
      ypoints(i,:,1) = (/binvalues(i,:)/)
      ypoints(i,:,2) = (/ypoints(i,:,1)/)
      ypoints(i,:,3) = (/1/)               ;;; Changed 0 to 1
      ypoints(i,:,4) = (/1/)               ;;; Changed 0 to 1

      xpoints(i,:,0) = (/barlocs(i,:)/)
      xpoints(i,:,1) = (/xpoints(i,:,0)/)
      xpoints(i,:,2) = (/barlocs(i,:) + barwidth/)
      xpoints(i,:,3) = (/xpoints(i,:,2)/)
      xpoints(i,:,4) = (/xpoints(i,:,0)/)
    end do
  else
    ypoints(:,0) = (/1/)               ;;; Changed 0 to 1
    ypoints(:,1) = (/binvalues/)
    ypoints(:,2) = (/ypoints(:,1)/)
    ypoints(:,3) = (/1/)               ;;; Changed 0 to 1
    ypoints(:,4) = (/1/)               ;;; Changed 0 to 1

    xpoints(:,0) = (/barlocs(0,:)/)
    xpoints(:,1) = (/xpoints(:,0)/)
    xpoints(:,2) = (/barlocs(0,:) + barwidth/)
    xpoints(:,3) = (/xpoints(:,2)/)
    xpoints(:,4) = (/xpoints(:,0)/)
  end if

  if(compare)
    fillindex = get_res_value(gsres,"gsFillIndex",(/0,6/))
  else
    if(multibars) then
      fillindex = get_res_value(gsres,"gsFillIndex",0)
      if(dimsizes(dimsizes(fillindex)).eq.1) then
         ftmp = new(nbarsinbin,typeof(fillindex))
         ftmp = fillindex
         delete(fillindex)
         fillindex = ftmp
         delete(ftmp)
      end if
    else
      fillindex = get_res_value(gsres,"gsFillIndex",0)
    end if
  end if
;
; Make sure fill indices are between 0 and 17.
; No more, because you have have added your own fill index.
;
;  fillindex = min((/max((/max(fillindex),0/)),17/))

  rgba_arr = convert_color_to_rgba(wks,colors)
  col_dims = dimsizes(rgba_arr)   ; better be 2D or 3D
  col_rank = dimsizes(col_dims)
  if(col_rank.gt.3.or.(col_rank.eq.3.and.col_dims(0).ne.nbarsinbin)) then
    print("Error: hist_columns: Invalid dimensionality of array of colors passed in.")
    print("   Was expecting a 1D or 2D array of colors")
    exit
  end if
  if(col_rank.eq.3) then
    ncolors = col_dims(1)
  else
    ncolors = col_dims(0)
  end if
  do i = 0, nbars - 1
    if(col_rank.eq.1) then
      gsres@gsFillColor = rgba_arr
    elseif(col_rank.eq.2) then
      gsres@gsFillColor = rgba_arr(i % ncolors,:)
    end if
    if(multibars) 
;
; Add the remaining histogram bars first, so that they will be drawn
; second. The first histogram bars will thus be drawn on top, if
; they are being stacked.
;
      do j=nbarsinbin-1,0,1
        if(col_rank.eq.3) then
          gsres@gsFillColor = rgba_arr(j,i % ncolors,:)
        end if
        gsres@gsFillIndex = fillindex(j)
        if(binvalues@horizontal) then
          polygons(j,i) = gsn_add_polygon(wks,xy,ypoints(j,i,:), \
                                         xpoints(j,i,:), gsres)
        else
          polygons(j,i) = gsn_add_polygon(wks,xy,xpoints(j,i,:), \
                                          ypoints(j,i,:), gsres)
        end if
      end do
    else
      gsres@gsFillIndex = fillindex
      if(binvalues@horizontal) then
        polygons(i) = gsn_add_polygon(wks,xy,ypoints(i,:),xpoints(i,:),gsres)
      else
        polygons(i) = gsn_add_polygon(wks,xy,xpoints(i,:),ypoints(i,:),gsres)
      end if
    end if
  end do
;
; Return the polygons created as attributes of the XY plot. This is
; necessary, b/c otherwise the polygons will go away when you exit the
; function.
;
  var_string      = unique_string("hpolygons")
  xy@$var_string$ = polygons
  return(xy)
end