很好玩么

Table[Table[
SphericalPlot3D[{Re[
SphericalHarmonicY[l, m, \[Theta], \[Phi]]]^2, -Re[
SphericalHarmonicY[l, m, \[Theta], \[Phi]]]^2}, {\[Theta], 0,
2 \[Pi]}, {\[Phi], 0, 2 \[Pi]}, PlotRange -> All ,
PlotPoints -> 70,
PlotStyle ->
Directive[Orange, Opacity[0.7], Specularity[White, 10]],
Mesh -> None, Boxed -> False, Axes -> False], {m, -l, l, 1}], {l,
0, 3, 1}]
mathematica有可能崩溃（当然不关我的事）

\[Psi][n_,l_,m_,r_,\[Theta]_,\[Phi]_]=E^(-((r a)/n)) LaguerreL[-1-l+n,1+2 l,r/(n a)] Sqrt[(2/(n a))^3 (-1-l+n)!/(2n (((l+n)!)^3) )] (r/(n a))^l Re[SphericalHarmonicY[l, m,\[Theta],\[Phi]]];
a=0.53;
Data={};
For[\[Theta]=0,\[Theta]<2\[Pi],\[Theta]=\[Theta]+0.02\[Pi],
For[\[Phi]=0,\[Phi]<2\[Pi],\[Phi]=\[Phi]+0.02\[Pi],
Do[r=RandomReal[{0,20}];
w=RandomReal[{0,0.003}];
If[w<=Abs[\[Psi][2,1,0,r,\[Theta],\[Phi]]]^2,AppendTo[Data,{r Sin[\[Theta]]Cos[\[Phi]],r Sin[\[Theta]]Sin[\[Phi]],r Cos[\[Theta]]}]]
,{10}]
]
]
Show[Graphics3D[{PointSize[Small],Hue[1/3,1,1,.2],Point[Data]}],ViewPoint->Front,Boxed->False]
\[Theta]=.;\[Phi]=.;r=.;w=.;