;$Id: pkt2_all.pro,v 1.2 2020/06/09 04:50:43 brandenb Exp $
if !d.name eq 'PS' then begin
  default,ysize,5
  default,thick,4
  device,xsize=18,ysize=ysize,yoffset=3
  !p.charthick=thick & !p.thick=thick & !x.thick=thick & !y.thick=thick
end
;
siz=1.7
yr1=[2e-8,5e-5]
yr2=[2e-6,5e-3]
yr3=[2e-5,5e-1]
;
@parameters
default,iread,0
default,ihydro,0
default,thick1,3
default,thick2,3
default,fourthird,1.
  iread=0
if iread eq 0 then begin
  power,'_mag','hel_mag',k=k,spec1=spec1,spec2=spec2,i=n,tt=t,/noplot
  power,'_Lor','hel_Lor',k=k,spec1=spec1h,spec2=spec2h,i=n,tt=t,/noplot
  power,'_2Lor','hel_2Lor',k=k,spec1=spec1c,spec2=spec2c,i=n,tt=t,/noplot
  iread=1
endif
nt=n_elements(t)
;
default,siz,1.7
!p.multi=[0,3,1]
half='!s!u 1!n!r!s-!r!d 2!n'
!x.title='!8k!6/!8k!6!d*!n'
!y.title='!8E!6!dM!n!6(!8k,t!6)'
!x.margin=[7.4,1.3]
!y.margin=[3.2,.4]
!p.charsize=siz
tilde='!9!s!aA!n!r!6'
;
default,w,.1
default,k0,60.
default,vA0,1.
s=1./(vA0^2*k0) & default,yr1,[2e-9,2e-2] & default,xr,[.9,576]/k0
tauA=1./(k0*vA0)
print,'tauA=',tauA
pc_read_param,/param2,o=param
eta=param.eta
E=t
D=t
II=t
IIK=t
;
iit=indgen(nt)
li1=indgen(nt)
li2=indgen(nt)
it0=iit[0]
nit=n_elements(iit)
;
for i=0,nit-1 do begin
  default,xstyle,1
  if xstyle then begin
    plot_oo,k/k0,s*spec2(*,0),xr=xr,yr=yr1,/nodata,xtickf='logticks_exp',ytickf='logticks_exp'
  endif else begin
    plot_oo,k/k0,s*spec2(*,0),xr=xr,yr=yr1,/nodata
  endelse
;
  it=iit(i)
  oplot,k/k0,abs(s*spec2(*,it)*k/2.),thick=thick1,col=55,li=1
  oplot,k/k0,s*spec2(*,it)*k/2.,thick=thick1,col=55
  oplot,k/k0,s*spec1(*,it),thick=thick2,col=122
  oplot,k/k0,k*0+1e-7
  ;
  ;  compute 1/k
  ;
  k1=k & k1(0)=1. & k1=1./k1 & k1(0)=0.
  EEM=total(spec2(*,it)) & xi=total(spec2(*,it)*k1)/EEM
;
;  plot II
;
default,xstyle,1
if xstyle then begin
  plot_oo,k/k0,s*spec2(*,0),xr=xr,yr=yr2,/nodata,xtickf='logticks_exp',ytickf='logticks_exp'
endif else begin
  plot_oo,k/k0,s*spec2(*,0),xr=xr,yr=yr2,/nodata
endelse
;
  it=iit(i)
  oplot,k/k0,abs(s*spec1h(*,it)),li=1,thick=thick2
  oplot,k/k0,    s*spec1h(*,it), li=0,thick=thick2
;
  it=iit(i)
  oplot,k/k0,abs(s*spec2h(*,it)),li=1,thick=thick2,co=188
  oplot,k/k0,    s*spec2h(*,it), li=0,thick=thick2,co=188
;
;  plot III
;
default,xstyle,1
if xstyle then begin
  plot_oo,k/k0,s*spec2(*,0),xr=xr,yr=yr3,/nodata,xtickf='logticks_exp',ytickf='logticks_exp'
endif else begin
  plot_oo,k/k0,s*spec2(*,0),xr=xr,yr=yr3,/nodata
endelse
;
  it=iit(i)
  ;oplot,k/k0,abs(s*spec1c(*,it)),li=1,thick=thick2
  ;oplot,k/k0,    s*spec1c(*,it), li=0,thick=thick2
  oplot,k/k0,abs(s*spec2c(*,it)*k/2.),thick=thick1,col=55,li=1
  oplot,k/k0,s*spec2c(*,it)*k/2.,thick=thick1,col=55
  oplot,k/k0,s*spec1c(*,it),thick=thick2,col=122
  print,'it,t=',it,t[it]
  wait,w
  stop
endfor
;
save,file='spec.sav',t,k,spec1,spec2,spec1h,spec2h
spawn,'cvs add -kb spec.sav'
;
if texdir eq 'inflation-hel' then ext='_all' else ext=''
;
cwd,run
print,'$mv idl.ps ~/tex/mhd/HallMHD/fig/pkt_'+run+ext+'.ps'
!p.multi=0
END