if !d.name eq 'PS' then begin
  device,xsize=15,ysize=26,yoffset=3
  !p.charthick=3 & !p.thick=3 & !x.thick=3 & !y.thick=3
end
if !d.name eq 'X' then begin
  window,0,xsize=480,ysize=480*2
end
;
; mv idl.ps scatter.ps
;
!p.charsize=1.6
!p.multi=[0,1,2]
nl=60 ; number of levels in countour plot
;ysize_cm=5.
;
; Read in data
;
file='scatterplot_t=2520.01.sav'
file='scatterplot_t=2763.01.sav'
;file='scatterplot_t=3111.99.sav'
restore,file=file
;
; Set some arrays etc.
;

nbins=100 ; Just for now (it is otherwise defined in the data file)

x=indgen(nbins)
xx = spread(x, 1, nbins)
yy = spread(x, 0, nbins)
;
; Plot v2_b2
;
;maxv=max(v2_b2)
;minv=0
;vect=indgen(nl)*(maxv-minv)/(nl-1)+minv
;xvar=xx*maxv2/nbins
;yvar=yy*maxb2/nbins
;contour,v2_b2,xvar,yvar,$
;  /fill,nlevels=nl,levels=vect,$
;  xtitle='!8v!6!u2!n',ytitle='!8B!6!u2!n',$
;  pos=[0.2, 0.7, 0.8, 0.9]
;
; Plot rho_b2
;
maxv=max(rho_b2)
minv=0
vect=indgen(nl)*(maxv-minv)/(nl-1)+minv
xvar=xx*(maxrho-minrho)/nbins+minrho
yvar=yy*maxb2/nbins
contour,rho_b2,xvar,yvar,$
  /fill,nlevels=nl,levels=vect,$
  xtitle='!4q!6',ytitle='!8B!6!u2!n',$
  pos=[0.25, 0.6, 0.85, 0.9]
;
; Plot rho_b2
;
maxv=max(rhov2_b2)
minv=0
vect=indgen(nl)*(maxv-minv)/(nl-1)+minv
xvar=xx*maxrhov2/nbins
yvar=yy*maxb2/nbins
contour,rhov2_b2,xvar,yvar,$
  /fill,nlevels=nl,levels=vect,$
  xtitle='!4q!8v!6!u2!n',ytitle='!8B!6!u2!n',$
  pos=[0.25, 0.1, 0.85, 0.4]


END