Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Log scale on radial contour plot with matplotlib

I have a sample script to make a radial contour plot:

import os
import math
import numpy as np
import matplotlib.pyplot as plt
import mpl_toolkits.axisartist.floating_axes as floating_axes
from matplotlib.projections import PolarAxes
from mpl_toolkits.axisartist.grid_finder import FixedLocator, MaxNLocator, DictFormatter
import random

# ------------------------------------ #

def setup_arc_radial_axes(fig, rect, angle_ticks, radius_ticks, min_rad, max_rad):

    tr = PolarAxes.PolarTransform()

    pi = np.pi

    grid_locator1 = FixedLocator([v for v, s in angle_ticks])
    tick_formatter1 = DictFormatter(dict(angle_ticks))

    grid_locator2 = FixedLocator([a for a, b in radius_ticks])
    tick_formatter2 = DictFormatter(dict(radius_ticks))

    grid_helper = floating_axes.GridHelperCurveLinear(tr,
                                extremes=((370.0*(pi/180.0)), (170.0*(pi/180.0)), max_rad, min_rad),
                                grid_locator1=grid_locator1,
                                grid_locator2=grid_locator2,
                                tick_formatter1=tick_formatter1,
                                tick_formatter2=tick_formatter2,
                                )

    ax1 = floating_axes.FloatingSubplot(fig, rect, grid_helper=grid_helper)
    fig.add_subplot(ax1)

    ax1.grid(True)

    # create a parasite axes whose transData in RA, cz
    aux_ax = ax1.get_aux_axes(tr)

    aux_ax.patch = ax1.patch
    ax1.patch.zorder=0.9

    #ax1.axis["left"].set_ticklabel_direction("+")

    return ax1, aux_ax

# ------------------------------------ #
# write angle values to the plotting array
angles = []
for mic_num in range(38):
    angle = float(mic_num)*(180.0/36.0)*(math.pi/180.0)+math.pi
    angles.append(angle)

# ------------------------------------ #
### these are merely the ticks that appear on the plot axis
### these don't actually get plotted

angle_ticks = range(0,190,10)
angle_ticks_rads = [a*math.pi/180.0 for a in angle_ticks]
angle_ticks_rads_plus_offset = [a+math.pi for a in angle_ticks_rads]
angle_ticks_for_plot = []
for i in range(len(angle_ticks)):
    angle_ticks_for_plot.append((angle_ticks_rads_plus_offset[i],r"$"+str(angle_ticks[i])+"$"))

# ------------------------------------ #

scale = 1.0
aspect = 1.50
height = 8.0
fig = plt.figure(1, figsize=(height*aspect*scale, height*scale))
fig.subplots_adjust(wspace=0.3, left=0.05, right=0.95, top=0.84)
fig.subplots_adjust()

plot_real_min = 30.0
plot_real_max = 100.0

plot_fake_min = 0.0
plot_fake_max = 5000.0

rad_tick_increment = 500.0

radius_ticks = []
for i in range(int(plot_fake_min),int(plot_fake_max)+int(rad_tick_increment),int(rad_tick_increment)):
    plot_fake_val = ((i-plot_fake_min)/(plot_fake_max-plot_fake_min))*(plot_real_max-plot_real_min)+plot_real_min
    radius_ticks.append((plot_fake_val, r"$"+str(i)+"$"))

ax2, aux_ax2 = setup_arc_radial_axes(fig, 111, angle_ticks_for_plot, radius_ticks, plot_real_min, plot_real_max)

azimuths = np.radians(np.linspace(0, 180, 91))
azimuths_adjusted = [ (x + math.pi) for x in azimuths ]
zeniths = np.arange(0, 5050, 50)
zeniths_adjusted = [((x-plot_fake_min)/(plot_fake_max-plot_fake_min))*(plot_real_max-plot_real_min)+plot_real_min for x in zeniths]

r, theta = np.meshgrid(zeniths_adjusted, azimuths_adjusted)
values = 90.0+5.0*np.random.random((len(azimuths), len(zeniths)))

aux_ax2.contourf(theta, r, values)

cbar = plt.colorbar(aux_ax2.contourf(theta, r, values), orientation='vertical')
cbar.ax.set_ylabel('Contour Value [Unit]', fontsize = 16)

plt.suptitle('Plot Title ', fontsize = 24, weight="bold")
plt.legend(loc=3,prop={'size':20})
plt.xlabel('Angle [deg]', fontsize=20, weight="bold")
plt.ylabel('Frequency [Hz]', fontsize=20, weight="bold")

# plt.show()
plt.savefig('plot.png', dpi=100)
plt.close()

...which gives me a plot that looks like:

enter image description here

However, I would like to have a logarithmic scale on the radius-axis. Does anyone know a convenient way to do this?

like image 719
HotDogCannon Avatar asked Jul 10 '15 10:07

HotDogCannon


2 Answers

Not elegant, alas, but you could alter the polar coordinate transform to do what you want. I got the code from here: https://github.com/matplotlib/matplotlib/blob/master/lib/matplotlib/projections/polar.py.

I changed the names to LogPolarTransform and InvertedLogPolarTransform, then altered the formulas to use a log scale. Basically, I changed these lines:

x[:] = np.where(mask, np.nan, r * np.cos(t))
y[:] = np.where(mask, np.nan, r * np.sin(t))

to these:

x[:] = np.where(mask, np.nan, np.log(r) * np.cos(t))
y[:] = np.where(mask, np.nan, np.log(r) * np.sin(t))

and this line:

r = np.sqrt(x*x + y*y)

to this:

r = np.exp(np.sqrt(x*x + y*y))

If you copy and paste the following code above what you already have, and change tr = PolarAxes.PolarTransform() to tr = LogPolarTransform(), you should end up with a log-scaled radial axis. Here's the resulting figure (I changed plot_real_min to 5.0 so things would show up better) :

Radial Contour Plot

from matplotlib.transforms import  Transform 

class LogPolarTransform(PolarAxes.PolarTransform):
    input_dims = 2
    output_dims = 2
    is_separable = False

    def __init__(self, axis=None, use_rmin=True):
        Transform.__init__(self)
        self._axis = axis
        self._use_rmin = use_rmin
    def transform_non_affine(self, tr):
        xy = np.empty(tr.shape, np.float_)
        if self._axis is not None:
            if self._use_rmin:
                rmin = self._axis.viewLim.ymin
            else:
                rmin = 0
            theta_offset = self._axis.get_theta_offset()
            theta_direction = self._axis.get_theta_direction()
        else:
            rmin = 0
            theta_offset = 0
            theta_direction = 1

        t = tr[:, 0:1]
        r = tr[:, 1:2]
        x = xy[:, 0:1]
        y = xy[:, 1:2]

        t *= theta_direction
        t += theta_offset

        r = r - rmin
        mask = r < 0
        x[:] = np.where(mask, np.nan, np.log(r) * np.cos(t))
        y[:] = np.where(mask, np.nan, np.log(r) * np.sin(t))

        return xy


    def inverted(self):
        return InvertedLogPolarTransform(self._axis, self._use_rmin)
    inverted.__doc__ = Transform.inverted.__doc__

class InvertedLogPolarTransform(Transform):
    """
    The inverse of the polar transform, mapping Cartesian
    coordinate space *x* and *y* back to *theta* and *r*.
    """
    input_dims = 2
    output_dims = 2
    is_separable = False

    def __init__(self, axis=None, use_rmin=True):
        Transform.__init__(self)
        self._axis = axis
        self._use_rmin = use_rmin

    def transform_non_affine(self, xy):
        if self._axis is not None:
            if self._use_rmin:
                rmin = self._axis.viewLim.ymin
            else:
                rmin = 0
            theta_offset = self._axis.get_theta_offset()
            theta_direction = self._axis.get_theta_direction()
        else:
            rmin = 0
            theta_offset = 0
            theta_direction = 1

        x = xy[:, 0:1]
        y = xy[:, 1:]
        r = np.exp(np.sqrt(x*x + y*y))
        with np.errstate(invalid='ignore'):
            # At x=y=r=0 this will raise an
            # invalid value warning when doing 0/0
            # Divide by zero warnings are only raised when
            # the numerator is different from 0. That
            # should not happen here.
            theta = np.arccos(x / r)
        theta = np.where(y < 0, 2 * np.pi - theta, theta)

        theta -= theta_offset
        theta *= theta_direction
        theta %= 2 * np.pi

        r += rmin

        return np.concatenate((theta, r), 1)

    def inverted(self):
        return PolarAxes.LogPolarTransform(self._axis, self._use_rmin)
like image 73
Amy Teegarden Avatar answered Nov 19 '22 01:11

Amy Teegarden


The class GridHelperCurveLinear is hardcoded as linear, so set_yscale('log') is not going to work. It would take a significant effort to use the logarithmic scale. The easiest thing to do is:

  • Apply a logarithmic filter to your data, so your graph will display correctly.
  • Use your own formatter for the labels.

You can find a similar use case here.

like image 26
enrico.bacis Avatar answered Nov 19 '22 02:11

enrico.bacis