pythonmatplotlibmeshgrid

Correct Visualization of data using np.meshgrid and ax.pcolormesh of matplotlib, error in visualization


Trying to make a visualization to understand the time-of-flight experiment

The Idea is simple provided 2 values for energy - calculate corresponding time of flight values and energies and visualize it on the xy plot where color stands for the corresponding neutron energy. But when I came to realization - I see that if I use pcolormesh - I have error in visualization - see attached Figure - in the upper right corner - right after the dashed vertical line.. meshgrid and pcolormesh

It seems that if I have provided a lower energy limit, then the border in time scale is strictly defined, and the corner of the triangle must be exactly where the dashed line is. But it's not. What is the error, and how can it be fixed?

Another thing to note here is that the lower the minimum value for energy is, the larger the error.

Code for Minimal reproducible example is shown below.

import pandas as pd
import numpy as np

import matplotlib.pyplot as plt
from matplotlib.colors import Normalize
import matplotlib.cm as cm

# diagram for one pulse case

def calculate_time_interval(FP,
                            energy_mev,
                            m_n=1.67492749804e-27,
                            k = 1.60218e-13,  # J/MeV
                            c=3e8):
    """ Calculates the time interval given the flight path (FP) and energy values. """

    energy_joules = energy_mev * k
    # Calculate velocity
    v = c * np.sqrt(1 - (m_n * c**2 / (energy_joules + m_n * c**2))**2)
    # Calculate time interval
    time_interval = FP / v
    return time_interval


min_E = 1e-2 #MeV
max_E = 1

FP = 100 # m

# calculating min and max TOF based on energy limits supplied
t_min = calculate_time_interval(FP=FP,
                                energy_mev=max_E)*1e3

t_max = calculate_time_interval(FP=FP,
                                energy_mev=min_E) * 1e3

print(t_min, t_max)

# Create a range for the flight path lengths up to 110 meters
flight_paths = np.linspace(0, FP, 500)

# Energy range 
energies = np.linspace(min_E, max_E, 500)

# Create meshgrid for flight paths and energies
FP_mesh, E_mesh = np.meshgrid(flight_paths, energies)

# Calculate velocities and corresponding ToF in ms
tof = calculate_time_interval(FP=FP_mesh, energy_mev=E_mesh) * 1e3


#  plot
fig, ax = plt.subplots(figsize=(8, 5))

c = ax.pcolormesh(tof, FP_mesh, E_mesh, cmap='coolwarm', shading='auto', norm=Normalize(vmin=min_E, vmax=max_E))

# Add colorbar
cbar = plt.colorbar(c, ax=ax, orientation='vertical')
cbar.set_label('Energy (MeV)')

# Add horizontal line at FP
ax.axhline(y=FP, color='black', linestyle='--')
ax.axvline(x=t_max, color='black', linestyle='--')

# Add title, labels
fig.suptitle('Time of Flight vs Flight Path Length with Energy (single pulse)')
ax.set_title(fr'(! $E_{{min}}$ = {np.min(E_mesh)} MeV => {np.round(t_min,6)} ms, $E_{{max}}$ = {np.max(E_mesh)} MeV => {np.round(t_max,6)} ms)')

ax.set_xlabel('Time of Flight (ms)')
ax.set_ylabel('Flight Path Length (m)')
ax.set_xlim(0, t_max * 2)
ax.set_ylim(0, 1.1*FP)

plt.grid(True)
plt.show()

Solution

  • For a triangular domain like that you should probably be using tripcolor rather than pcolormesh.

    ax.tripcolor(tof.flatten(), FP_mesh.flatten(), E_mesh.flatten(), cmap='coolwarm', norm=Normalize(vmin=min_E, vmax=max_E))
    

    Result: