AntSimi/py-eddy-tracker

Help in generating representative contour and save as .nc

hafez-ahmad opened this issue · 5 comments

Thank you so much for creating wonderfull Package. I am trying to reproduce following

https://py-eddy-tracker.readthedocs.io/en/latest/python_module/10_tracking_diagnostics/pet_pixel_used.html#sphx-glr-python-module-10-tracking-diagnostics-pet-pixel-used-py

Could you please help me on how to add some representative contour lines on top this plot? and how I can save as nc/ geotif file so I could further work on this.

Here my code

import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from matplotlib.colors import LogNorm

step = 0.125
bins = ((75, 100, step), (2, 23, step))
kwargs_pcolormesh = dict(
    cmap="terrain_r", vmin=0, vmax=0.75, factor=1 / a.nb_days, name="count"
)

fig, ((ax_a, ax_c), (ax_all, ax_ratio)) = plt.subplots(2, 2, figsize=(12, 12), subplot_kw={'projection': ccrs.PlateCarree()})

ax_a.set_title("Anticyclonic frequency")
ax_c.set_title("Cyclonic frequency")
ax_all.set_title("All eddies frequency")
ax_ratio.set_title("Ratio Cyclonic / Anticyclonic")

# Add labels
ax_a.text(-0.1, 1.1, 'A', transform=ax_a.transAxes, fontsize=16, fontweight='bold', va='top')
ax_c.text(-0.1, 1.1, 'B', transform=ax_c.transAxes, fontsize=16, fontweight='bold', va='top')
ax_all.text(-0.1, 1.1, 'C', transform=ax_all.transAxes, fontsize=16, fontweight='bold', va='top')
ax_ratio.text(-0.1, 1.1, 'D', transform=ax_ratio.transAxes, fontsize=16, fontweight='bold', va='top')

# Function to setup map features
def setup_map(ax):
    ax.set_extent([75, 100, 2, 23], crs=ccrs.PlateCarree())
    ax.coastlines()
    ax.add_feature(cfeature.LAND, color='gray')
    gl = ax.gridlines(draw_labels=True)
    gl.top_labels = False
    gl.right_labels = ax in [ax_c, ax_ratio]  # Show right labels only on the right plots
    gl.left_labels = ax in [ax_a, ax_all]  # Show left labels only on the left plots
    gl.bottom_labels = True

# Setup map features for each subplot
for ax in (ax_a, ax_c, ax_all, ax_ratio):
    setup_map(ax)

# Count pixel used for each contour
g_a = a.grid_count(bins, intern=True)
m_a = g_a.display(ax_a, **kwargs_pcolormesh)
g_c = c.grid_count(bins, intern=True)
m_c = g_c.display(ax_c, **kwargs_pcolormesh)

# Add a common color bar for the first two plots (ax_a and ax_c)
#cbar_ax = fig.add_axes([0.2, 0.53, 0.6, 0.02])  # Adjust the position and size of the color bar
#fig.colorbar(m_a, cax=cbar_ax, orientation='horizontal')

# Compute a ratio Cyclonic / Anticyclonic
ratio = g_c.vars["count"] / g_a.vars["count"]

# Mask manipulation to be able to sum the 2 grids
m_c = g_c.vars["count"].mask
m = m_c & g_a.vars["count"].mask
g_c.vars["count"][m_c] = 0
g_c.vars["count"] += g_a.vars["count"]
g_c.vars["count"].mask = m

m_all = g_c.display(ax_all, **kwargs_pcolormesh)
cbar_ax_all = fig.add_axes([0.92, 0.55, 0.02, 0.35])  # Adjust the position and size of the color bar for ax_all
fig.colorbar(m_all, cax=cbar_ax_all, orientation='vertical')

g_c.vars["count"] = ratio
m_ratio = g_c.display(
    ax_ratio, name="count", norm=LogNorm(vmin=0.1, vmax=10), cmap="coolwarm_r"
)
cbar_ax_ratio = fig.add_axes([0.92, 0.13, 0.02, 0.34])  # Adjust the position and size of the color bar for ax_ratio
fig.colorbar(m_ratio, cax=cbar_ax_ratio, orientation='vertical')

plt.tight_layout(rect=[0, 0, 0.9, 1])  # Adjust the layout to make space for the color bars

# Save figure
plt.savefig('Frequency.jpg', dpi=500)

To my knowledge, matplotlib does not allow you to create geotiffs.

Could you please help me on how to add some representative contour lines on top this plot?
Contour of which data?

Thank you so much for your response. @AntSimi for all four figures of Anticyclonic frequency, Cyclonic frequency
All eddies frequency, Ratio Cyclonic / Anticyclonic. I like to add contour on top of frequency.

Take a look at this method RegularGridDataset.contour

Could you please write demo like you created deom on frequency ?

Thank you

Use is similar to grid display method, keyword argument are listed on documentation and on matplotlib documentation linked in pet contour documentation.
Use could look like this if you want contour for some levels:

replace_by_my_grid.contour(replace_by_my_ax, 'replace_by_my_var', levels=[-3500,-2500,-1500])

Use could look like this if you want contour for a number of levels:

replace_by_my_grid.contour(replace_by_my_ax, 'replace_by_my_var', levels=5)