function expscl,a,vmn,vmx,cmn,cmx,tcol=tcol,nlevels=nl
s=size(a)
c=bytarr(s(1),s(2))
if (n_elements(nc)   EQ 0) then nc=18
if (n_elements(nl)   EQ 0) then nl=6
if (n_elements(tcol) EQ 0) then tcol=1
nc=cmx-cmn+1

; Set missing values to 0

index=where(a LT 0,cnt)
if (cnt GT 0) then c(index)=0

; Set values between 0.0 and 1.0*scl to 1

if (vmn EQ -vmx) then begin
  scl=vmx/2.0^(nl/2)
  c0=fix(nc/2)
  if (nl MOD 2 NE 0) then nl = nl + 1
endif else begin
  scl=vmx/2.0^nl
  lmin=0.0
  lmax=1.0*scl
  rind=where(a GE lmin AND a LT lmax,cnt)
  if (cnt GT 0) then c(rind)=tcol
endelse

; Apply scaling to remaining values

if (vmn EQ -vmx) then begin
  lmin=scl
  j = 1
  for i=c0-1,0,-1 do begin
    x=float(j*nl)/float(nc-1)
    lmax=2.0^x*scl
    print,i,j,lmin,lmax,lmax-lmin,cmn+i
    rind=where(a GE -lmax AND a LT -lmin,cnt)
    if (cnt GT 0) then c(rind)=cmn+i
    lmin=lmax
    j = j + 1;
  endfor
  x=float(nl)/float(nc-1)
  lmin=-scl
  lmax=scl
  print,i,j,lmin,lmax,lmax-lmin,cmn+c0
  rind=where(a GE lmin AND a LT lmax,cnt)
  if (cnt GT 0) then c(rind)=cmn+c0
  lmin=scl
  j = 1;
  for i=c0+1,nc-1 do begin
    x=float(j*nl)/float(nc-1)
    lmax=2.0^x*scl
    print,i,j,lmin,lmax,lmax-lmin,cmn+i
    rind=where(a GE lmin AND a LT lmax,cnt)
    if (cnt GT 0) then c(rind)=cmn+i
    lmin=lmax
    j = j + 1;
  endfor
endif else begin
  lmin=lmax
  for i=1,nc do begin
    x=float(i*nl)/float(nc)
    lmax=2.0^x*scl
    rind=where(a GE lmin AND a LT lmax,cnt)
    if (cnt GT 0) then c(rind)=cmn+i-1
    lmin=lmax
  endfor
endelse
  
; Values greater than vmax get set to cmx

if (vmn EQ -vmx) then begin
  rind=where(a GE vmx,cnt)
  if (cnt GT 0) then c(rind)=cmx
  rind=where(a LE vmn,cnt)
  if (cnt GT 0) then c(rind)=cmn
endif else begin
  rind=where(a GE vmx,cnt)
  if (cnt GT 0) then c(rind)=cmx
endelse

return,c
end