Pulsars and neutron stars/Pulsar properties/Ppdot figure script
Script for making p-pdot diagram
editThis figure was made using PSRCAT and Julia. First, to make the data files use the following PSRCAT commands
psrcat -c "p0 p1" -o short -l "type(radio) && exist(p0) && p1 > 0" -nohead -nonumber > all.dat
psrcat -c "p0 p1" -o short -l "type(radio) && exist(p0) && p1 > 0 && pb > 0" -nohead -nonumber > binary.dat
Then the following Julia commands:
#
using Winston;
all = readdlm("all.dat");
binary = readdlm("binary.dat");
all = sortrows(all,by=x->(x[1]));
binary = sortrows(binary,by=x->(x[1]));
# Use Chen, Ruderman 1993 death line
dline = 10^(2.0*78.0/7.0)*all[:,1].^(19.0/7.0)/3.2e19^2;
b10_14line = (1e14/3.2e19)^2./all[:,1];
b10_12line = (1e12/3.2e19)^2./all[:,1];
b10_10line = (1e10/3.2e19)^2./all[:,1];
b10_8line = (1e8/3.2e19)^2./all[:,1];
age_1000yr = all[:,1]/(1e3*365.25*86400.0*2);
age_100000yr = all[:,1]/(100e3*365.25*86400.0*2);
closefig()
# WARNING: color("red") is deprecated, use colorant"red" or parse(Colorant, "red")
# should be add(p, Curve(x, c, color=colorant"red")), but this doesn't work
#
p = FramedPlot(xlabel="Period (s)",ylabel="Period derivative (s/s)",xlog=true,ylog=true,yrange=(1e-22,1e-9));
pts_all = Points(all[:,1],all[:,2],"symbolsize",0.5,"color","black");
pts_binary = Points(binary[:,1],binary[:,2],"color","red");
setattr(pts_binary,label="Binary");
pts_dline = Curve(all[:,1],dline[:,1],"type","dotted","color","green");
plot_b10_14line = Curve(all[:,1],b10_14line[:,1],"type","dashed","color","blue");
plot_b10_12line = Curve(all[:,1],b10_12line[:,1],"type","dashed","color","blue");
plot_b10_10line = Curve(all[:,1],b10_10line[:,1],"type","dashed","color","blue");
plot_b10_8line = Curve(all[:,1],b10_8line[:,1],"type","dashed","color","blue");
plot_age_1000yr = Curve(all[:,1],age_1000yr[:,1],"type","dashed");
plot_age_100000yr = Curve(all[:,1],age_100000yr[:,1],"type","dashed");
l = Legend(.1, .9, {pts_binary})
text1 = PlotLabel(0.1,0.63,"1 kyr");
text2 = PlotLabel(0.1,0.47,"1 Myr");
text3 = PlotLabel(0.65,0.9,"10^{14} G","color","blue");
text4 = PlotLabel(0.15,0.77,"10^{12} G","color","blue");
text5 = PlotLabel(0.85,0.24,"10^{10} G","color","blue");
add(p,pts_all,pts_binary,pts_dline,plot_b10_14line,plot_b10_12line,plot_b10_10line,plot_b10_8line,plot_age_1000yr,plot_age_100000yr,l,text1,text2,text3,text4,text5)
savefig(p,"ppdot.png")