# ParametricPlot3D colorfunction for a specific region of an object

I have created a 3D parametric plot of an oblate spheroid, and I want to apply "TemperatureMap" as a colorfunction to its surface as a gradient in the Z direction. However, I only want to apply this coloring to part of the spheroid, and to color the rest of it gray. I am attempting this by first plotting a gray spheroid, and then overlaying that object with part of another spheroid, which should be colored. Here here is my attempt:

rx=2.53*10^9;
ry=2.375*10^9;
rz=8.65*10^8;
Req=1.46*10^9;
obl=0.19;
a=4.4879*10^9;

grayplot=ParametricPlot3D[Req/Sqrt[Sin[θ]^2+Cos[θ]^2/(1-obl)^2]*{Cos[ϕ]*Sin[θ],Sin[ϕ]*Sin[θ],Cos[θ]},{ϕ,0,2π},{θ,0,π},
Mesh->None,
ColorFunction->Function[{x,y,z},Glow[Darker[Gray, 0.4]]],
ColorFunctionScaling->False,
Boxed->False,
Axes->False,
PlotPoints->35,
Lighting->{{"Point",White,{0,0,0}}}
];

colorplot=ParametricPlot3D[(Req*1.02)/Sqrt[Sin[θ]^2+Cos[θ]^2/(1-obl)^2]*Boole[rx*Cos[ϕ]*Sin[θ]+ry*Sin[ϕ]*Sin[θ]+rz*Cos[θ]>=Req/Sqrt[Sin[θ]^2+Cos[\Theta]]^2/(1-obl)^2]]*{Cos[ϕ]*Sin[θ],Sin[\Phi]]*Sin[θ],Cos[θ]},{ϕ,0,2π},{θ,0,π},
ColorFunction->Function[{θ},If[θ<π/2,Glow@ColorData[{"TemperatureMap",{0,Req}}][#3]&},{Glow@       ColorData[{"TemperatureMap",{-Req,0}}]}]],
ColorFunctionScaling->False,
Mesh->None,
PlotPoints->24
];

planetgraphic=Graphics3D[{Blue,Sphere[{rx, ry, rz},7.15*10^9]}];

Show[{grayplot,colorplot,planetgraphic},
ViewPoint->{10,0,0},
ViewAngle->1/8,
ViewVertical->{1, 0, 0},
Boxed->False,
ImageSize->600,
Background->Black,
PlotRange->{{-1.1 a, 1.1 a}, {-1.1 a, 1.1 a},{-1.1 a, 1.1 a}}]


rx, ry, rz, and a, are all just constants based on my system. Req is the equatorial radius and obl is the spheroid's oblateness factor.

This code produces . What I want is something like this

which I colored by plotting thousands of individual points in Graphic3D (highly inefficient). Basically, I want to recreate that image but with a nice smooth parametric plot. I want the temperature gradient to scale from equator to pole for whatever region of the object I have colored in.

There is an obvious current problem with the lighting in the first image (which is why the colored region is black), but even beyond that the ColorFunction in the colorplot is not working properly. If I turn on ambient lighting, the "colored region" is just pure white.

• For clarity, can you give the specific ranges where your surface should be gray, and where it should have the color gradient? – J. M. is in limbo Feb 20 '16 at 19:36
• Yes, I forgot to explain. The colorplot is limited by the Boole function in the main argument. The (phi,theta) that satisfy the inequality inside of that Boole function describe the region that I want colored. – ahle6481 Feb 20 '16 at 19:41
• What do you mean by Z- direction. IN Oblate Spheroidal coordinate is a true 3D curvilinear system. You might be interested in this post. mathematica.stackexchange.com/questions/149401/… – Jose Enrique Calderon Dec 21 '18 at 21:13

### The case in hand

I would use MeshShading and create a single, unified graphic.

With[{ϕ0 = .5,
rx = 2.53*10^9,
ry = 2.375*10^9,
rz = 8.65*10^8,
Req = 1.46*10^9,
obl = 0.19},
ParametricPlot3D[
Req/Sqrt[Sin[θ]^2 +
Cos[θ]^2/(1 - obl)^2]*{Cos[ϕ]*Sin[θ],
Sin[ϕ]*Sin[θ], Cos[θ]}, {ϕ, 0,
2 π}, {θ, 0, π},
ColorFunction ->
Function[{x, y, z, ϕ, θ},
Glow@ColorData[{"TemperatureMap", Req {0, 1}}][Req - Abs@z]],
ColorFunctionScaling -> False, Boxed -> False, Axes -> False,
AxesLabel -> {x, y, z},
MeshFunctions -> {Function[{x, y,
z, ϕ, θ}, (rx*Cos[ϕ]*Sin[θ] +
ry*Sin[ϕ]*Sin[θ] + rz*Cos[θ]) - (Req/
Sqrt[Sin[θ]^2 + Cos[θ]^2/(1 - obl)^2])]},
Mesh -> {{0}},
PlotPoints -> 35, Lighting -> {{"Point", White, {0, 0, 0}}}]
]


One can coordinate the MeshFunctions and Mesh settings to create just about any combination of shaded and colored regions.

### Update: Further remarks on the use of Mesh etc.

Note that FEM and mesh regions are different and newer than Mesh and its related options of plotting. Conceptually they are related, but functionally within Mathematica they don't seem closely related.

Basic idea. For each mesh function mf, a set of values is prescribed or automatically chosen by Mathematica, according to option Mesh. (See MeshFunction for a description of the arguments, which are different for each plotter.) A line is drawn on the surface mf == c for each value c in the set. This will divide the surface up to regions, which will be colored according to the settings of MeshShading and other style options appropriate to the plotter. The value of MeshShading should be an array of depth equal to the number of mesh functions.

There are two special settings for a MeshShading entry: Automatic renders the surface according to the other styling directives. None means that no surface will be render (there will be a hole).

It is sometimes important to realize that the mesh lines are created by linearly interpolating along the boundaries of the polygons of the surface constructed by the plotter. If the numerical accuracy of the mesh lines is important, then it usually has to be addressed via PlotPoints. The option MaxRecursion may help, but only accidentally: The criteria for recursive subdivision depend on the curvature of the surface and not on the curvature of the mesh lines.

One mesh function. If a single mesh function is used, then the surface will be divided up into bands/stripes and possibly concentric annuli, depending on the function. The value of MeshShading should be of array-depth 1, that is, a simple, flat list. If $c_1, c_2,\dots$ form the set of mesh values, then regions for which the mesh function has a value less than $c_1$ are given the first color, those for which the value is between $c_1$ and $c_2$ are given the second color, and so on. The cycle around to the beginning when the end of the list is reached. If, say, dir is a list of 4 graphics directives and there are five mesh values, then MeshShading -> dir would color the surface according to the following scheme:

Example: Note that when the shading is Automatic, the PlotStyle texture colors the plot; when it is none, the white background shows through.

ParametricPlot[{u , v + (u - u^2)/4}, {u, 0, 1}, {v, 0, 1},
PlotStyle -> {Opacity[1], Texture[ExampleData[{"TestImage", "Mandrill"}]]},
Mesh -> 7, MeshShading -> {Automatic, Red, None},
MeshFunctions -> {Function[{x, y, u, v}, u]}]


Two mesh functions. The remarks about one mesh function still apply, except that the value of MeshShading should be of array depth 2 (i.e. a matrix), and the surface will usually be divided up into polygons with curved boundaries. If we have MeshShading -> dir, then when the first mesh function increases past a mesh value, the color will change from dir[[i, j]] to dir[[i + 1, j]]; when the second mesh function passes a mesh value, the color will change from dir[[i, j]] to dir[[i, j + 1]]. Of course the indices will wrap around when the corresponding edge of the matrix is reached.

Example. Suppose MeshShading is set to

MeshShading -> {{Automatic, Red}, {Blue, Pink}, {Lighter@Blue, None}}


Then we get, for the default mesh functions, u and v:

ParametricPlot[{u , v + (u - u^2)/4}, {u, 0, 1}, {v, 0, 1},
PlotStyle -> {Opacity[1], Texture[ExampleData[{"TestImage", "Mandrill"}]]},
Mesh -> {5, 4},
MeshShading -> {{Automatic, Red}, {Blue, Pink}, {Lighter@Blue, None}}]


Since v increases as you go up the image, the colors appear reversed (top to bottom) from how they appear in the matrix rows.

More mesh functions. As the number of functions and depth of the array dir of directives in MeshShading increases, things get more complicated, but they work in the same way. When a mesh function crosses a boundary value, the corresponding index in

dir[[i1, i2, i3, ...]]


changes. It is rare that one needs to consider this case though.

Complicated mesh functions. The complicated pattern produced in the plot on the left below by the mesh functions

mf1 = Function[{x, y, u, v}, Cos[5 x] - y];
mf2 = Function[{x, y, u, v}, x - Cos[5 y]];


is hard to describe. However by plotting each mesh function, we can see where it crosses the different mesh values by Mathemetica automatically for the settings of Mesh to 5 and 4 respectively.

One can also observe in the lighter blue region below the eye (around the vertex of the cosine minimum), the effect of linear interpolation. To get a less jagged appearance, increase PlotPoints.

ParametricPlot[{u , v + (u - u^2)/4}, {u, 0, 1}, {v, 0, 1},
PlotStyle -> {Opacity[1], Texture[ExampleData[{"TestImage", "Mandrill"}]]},
Mesh -> 5,
MeshFunctions -> {mf1, mf2},
MeshShading -> {{Automatic, Red}, {Blue, Pink}, {Lighter@Blue, None}}]

ParametricPlot3D[{u , v + (u - u^2)/4, mf1[u , v + (u - u^2)/4, u, v]},
{u, 0, 1}, {v, 0, 1},
PlotStyle -> {Opacity[1], Texture[ExampleData[{"TestImage", "Mandrill"}]]},
Mesh -> 5,
MeshFunctions -> {mf1},
MeshShading -> {{Automatic, Red}, {Blue, Pink}, {Lighter@Blue, None}}[[2, All]],
BoxRatios -> {1, 1, 0.6}]

ParametricPlot3D[{u , v + (u - u^2)/4, mf2[u , v + (u - u^2)/4, u, v]},
{u, 0, 1}, {v, 0, 1},
PlotStyle -> {Opacity[1], Texture[ExampleData[{"TestImage", "Mandrill"}]]},
Mesh -> 5,
MeshFunctions -> {mf2},
MeshShading -> {{Automatic, Red}, {Blue, Pink}, {Lighter@Blue, None}}[[All, 2]],
BoxRatios -> {1, 1, 0.6}]

• This didn't come to me as I don't think about mesh, not intimately familiar with it. Here you seem to use it like a second color function but I normally just think of it as those grid lines on top of a 3d graphic. I know it is intimately connected with the new region functions as well. Do you know of a resource where I could learn more about what a mesh is? Sorry for the vague question lol – Jason B. Feb 24 '16 at 20:57
• @JasonB I mean to expand my answer at some point to explain what's going on here. I just don't have time today. I figured it out the hard way, by playing with the examples in the docs. Basically, each time a Mesh value is crossed, the position in the MeshShading array shifts one (and cycles around to the beginning if the end is reached). – Michael E2 Feb 24 '16 at 21:57
• @JasonB I added some :) explanation. AFAIK, it's not explained in the documentation, except from what you can infer from examples. – Michael E2 Feb 28 '16 at 2:12
• Thanks again for the lesson on MeshShading, as soon as I saw this question, I knew this was the method to use. – Jason B. May 4 '16 at 10:24

You were mixing your pure function definition types in your ColorFunction. The safest way to make a color function for a ParametricPlot3D is of the form Function[{x,y,z,phi,theta},......] - you had a Function[{theta}, ....#3 &...] in which you wanted the #3 to apply to the z coordinate.

rx = 2.53*10^9;
ry = 2.375*10^9;
rz = 8.65*10^8;
Req = 1.46*10^9;
obl = 0.19;
a = 4.4879*10^9;

grayplot =
ParametricPlot3D[
Req/Sqrt[
Sin[θ]^2 + Cos[θ]^2/(1 - obl)^2]*{Cos[ϕ]*
Sin[θ], Sin[ϕ]*Sin[θ],
Cos[θ]}, {ϕ, 0, 2 π}, {θ, 0, π},
Mesh -> None,
ColorFunction -> Function[{x, y, z}, Glow[Darker[Gray, 0.4]]],
ColorFunctionScaling -> False,
Boxed -> False,
Axes -> False,
PlotPoints -> 35,
Lighting -> {{"Point", White, {0, 0, 0}}}
];

colorplot =
ParametricPlot3D[(Req*1.02)/
Sqrt[Sin[θ]^2 + Cos[θ]^2/(1 - obl)^2]*
Boole[rx*Cos[ϕ]*Sin[θ] +
ry*Sin[ϕ]*Sin[θ] + rz*Cos[θ] >=
Req/Sqrt[
Sin[θ]^2 + Cos[θ]^2/(1 - obl)^2]]*{Cos[ϕ]*
Sin[θ], Sin[ϕ]*Sin[θ],
Cos[θ]}, {ϕ, 0, 2 π}, {θ, 0, π},
PlotRange -> Chartingget3DPlotRange@grayplot,
ColorFunction -> Function[{x, y, z, ϕ, θ},
Glow@ColorData[{"TemperatureMap", {Req, 0}}][Abs@z]],
ColorFunctionScaling -> False,
Mesh -> None,
PlotPoints -> 100
];

Show[{grayplot, colorplot},
ViewPoint -> {10, 0, 0},
ViewAngle -> 1/8,
ViewVertical -> {1, 0, 0},
Boxed -> False,
ImageSize -> 600,
Background -> Black,
PlotRange -> {{-1.1 a, 1.1 a}, {-1.1 a, 1.1 a}, {-1.1 a, 1.1 a}}]


I left out planetgraphics as it isn't relevant to the question at hand, and I increased the number of PlotPoints in the colorplot to get a less jagged edge.

In fixing the color function, I may have changed it to something different than what you were looking for. The syntax you had in the third argument to the If is confusing, so let me know if you wanted something different.

I think you could get all this out of one plot, if you just make the color function a little more complicated so as to incorporate your Boole conditions,

ParametricPlot3D[
Req/Sqrt[Sin[θ]^2 + Cos[θ]^2/(1 - obl)^2]*{Cos[ϕ]*
Sin[θ], Sin[ϕ]*Sin[θ],
Cos[θ]}, {ϕ,
0, 2 π}, {θ, 0, π},
Mesh -> None,
ColorFunction -> Function[{x, y, z, ϕ, θ},
If[rx*Cos[ϕ]*Sin[θ] + ry*Sin[ϕ]*Sin[θ] +
rz*Cos[θ] >=
Req/Sqrt[Sin[θ]^2 + Cos[θ]^2/(1 - obl)^2],
Glow@ColorData[{"TemperatureMap", {Req, 0}}][Abs@z],
Glow[Darker[Gray, 0.4]]
]
], PlotPoints -> 300,
ColorFunctionScaling -> False,
Boxed -> False,
Axes -> False,
Lighting -> {{"Point", White, {0, 0, 0}}},
ViewPoint -> {10, 0, 0},
ViewAngle -> 1/8,
ViewVertical -> {1, 0, 0},
Boxed -> False,
ImageSize -> 600,
Background -> Black,
PlotRange -> {{-1.1 a, 1.1 a}, {-1.1 a, 1.1 a}, {-1.1 a,
1.1 a}}]


But I had to set the PlotPoints` to a very high level to get a decent output, and it still isn't as good as the plot above,

• Thank you so much for your help! I wish I could divide the bounty between you and the answer below. However, that other answer is the method that I will actually use. This was a great answer though! – ahle6481 Feb 24 '16 at 20:34