Python Plots
The title of this page is far too grandiose: we will discuss only some simple aspects of the pylab package for python. Pylab is a collection of functions designed to mimic much of the function of Matlab in python, and in particular, to mimic Matlab's powerful plotting capabilities. Pylab imports and capitalizes on numpy and another python package called matplotlib. Indeed, you will probably find henceforth that you simply load pylab at the beginning of every session and leave it at that.
In order to illustrate the use of the plotting facilites in pylab,
let's write a python function to make a classic picture. We'll plot
a function passed as an argument from a pylab session, over
an interval specified as arguments to the python function, and then
plot the rectangles of a Riemann left sum to approximate the
integral of the given curve. We start with just a function
to plot some curve. In our file "myfunctions.py" we change
the import line at the top to
from pylab import *
This replaces the line we had before that imported numpy.
We type the following code e.g. at the end of the file.
def riemannplot(fct,a,b,n): if n<=0: print "Error: n must be positive" return False smoothh= (b-a)/100.0 x = arange(a,b+smoothh,smoothh) plot(x,fct(x)) show()
We recycled code to make sure that n is positive but then computed a reasonably fine array of points where we want to plot the function given as fct. Then we just call the plot function. The plot command requires two arguments: an array of abscissae (x-coordinates) for the points to be plotted, and an array of ordinates (y-coordinates). The plot command creates the figure, but does not display it on the screen. To do that you must call the show command. This allows us to put more curves on the plot being created, or alter it in other ways.
We call our new function as follows.
>>> import myfunctions as mf >>> def f(x): ... return x**3-2*x**2+0.5*x+0.5 ... >>> mf.riemannplot(f,-1,2,8)
At this point a new window appears on our screen, containing the plot we made. The plot is shown in Figure 1.
Figure 1
Next, we want to put in all the boxes showing the elements of the Riemann sum. We can use the bar command for this. The bar command again takes two arguments. The first is an array of x-coordinates for the left corner of each bar, and the second is a height of the bar, which can be negative. Thus, we make a new variable h that represents the width of each bar, and then use that to create an array riemannx of x-coordinates for the left corners of the bars. The array riemanny is just the value of the input function fct at those points. Thus, the picture we are trying to draw corresponds to a Riemann left sum. The picture we actually drew is ugly - it appears in Figure 2. Our code appears below.
def riemannplot(fct,a,b,n): if n<=0: print "Error: n must be positive" return False smoothh= (b-a)/100.0 x = arange(a,b+smoothh,smoothh) plot(x,fct(x)) h = (double(b)-double(a))/double(n) riemannx = arange(a,b,h) riemanny = fct(riemannx) bar(riemannx,riemanny) show()
Figure 2
Figure 3
Figure 4
The problems with this new picture include our inability to
see the original curve through the bars; the bars are too wide; and
we would really like the bars to be a different color to
contrast with the curve. We can fix all these things by adding
attributes to our bar command. This means that we
add arguments to the command, each of which has the form
name=value. In our case, the three things we want to
change include the color: facecolor='orange'. The width
of the bars should be h, so we can add an
argument width=h. Finally, we can make the
bars 50% transparent by setting alpha=0.5. Our
bar command now looks like
bar(riemannx,riemanny,width=h,alpha=0.5,facecolor='orange')
The new plot appears in Figure 3.
We should probably also label the axes, maybe provide a title, and maybe restrict the window to the part we care about. In particular, the default plot ranged out to 2.5 for reasons we do not fully understand. We'll restrict it again. Here is the code.
def riemannplot(fct,a,b,n): if n<=0: print "Error: n must be positive" return False smoothh= (b-a)/100.0 x = arange(a,b+smoothh,smoothh) plot(x,fct(x)) h = (double(b)-double(a))/double(n) riemannx = arange(a,b,h) riemanny = fct(riemannx) bar(riemannx,riemanny,width=h,alpha=0.5,facecolor='orange') xlabel('x') ylabel('f(x)') title('Riemann Left Sum for f(x)') xlim(-1,2) show()
The commands we just added are pretty intuitive. Note in particular the xlim command, which takes two arguments: the left limit; and the right limit. There is a similar ylim command that we could have used, had we wanted to. The picture appears in Figure 4.
At this point we have a function that has some wide applicability.
For example, we can use it to see how the Riemann sum improves as
we refined the partition by increasing the number of points.
In Figure 5 we see the results of a call to the function with n=20.
>>> mf.riemannplot(f,-1,2,20)
Figure 5