# Intensity pattern [closed]

I'm trying to get the following intensity pattern from this paper

In the Newtonian case, the equations for the intensity (eq. 9 and 10 in the article), I put them in Mathematica with the change $$b^2 = x^2 + y^2$$ because the intensity is circularly symmetric with the parameter impact $$b$$ as radius (page 8 in this article)

H[x_, y_] := 10^3 1/(32 \[Pi] (Sqrt[x^2 + y^2])^3);
G[x_, y_] := 10^3 (1/(32 \[Pi]^2 (Sqrt[x^2 + y^2])^3)) ArcTan[Sqrt[x^2 + y^2]/Sqrt[4 - (Sqrt[x^2 + y^2])^2]] -10^3 Sqrt[4 - (Sqrt[x^2 + y^2])^2]/(128 \[Pi]^2 (Sqrt[x^2 + y^2])^2);
combinedFunction[x_, y_] := Piecewise[{{G[x, y], x^2 + y^2 < 2}, {H[x, y], x^2 + y^2 > 2}}]
DensityPlot[combinedFunction[x, y], {x, -10, 10}, {y, -10, 10}, ColorFunction -> "SunsetColors", Exclusions -> All,PlotPoints -> 100, PlotLegends -> Automatic]


and I get

However, it does not look like the image in the article. Could anyone advise me or give me a suggestion? You should obtain a circumference of radius $$2$$ as in the paper. I'm also trying to get the intensity pattern for a more general case where numerical integration was performed (I actually asked the question here) but I can't find how to make the change $$b^2 = x^2 + y^2$$. I appreciate suggestions or observations

• Your plot is only showing a range up to ~0.2, and showing everything above that as white. Add PlotRange -> {0, 1.2}. Nov 3, 2023 at 22:02
• Also, shouldn't the combinedFunction be Piecewise[{{G[x, y], Sqrt[x^2 + y^2] < 2}, {H[x, y], Sqrt[x^2 + y^2] > 2}}]? Nov 3, 2023 at 22:03
• @MelaGo You're right. I honestly didn't notice that detail. Thank you so much! Nov 3, 2023 at 22:26