I am not sure exactly what you mean by "discrete 3d set of plots for magnetic field lines".
Is this the correct WL definition of the equation?
r[?_, ?_] := Module[{rn = 1, ?0 = 0, ? = -0.18}, rn Cos[? + ? Cos[? + ?0]]^2]
If so, you can use Table to generate a list of values of the function for a range of ? and ? that you can then plot. e.g.
Table[{?, ?, r[?, ?]}, {?, 0, 2 Pi, Pi/12}, {?, 0, 0.4, 0.1}];
