;$Id: pspecm_hel_delk.pro,v 1.4 2021/03/25 14:19:36 brandenb Exp $ if !d.name eq 'PS' then begin ;device,xsize=18,ysize=8,yoffset=3 device,xsize=18,ysize=6,yoffset=3 !p.charthick=1.6 & !p.thick=1.6 & !x.thick=1.6 & !y.thick=1.6 end ; dir1='M512sig1_k6_wav03_1p9em1b' dir2='M512sig1_k6_wav03_1p9em1b_delk03' dir3='M512sig1_k6_wav03_1p9em1b_delk1' dir4='M512sig1_k6_wav03_1p9em1b_delk3' dir5='M512sig1_k6_wav03_1p9em1b_delk10' ; xs=.34 ys1=.06 ys2=.6 ; yr_pspecm_scl3=[1d-20,3d-7] yr_pspechm_scl3=ys2*[1d-20,2.2d-13] ; thick4=5 siz=1.0 si2=1.0 !p.charsize=1.2 !x.margin=[7.8,0.5] !y.margin=[3.2,0.2] !p.multi=[0,2,1] ;fac_h=1.263e-18 fac2=2. half='!s!u 1!n!r!s-!r!d 2!n' specm='/specm.sav' wspecm='/specm.txt' label='' ; ; determine normalization ; default,maglabel,'mag' default,green_lines,0 default,plot_polarization,1 ytickf='logticks_exp' yr=yr_pspecm_scl3 ; ; new units ; T_in100GeV=1. & gstar_in100=1. T_in100GeV=1.5d-3 & gstar_in100=.15 ;gS_in100=1. ;gS_in100=.0391 tend=1. Hstar=2.066d10*T_in100GeV^2*gstar_in100^.5 H0=3.241d-18 ;a0=1.254d15*T_in100GeV*gS_in100^(1./3.) a0=1d12 fact=(Hstar/H0)^2*(tend/a0)^4 print,'fact=',fact ; ;Hstar=2.066d10 ;a0=1.2545d15 freq=Hstar/(2.*!pi*a0) H0=3.24d-18 ;changed from H0 = 100 km/s/Mpc to 70 km/s/Mpc ;H0=2.268545d-18 rhocrit_fac=3./(8.*!pi) ;fac_h=1.263e-18 ;changed H0 from 100 to 70 km/s/Mpc, and computed h_c ;from formula fac_h=sqrt(3./2.)*H0/!pi fac=(Hstar/H0)^2/rhocrit_fac*(1./a0)^4 facGW=1./(16.*!pi) ;xr=xs*[2e-9,1.6e-5] xr=[7e-10,3e-7] circ_sym,0.7,1 mixed=0. ; ; 1st panel ; !p.title='!6';+run !x.title='!8f!6 [Hz]' !y.title='!6!8h!6!s!d0!n!r!u2!n!7X!6!dGW!n(!8f!6)' ; loadct,0 plot_oo,xr,yr,ytickf=ytickf,/nodata up=rtable('data/GreenUpper2021.csv',2) dn=rtable('data/GreenLower2021.csv',2) xxx=10.^[reform(dn[0,*]),reverse(reform(up[0,*])),dn[0,0]] yyy=10.^[reform(dn[1,*]),reverse(reform(up[1,*])),dn[1,0]] polyfill,xxx,yyy<(.9*yr[1]),col=222 loadct,5 ; restore,'../'+dir1+specm grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot ;plot_oo,freq*k,OmegaGW,xr=xr,yr=yr,ytickf=ytickf,/nodata oplot,freq*k,OmegaGW ; restore,'../'+dir2+specm grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot oplot,freq*k,OmegaGW,col=55 ; loadct,6 restore,'../'+dir3+specm grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot oplot,freq*k,OmegaGW,col=122 loadct,5 ; restore,'../'+dir4+specm grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot oplot,freq*k,OmegaGW,col=155 ; restore,'../'+dir5+specm grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot oplot,freq*k,OmegaGW,col=122 ; ;xyouts,xs*5.5e-8,ys1*6.0e-14,'!8k!6!df!n=2',siz=siz,col=155 ;xyouts,xs*2.0e-7,ys1*1.0e-13,'!66',siz=siz ;xyouts,xs*7.0e-7,ys1*7.3e-13,'!620',siz=siz,col=55 ;xyouts,xs*2.1e-6,ys1*4.8e-12,'!660',siz=siz,col=122 ; ;xx=xs*[3e-9,1.0e-7] & oplot,li=3,xx,ys1*3e-2*xx ;xyouts,xs*1.4e-8,ys1*2.3e-11,'!9A!8f!6',siz=siz xout=1.6e-7 xyouts,xout,2.0e-8,'!6(a)',siz=si2 ; ; 2nd panel ; !y.title='!8h!6!dc!n(!8f!6)' default,kin,0 default,hdone,1 ytickf='logticks_exp' yr=yr_pspechm_scl3 circ_sym,0.7,1 ; restore,'../'+dir1+specm fach=1.263e-18/(freq*k[1:*]) grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot plot_oo,freq*k,sqrt(k*grah1m)/a0,xr=xr,yr=yr,ytickf=ytickf,/nodata oplot,freq*k,sqrt(k*grah1m)/a0 ; ; 2nd line ; restore,'../'+dir2+specm fach=1.263e-18/(freq*k[1:*]) grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot oplot,freq*k,sqrt(k*grah1m)/a0,col=55 ; ; 3rd line ; loadct,6 restore,'../'+dir3+specm fach=1.263e-18/(freq*k[1:*]) grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot oplot,freq*k,sqrt(k*grah1m)/a0,col=122 loadct,5 ; ; 4th line ; restore,'../'+dir4+specm fach=1.263e-18/(freq*k[1:*]) grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot oplot,freq*k,sqrt(k*grah1m)/a0,col=155 ; ; 5th line ; restore,'../'+dir5+specm fach=1.263e-18/(freq*k[1:*]) grav_tot=grav1m+grah1m/tm[0]^2+mixed OmegaGW=facGW*fac*k*grav_tot oplot,freq*k,sqrt(k*grah1m)/a0,col=122 ; ;xyouts,xs*3.0e-8,ys2*2.0e-17,'!8k!6!df!n=2',siz=siz,col=155 ;xyouts,xs*1.5e-7,ys2*2.3e-17,'!66',siz=siz ;xyouts,xs*7.0e-7,ys2*1.8e-18,'!620',siz=siz,col=55 ;xyouts,xs*1.4e-6,ys2*1.6e-17,'!660',siz=siz,col=122 ; ;xx=xs*[2.6e-9,1.3e-7] & oplot,li=3,xx,ys2*6e-18/xx^.5 ;xyouts,xs*1.9e-7,ys2*1.0e-14,'!9A!8f!6!u-1/2!n',siz=siz xyouts,xout,3.3e-14,'!6(b)',siz=si2 ; print,"$mv idl.ps ~/GitHub/Yutong/MGW-NANOGrav/Figures/pspecm_hel_delk.eps" END