;$Id: ppower_all2.pro,v 1.8 2026/06/07 08:10:07 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','_kin',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)
xiMag=fltarr(nt)
Lu=fltarr(nt)
IA=fltarr(nt)
tA=fltarr(nt)
;
;  This routine also computes Cxi and CEE.
;  We take the value at t=1e5.
;
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          IA          Cxi         CEE"
fo="(i4,e10.2,7e12.4)"
for i=i1,i2,istride do begin
  plot_oo,k,spec1(*,i),xr=xr,yr=yr,tit='t='+str(t(i))
  oplot,k,spec2(*,i),li=1
  ;
  EMag[i]=total(spec1(*,i))
  xiMag[i]=total(k1*spec1(*,i))/EMag[i]
  CM=t[i]*sqrt(2.*EMag[i])/xiMag[i]
  Lu[i]=sqrt(2.*EMag[i])*xiMag[i]/eta
  its=findex(t[i],ts.t)
  IA[i]=ts.arms[its]^2
  Cxi=xiMag[i]/(IA[i]^.25*t[i]^.5)
  CEE=EMag[i]*t[i]/IA[i]^.5
  tA[i]=xiMag[i]/sqrt(2.*EMag[i])
  ;
  print,i,t[i],Emag[i],xiMag[i],CM,Lu[i],IA[i],Cxi,CEE,fo=fo
  printf,1,i,t[i],Emag[i],xiMag[i],CM,Lu[i],IA[i],Cxi,CEE,fo=fo
  oplot,[1,1]/xiMag[i],yr,li=3
  oplot,xr,[1,1]*EMag[i]/xiMag[i],li=3
  ;xx=[.012,.3] & oplot,xx,1e-4*xx^3
  xx=[.012,.3] & oplot,xx,1e-4*xx^5
  wait,w
endfor
close,1
;
spawn,'cvs add kpm.txt'
;save,file='grav1_last.sav',k,grav1_last,grav1_next
END
