@parameters
default,kin,0
default,hdone,1
get_lhalf_factor_in_GW,lhalf_factor_in_GW
pc_read_param,/param2,obj=param2
cstress_prefactor=param2.cstress_prefactor
;
power,'_GWs','hel_GWs',k=k,spec1=grav1,spec2=grav2,i=n,tt=t,/noplot
if hdone then $
power,'_GWh','hel_GWh',k=k,spec1=grah1,spec2=grah2,i=n,tt=t,/noplot
;power,'_Str','hel_Str',k=k,spec1=stre1,spec2=stre2,i=n,tt=t,/noplot
;power,'_Str','hel_Str',k=k,spec1=grav1,spec2=grav2,i=n,tt=t,/noplot
;
;  $Id: polarization_degree.pro,v 1.1 2018/03/22 18:40:17 roper Exp $
;
;  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
;
pc_read_ts,o=ts

int_s = total(grav1, 1)
int_as = total(grav2, 1)
pol_energy = int_as/int_s
;print,'polarization energy', pol_energy

int_s = total(grah1, 1)
int_as = total(grah2, 1)
pol_hrms = int_as/int_s
;print,'polarization hrms', pol_hrms

t1=5.
tmin=1e-3
good=where(t gt t1)
!p.multi=[0,3,1]
;
plot,t(good),pol_energy(good)
oplot,t(good),accum(pol_energy(good))
;

plot,t(good),pol_hrms(good)
oplot,t(good),accum(pol_hrms(good))
;

good=where(ts.t gt t1)
a = ts.oum/sqrt(ts.o2m*ts.u2m)
plot,ts.t(good), a(good)
oplot,ts.t(good), accum(a(good))

!p.multi=0




END