Producing polar contour plots with matplotlib
In my field I often need to plot polar contour plots, and generally plotting tools don’t make this easy. In fact, I think I could rate every single graphing/plotting package in the world by the ease of producing a polar contour plot – and most would fail entirely! Still, I have managed to find a fairly nice way of doing this using my plotting package of choice: matplotlib.
I must warn you first – a Google search for matplotlib polar contour or a similar search term will produce a lot of completely out-dated answers. The most commonly found answers are those such as this StackOverflow question and this forum post. In fact, the first question was asked by me last year – and got an answer which is now completely out of date. Basically, all of these answers tell you that you can’t do a polar contour plot directly in matplotlib, and you must convert your points from polar co-ordinates to cartesian co-ordinates first. This isn’t difficult, but is a pain to do, and of course you then end up with cartesian axes which doesn’t look great. The great news is that you can now do polar contour plots directly with matplotlib!
So, how do you do them? Simple really, you just create some polar axes and plot a contour plot on them:
fig, ax = subplots(subplot_kw=dict(projection='polar')) cax = ax.contourf(theta, r, values, nlevels)
This produces a filled contour plot, as it uses the contourf function, using the contour function would give simple contour lines. The first three parameters which must be given to this function are all two-dimensional arrays containing: the radii, the angles (theta) and the actual values to contour. The final parameter is the number of contour levels to plot – you tend to want lower numbers for line contours and higher numbers for filled contour plots (to get a smooth look).
I never quite understood these two-dimensional arrays, and why they were needed. I normally had my data in the form of three lists that were basically columns of a table, where each row of the table defined a point and value. For example:
Radius | Theta | Value |
10 | 0 | 0.7 |
10 | 90 | 0.45 |
10 | 180 | 0.9 |
10 | 270 | 0.23 |
20 | 0 | 0.5 |
20 | 90 | 0.13 |
20 | 180 | 0.52 |
20 | 270 | 0.98 |
Each of these rows define a point – for example, the first row defines a point with a radius of 10, an angle of 0 degrees and a value of 0.7. I could never understand why the contour function didn’t just take these three lists and plot me a contour plot. In fact, I’ve written a function that will do just that, which I will describe below, but first let me explain how those values are converted to two-dimensional arrays.
First of all, lets think of the dimensions: we obviously have dimensions here in our data, the radius and the angle. In fact, we could re-shape our values array so that it is two-dimensional fairly easily. We can see from the table above that we are doing all the azimuth angles for a radius of 10 degrees, then the same azimuth angles for a radius of 20 degrees, etc. Thus, rather than our values being stored in a one-dimensional list, we could put them in a two-dimensional table where the columns are the azimuth angles, and the rows are the radii:
0 | 90 | 180 | 270 | |
10 | 0.7 | 0.45 | 0.9 | 0.23 |
20 | 0.5 | 0.13 | 0.52 | 0.98 |
This is exactly the sort of two dimensional array that we need to give to the contourf function. That’s not too hard to understand – but why on earth do the radii and angle arrays have to be two-dimensional too. Well, basically we just need two arrays like the one above, but with the relevant radii and angles in the cells, rather than the values. So, for the angles, we’d have:
0 | 90 | 180 | 270 | |
10 | 0 | 90 | 180 | 270 |
20 | 0 | 90 | 180 | 270 |
And for the radii we’d have:
0 | 90 | 180 | 270 | |
10 | 10 | 10 | 10 | 10 |
20 | 0 | 20 | 20 | 20 |
Then, when we take all three arrays together, each cell will define the three bits of information we need. So, the top left cell gives us an angle of 0, a radius of 10 and a value of 0.7. Luckily, you don’t have to make these arrays by hand – a handy NumPy function called meshgrid will do it for you:
>>> radii = np.arange(0, 60, 10) >>> print radii [ 0 10 20 30 40 50] >>> angles = np.arange(0, 360, 90) >>> print angles [ 0 90 180 270] >>> np.meshgrid(angles, radii) (array([[ 0, 90, 180, 270], [ 0, 90, 180, 270], [ 0, 90, 180, 270], [ 0, 90, 180, 270], [ 0, 90, 180, 270], [ 0, 90, 180, 270]]),
array([[ 0, 0, 0, 0], [10, 10, 10, 10], [20, 20, 20, 20], [30, 30, 30, 30], [40, 40, 40, 40], [50, 50, 50, 50]]))
One thing to remember is that the plotting function requires the angle (theta) in radians, not degrees, so if your data is in degrees (as it often is) then you’ll need to convert it to radians using the NumPy radians function.
After doing all of this you can get your data into the contour plotting function correctly, and you can get some polar axes for it to be plotted on. However, if you do this, you’ll find that your axes look something like this:
You can see that zero degrees isn’t at the top, it’s at the ‘East’ or ‘3 o’clock’ position, and the angles go round the wrong way! Apparently that’s how these things are often done in maths – but in my field particularly people want to have a polar plot like a compass, with zero at the top!
If you try and find how to do this, you’ll find a StackOverflow answer with a brilliant subclass of PolarAxes which does this for you. It’s brilliant that matplotlib allows you do this sort of customisation, but if you look below the accepted answer you’ll find a link to the matplotlib documentation for a function called set_theta_zero_location. This function very nicely takes a compass direction (“N” or “E” or “NE” etc) for where zero should be, and puts it there! Similarly, the function set_theta_direction sets the direction in which the angles will increase. All you need to do to use these is call them from the axes object:
ax.set_theta_zero_location("N") ax.set_theta_direction(-1)
The example above will set up the plot for a ‘normal’ compass-style plot with zero degrees at the north, and the angles increasing clockwise. If you find that these lines of code give an error you need to update your matplotlib version – these methods were only added in the latest version (v1.1.0).
So, now we’ve covered everything that I’ve gradually learnt about doing this, we can put it all together in a function. I use the function below whenever I want to plot a polar contour plot, and it works fine for me. It is documented through the docstring shown in the code below.
I can’t guarantee the code will work for you, but hopefully this post has been helpful and you’ll now be able to go away and create polar contour plots in matplotlib.
import numpy as np from matplotlib.pyplot import * def plot_polar_contour(values, azimuths, zeniths): """Plot a polar contour plot, with 0 degrees at the North. Arguments: * `values` -- A list (or other iterable - eg. a NumPy array) of the values to plot on the contour plot (the `z` values) * `azimuths` -- A list of azimuths (in degrees) * `zeniths` -- A list of zeniths (that is, radii) The shapes of these lists are important, and are designed for a particular use case (but should be more generally useful). The values list should be `len(azimuths) * len(zeniths)` long with data for the first azimuth for all the zeniths, then the second azimuth for all the zeniths etc. This is designed to work nicely with data that is produced using a loop as follows: values = [] for azimuth in azimuths: for zenith in zeniths: # Do something and get a result values.append(result) After that code the azimuths, zeniths and values lists will be ready to be passed into this function. """ theta = np.radians(azimuths) zeniths = np.array(zeniths) values = np.array(values) values = values.reshape(len(azimuths), len(zeniths)) r, theta = np.meshgrid(zeniths, np.radians(azimuths)) fig, ax = subplots(subplot_kw=dict(projection='polar')) ax.set_theta_zero_location("N") ax.set_theta_direction(-1) autumn() cax = ax.contourf(theta, r, values, 30) autumn() cb = fig.colorbar(cax) cb.set_label("Pixel reflectance") return fig, ax, cax
Categorised as: Programming, Python, Remote Sensing
Thanks, this is a good write up. I’ve been struggling with plotting spectral wave data in a polar plot in matplotlib.
I found your post after seeing your question on stackoverflow, might be worth linking them together.
Thanks a lot for posting this. I can’t use it yet because my work computer won’t allow me to download matplotlib v1.1.0 but I will hopefully soon.
To James Morrison: I’m working with wave spectral data too. If you haven’t heard of it, the WAFO toolbox could be good for you – http://www.maths.lth.se/matstat/wafo/. It has a python version too.
I tried this out but it’s still quite a bit slower than what my previous method was: for each point along a ray, add an x,y,intensity to a scatter plot. This is still roughly 4x faster than a polar contour as above, but the output ends up looking pretty similar.
Thanks! This was super useful.
Thank you very much for sharing the useful information! I have one question. Is there any chance to bring the polar contour plot to 3D (to make cylindrical plot)? I work with computational chemistry and performed a potential energy surface scan of a small 3 atomic molecule, where I was stretching one bond at different angles. The data which I produced would form sort of Mexican hat. I can plot them 3D in Cartesian coordinates, but the plot is not as easily readable as if it was plotted in cylindrical coordinates. Thanks!
Hi there
I have tried using this but it didn’t work. I have a (9,6) array that contains the data I would like to analyse (MISR reflectance values, so 9 for the 9 MISR sensors/zeniths ±70, ± 60, ±45, ±25, 0; 6 for View Zenith, Sun Zenith, View Azimuth, Sun Azimuth, Relative Azimuth and Reflectance). So I set the arguments:
values = array[:,5]
azimuths = array[:,2]
zeniths = array[:,0]
then I run plot_polar_contour(values,azimuths,zeniths) but I get the following error message:ValueError Traceback (most recent call last)
in ()
—-> 1 plot_polar_contour(values,azimuths,zeniths)
in plot_polar_contour(values, azimuths, zeniths)
28
29 values = np.array(values)
—> 30 values = values.reshape(len(azimuths), len(zeniths))
31
32 r, theta = np.meshgrid(zeniths, np.radians(azimuths))
ValueError: total size of new array must be unchanged
len(values), len(azimuths) and len(zeniths) is 9 in all cases, so not sure what the problem is!
Also I’m not quite sure what you mean with this bit:
values = []
for azimuth in azimuths:
for zenith in zeniths:
# Do something and get a result
values.append(result)
so I skipped it, can you explain what you mean?
Cheers
Alex
Hi,
I think the issue is that all of your arrays are the same length. I don’t quite understand exactly what you’re trying to do: you say that you have a 9×6 array, but it seems like your data itself is only two dimensional, not three dimensional. To plot a polar contour plot you need to have three pieces of information for each point on the graph: the value (which will be shown as the colour), the zenith location (or radius) and the azimuth location (or angle). You have the MISR reflectance (your value), and the zenith angle for each of these reflectances (your zenith/radius). You also have an azimuth location for each of your measurements – but combined together these don’t form a complete set of measurements. That is, you may know the value at azimuth=100, zenith=70, and also the value at azimuth=140 for and zenith=60 – but actually you need to have values at all combinations of azimuth and zenith angles. That is, if you have 9 zenith angles and 6 azimuth angles you need 9*6 = 54 measurements. All you have is 9 measurements.
Does that make sense? If you look at the examples in the blog post you’ll see that I always have values for all combinations of zeniths and azimuths. I don’t know how to go about producing a graph when you don’t have all of the data – you may need to interpolate somehow, or it may be that a polar plot isn’t the best way to visualise the data.
In terms of the commented for loop that you didn’t understand: that was just an example of how you could generate the required data. You already have the data, so you don’t need to worry about that. However, it does illustrate how the data needs to be formatted – you can see that you are producing a value for every combination of zeniths and azimuths by using a nested for loop.
Hope that helps,
Robin
Hi, thank you for this!!
I’m having a problem in that it will only plot my data up until the last angle before 2pi. (my data goes up to an increment lower than 2pi). I’d like it to fill in the gap between this last index and 2pi so that it doesn’t look like someone took a slice of pie out of my plot.
Thanks,
Naomi
Hi Naomi,
You should be able to do this by repeating the data for angle 0 at angle 2pi: that is, just repeat the first element of the array at the end, and it should fill in the gap.
Hope that makes sense,
Robin
Great 🙂
[…] Producing polar contour plots with matplotlib […]
Amazing help! Thank you! Worked like a charm!
Hi and thank you for this useful tool!
I am not very familiar with Python and would need some help.
I am trying to create this kind of plot, excluding all values equal to zero (i.e. I would like them to appear in white). Could you indicate me a tip please?
Unfortunately I’m not sure how to best to do this – maybe try asking on StackOverflow?
Hello, I’m having quite a bit of trouble as I’m new to matplotlib and simply do not understand the variables that you have used there. My goal is to simply plot a limacon equation determined by the user when asked for a, b, and their choice of a sin or cos function. I have used this to generate a list of radii and thetas that I would like to plot. However, this is where I’m lost and confused. Here is my code? What should I do? This is for my computer science class project.
def genpolcoordlist(a, b, trig):
rlist = []
thetalist = []
count = 0
theta = 0
while count <= 12:
if trig == "Cos" or trig == "cos" or trig == "Cosine" or trig == "cosine":
r = ((a) + ((b)(math.cos(math.radians(theta)))))
thetalist = thetalist + [theta]
rlist = rlist + [r]
count = count + 1
theta = theta + 30
elif trig == "Sin" or trig == "sin" or trig == "Sine" or trig == "sine":
r = ((a) + ((b)(math.sin(math.radians(theta)))))
thetalist = thetalist + [theta]
rlist = rlist + [r]
count = count + 1
theta = theta + 30
return rlist, thetalist
rlist, thetalist =gencartcoordlist(a, b, trig)
I just dont know what to do with these lists.
Unfortunately I haven’t got the time at the moment to try and debug this code – try posting on StackOverflow and see if someone there can help.
I think you will need to .append each theta to thetalist, surely? Similarly with r and rlist.