Wednesday 29 July 2009

Adding shadow to the key

I think, it changes the impression a graph can make by quite a lot, if it looks sort of three-dimensional. Even a tiny twist to the boring 2D lookout can make a difference, though it doesn't give any new information to the reader. However, in a presentation, it is almost expected that one produces "exciting" graphs. Today, we will discuss an easy way to add a slight shadow to the key on a graph, as if it were lifted from the plane of the curves. At the end of our exercise, we will have the following figure




I believe the steps are quite straightforward, so I won't spend too much time on explaining every detail. Here is the script that we need

reset

xl=0; xh=1; yl=-1; yh=1;
eps=0.01;
rx=0.6; ry=0.8; kw=0.35; kh=0.15
lh=0.06; al=0.1

key1="First function"
key2="Second function"

set table 'shadowkey.dat'
splot [xl:xh] [yl:yh] x/(xh-xl)
unset table

set object 1 rect from graph rx,ry rto kw,kh fc rgb "#aaaaaa" fs solid 1.0 front lw 0
set object 2 rect from graph rx-eps,ry+eps rto kw,kh front fs empty
set label 1 at graph 1.1*al+rx, ry+2*lh key1 front
set label 2 at graph 1.1*al+rx, ry+lh key2 front
set arrow from graph rx, ry+2*lh rto al, 0 lt 1 lw 1.5 nohead front
set arrow from graph rx, ry+lh rto al, 0 lt 3 lw 1.5 nohead front

unset colorbox
unset key
set palette defined (0 "#8888ff", 1 "#ffffff")
plot [xl:xh] [yl:yh] 'shadowkey.dat' w ima, \
x*x*exp(-x) lw 1.5, cos(13*x)*exp(-3*x) lt 3 lw 1.5

I collected all variables at the beginning of the script, so that it will be easier to adapt it to any situations. The first four numbers define the xrange and the yrange. The next four numbers specify the position and the size of the key. We will have to make our own key, and the size of the key will depend on the particular terminal one uses. This is why it is handy to define them at the beginning of the script, so if you are not satisfied with the results, you can easily adjust both the position and the size. The next number (lh) gives the hight of one key line, while 'al' will be the length of the line which represents the curve. 'key1' and 'key2' are just two arbitrary strings, holding the text of the key.

The next three lines are necessary only if you want to draw the background gradient, but if you are happy with a white plot, you can skip this. I should also point out that the first four numbers that are needed for the xrange and yrange, are needed only because we have to "synchronise" the plot with its background. If you don't want the gradient, you can drop these definitions, and you haven't got to specify the range in the plot either. This might be useful, when you don't know the plotting range beforehand, and want to let gnuplot calculate it for you.

The next step is to draw two rectangles in the front of the graph: one gray, and one white. Note that the rectangle in gray has no boundary, i.e., that is set to zero width. Also note that we specify the coordinates in terms of the numbers we defined at the beginning, so both rectangles will change accordingly, if you change those numbers. This means that the white rectangle and its shadow will always be linked, given by the variable 'eps'. When we are done with the white rectangle, we can place the text of the keys, and the two lines representing our plots. For this purpose we draw two arrows without heads, and with the linetypes that we will use for plotting the curves.

The final step is to actually draw the curves. Before that, we plot our file, 'shadowkey.dat', so the graph will have a background gradient. If you skipped those three lines writing the data file to disc, you should replace the plotting command by

plot [xl:xh] [yl:yh] x*x*exp(-x) lw 1.5, cos(13*x)*exp(-3*x) lt 3 lw 1.5

It should go without saying that in this case, you can also drop the palette definition, for it will not be used. Well, this is it for today. I know that this was a simple trick, but we can't have something complicated every day!

Sunday 26 July 2009

Maps - Contour plots with labels

So, you have always wondered how on Earth one can make a real map with gnuplot. Well, there is a simple and a not-so-simple way to this. First, let us see the simple way. Since it is simple, this method won't have one of the features, isoline labelling, that the second one has. As we go along, the map will become more and more complicated, but I hope that I set the right pace, and it will be easy to follow.

A map is nothing but a colour-coded 3D plot, with the isolines attached to it. We will produce a simple map using the function
sin(1.3*x)*cos(0.9*y)+cos(.8*x)*sin(1.9*y)+cos(y*.2*x)
I don't think that this function has any particular meaning, but it looked quite all right, at least, as far as maps and isolines are concerned. If you have a matrix to plot, you can replace this function with that file. In principle, we could plot both the colour-coded map and the isolines from gnuplot, but we will have much greater flexibility, if we first direct the plots to a file, and then call the data from those files. One of the advantages of doing this is that in this way, we can make sure that the various plots have the same size. This requires a small overhead in terms of scripting, namely, we have got to issue the command
set table 'something.dat'
splot something
unset table
While this might appear superfluous, there are good reasons to do this. If we plot our function through the file 'something.dat', then we can use the 'with image' modifier of the plot command, and this means that the plot will be of the same size as would a normal 2D plot. Otherwise, if we use
set pm3d map
splot something with pm3d
the plot will actually be a bit smaller. The reason behind this is that in a 3D plot, we have to have some space for the 'z' axis, and even if we drop it in the map view, the space-holder for the 'z' axis is still there, therefore, gnuplot makes the whole plot a bit smaller. If, however, we plot the map through a file, we can use 'plot', in splot's stead, therefore, the 'z' axis will not appear anywhere in the processing of the plot.

After this interlude, let us see our first version of the map.
reset
f(x,y)=sin(1.3*x)*cos(.9*y)+cos(.8*x)*sin(1.9*y)+cos(y*.2*x)
set xrange [-5:5]
set yrange [-5:5]
set isosample 250, 250
set table 'test.dat'
splot f(x,y)
unset table

set contour base
set cntrparam level incremental -3, 0.5, 3
unset surface
set table 'cont.dat'
splot f(x,y)
unset table

reset
set xrange [-5:5]
set yrange [-5:5]
unset key
set palette rgbformulae 33,13,10
p 'test.dat' with image, 'cont.dat' w l lt -1 lw 1.5
In the first couple of lines, we plot the actual function into a file, called 'test.dat'. There will be some 25000 data points in this file, for there are 100 samples (default), and 250 isosamples. Having written these points into a file, we plot the contour. We just happen to know that the value of f(x,y) is between -3 and 3, so we set the contour levels at -3, -2.5...2.5, 3. When we are done with all this, we simply reset our plot, set the x and yrange, specify the palette that we want to use with the colouring, and call the plot command. For the actual data points, we use the 'with image' modifier, while for the contours, we increase the linewidth to 1.5, instead of using the default value of 1. We, thus, have just produced the following image.



Now, this was sort of standard, but the question inevitably comes up, whether we could do something more with this. When one has contour lines, one wants to label them, I presume. Petr Mikulik has a gawk script that will do this, in particular, it will place the appropriate labels next to the contour lines. I will discuss another approach here. The main difference with Petr Mikulik's method is that here we will put the labels on the contour lines, with some white space, of course.

First, it is quite instructive to look at the file containing the contour lines. In order not to cobble the plot, we will restrict it to the range [-5:0] (x) and [2:5], but it is really just for the sake of simplicity. With this modification, the file should read as
#Surface 0 of 1 surfaces

# Contour 0, label:        2
-0.482969  3.59036  2
-0.505051  3.58024  2
-0.509852  3.57831  2
-0.539895  3.56627  2
-0.555556  3.55994  2
-0.572136  3.55422  2
...

# Contour 1, label:      1.5
-1.11111  4.57693  1.5
-1.10874  4.57831  1.5
-1.0877  4.59036  1.5
-1.06617  4.60241  1.5
-1.06061  4.60546  1.5
-1.04273  4.61446  1.5
...
What this means is that the 0th contour line is drawn at z=2, and then the (x,y,z) values are listed. Of course, since we know that this particular line is at z=2, the last coordinate doesn't play any role, but it is listed all the same. When we plot this file, consecutive (x,y) points will be connected by a straight line segment. The next contour is at 1.5 etc. It might happen that one contour line is made up of several blocks, both for technical, and fundamental reasons. The technical reason is gnuplot's inner procedure of computing the an isoline, while the fundamental is simply that there might be several disconnected regions with the same isolevel. Anyway, if there are several blocks, they are separated by a blank line within the same contour.

In order to create the labels and the white spaces for them, we will use the following short script:
#!/bin/bash

gawk -v d=$2 -v w=$3 -v os=$4 'function abs(x) { return (x>=0?x:-x) }
    {
            if($0~/# Contour/) nr=0
            if(nr==int(os+w/2) && d==0) {i++; a[i]=$1; b[i]=$2; c[i]=$3;}
            if(abs(nr-os-w/2)>w/2 && d==1) print $0
            nr++
    }
    END {   if(d==0) {
                    for(j=1;j<=i;j++)
                    printf "set label %d \"%g\" at %g, %g centre front\n", j, c[j], a[j], b[j]
            }
    }' $1
This script operates on the data file of contour lines that we had before, and simply takes out a couple of points of the contours, while keeping one of the coordinates as the position of the labels. You can change the length of the white space by setting the 'w' argument, while the position of the white space can be modified by changing the 'os' argument. The first argument, 'd', determines what we want to do with the script: make the labels, or tweak the input file.
Now, we can call the script from our gnuplot script as follows:
reset
set xrange [-5:0]
set yrange [2:5]
set isosample 150, 150
set table 'test.dat'
splot sin(1.3*x)*cos(.9*y)+cos(.8*x)*sin(1.9*y)+cos(y*.2*x)
unset table
set cont base
set cntrparam level incremental -3, 0.5, 3
unset surf
set table 'cont.dat'
splot sin(1.3*x)*cos(0.9*y)+cos(.8*x)*sin(1.9*y)+cos(y*.2*x)
unset table
reset
set xrange [-5:0]
set yrange [2:5]
unset key
set palette rgbformulae 33,13,10
l '<./cont.sh cont.dat 0 15 0'
p 'test.dat' w ima, '<./cont.sh cont.dat 1 15 0' w l lt -1 lw 1.5



There are two more twists that we could add to the figures above. One is that the labels can be rotated in such a way that they are parallel to the curves. This we can easily achieve by adding two more lines to our script above: we take the first and last dropped coordinate, and calculate the angle that that line segment would make with the horizontal. Using this angle, we rotate the labels accordingly. With this modification, the script now becomes
#!/bin/bash

gawk -v d=$2 -v w=$3 -v os=$4 'function abs(x) { return (x>=0?x:-x) }
    {
            if($0~/# Contour/) nr=0
            if(nr==int(os+w/2) && d==0) {a[i]=$1; b[i]=$2; c[i]=$3;}
            if(nr==int(os+w/2)-1 && d==0) {i++; x = $1; y = $2;}
            if(nr==int(os+w/2)+1 && d==0) r[i]= 180.0*atan2(y-$2, x-$1)/3.14
            if(abs(nr-os-w/2)>w/2 && d==1) print $0
            nr++
    }
    END {   if(d==0) {
                    for(j=1;j<=i;j++)
                    printf "set label %d \"%g\" at %g, %g centre front rotate by %d\n", j, c[j], a[j], b[j], r[j]
            }
    }' $1


and the resulting graph is shown here:




Note that the only difference between this and the previous script is the following two lines
if(nr==int(os+w/2)-1 && d==0) {i++; x = $1; y = $2;}
if(nr==int(os+w/2)+1 && d==0) r[i]= 180.0*atan2(y-$2, x-$1)/3.14
where we calculate the appropriate angles.

The second twist is that it might make interpretation of the figure easier, if the contours are coloured, as are the corresponding labels. We can use the script above with some small modification. The key is to use the 'index' modifier. Recall that the contours are rendered in blocks, and with the help of the 'index' keyword, we can specify which block we want to plot. We, therefore, modify our script as follows:
#!/bin/bash

gawk -v d=$2 -v w=$3 -v os=$4 'function abs(x) { return (x>=0?x:-x) }
   {
           if($0~/# Contour/) nr=0
           if(nr==int(os+w/2) && (d%2)==0) {a[i]=$1; b[i]=$2; c[i]=$3;}
           if(nr==int(os+w/2)-1 && (d%2)==0) {i++; x = $1; y = $2;}
           if(nr==int(os+w/2)+1 && (d%2)==0) r[i]= 180.0*atan2(y-$2, x-$1)
           if(abs(nr-os-w/2)>w/2 && (d%2)==1) print $0
           nr++
   }
   END {   if(d==0) {
                   for(j=1;j<=i;j++)
                   printf "set label %d \"%g\" at %g, %g centre front rotate by %d tc lt %d \n", j, c[j], a[j], b[j], r[j], j

           }
           if(d==2) {
                   printf "plot \"test.dat\" w ima, \"cont.plt\" index 0 w l lt 1,\\\n"
                   for(j=2;j<i;j++) printf "\"\" index %d w l lt %d,\\\n", j-1, j
                   printf "\"\" index %d w l lt %d\n", i-1, i

           }
   }' $1

Calling the script produces the following graph:



Note that we had to call plot through a new temporary file, cont.plt, which should be produced by redirecting the output of
cont5.sh 2 15 0 > cont.plt

Wednesday 15 July 2009

Grid with multiple colours

A couple of days ago, I very badly needed a graph on which all four axes were used, and in addition to that, I needed a grid on the graph. Now, the problem is, if you just plot a grid, it will be hard to follow which grid line corresponds to which axis. So, I decided to colour both the grid, and the axes, and use the same colour for corresponding grid lines and axes. Well, this is what I decided, but then it turned out to be a rather hard nut to crack. (In gnuplot 4.3, the situation is a bit better, I will discuss that tomorrow.) I have found a work-around, however, so if you are interested, you can find the details below. But let us see the graph first!



The standard way of defining a grid is a line similar to this:
set grid ytics lt 0 lw 1 lc rgb "#880000"


First, we define where we want to have the grid (everywhere, where there is an ytics), then the linetype, the line width, perhaps, and the colour. The problem is that when you issue this command multiple times, for some funny reason, only the last colour setting will take effect, i.e., if after the line above, you happened to have
set grid y2tics lt 0 lw 1 lc rgb "#008800"

then gnuplot would place the grid lines at y2tics all right, but the grid lines at ytics would also be green, and this is not what we want. The way out of this difficulty is to trick gnuplot into thinking that we have more than one plot, and re-set the grid before each new plot. Of course, we don't actually want to plot anything after the first one, therefore, we will just plot 1/0. It is a vicious way to deceive gnuplot, but it does the trick. Now, here is the script

reset
h=4.136e-15 # eV s
hbar=6.57e-16 # eV s
c=3e17  # nm/s

set multiplot
unset key
set xlabel 'Wavelength [nm]' tc rgb "#0000ff"
set x2label 'Temperature [K]'
set ylabel 'Energy [eV]' tc rgb "#880000"
set y2label 'dE/d{/Symbol l} [meV/nm]' tc rgb "#008800"
set xrange [200:2000]
set x2range [0:10]
set yrange [0:6]
set y2range [1:-7]
set y2tics 1,-1,-7
set ytics nomirror tc rgb "#880000"
set y2tics tc rgb "#008800"
set xtics nomirror tc rgb "#0000ff"
set x2tics
set mxtics 2
set mytics 2
set grid xtics lt 0 lw 1 lc rgb "#0000ff"

plot h*c/x w l lc rgb "#880000", -1000*h*c/x/x axes x1y2 w l lc rgb "#008800"

unset grid
set xlabel " "
set xtics format " "
set x2tics format " "
set x2label " "
set ylabel " "
set y2label " "
set grid x2tics lt 0 lw 1 lc rgb "#000000"
plot 1/0

unset grid
set grid ytics lt 0 lw 1 lc rgb "#880000"
plot 1/0

unset grid
set grid y2tics lt 0 lw 1 lc rgb "#008800"
plot 1/0

unset multiplot


First, we define a couple of constants, and then set up the axis labels. Note that when doing so, we assign a colour to each axis, so that we will know which axis the plot belongs to. Of course, we could do this with arrows, but since we want to have coloured grids, it is better to do it in this way. Next, we set the various ranges, and then set up the first grid line, which will be drawn at every xtics, and in blue. (Blue was the label colour of x.) After plotting the curves, we unset the grid, and we also take off the labels and ticlabels. Note that when doing so, we, in fact, re-set the labels, but with an empty string. The reason for this is that, if we simply take off the labels, then the margins of the graph will change, and when we plot our next graph (which is empty, but we still need the grid), it would be of a different size, and would not align with our original plot. Having re-labelled the graph, we re-set the grid, this time putting the lines on x2tics in black, then on ytics in a dark shade of red, and finally, on y2tics, in dark green. At the very end, we unset the multiplot, i.e., escape from the whole plot.
This is a bit dirty, but if you are pressed for it, it should still be OK. Next time, I will show the easy route, provided you have gnuplot 4.3.