;$Id: ppower_all0.pro,v 1.2 2026/10/24 15:00:58 brandenb Exp $
@parameters
print,"$sed.csh data/param2.nml"
pc_read_param,/param2,obj=param2
pc_read_param,obj=param
pc_read_ts,obj=ts
eta=param2.eta
;
;power,'_mag','_ele',k=k,spec1=spec1,spec2=spec2,i=n,tt=t,/noplot,/lks
power,'_mag','hel_mag',k=k,spec1=spec1,spec2=spec2,i=n,tt=t,/noplot,/lks
dk=k[1]-k[0]
;
;  This routine is used to view all power spectra
;  that were written during the run (controlled by dspec)
;  The plot range yr is set to a default value.
;  The time interval of plotting is "w" (default value)
;
;   5-oct-02/axel: written
;
spawn,'if (! -f param.pro) touch parameters.pro'
nt=n_elements(t)
xr=[k[1],max(k)]
k1=k & k1[0]=1. & k1=1./k1 & k1[0]=0.
EMag=fltarr(nt)
HMag=fltarr(nt)
xiMag=fltarr(nt)
k2Mag=fltarr(nt)
dk=k[1]-k[0]
;
default,w,.01
default,i1,0
default,i2,nt-1
i2=nt-1
default,istride,1
;
close,1
openw,1,'kpm.txt'
printf,1,"   i      t        EMag          xiMag        CM          Lu"
fo="(i4,e10.2,13e13.4)"
for i=i1,i2,istride do begin
  plot_oo,k,spec1(*,i),xr=xr,yr=yr,tit='t='+str(t(i))
  xx=[10.,300.] & oplot,xx,1./xx^(7./3.)
  xx=[1.2,10.] & oplot,xx,1e-7*xx^4
  oplot,k,.5*k*abs(spec2(*,i)),li=1
  oplot,k,.5*k*spec2(*,i),ps=1,col=122
  oplot,k,-.5*k*spec2(*,i),ps=1,col=55
  ;
  its=findex(t[i],ts.t)
  ;print,t[i],ts.tph[its],ts.ascale[its],ts.t[its]
  ;
  ;EMag[i]=dk*total(spec1(*,i))
  EMag[i]=total(spec1(*,i))
  HMag[i]=total(spec2(*,i))
  xiMag[i]=total(k1*spec1(*,i))/EMag[i]   ;(comoving)
  k2Mag[i]=total(k^2*spec1(*,i))/EMag[i]
  ;
  ; (same for electric field)
  ;EMag[i]=total(spec2(*,i))
  ;xiMag[i]=total(k1*spec2(*,i))/EMag[i]
  ;
  CM=t[i]*sqrt(2.*EMag[i])/xiMag[i]
  Lu=sqrt(2.*EMag[i])*xiMag[i]/eta
  vA=sqrt(2.*EMag[i])
  tau=xiMag[i]/vA
  print,   i,t[i],Emag[i],xiMag[i],CM,Lu,fo=fo
  printf,1,i,t[i],Emag[i],xiMag[i],CM,Lu,fo=fo
  oplot,[1,1]/xiMag[i],yr,li=3
  oplot,xr,[1,1]*EMag[i]/xiMag[i],li=3
  wait,w
endfor
close,1
;
spawn,'cvs add kpm.txt'
;save,file='grav1_last.sav',k,grav1_last,grav1_next
END