scala logaritmica su terreno radiale profilo con matplotlib

Question

Ho uno script di esempio per fare una trama di contorno radiale:

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() 

... che mi dà una trama che assomiglia a:

Tuttavia, mi piacerebbe avere una scala logaritmica sull'asse del raggio. Qualcuno sa un modo conveniente per farlo?

Answer 1

Non elegante, purtroppo, ma è possibile modificare la trasformazione delle coordinate polari per fare ciò che si desidera. Ho ricevuto il codice da qui: https://github.com/matplotlib/matplotlib/blob/master/lib/matplotlib/projections/polar.py.

Ho modificato i nomi in LogPolarTransform e InvertedLogPolarTransform, quindi ho modificato le formule per utilizzare una scala di registro. Fondamentalmente, ho cambiato queste righe:

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

a questi:

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

e questa linea:

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

a questo:

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

Se si copia e incolla il seguente codice sopra quello che hai già, e c affiancare tr = PolarAxes.PolarTransform() a tr = LogPolarTransform(), si dovrebbe finire con un asse radiale scalato nel registro. Ecco la cifra risultante (ho cambiato plot_real_min a 5,0 così le cose sarebbero presentarsi meglio):

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) 

Answer 2

La classe GridHelperCurveLinear è hardcoded come lineare, quindi set_yscale('log') non funzionerà. Ci vorrebbe uno sforzo significativo per usare la scala logaritmica. La cosa più semplice da fare è:

• applicare un filtro logaritmica ai vostri dati , quindi il grafico verrà visualizzato correttamente.
• Utilizzare il proprio formattatore per le etichette.

È possibile trovare un caso simile here.