Let me give a model solution which can easily be adapted.
2 dimensions
Consider a photon which moves in the x-y-plane, starting at time t = 0 in the origin {0,0} and moving towards the positive x-axis. At each tick of the clock, corresponding to a constant distance 1 travelled by the photon, the photon will experience a scattering event which leads to a deviation of the original direction by an angle "a". The angle "a" is a random variable with a certain given distibution dist.
We ask for the random trajectory of the photon after a given number n of scattering events. We would like to see several of those trajectories in one picture.
For simplicity, I assume that the angle distribution is a symmetric NormaDistribution[0,s]
with a given width s.
dist[s_] = NormalDistribution[0, s];
The scattering angle is then
a[s_] := RandomVariate[dist[s]]
A sample of angles is
Array[a[1] &, 10];
We need to cumulate these angeles because the new one gives the deviation with respect to the current one:
Accumulate[%];
Now the x-y-coordinates of each subtrajectory of one tick of the clock are
{Cos[#], Sin[#]} & /@ %;
This again needs to be cumulated in order to get the coordinates of the trajectory:
Accumulate[%];
Now we have one trajectory an plot it
ListLinePlot[%]
(* 150304 plot_1_trajectory.jpg *)

Putting everything together gives us a trajectory by
tr[s_, n_] :=
Accumulate[{Cos[#], Sin[#]} & /@ (Accumulate[Array[a[s] &, n]])]
Plotting a bunch of k = 11 trajectories, each of length n = 100 with a width s of the normal Distribution
With[{s = 0.15, n = 100, k = 11}, ListLinePlot[Table[tr[s, n], {k}]]]
(* 150304 plot_k_trajectories.jpg *)

Varying s from small values (s = .01, prefers forward directions) to larger ones (s = 1, gives strong scattering) shows interesting outcomes.
3 dimensions
dist[s_] =
NormalDistribution[0, s]; (* distribution of scattering angle *)
a[s_] := RandomVariate[dist[s]]; (* scattering angle *)
f := 2 \[Pi] Random []; (* rotation angle *)
k = 30; (* length of a trajectory *)
m = 8; (* number of trajectories *)
s = 0.3; (* width of distribution of scattering angle *)
(* run starts here *)
tbtrk = Table[af = Accumulate[Array[{a[s], f} &, k]];
{Hue[j/m],
trk = Line @@ {Accumulate[
Prepend[Table[{Sin[#1] Cos[#2], Sin[#1] Sin[#2],
Cos[#1]} &@(Sequence @@ af[[i]]), {i, 1,
Length[af]}], {0, 0, 0}]]}}, {j, 0, m}];
d = 5; (* size of the display box *)
Graphics3D[{{PointSize[0.02], Point[{0, 0, 0}]}, tbtrk},
PlotRange -> {{-d, d}, {-d, d}, {-5, k}}, Boxed -> True]
(* 150304_plot _ 3D_k _trajectories.jpg *)
