trento/scripts/quick_view.py

0001 #!/usr/bin/env python3
0002 import numpy as np
0003 import matplotlib.pyplot as plt
0004 import h5py, sys
0005
0006 def help():
0007
0008     """
0009   This script plots
0010     (1) transverse intergated entropy/density as function of rapidity
0011     (2) Projection of entropy/density onto x-y (eta=0) plane
0012     (3) ............................. onto y-eta (x=0) plane.
0013     (4) ............................. onto x-eta (y=0) plane.
0014
0015   Usage:
0016     {:s} trento3d-hdf5-output [list-of-event-id-to-convert]
0017
0018   For example, to view all events:
0019     {:s} ic.hdf5
0020   To view only events #2 and #3:
0021     {:s} ic.hdf5 2 3
0022 """
0023     print(help.__doc__.format(__file__, __file__, __file__))
0024
0025 def plot(dataset):
0026     fig, axes = plt.subplots(
0027         nrows=2, ncols=2,
0028         figsize=(8, 7)
0029     )
0030     ax = axes[0,0]
0031     field = dataset.value
0032     Nx, Ny, Neta = field.shape
0033     deta = dataset.attrs['deta']
0034     dxy = dataset.attrs['dxy']
0035     Lx, Ly, Leta = Nx*dxy/2., Ny*dxy/2., (Neta-1)*deta/2.
0036     x = np.linspace(-Lx, Lx, Nx)
0037     y = np.linspace(-Ly, Ly, Ny)
0038     eta = np.linspace(-Leta, Leta, Neta)
0039     ax.plot(eta, np.sum(field, axis=(0,1))*dxy*dxy)
0040     ax.set_xlabel(r'$\eta$')
0041     ax.set_ylabel(r'$dS/d\eta$')
0042     ax.set_title('Tranverse integrated entropy')
0043     Nxmid, Nymid, Netamid = int((Nx-1)/2), int((Ny-1)/2), int((Neta-1)/2)
0044     for ax, projection, ranges, xlabel, ylabel, title in zip(
0045         axes.flatten()[1:],
0046         [field[:,:,Netamid], field[Nxmid,:,:], field[:,Nymid,:]],
0047         [[-Lx, Lx, -Ly, Ly], [-Ly, Ly, -Leta, Leta], [-Lx, Lx, -Leta, Leta]],
0048         [r'$y$ [fm]', r'$\eta$', r'$\eta$'],
0049         [r'$x$ [fm]', r'$y$ [fm]', r'$x$ [fm]'],
0050         [r'$x$-$y$ projection, $\eta=0$',r'$y$-$\eta$ projection, $x=0$',
0051          r'$x$-$\eta$ projection, $y=0$']
0052         ):
0053         ax.contourf(projection, extent=ranges)
0054         ax.set_xlabel(xlabel)
0055         ax.set_ylabel(ylabel)
0056         ax.set_title(title)
0057         ax.axis('equal')
0058     plt.subplots_adjust(wspace=0.4, hspace=0.4)
0059     plt.show()
0060
0061 def main():
0062     if len(sys.argv) <= 1:
0063         help()
0064         exit()
0065     f = h5py.File(sys.argv[1], 'r')
0066     elist = ['event_{}'.format(index) for index in sys.argv[2:]] \
0067             if len(sys.argv)>=3 else list(f.keys())
0068     for eid in elist:
0069         plot(f[eid])
0070
0071 if __name__ == '__main__':
0072     main()