;$Id: pkpm8192a.pro,v 1.10 2026/05/24 17:20:47 brandenb Exp $
if !d.name eq 'PS' then begin
  device,xsize=18,ysize=12,yoffset=3
  !p.charthick=3 & !p.thick=3 & !x.thick=3 & !y.thick=3
end
;
!p.charsize=1.7
!x.margin=[9.4,.5]
!y.margin=[3.2,.3]
;
siz=1.7
si2=1.5
xr=[1.,100]
xr=[.4,240]
yr=[1e-4,2.2e-2]
!x.title='!7n!6!dM!n(!8t!6)'
!y.title='!8v!6!dA!n(!8t!6)'
lstop=1
lstop=0
;
file='/kpm.txt'
ss1=.5
ss2=1.3
fo='(a,2x,9e11.1)'
;
;-----------------------------------------------------------------------------
dir1='1024_k1em2_k1em0_vA1em1_nu1em4b' & istart1=2+6*2  ;black
dir2='1024_k1em3_k1em1_vA7em3_nu1em4b' & istart2=2+6*3  ;red
dir3='1024_k1em3_k1em1_vA1em1_nu1em3b' & istart3=2+6*3  ;blue
dir4='8192_k1em2_k1em0_vA1em1_nu1em5' & istart4=2+6 ;green, running, CB=2.9-2.95
dir5='8192_k1em3_k1em1_vA9em3_nu1em5' & istart5=0+6*3  ;Orange
dir5='8192_k1em3_k1em1_vA9em3_nu1em5u2' & istart5=0+6*0 & i1_dir2=21 & i2_dir2=29  ;Orange
;
;  i   t      EMag       xiMag        CM          Lu
;  & CM=a[4,*] & Lu=a[5,*]
;
default,CM1,10.
default,CM4,22.
dir=dir1 & i1=istart1
a=rtable('../'+dir+file,head=1,6)
t=a[1,*] & vA=sqrt(2.*a[2,*]) & xiM=a[3,*]
tA=xiM/vA
print,'black: tproper(0)=',CM1*tA[0]
!psym=8 & circ_sym,ss1,0
plot_oo,xiM,vA,ps=8,xr=xr,yr=yr
tt=every(t[0,i1:*],6) & vA=every(vA[0,i1:*],6) & xiM=every(xiM[0,i1:*],6)
!psym=8 & circ_sym,ss2,1
oplot,xiM,vA,ps=8
print,'black: ',tt,fo=fo
print,'black: ',t,fo=fo
xyouts,0.42,1.22e-2,siz=si2,'!8t!6=0'
xyouts,0.65,5.60e-3,siz=si2,'!8t!6=10!u3!n'
xyouts,1.60,2.01e-3,siz=si2,'!8t!6=10!u4!n'
if lstop then stop
;
dir=dir2 & i1=istart2
a=rtable('../'+dir+file,head=1,6)
t=a[1,*] & vA=sqrt(2.*a[2,*]) & xiM=a[3,*]
tA=xiM/vA
print,'red: tproper(0)=',CM1*tA[0]
;print,'red:   ',reform(t)
!psym=8 & circ_sym,ss1,0
oplot,xiM,vA,ps=8,col=122
tt=every(t[0,i1:*],6) & vA=every(vA[0,i1:*],6) & xiM=every(xiM[0,i1:*],6)
!psym=8 & circ_sym,ss2,1
oplot,xiM,vA,ps=8,col=122
print,'red:   ',tt,fo=fo
xyouts,30.,2.3e-4,siz=siz,'!8t!6=10!u6!n',col=122
if lstop then stop
;stop
;
dir=dir3 & i1=istart3
a=rtable('../'+dir+file,head=1,6)
t=a[1,*] & vA=sqrt(2.*a[2,*]) & xiM=a[3,*]
tA=xiM/vA
print,'blue: tproper(0)=',CM1*tA[0]
;print,'blue:  ',reform(t)
!psym=8 & circ_sym,ss1,0
oplot,xiM,vA,ps=8,col=55,th=4
tt=every(t[0,i1:*],6) & vA=every(vA[0,i1:*],6) & xiM=every(xiM[0,i1:*],6)
!psym=8 & circ_sym,ss2,1
oplot,xiM,vA,ps=8,col=55
print,'blue:  ',tt,fo=fo
xyouts,10.5,1.34e-2,siz=siz,'!8t!6=10!u3!n',col=55
xyouts,15.,9.1e-3,siz=siz,'!8t!6=10!u4!n',col=55
xyouts,35.4,2.8e-3,siz=siz,'!8t!6=10!u5!n',col=55
xyouts,100.,8e-4,siz=si2,'!8t!6=10!u6!n',col=55
if lstop then stop
;
dir=dir4 & i1=istart4
a=rtable('../'+dir+file,head=1,6)
t=a[1,*] & vA=sqrt(2.*a[2,*]) & xiM=a[3,*]
tA=xiM/vA
print,'green: tproper(0)=',CM4*tA[0]
loadct,6
!psym=8 & circ_sym,ss1,0
oplot,xiM,vA,ps=8,col=144
tt=every(t[0,i1:*],6) & vA=every(vA[0,i1:*],6) & xiM=every(xiM[0,i1:*],6)
!psym=8 & circ_sym,ss2,1
oplot,xiM,vA,ps=8,col=144
print,'green: ',tt,fo=fo
loadct,5
;
dir=dir5 & i1=istart5
a=rtable('../'+dir+file,head=1,6)
t=a[1,*] & vA=sqrt(2.*a[2,*]) & xiM=a[3,*]
tA=xiM/vA
print,'orange: tproper(0)=',CM4*tA[0]
!psym=8 & circ_sym,ss1,0
oplot,xiM,vA,ps=8,col=155
tt=every(t[0,i1:*],6) & vA=every(vA[0,i1:*],6) & xiM=every(xiM[0,i1:*],6)
!psym=8 & circ_sym,ss2,1
oplot,xiM,vA,ps=8,col=144
xyouts,24.5,4.8e-4,siz=siz,'!8t!6=10!u6!n',col=155
print,'orange: ',tt,fo=fo
;
;-----------------------------------------------------------------------------
!psym=0
xx=[0.53,9.] & oplot,xx,xx^.0/.60e2,li=1
xx=[1.03,10.] & oplot,xx,xx^.25/1.24e2,li=1
xx=[2.4,15.] & oplot,xx,xx^.7/7e2,li=1
;
xx=[1.3,13.] & oplot,xx,xx/450,li=2
;
xx=[5.,50.] & oplot,xx,xx/1e4,li=3
xx=[17.,140.] & oplot,xx,xx/1e5,li=3
xx=[.4,80.] & oplot,xx,1.1e-2/xx,li=3
;
xyouts,.5,.0021,siz=siz,'!7n!6!dM!n/450'
xyouts,3.,.00034,siz=siz,'!7n!6!dM!n/10!u4!n'
xyouts,9.,.00012,siz=siz,'!7n!6!dM!n/10!u5!n'
;
print,"$mv idl.ps ~/Overleaf/Mattia/Isochrones/fig/pkpm8192a.eps"
END