you will immediately notice that apart from the phong, not too much can be seen, for the simple reason that the surface is not coloured, so it is very hard to get a feeling of the heights and depths. Using pm3d, we can colour the surface, but we can't add phong. So, the question is, how do we combine the figure above with this one:
While it appears to be trivial, it turns out to be rather hard. But will manage, don't worry, and at the end of the day, we will have figures that look like these:
Understanding the colouring scheme
To begin, we have got to understand how the colouring in pm3d is done. When you plot without specifying the colour scheme, you are given the default, which looks like this (look at the colour box at the bottom, don't worry about the curves)
If you dig a bit deeper, e.g., by issuing the help command
you will be presented with cryptic messages, like
Some nice schemes in RGB color space 7,5,15 ... traditional pm3d (black-blue-red-yellow) 3,11,6 ... green-red-violet 23,28,3 ... ocean (green-blue-white); try also all other permutations 21,22,23 ... hot (black-red-yellow-white) 30,31,32 ... color printable on gray (black-blue-violet-yellow-white) 33,13,10 ... rainbow (blue-green-yellow-red) 34,35,36 ... AFM hot (black-red-yellow-white)
So, what does this all mean?
When invoking pm3d, first, the colour range (cbrange) is calculated (or specified), and then, using the colour range, all values of the function z(x,y) are mapped to the [0:1] interval. Using this value, which is called
set palette rgbformulae 7, 5, 15
(This is the default, so you actually haven't got to invoke this command...), you assign function 7 to red(gray), function 5 to green(gray) and function 15 to blue(gray). But then, how on Earth does one figure out which function to use? Well, the easy solution is to issue
show palette rgbformulae
upon which, this list is shown:
* there are 37 available rgb color mapping formulae: 0: 0 1: 0.5 2: 1 3: x 4: x^2 5: x^3 6: x^4 7: sqrt(x) 8: sqrt(sqrt(x)) 9: sin(90x) 10: cos(90x) 11: |x-0.5| 12: (2x-1)^2 13: sin(180x) 14: |cos(180x)| 15: sin(360x) 16: cos(360x) 17: |sin(360x)| 18: |cos(360x)| 19: |sin(720x)| 20: |cos(720x)| 21: 3x 22: 3x-1 23: 3x-2 24: |3x-1| 25: |3x-2| 26: (3x-1)/2 27: (3x-2)/2 28: |(3x-1)/2| 29: |(3x-2)/2| 30: x/0.32-0.78125 31: 2*x-0.84 32: 4x;1;-2x+1.84;x/0.08-11.5 33: |2*x - 0.5| 34: 2*x 35: 2*x - 0.5 36: 2*x - 1 * negative numbers mean inverted=negative colour component * thus the ranges in `set pm3d rgbformulae' are -36..36
So, with these in mind, the default pm3d colouring scheme is sqrt(gray) for red, x^3 for green and sin(360x) for blue. These three functions are shown in the previous figure.
So far, so good, but how do we turn all this into a phong on a general surface? If you recall from the post on phonged surfaces, what we did was to make the colour whiter when the value of the Gaussian function was high, and redder, when it was small. This worked very well for a single colour, because in that case, we simply fixed the red function at 1, and changed the value of the green and blue components according to the value of the Gaussian. But if we want to retain the pm3d colouring and apply phong at the same time, we have a bit of a problem here, because all of a sudden, the colour of a point on the surface will depend on two quantities: the value if z(x,y), and the value of f(x,y,z), which is our Gaussian function. However, the colour functions take only one argument, and nothing more. How will we parametrize a 2D space with only one variable? Well, we will discretise it in a smart way. Let us suppose, that we want to have only 3 white levels, and let us suppose, for a moment, that we are not restricted to the [0:1] interval for gray, but it can run between 0 and 4. Then, using only the red function for now, we will take the following scheme:
red(x) = sqrt(x)*3/3 + 0, if x is in [0:1] red(x) = sqrt(x-1)*2/3 + 1/3, if x is in [1:2] red(x) = sqrt(x-2)/3 + 2/3, if x is in [2:3] red(x) = sqrt(x-3)*0 + 3/3, if x is in [3:4]
So, what happens here is that the red component has the same functional behaviour (namely, sqrt(x)) in all four intervals, but the base value is shifted upwards (to the white), while the amplitude is decreased in such a way that red(x) is always in the [0:1] interval. Obviously, we can play the same trick with the green and blue components. What we should notice is that in each interval, the argument of the original function is the fractional part of x, and that the coefficient of sqrt(x) plus the base line adds up to one. Therefore, we can conveniently write
red(x) = RED(frac(x))*(1-floor(x)/N) + floor(x)/N
where frac(x) is just the fractional part of x,
frac(x) = x - floor(x)
and RED(x) is the original colouring function (sqrt(x) in the example above). N was the number of white levels (3 above) that we wanted to have. Now we have almost everything at our disposal in order to add phong to an arbitrary surface. The last thing to do is to make sure that the argument of red(x) etc. falls into the [0:1] interval. We achieve this buy expanding the colour range with respect to the zrange (which is the interval of possible z values). With these in mind, the gnu file would read as follows:
reset unset colorbox unset key unset border unset tics f(x,y) = sin(x*x+y*y)/(x*x+y*y) A = 100.0 zi = -0.4 za = 1.0 ci = zi ca = zi + (za-zi)*A set isosample 150, 150 set cbrange [ci:ca] set xrange [-3:3] set yrange [-3:3] set zrange [zi:za] set table 'phong2.dat' splot f(x,y) unset table # These are the functions of the original colour palette r(x) = sqrt(x) g(x) = x**3.0 b(x) = q(x<0.5?sin(x*2.0*pi):0.0) # The Gaussian phong h(x,a,b) = exp(-(x-a)*(x-a)/b/b) # These are helper functions R(x) = A*x - A*floor(x) # Fractional part of x L(x) = floor(x) # This is just a short-hand for the floor function F(x) = L(A*x)/(ca-ci) # The floor function divided by the colour range # This is where we want to place the white spot x0 = -0.3 y0 = -1.0 z0 = f(x0, y0) # And these are the new colour scheme red(gray) = r(R(gray))*(1.0-F(gray))+F(gray) green(gray) = green(R(gray))*(1.0-F(gray))+F(gray) blue(gray) = blue(R(gray))*(1.0-F(gray))+F(gray) set palette model RGB functions red(gray), green(grey), blue(gray) splot 'phong2.dat' u 1:2:3:((A-1)*(za-zi)*h($3,z0,0.08)*h($1,x0,0.3)*h($2,y0,0.3)+$3)) w pm3d
That is, first we decide that we want to have 100 white levels, this is A = 100. Then we fix the zrange according to the function that we want to plot, and after that, we specify the colour range. Here we shouldn't forget that have got to stretch the colour range to 100 times the zrange. Next, we identify the functions that we will need. For a start, we begin with the default pm3d scheme, whose functions we have already discussed above. Then a couple of helper functions follows, with the definition of the new colour scheme. Immediately before plotting the surface, we set the palette, and then we call the file that we produced a couple of lines earlier. We plot the surface using the x, y, z values, and for the colouring, we stretch the range of gray as we discussed above: the value of z(x,y) (this is the 3rd column, $3) will determine the base colour, and the value of the Gaussian will determine how much white tinge this base colour will have.
Well, there is not too much more to it. Once you decided your colour scheme and the number of white levels, the rest should be plain sailing.