Relativistic Aberration and the Appearance of Rapidly Moving Objects

Peter Cramer, TEG

(RowDefault[, (Spacer1), , (Spacer1), ])  #WolframTechConf

Outline

•   Introduction

•   Mathematical Model

•   Interactive Simulation

•   Conclusions

•   Discussion

ABSTRACT

I present a simulation of how the appearance of rapidly moving objects is distorted due to relativistic effects.  It is based on a computationally efficient model of aberration that I will derive.  I will use the Manipulate[] function to demonstrate a number of interesting visual effects.

Introduction                 (Mr. Tompkins, Gamow, Cambridge University Press)

•  Relative motion affects our perception of objects.

•  The frequency of detected light changes due to the Doppler effect.

•  Objects appear to shrink the direction of motion (Lorentz contraction).

•  Shapes are distorted by aberration.

•  Up to now, no general solution has been found for predicting the apparent shape of objects moving at relativistic speeds.

A Different Approach

•   Treat objects as collections of point sources of light

•   Treat aberration as a latency effect

•   Perform analysis in a single inertial frame of reference

•   Assume constant velocity

•   Use the Manipulate[] function

Aberration as Angular Effect

Aberration was first discovered by astronomer James Bradley in 1727.  He noted that the angular positions of stars varies with the time of year.  The effect is of order  v/c = 20 arcsec.   This is angular aberration.

   (from Wikipedia)

But this is not the whole picture.

Aberration as Latency Effect

Aberration:

Experiment

•  Construct object of interest.
•  Accelerate object to desired velocity.
•  Observe the object.

Simulation

•   Construct the desired collection of point sources.
•   Boost the object, i.e., apply Lorentz contraction.
•   Use aberration model to calculate the apparent position of each source.
•   Use optics model to simulate what a telescope would see.
•   Render the apparent positions and the telescopic view.

Latency Model

The ray is detected when it reaches the observer at the origin:

(Subscript[(OverVector[r]), o]) + (OverVector[v]) (Subscript[t, e  ]) + c (t-(Subscript[t, e])) (OverHat[u]) = 0                     

(Subscript[(OverVector[r]), o]) + (OverVector[v]) (Subscript[t, e  ]) =  c (t-(Subscript[t, e])) (OverHat[u])                                   

((Subscript[r, o])^2) + (((v (Subscript[t, e])))^2) + 2 (Subscript[t, e]) ((Subscript[(OverVector[r]), o]) · (OverVector[v]))  =  (c^2) (((t - (Subscript[t, e])))^2)           

Solving the quadratic gives two solutions:

	(Subscript[t, e])  =  (γ^2) (t + d ± (Sqrt[(d^2)+2 d t+(((v t))^2)+((((Subscript[r, o]) /γ))^2)]))       

	where       c->1,       γ = (1/(Sqrt[1-(v^2)])),      d = (Subscript[(OverVector[r]), 0]) · (OverVector[v])

Which Solution?

Retain the solution that is real and satisfies causality  (Subscript[t, e]) < t.  

Once (Subscript[t, e]) is known, one can compute the emit position  (Subscript[(OverVector[r]), e]) , which is also the apparent position, from 

	(Subscript[(OverVector[r]), e]) = (Subscript[(OverVector[r]), 0]) + (OverVector[v]) (Subscript[t, e])(t)

Summary

•  Specify (Subscript[(OverVector[r]), 0]), (OverVector[v]), t
•  Compute (Subscript[t, e])(t)
•  Compute x((Subscript[t, e]))
•  Repeat for all sources

Lorentz Contraction Model

We apply Lorentz contraction immediately after specifying the initial positions of all the point sources.   In scalar form it is given by

L = (Subscript[L, o])/γ

We will use the vector form:


   
where  is a unit vector in the direction of the velocity.

Optics Model

We model the observer’s telescope as an ideal camera (pinhole camera).

Simulation

Manipulate[
Module[
{sources,tes,te,xe,ye,ze,uC,uCe,aberration,colors,apparentPoints,rotatedPositions,projectedPositions,te1Desired,tDesired,nSources,initialDistances,t,generateSources,vy,y0=0.},

(* The user inputs a magnitude, negative sign gives vector the correct sense. *)
(* Absolute value prevents user from seeing zero stars. *)
vy=-Abs[vyin];

(* construct the objectof interest which consists of a set of point sources. *)
Which[
shape=="cube",
sources=generateCubeEdges[x0,y0,0,spread,nSourcesDesired],

shape=="cluster",
SeedRandom[1];
sources=RandomPlummerCluster[nSourcesDesired,spread];
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="galaxy_x",
SeedRandom[1];
sources=RandomSpiralGalaxy[nSourcesDesired,1.5,0.2,spread,0.23*spread];
sources=Transpose[RotateLeft[Transpose[sources],2]];
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="galaxy_y",
SeedRandom[1];
sources=RandomSpiralGalaxy[nSourcesDesired,1.5,0.2,spread,0.23*spread];
sources=Transpose[RotateLeft[Transpose[sources],1]];
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="galaxy_z",
SeedRandom[1];
sources=RandomSpiralGalaxy[nSourcesDesired,1.5,0.2,spread,0.23*spread];
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="line_x",
sources=2*{spread*(Range[0,nSourcesDesired-1]/(nSourcesDesired-1 )-.5),ConstantArray[0,nSourcesDesired],ConstantArray[0,nSourcesDesired]}//Transpose;
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="line_z",
sources=2*{ConstantArray[0,nSourcesDesired],ConstantArray[0,nSourcesDesired],spread*(Range[0,nSourcesDesired-1]/(nSourcesDesired-1 )-.5)}//Transpose;
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="line_y",
sources=2*{ConstantArray[0,nSourcesDesired],spread*(Range[0,nSourcesDesired-1]/(nSourcesDesired-1 )-.5),ConstantArray[0,nSourcesDesired]}//Transpose;
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="circle_y",
sources=generateDiskPoints[0,0,0,spread,nSourcesDesired];
sources=Transpose[RotateLeft[Transpose[sources],1]];
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="circle_x",
sources=generateDiskPoints[0,0,0,spread,nSourcesDesired];
sources=Transpose[RotateLeft[Transpose[sources],2]];
sources+=ConstantArray[{x0,y0,0},Length[sources]],

shape=="circle_z",
sources=generateDiskPoints[x0,y0,0,spread,nSourcesDesired],

shape=="sphere",
sources=generateSpherePoints[x0,y0,0,spread,nSourcesDesired]
];
sources=Transpose[Join[Transpose[sources],Transpose[ConstantArray[{0,vy,0},Length[sources]]]]];

(* place a point source at the selected shape's centroid for visual reference *)
(* The simulated telescope tracks this point. *)
PrependTo[sources,{x0,y0+.001*spread,0,0,vy,0}];  

(* This point is used by the velocity vector display *)
PrependTo[sources,{x0,y0,0,0,vy,0}];

(* Lorentz-contract the object *)
nSources=Length[sources];
If[applyLorentzContraction,
sources=lorentzContract[nSources,sources,vy,c];
];

(* Each point source has the same color in both displays. *)
colors=ColorData["Rainbow"]/@Range[0,1,1/(nSources-1)];

(* apparent position of rest-centroid is set by user *)
te1Desired=(yeControl-y0)/vy;
tDesired=te1Desired+Norm[{x0,yeControl,0}]/c;
t=tDesired;

(*Calculate apparent positions*)
apparentPoints={};
Table[Module[{x0s,y0s,z0s,vxs,vys,vzs,color},
{x0s,y0s,z0s,vxs,vys,vzs}=sources[[i]];
color=colors[[i]];
tes=emissionTimes[x0s,y0s,z0s,vxs,vys,vzs,c,t];
tes=filterEmitTimes[tes,t];
te=If[Length[tes]==1,tes[[1]],Null];
xe=If[NumberQ[te],x0s+vxs*te,Null];
ye=If[NumberQ[te],y0s+vys*te,Null];
ze=If[NumberQ[te],z0s+vzs*te,Null];
If[NumberQ[xe]&&NumberQ[ye]&&NumberQ[ze],AppendTo[apparentPoints,{xe,ye,ze}]];
{color,PointSize[Medium],
If[NumberQ[xe]&&NumberQ[ye]&&NumberQ[ze],{color,PointSize[Medium],Point[{xe,ye,ze}]}]
}],{i,nSources}];




(* Aberration is angle between actual and apparent position vectors of centroid source. *)
(*uCe=apparentPoints[[1]]/Norm[apparentPoints[[1]]];*)       (* normalized apparent centroid position vector *)
(*uC={x0,y0+vy*t,0.0};*)
(*uC /= Norm[uC]; *)                                                                                             (* normalized actual centroid position vector *)
(*aberration=ArcCos[uCe.uC];*)



(* "Aim the telescope" by rotating the apparent positions of the individual stars. *)
rotatedPositions=rotatePositions[apparentPoints];

(* Simulate imaging by telescope. Quiet suppresses divide by zero warnings. *)
projectedPositions=pinholeProjection[rotatedPositions,focalDistance]//Quiet;

projectedPositions*=-1;

(* Create the display *)
Grid[{{
(* Left hand display (observer frame of reference) *)
Graphics3D[{Table[Module[{x0s,y0s,z0s,vxs,vys,vzs,color},
{x0s,y0s,z0s,vxs,vys,vzs}=sources[[i]];
color=colors[[i]];
{
{color,PointSize[Medium],Point[{apparentPoints[[i,1]],apparentPoints[[i,2]],apparentPoints[[i,3]]}]}
}],{i,nSources}
],
Black,PointSize[Large],Point[{0,0,0}],  (* observer position *)
If[showVelocityVector,{Black,Arrowheads[{.02}],Arrow[{apparentPoints[[1]],apparentPoints[[1]]+{0,-10,0}}]}],
If[showAim,{Purple,Thickness[Medium],Arrowheads[{.02}],Arrow[{{0,0,0},apparentPoints[[1]]}]}],
If[showLabels,{Text[Style["telescope position",Medium],{4,-2,0.5}]}],
If[showLabels&&NumberQ[te],{Text[Style["apparent position",Red],{xe,ye,ze}]}]
},
Boxed->False,PlotRange->{{-size3D,size3D},{-size3D,size3D},{-size3D,size3D}},Axes->showAxes,AxesOrigin->{0,0,0},AxesLabel->{"X","Y","Z"},
ImageSize->600,PlotLabel->"apparent positions",Background->Lighter[LightGray],ImageMargins->0
],

(* Right hand display (telescope view)  *)
Graphics[{
Table[Module[{x0s,y0s,z0s,vxs,vys,vzs,color},
{x0s,y0s,z0s,vxs,vys,vzs}=sources[[i]];
color=colors[[i]];
{
{color,PointSize[Medium],Point[{projectedPositions[[i,1]],projectedPositions[[i,2]]}]}
}],{i,nSources}
](*Text[Style["aberration = "<>ToString[aberration],Medium,Red],Scaled[{.2,.1}]]*)
},
Axes->True,
PlotRange->{{-10,10},{-10,10}},
LabelStyle->Directive[Red],
PlotLabel->"telescope view",
Background->Nest[Darker,Blue,4],
ImageSize->600]}},Alignment->{Bottom,Bottom}
]
],

(* controls *)
{{x0,-10,"x coordinate (ly)"},-20,20,Appearance->"Labeled"},
{{yeControl,0,"y coordinate (ly)"},-20,20,Appearance->"Labeled"},
Delimiter,
{{vyin,0.001," speed (ly/y)"},0.001,0.999,Appearance->"Labeled"},
{{c,1,"speed of light (ly/yr)"},1,3,DefaultValue->1,Appearance->"Labeled"},
Delimiter,
{{spread,5," size (ly)"},1,10,Appearance->"Labeled"},
{{nSourcesDesired,100,"number of sources"},2,2000,1},
{{shape,"cube"," shape"},{"sphere","line_x","line_y","line_z","cube","circle_x","circle_y","circle_z","cluster","galaxy_x","galaxy_y","galaxy_z"}},
Delimiter,
{{applyLorentzContraction,True,"apply Lorentz contraction"},{True,False}},
Delimiter,
{{showLabels,True,"show labels"},{True,False}},
{{showAxes,True,"show axes"},{True,False}},
{{showAim,True,"show position vector"},{True,False}},
{{showVelocityVector,True,"show velocity direction"},{True,False}},
{{size3D,20,"observer frame scale"},1,50,Appearance->"Labeled"},
{{focalDistance,10,"telescope magnification"},1,20,Appearance->"Labeled"},

(* functions *)
Initialization:>(

(* select the phyiscally meaningful of the two quadratic solutions for the time of emission. *)
filterEmitTimes[tesi_,t_]:=Module[{emitTimes=tesi},

emitTimes=Select[emitTimes,Im[#1]==0&];      (* emission time must be real! *)
emitTimes=Select[emitTimes,#1≤t&];               (* emission must precede detection! *)
emitTimes=Union[emitTimes];                                      (* remove duplicates *)

(* catch pathological cases *)
If[Length[emitTimes]!=1,Beep[]; Print["emit times after filtering = ",emitTimes]];

emitTimes
];


(* apply Lorentz contraction to the "galaxy". *)
lorentzContract[nSources_,sourcesIn_,vy_,c_]:=Module[{scs=sourcesIn,lorentzFactor,yCentroidArray,deltaYs,contractedYs},

(* Apply Lorentz contraction, holding centroid fixed. Assumes vx=vz=0.  *)
lorentzFactor=Sqrt[1-vy^2/c^2];

yCentroidArray=ConstantArray[scs[[1,2]],nSources];

deltaYs=scs[[;;,2]]-yCentroidArray;
deltaYs*=lorentzFactor;

contractedYs=yCentroidArray+deltaYs;
scs[[;;,2]]=contractedYs;

scs
];


(* Simulate the observation of the object by a telescope. *)
(* Model the telescope as an ideal, aberration-free camera. *)
pinholeProjection[points_,distance_]:=Module[
{projectedPoints,desiredYaxis,myRotationMatrix,rotationAngle,rotatedPoints},

If[Length[points]<2,Return[points]];

projectedPoints=Table[{((points[[i,1]])/(points[[i,3]])),((points[[i,2]])/(points[[i,3]]))},{i,Length[points]}];
projectedPoints*distance
];



(* "Aim" the simulated telscope at "centroid" of object *)
rotatePositions[array_]:=Module[{xc,yc,thetaCentroid,mr,mrz,mrx,mrzpp},

If[Length[array]<2,Return[array]];

xc=array[[1,1]];
yc=array[[1,2]];
thetaCentroid=ArcTan[yc,xc];
mrz=RotationMatrix[thetaCentroid,{0,0,1}];
mrx=RotationMatrix[-π/2,{1,0,0}];
(*mrzpp=RotationMatrix[π/2,{0,0,1}];*)

mr=mrx.mrz;
(*mr=mrzpp.mrx.mrz;*)
(mr.#1&)/@array
];


(* solve quadratic equation for the possible emission times corresponding to a detection at time t *)
emissionTimes[x0_,y0_,z0_,vx_,vy_,vz_,c_,t_]:=Module[
{r0v={x0,y0,z0},vv={vx,vy,vz},r0,v},

r0=Norm[r0v];
v=Norm[vv];

{(1/((c^2)-(v^2)))((c^2) t+vv.r0v-√((c^2) (r0^2)-(r0^2) (v^2)+(c^2) (t^2) (v^2)+2 (c^2) t vv.r0v+(((vv.r0v))^2))),(1/((c^2)-(v^2)))((c^2) t+vv.r0v+√((c^2) (r0^2)-(r0^2) (v^2)+(c^2) (t^2) (v^2)+2 (c^2) t vv.r0v+(((vv.r0v))^2)))}
];


(* Create a spherical shell of uniformly distributed stars using Fibonacci lattice. *)
generateSpherePoints[x0_, y0_, z0_, radius_, nPoints_] := Module[{ theta, z, r, points},
 
(* Generate points on the sphere *)
points = Table[
z = 1. - 2*(i - 1)/(nPoints - 1); (* z-coordinate *)
r = Sqrt[1. - z^2]; (* Radius of the circle at height z *)
theta = 2.*Pi*i/GoldenRatio; (* Angle around the circle *)
{x0 + radius*r*Cos[theta], y0 + radius*r*Sin[theta], z0 + radius*z},
{i, 1, nPoints}
];

points
];



(* Create a cubic shell of uniformly distributed stars. *)
generateCubeEdges[x0_,y0_,z0_,edgeLength_,nPointsDesired_]:=Module[
{nPointsPerEdge,corners,edges,nPoints,pointsOnEdges,edgeVectors,pointsAlongEdges},

nPointsPerEdge=Round[nPointsDesired/12];

(*Define the 8 corners of the cube centered at {x0,y0,z0}*)
corners=edgeLength/2*Tuples[{1,-1},3];
corners+=ConstantArray[{x0,y0,z0},8];

(*Define the 12 edges of the cube as pairs of corners*)
edges={{1,2},{1,3},{1,5},{2,4},{2,6},{3,4},{3,7},{4,8},{5,6},{5,7},{6,8},{7,8}};

(*Number of points along each edge*)
nPoints=Max[nPointsPerEdge,2];

(*At least two points per edge*)(*Generate points along each edge*)pointsAlongEdges=Flatten[Table[Table[corners[[edges[[i,1]]]]+(corners[[edges[[i,2]]]]-corners[[edges[[i,1]]]])*(j-1)/(nPoints-1),{j,1,nPoints}],{i,1,Length[edges]}],1];

(*Return unique points*)
DeleteDuplicates[pointsAlongEdges]
];



RandomPlummerCluster[n_,scale_]:=Module[{r,theta,phi,x,y,z,ret},

ret=Table[
r=scale Sqrt[RandomReal[]/(1-RandomReal[])];(*Radial distance*)
theta=RandomReal[{0,Pi}];(*Polar angle*)
phi=RandomReal[{0,2 Pi}];(*Azimuthal angle*)
x=r Sin[theta] Cos[phi];
y=r Sin[theta] Sin[phi];
z=r Cos[theta];
{x,y,z},{n}
];

ret-ConstantArray[Mean[ret],Length[ret]]
];


RandomSpiralGalaxy[n_,a_,b_,diskScale_,bulgeScale_]:=Module[
{r,theta,x,y,z,diskStars,bulgeStars,phi},

(*Generate disk stars along spiral arms*)
diskStars=Table[

r=diskScale RandomReal[];

(*Exponential radial falloff*)
theta=RandomReal[{0,2 Pi}];

(*Random angular position*)
phi=theta+a Log[r];

(*Spiral arm equation*)
{r Cos[phi],r Sin[phi],RandomReal[{-0.05,0.05}]},

(*Thin disk*)
{n}
];

(*Generate bulge stars as a spherical distribution*)
bulgeStars=Table[
r=bulgeScale RandomReal[]^(1/3);
(*Spherical bulge,scaled*)
theta=RandomReal[{0,Pi}];
phi=RandomReal[{0,2 Pi}];
{r Sin[theta] Cos[phi],r Sin[theta] Sin[phi],r Cos[theta]},{n/5}];
(*Combine the disk and bulge stars*)
Join[diskStars,bulgeStars]
];


(* Create a set of uniformly distributed stars that forms a circle. *)
generateDiskPoints[x0_, y0_, z0_, radius_, nBoundaryPoints_, includeInterior_: False, nInteriorRings_: 10] :=
 Module[{boundaryPoints, interiorPoints, theta, r, ringPoints},
  
  (* Generate points along the boundary of the disk *)
  boundaryPoints = 
   Table[
    {x0 + radius*Cos[theta], y0 + radius*Sin[theta], z0},
    {theta, 0, 2*Pi - 2*Pi/nBoundaryPoints, 2*Pi/nBoundaryPoints}
    ];
  
  (* Generate uniformly distributed points inside the disk if includeInterior is True *)
  interiorPoints = 
   If[includeInterior,
    Flatten[
     Table[
      r = radius*(i/nInteriorRings);
      ringPoints = Round[nBoundaryPoints*(i/nInteriorRings)];
      Table[
       {x0 + r*Cos[theta], y0 + r*Sin[theta], z0},
       {theta, 0, 2*Pi - 2*Pi/ringPoints, 2*Pi/ringPoints}
       ],
      {i, 1, nInteriorRings}
      ], 1],
    {}
    ];
  
  (* Combine boundary and interior points *)
  Join[boundaryPoints, interiorPoints]
];

)]

Conclusions

•  I have presented a method for modeling the appearance of objects moving at relativistic speeds, accounting for Lorentz contraction and relativistic aberration.  

•  The method applies to any object that can be represented by a collection of point sources of light.  

•  Apparent rotation predicted by Terrell, Penrose and other authors is clearly visible in the simulation.  

•  Spherical objects always appear to have a circular outline to the observer, as predicted by Penrose.

•  Lorentz contraction is in principle directly observable in certain geometries where aberration effects cancel.

•  Relativistic aberration effects may be useful in deducing the velocity vectors of non-spherically symmetric objects such as galaxies.

•  Aberration effects for approaching and receding objects are very different.  Approaching objects appear stretched while receding objects appear compressed by aberration.   The stretching is in principle observable, but there are no known instances of objects approaching Earth with relativistic velocities.

•  Aberration produces apparent elongation of approaching objects and apparent compression of receding objects.

For approaching objects, aberration and Lorentz contraction effects tend to oppose one another.   For receding objects they both produce apparent compression.   

•  Traditional approaches to aberration, which treat it as an angular effect, entirely miss the most interesting case, that of purely radial motion.

•  Manipulate[] is a tour de force.

Discussion