Application of H5 Files with PENMSH

In [1]:

    

%pylab inline


    

    




Welcome to pylab, a matplotlib-based Python environment [backend: module://IPython.zmq.pylab.backend_inline].
For more information, type 'help(pylab)'.

In [2]:

    

import numpy as np
import h5py
import matplotlib.pyplot as plt
import matplotlib.cm as cm

def boundaries_xy(h5_file_handle):
    """
    Return a list of tuples (xmin, xmax, ymin, ymax) for each CM, as 
    well as a global tuple defining the XY view

    Args: h5_file_handle

    Returns: pm_data, pm_global
    pm_data is a list of tuples and pm_global is a tuple list
    """

    xdim = h5_file_handle.get('xdim')[:]
    ydim = h5_file_handle.get('ydim')[:]
    zdim = h5_file_handle.get('zdim')[:]
    
    pm_data = []
    for k, cbz in enumerate(zdim[:-1]):
        for j, cby in enumerate(ydim[:-1]):
            for i, cbz in enumerate(xdim[:-1]):
               pm_data.append((xdim[i], xdim[i+1], ydim[j], ydim[j+1]))
    
    pm_global = (np.min(xdim), np.max(xdim), np.min(ydim), np.max(ydim))
    return pm_data, pm_global

def grab_macro_cm_dims(h5_file_handle):

    xdim = h5_file_handle.get('xdim')[:]
    ydim = h5_file_handle.get('ydim')[:]
    zdim = h5_file_handle.get('zdim')[:]
    cm_span = dict()
    cm_span['x'] = np.shape(xdim[:])[0] - 1
    cm_span['y'] = np.shape(ydim[:])[0] - 1
    cm_span['z'] = np.shape(zdim[:])[0] - 1

    return cm_span

def grab_coarse_meshes(h5_file_handle):
    
    cm_data = []
    
    cm_span = grab_macro_cm_dims(h5_file_handle)
    for i in range(cm_span['x']*cm_span['y']*cm_span['z']):
        cm_name = 'cm' + str(i+1).zfill(3)
        data = h5_file_handle.get(cm_name)
        cm_data.append(data)
        
    return cm_data, cm_span

def cm_view(cm_data, i):
    """
    Return a compatible numpy array providing 2D plane view of CM
    (e.g. 'XY' for z-cut).

    For now, the view is 'XY' and taking the first z fine mesh.
    [i-1] beacuse cm_data is zero-based based on Python zero-index
    [0] for the first z fine mesh
    [::-1] reverses the array as PENMSH data is in 'YX' form, flip back to 'XY'
    """
    return np.array(cm_data[i-1][0][::-1])

# Not sure if the pixel is proportioned to appropriate dx/dy, will need to investigate
def plot_cm(cm_data, cm_idx):
    """
    Return an image of the coarse mesh #cm_idx in 2D plane view using nearest intepolation and rainbow colormap
    """
    return plt.imshow(cm_view(cm_data, cm_idx),interpolation='nearest', cmap=cm.gist_rainbow_r, extent=pm_data[cm_idx-1])

def plot_xy(cm_data, cm_span, z_level):
    """
    Plot a 2D plane view, given the coarse-mesh data, the coarse-mesh span, and z-level[int]
    """
    assert 0 < z_level <= cm_span['z']
    cm_min = cm_span['x'] * cm_span['y'] * (z_level - 1)
    cm_max = cm_span['x'] * cm_span['y'] * z_level
    for cm_idx in range(cm_min, cm_max):
        plot_cm(cm_data, cm_idx)


    

In [3]:

    

# This is the 'h5' file that is exported by '-h5' using PENMSH for case 1 hexagon
myfile = 'he3demo.h5'
h5_file_handle = h5py.File(myfile,'r')

# Collect the coarse mesh data from the h5 file
cm_data, cm_span = grab_coarse_meshes(h5_file_handle)    

# Collect the physical mesh boundary data
pm_data, pm_global = boundaries_xy(h5_file_handle)

# Size
plt.rcParams['figure.figsize'] = 10,10
# Assemble the XY mesh data (right first zfm in each z-level)
plot_xy(cm_data, cm_span, 1)

# Set XY View, and Show
plt.axis(pm_global)

plt.show()