Lightbending ptrace (xfig example)
Jump to navigation
Jump to search
Relativistic photon trajectory
require("/home/falkner/Public/lightbending.git/lbscripts"); define lb_plot_ptracesketch( R, b ){ variable size = qualifier("size",16); variable rmax = qualifier("rmax",10.); variable index = qualifier("index",""); if( index !="" ) index = string( "\\ensuremath{_"+index+"}" ); variable dx = qualifier("dx",1.); variable a1 = lb_a_ub( 1./R, b ); variable a2 = PI-a1; variable p1 = lb_psi_ua_belob( 1./R, a1 ); variable p2 = lb_psi_ua_belob( 1./R, a2 ; b=b ); variable p = lb_p_b( b ); variable psip = lb_psi_ua_belob( 1./p, PI/2. ; b=b ); variable pon = 0; if( qualifier_exists("periastron") ) pon = 1; variable pmax = p1; if(pon) pmax = p2; variable trace = lb_ptrace_belob( b; pmax=pmax, rmax=rmax); variable xmax = abs(max(trace.r)*cos(min(trace.p))); variable ymax = abs(max(trace.r)*sin(min(trace.p))); variable tmp,dr; variable xf; if( qualifier_exists("addto") ) xf=qualifier("addto"); else{ xf = xfig_plot_new(size,size); xf.world(-1.2*xmax,1.2*xmax,-1.2*xmax,1.2*xmax;xpad,ypad); xf.axis(;off); } %Circle at R tmp = [0:2*PI:#100]; xf.plot( R*cos(tmp), R*sin(tmp) ; line=2 ); tmp = 2.5/R*7./4.*PI; xf.plot( [0,R*cos(tmp)], [0,R*sin(tmp)] ; line=2 ); xf.xylabel(.5*R*cos(0.93*tmp), .5*R*sin(0.93*tmp) , "$R$"+index ); %trace xf.plot( trace.r*cos(trace.p), trace.r*sin(trace.p) ; width=4, color="red", depth=20); xf.xylabel(dx*xmax,b,"photon"+index,.5,-.5); % impact parameter xf.plot( dx*[0.8*xmax, 0.8*xmax],[0,b] ); xf.xylabel( dx*0.8*xmax, b/2., "$b$"+index, .5,0); ifnot( qualifier_exists("addto") ){ % observer sky xf.plot( [1.1*xmax,1.1*xmax],[-1.2*ymax,1.2*ymax] ; line=3 ); xf.xylabel( 1.12*xmax, 0, "{\\tiny observer sky at infinity}",-.5,0 ; rotate=90 ); % line of sight xf.plot( [0,xmax],[0,0] ); xf.xylabel( .5*(xmax+R),-.05*R,"\\tiny line of sight",.0,.5); } ifnot( qualifier_exists("nopsi") ){ % Psi R xf.plot( [0,1.5*R*cos(p1)],[0,1.5*R*sin(p1)] ; depth=19 ); tmp = [0:p1:#100]; dr = .45; xf.plot( dr*R*cos(tmp),dr*R*sin(tmp) ); dr += .2; xf.xylabel( dr*R*cos(p1/2.), dr*R*sin(p1/2.), "$\\Psi_{R"+index+"}$"); } % alpha if( qualifier_exists("alpha") ){ tmp = [p1-a1:p1:#100]; dr=.4; xf.plot( R*cos(p1)+dr*R*cos(tmp), R*sin(p1)+dr*R*sin(tmp) ;depth=19); xf.xylabel( R*cos(p1)+.7*dr*R*cos(p1-a1/2.), R*sin(p1)+.7*dr*R*sin(p1-a1/2.), "$\\alpha$"+index ); dr = 1.; xf.plot( R*cos(p1)+[0,dr*R*cos(p1-a1)],R*sin(p1)+[0,dr*R*sin(p1-a1)] ; line=1,depth=19); } if(pon){ ifnot( qualifier_exists("nopsi") ){ % Psi R prime xf.plot( [0,1.5*R*cos(p2)],[0,1.5*R*sin(p2)] ; depth=19); tmp = [0:p2:#100]; dr = .15; xf.plot( dr*R*cos(tmp),dr*R*sin(tmp) ); dr += .3; xf.xylabel( dr*R*cos(.5*(psip+p2)), dr*R*sin(.5*(psip+p2)), "$\\Psi^*_{R"+index+"}$"); % periastron xf.plot( [0,p*cos(psip)],[0,p*sin(psip)] ); tmp = [0:psip:#100]; dr = .3; xf.plot( dr*R*cos(tmp), dr*R*sin(tmp) ); dr += .2; xf.xylabel( dr*R*cos(.5*(psip+p1)), dr*R*sin(.5*(psip+p1)), "$\\Psi_{p"+index+"}$"); } if( qualifier_exists("alpha") ){ % alpha prime tmp = [p2-a2:p2:#100]; dr = .4; xf.xylabel( R*cos(p2)+.6*dr*R*cos(p2-a2/2.), R*sin(p2)+.6*dr*R*sin(p2-a2/2.), "$\\alpha^*$"+index ); xf.plot( R*cos(p2)+dr*R*cos(tmp), R*sin(p2)+dr*R*sin(tmp) ;depth=19); dr = 1.5; xf.plot( R*cos(p2)+[0,dr*R*cos(p2-a2)],R*sin(p2)+[0,dr*R*sin(p2-a2)] ; line=1,depth=19); if(qualifier_exists("symm")){ % alpha prime at alpha tmp = [p1:p1+a2:#100]; dr = .4; xf.xylabel( R*cos(p1)+.6*dr*R*cos(p1+a2/2.), R*sin(p1)+.6*dr*R*sin(p1+a2/2.), "$\\alpha^*$"+index ); xf.plot( R*cos(p1)+dr*R*cos(tmp), R*sin(p1)+dr*R*sin(tmp);depth=19 ); dr = 1.; xf.plot( R*cos(p1)+[0,dr*R*cos(p1+a2)],R*sin(p1)+[0,dr*R*sin(p1+a2)] ; line=1,depth=19); } } } return xf; } variable R = 2.5; variable b = 2.9; variable xf; xf = lb_plot_ptracesketch( R, b ); xf.render("ptrace1.pdf"); xf = lb_plot_ptracesketch( R, b ; periastron ); xf.render("ptrace2.pdf"); xf = lb_plot_ptracesketch( R, b ; alpha,periastron,symm ); xf.render("ptrace3.pdf"); variable R2 = 1.5*R; variable b2 = 1.3*b; variable p1 = lb_psi_ua_belob( 1./R, lb_a_ub( 1./R, b ) ); variable p2 = lb_psi_ua_belob( 1./R2, lb_a_ub( 1./R2, b2 ) ); variable p3 = lb_psi_ua_belob( 1./R2, PI-lb_a_ub( 1./R2, b2 )); xf = lb_plot_ptracesketch( R2, b2 ; nopsi,index="2",periastron,dx=0.9); xf = lb_plot_ptracesketch( R, b ; nopsi, addto=@xf, index="1" ); xf.plot( R,0; sym=1, size=1,fill ); xf.plot( R*cos(p1),R*sin(p1); sym=1, size=1,fill,fillcolor="red",depth=18 ); xf.plot( R2,0; sym=0, size=4,fill ); xf.plot( R2*cos(p2),R2*sin(p2); sym=0, size=4,fill,fillcolor="red",depth=18 ); xf.plot( R2*cos(p3),R2*sin(p3); sym=0, size=4,fill,fillcolor="white",depth=18 ); xf.render("ptrace4.pdf"); exit;