I would like to make the surface of my 3D plot have reflection and shadows, like in this plot:

(from here).

Is there any way to do this with Mathematica? If not, what software can make plots like these?

• You will want to use a combination of Specularity and Lighting. Perhaps this answer 221946 will give some ideas. – Tim Laska May 16 '20 at 0:48

Here is an example based on the answer alluded to in my above comment. I will create some fake data of Gaussian's for illustration purpose and apply Textures, Specularity, and Point Lighting.

f[x_, y_] := E^(-x^2 - y^2)
Plot3D[1.5 f[x, y] + f[x + 2, y + 2] + f[x - 2, y - 2] +
1/2 (f[x + 2, y - 2] + f[x - 2, y + 2]), {x, -4, 4}, {y, -4, 4},
Mesh -> 10, ImageSize -> 400, Ticks -> None, PlotRange -> All,
TextureCoordinateFunction -> ({#3, #3} &),
PlotStyle ->
Directive[Specularity[White, 50],
FaceForm[Texture[ColorData["SouthwestColors", "Image"]]]],
Axes -> False, Lighting -> {{"Point", White, {4, 0, 8}}},
MeshFunctions -> {#3 &}, PlotPoints -> 50, Boxed -> False,
ViewPoint -> Top, PerformanceGoal -> "Quality"]


To investigate directional lighting we can borrow and adapt the code from the Wolfram U course Advanced 3D Graphics in the Wolfram Language to create a nice dynamic visualization.

Manipulate[
With[{pt =
6 {Cos[θ] Cos[ρ], Sin[θ] Cos[ρ],
Sin[ρ]}, center = {0, 0, 0}},
Show[
Plot3D[
1.5 f[x, y] + f[x + 2, y + 2] + f[x - 2, y - 2] +
1/2 (f[x + 2, y - 2] + f[x - 2, y + 2]), {x, -4, 4}, {y, -4, 4},
Mesh -> 10, ImageSize -> 400, Ticks -> None, PlotRange -> All,
TextureCoordinateFunction -> ({#3, #3} &),
PlotStyle ->
Directive[Specularity[White, 50],
FaceForm[Texture[ColorData["SouthwestColors", "Image"]]]],
Axes -> False, Lighting -> {{"Directional", color, {pt, center}}},
MeshFunctions -> {#3 &}, PlotPoints -> 50, Boxed -> False],
Graphics3D[{
Specularity[White, 50],
Dynamic[
If[show, {color,
Style[Cylinder[{pt, 1.05 pt}, .5],
Lighting -> "Neutral"]}, {}]],
Opacity[.1], Sphere[center, 6]},
Lighting -> {{"Directional", color, {pt, center}}}],

dirplane,

Graphics3D[
Dynamic[If[
show, {Gray, Dashed,
Line[N@Table[

6 {Cos[t] Cos[ρ], Sin[t] Cos[ρ], Sin[ρ]}, {t,
0, 2 π, (2 π)/50}]]}, {}]]],

Graphics3D[
Dynamic[If[
show, {Gray, Dashed,
Line[N@Table[
6 {Cos[θ] Cos[s], Sin[θ] Cos[s], Sin[s]}, {s,
0, Pi}, {s, 0, 2 π, (2 π)/50}]]}, {}]]],

PlotRange -> {{-6, 6}, {-6, 6}, {.5, 6}}, PlotRangePadding -> .5,
Boxed -> False,
ViewAngle -> .283, ViewCenter -> {{.5, .5, .5}, {.5, .6}}
]
],
{color, LightBlue},
Delimiter,
"Direction",
{{θ, 0}, 0, 2 π}, {{ρ, π/4}, 0, π},
Delimiter,
{{show, True, "show light"}, {False, True}}, ControlPlacement -> Left,
Initialization :> {
Clear[dirplane];
dirplane =
ParametricPlot3D[{x, y, 0}, {x, -4, 4}, {y, -4, 4}, Mesh -> None,
BoundaryStyle -> Gray, PerformanceGoal -> "Quality"]
}]


# Attempt at a Banded Custom Color Map

The sample image provided in the OP has a color map that as far as I can tell is not in the ColorData collection. The colors in the image seem to go from Purple->Green->Yellow similar to the perceptually uniform Viridis color map. The following workflow shows how to create a custom ColorData["xxxx","Image"] function using a color map you like from the web and create a banding effect by brightening the colors. By using MeshShading, one can alternate between bright and normal colors creating a banding effect.

(* Create Viridis ColorFunction *)
img = ImageCrop@
Import["http://makie.juliaplots.org/stable/assets/viridis.png"];
dims = ImageDimensions[img];
colorrgb1 = {ImageData[img][[IntegerPart@(dims[[2]]/2), All]]};
(* Brighten Slightly *)
colorrgb2 = 1.1*{ImageData[img][[IntegerPart@(dims[[2]]/2), All]]};
(* Look at InputForm of ColorData["SouthwestColors","Image"]*)
cmfn = Graphics[{Raster[#, {{0, 0}, {1, 1}}]},
{ImageSize -> 250, ContentSelectable -> False,
AspectRatio -> 1/8,
PlotRange -> {{0, 1}, {0, 1}}}] &;
(* Create ColorData suitable for textures *)
cmimg1 = cmfn@colorrgb1;
cmimg2 = cmfn@colorrgb2;
(* Create Test Function *)
f[x_, y_] := E^(-x^2 - y^2)
fs[x_, y_] := (6 f[x, y] + 4 (f[x + 2, y + 2] + f[x - 2, y - 2]) +
9 (f[x + 2, y - 2] + f[x - 2, y + 2]))/9
(* Alternating FaceForms for Mesh *)
ff1 = FaceForm[Texture[cmimg1]];
ff2 = FaceForm[Texture[cmimg2]];
(* Plot the function *)
ParametricPlot3D[{x, y, fs[x, y]}, {x, -4, 4}, {y, -4, 4}, Mesh -> 25,
ImageSize -> 400, Ticks -> None, PlotRange -> All,
TextureCoordinateFunction -> ({#3, #3} &), Axes -> False,
Lighting -> {{"Directional", White, {{10, 5, 10}, {0, 0, 0}}}},
MeshFunctions -> {#3 &}, PlotPoints -> 75, Boxed -> False,
ViewPoint -> Top, PerformanceGoal -> "Quality",
Directive[Specularity[White, 80], ff1]}]


# Update Using Custom Color Maps Following @JM's Advice

This Thread contains a number of custom color maps including those of the Vidiris family. Here is a workflow using color map that goes Red->Blue->Green->Yellow that appears to match the colors better.

(* Create Viridis ColorFunction *)
ClearAll[MPLColorMap]
<< "http://pastebin.com/raw/pFsb4ZBS";
colorrgb1 = {MPLColorMap["EricsRdBuGnYl2"] /@ Subdivide[1, 50]} /.
RGBColor[x__] :> List[x];
(* Darken Slightly *)
colorrgb2 = 0.85*colorrgb1;
(* Look at InputForm of ColorData["SouthwestColors","Image"]*)
cmfn = Graphics[{Raster[#, {{0, 0}, {1, 1}}]},
{ImageSize -> 250, ContentSelectable -> False,
AspectRatio -> 1/8,
PlotRange -> {{0, 1}, {0, 1}}}] &;
(* Create ColorData suitable for textures *)
cmimg1 = cmfn@colorrgb1;
cmimg2 = cmfn@colorrgb2;
(* Create Test Function *)
f[x_, y_] := E^(-x^2 - y^2)
fs[x_, y_] := (6 f[x, y] + 4 (f[x + 2, y + 2] + f[x - 2, y - 2]) +
9 (f[x + 2, y - 2] + f[x - 2, y + 2]))/9
(* Alternating FaceForms for Mesh *)
ff1 = FaceForm[Texture[cmimg1]];
ff2 = FaceForm[Texture[cmimg2]];
(* Plot the function *)
ParametricPlot3D[{x, y, fs[x, y]}, {x, -4, 4}, {y, -4, 4}, Mesh -> 25,
ImageSize -> 400, Ticks -> None, PlotRange -> All,
TextureCoordinateFunction -> ({#3, #3} &), Axes -> False,
Lighting -> {{"Directional", White, {{10, 5, 20}, {0, 0, 0}}}},
MeshFunctions -> {#3 &}, PlotPoints -> 75, Boxed -> False,
ViewPoint -> Top, PerformanceGoal -> "Quality",