ensure closed file
This commit is contained in:
parent
36d2ae1c2a
commit
4a943ff844
|
@ -111,7 +111,7 @@ class ConfigMaterial(Config):
|
||||||
|
|
||||||
Parameters
|
Parameters
|
||||||
----------
|
----------
|
||||||
fname : str
|
fname : str or pathlib.Path
|
||||||
Filename of the DREAM.3D (HDF5) file.
|
Filename of the DREAM.3D (HDF5) file.
|
||||||
grain_data : str
|
grain_data : str
|
||||||
Name of the group (folder) containing grain-wise data. Defaults
|
Name of the group (folder) containing grain-wise data. Defaults
|
||||||
|
@ -154,26 +154,26 @@ class ConfigMaterial(Config):
|
||||||
defined separately.
|
defined separately.
|
||||||
|
|
||||||
"""
|
"""
|
||||||
b = util.DREAM3D_base_group(fname) if base_group is None else base_group
|
with h5py.File(fname, 'r') as f:
|
||||||
c = util.DREAM3D_cell_data_group(fname) if cell_data is None else cell_data
|
b = util.DREAM3D_base_group(f) if base_group is None else base_group
|
||||||
f = h5py.File(fname,'r')
|
c = util.DREAM3D_cell_data_group(f) if cell_data is None else cell_data
|
||||||
|
|
||||||
if grain_data is None:
|
if grain_data is None:
|
||||||
phase = f['/'.join([b,c,phases])][()].flatten()
|
phase = f['/'.join([b,c,phases])][()].flatten()
|
||||||
O = Rotation.from_Euler_angles(f['/'.join([b,c,Euler_angles])]).as_quaternion().reshape(-1,4) # noqa
|
O = Rotation.from_Euler_angles(f['/'.join([b,c,Euler_angles])]).as_quaternion().reshape(-1,4) # noqa
|
||||||
_,idx = np.unique(np.hstack([O,phase.reshape(-1,1)]),return_index=True,axis=0)
|
_,idx = np.unique(np.hstack([O,phase.reshape(-1,1)]),return_index=True,axis=0)
|
||||||
idx = np.sort(idx)
|
idx = np.sort(idx)
|
||||||
else:
|
else:
|
||||||
phase = f['/'.join([b,grain_data,phases])][()]
|
phase = f['/'.join([b,grain_data,phases])][()]
|
||||||
O = Rotation.from_Euler_angles(f['/'.join([b,grain_data,Euler_angles])]).as_quaternion() # noqa
|
O = Rotation.from_Euler_angles(f['/'.join([b,grain_data,Euler_angles])]).as_quaternion() # noqa
|
||||||
idx = np.arange(phase.size)
|
idx = np.arange(phase.size)
|
||||||
|
|
||||||
if cell_ensemble_data is not None and phase_names is not None:
|
if cell_ensemble_data is not None and phase_names is not None:
|
||||||
try:
|
try:
|
||||||
names = np.array([s.decode() for s in f['/'.join([b,cell_ensemble_data,phase_names])]])
|
names = np.array([s.decode() for s in f['/'.join([b,cell_ensemble_data,phase_names])]])
|
||||||
phase = names[phase]
|
phase = names[phase]
|
||||||
except KeyError:
|
except KeyError:
|
||||||
pass
|
pass
|
||||||
|
|
||||||
|
|
||||||
base_config = ConfigMaterial({'phase':{k if isinstance(k,int) else str(k): None for k in np.unique(phase)},
|
base_config = ConfigMaterial({'phase':{k if isinstance(k,int) else str(k): None for k in np.unique(phase)},
|
||||||
|
|
|
@ -365,7 +365,7 @@ class Grid:
|
||||||
|
|
||||||
Parameters
|
Parameters
|
||||||
----------
|
----------
|
||||||
fname : str or or pathlib.Path
|
fname : str or pathlib.Path
|
||||||
Filename of the DREAM.3D (HDF5) file.
|
Filename of the DREAM.3D (HDF5) file.
|
||||||
feature_IDs : str, optional
|
feature_IDs : str, optional
|
||||||
Name of the dataset containing the mapping between cells and
|
Name of the dataset containing the mapping between cells and
|
||||||
|
@ -401,22 +401,22 @@ class Grid:
|
||||||
orientation and phase are considered.
|
orientation and phase are considered.
|
||||||
|
|
||||||
"""
|
"""
|
||||||
b = util.DREAM3D_base_group(fname) if base_group is None else base_group
|
with h5py.File(fname, 'r') as f:
|
||||||
c = util.DREAM3D_cell_data_group(fname) if cell_data is None else cell_data
|
b = util.DREAM3D_base_group(f) if base_group is None else base_group
|
||||||
f = h5py.File(fname, 'r')
|
c = util.DREAM3D_cell_data_group(f) if cell_data is None else cell_data
|
||||||
|
|
||||||
cells = f['/'.join([b,'_SIMPL_GEOMETRY','DIMENSIONS'])][()]
|
cells = f['/'.join([b,'_SIMPL_GEOMETRY','DIMENSIONS'])][()]
|
||||||
size = f['/'.join([b,'_SIMPL_GEOMETRY','SPACING'])] * cells
|
size = f['/'.join([b,'_SIMPL_GEOMETRY','SPACING'])] * cells
|
||||||
origin = f['/'.join([b,'_SIMPL_GEOMETRY','ORIGIN'])][()]
|
origin = f['/'.join([b,'_SIMPL_GEOMETRY','ORIGIN'])][()]
|
||||||
|
|
||||||
if feature_IDs is None:
|
if feature_IDs is None:
|
||||||
phase = f['/'.join([b,c,phases])][()].reshape(-1,1)
|
phase = f['/'.join([b,c,phases])][()].reshape(-1,1)
|
||||||
O = Rotation.from_Euler_angles(f['/'.join([b,c,Euler_angles])]).as_quaternion().reshape(-1,4) # noqa
|
O = Rotation.from_Euler_angles(f['/'.join([b,c,Euler_angles])]).as_quaternion().reshape(-1,4) # noqa
|
||||||
unique,unique_inverse = np.unique(np.hstack([O,phase]),return_inverse=True,axis=0)
|
unique,unique_inverse = np.unique(np.hstack([O,phase]),return_inverse=True,axis=0)
|
||||||
ma = np.arange(cells.prod()) if len(unique) == cells.prod() else \
|
ma = np.arange(cells.prod()) if len(unique) == cells.prod() else \
|
||||||
np.arange(unique.size)[np.argsort(pd.unique(unique_inverse))][unique_inverse]
|
np.arange(unique.size)[np.argsort(pd.unique(unique_inverse))][unique_inverse]
|
||||||
else:
|
else:
|
||||||
ma = f['/'.join([b,c,feature_IDs])][()].flatten()
|
ma = f['/'.join([b,c,feature_IDs])][()].flatten()
|
||||||
|
|
||||||
return Grid(material = ma.reshape(cells,order='F'),
|
return Grid(material = ma.reshape(cells,order='F'),
|
||||||
size = size,
|
size = size,
|
||||||
|
|
Loading…
Reference in New Issue