pyerrors.input.hadrons
View Source
0import os 1import warnings 2from collections import Counter 3import h5py 4import numpy as np 5from ..obs import Obs, CObs 6from ..correlators import Corr 7from ..dirac import epsilon_tensor_rank4 8 9 10def _get_files(path, filestem, idl): 11 ls = os.listdir(path) 12 13 # Clean up file list 14 files = list(filter(lambda x: x.startswith(filestem), ls)) 15 16 if not files: 17 raise Exception('No files starting with', filestem, 'in folder', path) 18 19 def get_cnfg_number(n): 20 return int(n.replace(".h5", "")[len(filestem) + 1:]) # From python 3.9 onward the safer 'removesuffix' method can be used. 21 22 # Sort according to configuration number 23 files.sort(key=get_cnfg_number) 24 25 cnfg_numbers = [] 26 filtered_files = [] 27 for line in files: 28 no = get_cnfg_number(line) 29 if idl: 30 if no in list(idl): 31 filtered_files.append(line) 32 cnfg_numbers.append(no) 33 else: 34 filtered_files.append(line) 35 cnfg_numbers.append(no) 36 37 if idl: 38 if Counter(list(idl)) != Counter(cnfg_numbers): 39 raise Exception("Not all configurations specified in idl found, configurations " + str(list(Counter(list(idl)) - Counter(cnfg_numbers))) + " are missing.") 40 41 # Check that configurations are evenly spaced 42 dc = np.unique(np.diff(cnfg_numbers)) 43 if np.any(dc < 0): 44 raise Exception("Unsorted files") 45 if len(dc) == 1: 46 idx = range(cnfg_numbers[0], cnfg_numbers[-1] + dc[0], dc[0]) 47 elif idl: 48 idx = idl 49 warnings.warn("Configurations are not evenly spaced.", RuntimeWarning) 50 else: 51 raise Exception("Configurations are not evenly spaced.") 52 53 return filtered_files, idx 54 55 56def read_meson_hd5(path, filestem, ens_id, meson='meson_0', idl=None): 57 """Read hadrons meson hdf5 file and extract the meson labeled 'meson' 58 59 Parameters 60 ----------------- 61 path : str 62 path to the files to read 63 filestem : str 64 namestem of the files to read 65 ens_id : str 66 name of the ensemble, required for internal bookkeeping 67 meson : str 68 label of the meson to be extracted, standard value meson_0 which 69 corresponds to the pseudoscalar pseudoscalar two-point function. 70 idl : range 71 If specified only configurations in the given range are read in. 72 """ 73 74 files, idx = _get_files(path, filestem, idl) 75 76 tree = meson.rsplit('_')[0] 77 corr_data = [] 78 infos = [] 79 for hd5_file in files: 80 h5file = h5py.File(path + '/' + hd5_file, "r") 81 if not tree + '/' + meson in h5file: 82 raise Exception("Entry '" + meson + "' not contained in the files.") 83 raw_data = list(h5file[tree + '/' + meson + '/corr']) 84 real_data = [o[0] for o in raw_data] 85 corr_data.append(real_data) 86 if not infos: 87 for k, i in h5file[tree + '/' + meson].attrs.items(): 88 infos.append(k + ': ' + i[0].decode()) 89 h5file.close() 90 corr_data = np.array(corr_data) 91 92 l_obs = [] 93 for c in corr_data.T: 94 l_obs.append(Obs([c], [ens_id], idl=[idx])) 95 96 corr = Corr(l_obs) 97 corr.tag = r", ".join(infos) 98 return corr 99 100 101class Npr_matrix(np.ndarray): 102 103 def __new__(cls, input_array, mom_in=None, mom_out=None): 104 obj = np.asarray(input_array).view(cls) 105 obj.mom_in = mom_in 106 obj.mom_out = mom_out 107 return obj 108 109 @property 110 def g5H(self): 111 """Gamma_5 hermitean conjugate 112 113 Uses the fact that the propagator is gamma5 hermitean, so just the 114 in and out momenta of the propagator are exchanged. 115 """ 116 return Npr_matrix(self, 117 mom_in=self.mom_out, 118 mom_out=self.mom_in) 119 120 def _propagate_mom(self, other, name): 121 s_mom = getattr(self, name, None) 122 o_mom = getattr(other, name, None) 123 if s_mom is not None and o_mom is not None: 124 if not np.allclose(s_mom, o_mom): 125 raise Exception(name + ' does not match.') 126 return o_mom if o_mom is not None else s_mom 127 128 def __matmul__(self, other): 129 return self.__new__(Npr_matrix, 130 super().__matmul__(other), 131 self._propagate_mom(other, 'mom_in'), 132 self._propagate_mom(other, 'mom_out')) 133 134 def __array_finalize__(self, obj): 135 if obj is None: 136 return 137 self.mom_in = getattr(obj, 'mom_in', None) 138 self.mom_out = getattr(obj, 'mom_out', None) 139 140 141def read_ExternalLeg_hd5(path, filestem, ens_id, idl=None): 142 """Read hadrons ExternalLeg hdf5 file and output an array of CObs 143 144 Parameters 145 ---------- 146 path : str 147 path to the files to read 148 filestem : str 149 namestem of the files to read 150 ens_id : str 151 name of the ensemble, required for internal bookkeeping 152 idl : range 153 If specified only configurations in the given range are read in. 154 """ 155 156 files, idx = _get_files(path, filestem, idl) 157 158 mom = None 159 160 corr_data = [] 161 for hd5_file in files: 162 file = h5py.File(path + '/' + hd5_file, "r") 163 raw_data = file['ExternalLeg/corr'][0][0].view('complex') 164 corr_data.append(raw_data) 165 if mom is None: 166 mom = np.array(str(file['ExternalLeg/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 167 file.close() 168 corr_data = np.array(corr_data) 169 170 rolled_array = np.rollaxis(corr_data, 0, 5) 171 172 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 173 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]): 174 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx]) 175 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx]) 176 matrix[si, sj, ci, cj] = CObs(real, imag) 177 178 return Npr_matrix(matrix, mom_in=mom) 179 180 181def read_Bilinear_hd5(path, filestem, ens_id, idl=None): 182 """Read hadrons Bilinear hdf5 file and output an array of CObs 183 184 Parameters 185 ---------- 186 path : str 187 path to the files to read 188 filestem : str 189 namestem of the files to read 190 ens_id : str 191 name of the ensemble, required for internal bookkeeping 192 idl : range 193 If specified only configurations in the given range are read in. 194 """ 195 196 files, idx = _get_files(path, filestem, idl) 197 198 mom_in = None 199 mom_out = None 200 201 corr_data = {} 202 for hd5_file in files: 203 file = h5py.File(path + '/' + hd5_file, "r") 204 for i in range(16): 205 name = file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['gamma'][0].decode('UTF-8') 206 if name not in corr_data: 207 corr_data[name] = [] 208 raw_data = file['Bilinear/Bilinear_' + str(i) + '/corr'][0][0].view('complex') 209 corr_data[name].append(raw_data) 210 if mom_in is None: 211 mom_in = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 212 if mom_out is None: 213 mom_out = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float) 214 215 file.close() 216 217 result_dict = {} 218 219 for key, data in corr_data.items(): 220 local_data = np.array(data) 221 222 rolled_array = np.rollaxis(local_data, 0, 5) 223 224 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 225 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]): 226 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx]) 227 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx]) 228 matrix[si, sj, ci, cj] = CObs(real, imag) 229 230 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out) 231 232 return result_dict 233 234 235def read_Fourquark_hd5(path, filestem, ens_id, idl=None, vertices=["VA", "AV"]): 236 """Read hadrons FourquarkFullyConnected hdf5 file and output an array of CObs 237 238 Parameters 239 ---------- 240 path : str 241 path to the files to read 242 filestem : str 243 namestem of the files to read 244 ens_id : str 245 name of the ensemble, required for internal bookkeeping 246 idl : range 247 If specified only configurations in the given range are read in. 248 vertices : list 249 Vertex functions to be extracted. 250 """ 251 252 files, idx = _get_files(path, filestem, idl) 253 254 mom_in = None 255 mom_out = None 256 257 vertex_names = [] 258 for vertex in vertices: 259 vertex_names += _get_lorentz_names(vertex) 260 261 corr_data = {} 262 263 tree = 'FourQuarkFullyConnected/FourQuarkFullyConnected_' 264 265 for hd5_file in files: 266 file = h5py.File(path + '/' + hd5_file, "r") 267 268 for i in range(32): 269 name = (file[tree + str(i) + '/info'].attrs['gammaA'][0].decode('UTF-8'), file[tree + str(i) + '/info'].attrs['gammaB'][0].decode('UTF-8')) 270 if name in vertex_names: 271 if name not in corr_data: 272 corr_data[name] = [] 273 raw_data = file[tree + str(i) + '/corr'][0][0].view('complex') 274 corr_data[name].append(raw_data) 275 if mom_in is None: 276 mom_in = np.array(str(file[tree + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 277 if mom_out is None: 278 mom_out = np.array(str(file[tree + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float) 279 280 file.close() 281 282 intermediate_dict = {} 283 284 for vertex in vertices: 285 lorentz_names = _get_lorentz_names(vertex) 286 for v_name in lorentz_names: 287 if v_name in [('SigmaXY', 'SigmaZT'), 288 ('SigmaXT', 'SigmaYZ'), 289 ('SigmaYZ', 'SigmaXT'), 290 ('SigmaZT', 'SigmaXY')]: 291 sign = -1 292 else: 293 sign = 1 294 if vertex not in intermediate_dict: 295 intermediate_dict[vertex] = sign * np.array(corr_data[v_name]) 296 else: 297 intermediate_dict[vertex] += sign * np.array(corr_data[v_name]) 298 299 result_dict = {} 300 301 for key, data in intermediate_dict.items(): 302 303 rolled_array = np.moveaxis(data, 0, 8) 304 305 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 306 for index in np.ndindex(rolled_array.shape[:-1]): 307 real = Obs([rolled_array[index].real], [ens_id], idl=[idx]) 308 imag = Obs([rolled_array[index].imag], [ens_id], idl=[idx]) 309 matrix[index] = CObs(real, imag) 310 311 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out) 312 313 return result_dict 314 315 316def _get_lorentz_names(name): 317 lorentz_index = ['X', 'Y', 'Z', 'T'] 318 319 res = [] 320 321 if name == "TT": 322 for i in range(4): 323 for j in range(i + 1, 4): 324 res.append(("Sigma" + lorentz_index[i] + lorentz_index[j], "Sigma" + lorentz_index[i] + lorentz_index[j])) 325 return res 326 327 if name == "TTtilde": 328 for i in range(4): 329 for j in range(i + 1, 4): 330 for k in range(4): 331 for o in range(k + 1, 4): 332 fac = epsilon_tensor_rank4(i, j, k, o) 333 if not np.isclose(fac, 0.0): 334 res.append(("Sigma" + lorentz_index[i] + lorentz_index[j], "Sigma" + lorentz_index[k] + lorentz_index[o])) 335 return res 336 337 assert len(name) == 2 338 339 if 'S' in name or 'P' in name: 340 if not set(name) <= set(['S', 'P']): 341 raise Exception("'" + name + "' is not a Lorentz scalar") 342 343 g_names = {'S': 'Identity', 344 'P': 'Gamma5'} 345 346 res.append((g_names[name[0]], g_names[name[1]])) 347 348 else: 349 if not set(name) <= set(['V', 'A']): 350 raise Exception("'" + name + "' is not a Lorentz scalar") 351 352 for ind in lorentz_index: 353 res.append(('Gamma' + ind + (name[0] == 'A') * 'Gamma5', 354 'Gamma' + ind + (name[1] == 'A') * 'Gamma5')) 355 356 return res
View Source
57def read_meson_hd5(path, filestem, ens_id, meson='meson_0', idl=None): 58 """Read hadrons meson hdf5 file and extract the meson labeled 'meson' 59 60 Parameters 61 ----------------- 62 path : str 63 path to the files to read 64 filestem : str 65 namestem of the files to read 66 ens_id : str 67 name of the ensemble, required for internal bookkeeping 68 meson : str 69 label of the meson to be extracted, standard value meson_0 which 70 corresponds to the pseudoscalar pseudoscalar two-point function. 71 idl : range 72 If specified only configurations in the given range are read in. 73 """ 74 75 files, idx = _get_files(path, filestem, idl) 76 77 tree = meson.rsplit('_')[0] 78 corr_data = [] 79 infos = [] 80 for hd5_file in files: 81 h5file = h5py.File(path + '/' + hd5_file, "r") 82 if not tree + '/' + meson in h5file: 83 raise Exception("Entry '" + meson + "' not contained in the files.") 84 raw_data = list(h5file[tree + '/' + meson + '/corr']) 85 real_data = [o[0] for o in raw_data] 86 corr_data.append(real_data) 87 if not infos: 88 for k, i in h5file[tree + '/' + meson].attrs.items(): 89 infos.append(k + ': ' + i[0].decode()) 90 h5file.close() 91 corr_data = np.array(corr_data) 92 93 l_obs = [] 94 for c in corr_data.T: 95 l_obs.append(Obs([c], [ens_id], idl=[idx])) 96 97 corr = Corr(l_obs) 98 corr.tag = r", ".join(infos) 99 return corr
Read hadrons meson hdf5 file and extract the meson labeled 'meson'
Parameters
- path (str): path to the files to read
- filestem (str): namestem of the files to read
- ens_id (str): name of the ensemble, required for internal bookkeeping
- meson (str): label of the meson to be extracted, standard value meson_0 which corresponds to the pseudoscalar pseudoscalar two-point function.
- idl (range): If specified only configurations in the given range are read in.
View Source
102class Npr_matrix(np.ndarray): 103 104 def __new__(cls, input_array, mom_in=None, mom_out=None): 105 obj = np.asarray(input_array).view(cls) 106 obj.mom_in = mom_in 107 obj.mom_out = mom_out 108 return obj 109 110 @property 111 def g5H(self): 112 """Gamma_5 hermitean conjugate 113 114 Uses the fact that the propagator is gamma5 hermitean, so just the 115 in and out momenta of the propagator are exchanged. 116 """ 117 return Npr_matrix(self, 118 mom_in=self.mom_out, 119 mom_out=self.mom_in) 120 121 def _propagate_mom(self, other, name): 122 s_mom = getattr(self, name, None) 123 o_mom = getattr(other, name, None) 124 if s_mom is not None and o_mom is not None: 125 if not np.allclose(s_mom, o_mom): 126 raise Exception(name + ' does not match.') 127 return o_mom if o_mom is not None else s_mom 128 129 def __matmul__(self, other): 130 return self.__new__(Npr_matrix, 131 super().__matmul__(other), 132 self._propagate_mom(other, 'mom_in'), 133 self._propagate_mom(other, 'mom_out')) 134 135 def __array_finalize__(self, obj): 136 if obj is None: 137 return 138 self.mom_in = getattr(obj, 'mom_in', None) 139 self.mom_out = getattr(obj, 'mom_out', None)
ndarray(shape, dtype=float, buffer=None, offset=0, strides=None, order=None)
An array object represents a multidimensional, homogeneous array of fixed-size items. An associated data-type object describes the format of each element in the array (its byte-order, how many bytes it occupies in memory, whether it is an integer, a floating point number, or something else, etc.)
Arrays should be constructed using array
, zeros
or empty
(refer
to the See Also section below). The parameters given here refer to
a low-level method (ndarray(...)
) for instantiating an array.
For more information, refer to the numpy
module and examine the
methods and attributes of an array.
Parameters
- (for the __new__ method; see Notes below)
- shape (tuple of ints): Shape of created array.
- dtype (data-type, optional): Any object that can be interpreted as a numpy data type.
- buffer (object exposing buffer interface, optional): Used to fill the array with data.
- offset (int, optional): Offset of array data in buffer.
- strides (tuple of ints, optional): Strides of data in memory.
- order ({'C', 'F'}, optional): Row-major (C-style) or column-major (Fortran-style) order.
Attributes
- T (ndarray): Transpose of the array.
- data (buffer): The array's elements, in memory.
- dtype (dtype object): Describes the format of the elements in the array.
- flags (dict): Dictionary containing information related to memory use, e.g., 'C_CONTIGUOUS', 'OWNDATA', 'WRITEABLE', etc.
- flat (numpy.flatiter object):
Flattened version of the array as an iterator. The iterator
allows assignments, e.g.,
x.flat = 3
(Seendarray.flat
for assignment examples; TODO). - imag (ndarray): Imaginary part of the array.
- real (ndarray): Real part of the array.
- size (int): Number of elements in the array.
- itemsize (int): The memory use of each array element in bytes.
- nbytes (int):
The total number of bytes required to store the array data,
i.e.,
itemsize * size
. - ndim (int): The array's number of dimensions.
- shape (tuple of ints): Shape of the array.
- strides (tuple of ints):
The step-size required to move from one element to the next in
memory. For example, a contiguous
(3, 4)
array of typeint16
in C-order has strides(8, 2)
. This implies that to move from element to element in memory requires jumps of 2 bytes. To move from row-to-row, one needs to jump 8 bytes at a time (2 * 4
). - ctypes (ctypes object): Class containing properties of the array needed for interaction with ctypes.
- base (ndarray):
If the array is a view into another array, that array is its
base
(unless that array is also a view). Thebase
array is where the array data is actually stored.
See Also
array
: Construct an array.
zeros
: Create an array, each element of which is zero.
empty
: Create an array, but leave its allocated memory unchanged (i.e.,
it contains "garbage").
dtype
: Create a data-type.
numpy.typing.NDArray
: An ndarray alias :term:generic <generic type>
w.r.t. its dtype.type <numpy.dtype.type>
.
Notes
There are two modes of creating an array using __new__
:
- If
buffer
is None, then onlyshape
,dtype
, andorder
are used. - If
buffer
is an object exposing the buffer interface, then all keywords are interpreted.
No __init__
method is needed because the array is fully initialized
after the __new__
method.
Examples
These examples illustrate the low-level ndarray
constructor. Refer
to the See Also
section above for easier ways of constructing an
ndarray.
First mode, buffer
is None:
>>> np.ndarray(shape=(2,2), dtype=float, order='F')
array([[0.0e+000, 0.0e+000], # random
[ nan, 2.5e-323]])
Second mode:
>>> np.ndarray((2,), buffer=np.array([1,2,3]),
... offset=np.int_().itemsize,
... dtype=int) # offset = 1*itemsize, i.e. skip first element
array([2, 3])
Gamma_5 hermitean conjugate
Uses the fact that the propagator is gamma5 hermitean, so just the in and out momenta of the propagator are exchanged.
Inherited Members
- numpy.ndarray
- dumps
- dump
- all
- any
- argmax
- argmin
- argpartition
- argsort
- astype
- byteswap
- choose
- clip
- compress
- conj
- conjugate
- copy
- cumprod
- cumsum
- diagonal
- dot
- fill
- flatten
- getfield
- item
- itemset
- max
- mean
- min
- newbyteorder
- nonzero
- partition
- prod
- ptp
- put
- ravel
- repeat
- reshape
- resize
- round
- searchsorted
- setfield
- setflags
- sort
- squeeze
- std
- sum
- swapaxes
- take
- tobytes
- tofile
- tolist
- tostring
- trace
- transpose
- var
- view
- ndim
- flags
- shape
- strides
- data
- itemsize
- size
- nbytes
- base
- dtype
- real
- imag
- flat
- ctypes
- T
View Source
142def read_ExternalLeg_hd5(path, filestem, ens_id, idl=None): 143 """Read hadrons ExternalLeg hdf5 file and output an array of CObs 144 145 Parameters 146 ---------- 147 path : str 148 path to the files to read 149 filestem : str 150 namestem of the files to read 151 ens_id : str 152 name of the ensemble, required for internal bookkeeping 153 idl : range 154 If specified only configurations in the given range are read in. 155 """ 156 157 files, idx = _get_files(path, filestem, idl) 158 159 mom = None 160 161 corr_data = [] 162 for hd5_file in files: 163 file = h5py.File(path + '/' + hd5_file, "r") 164 raw_data = file['ExternalLeg/corr'][0][0].view('complex') 165 corr_data.append(raw_data) 166 if mom is None: 167 mom = np.array(str(file['ExternalLeg/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 168 file.close() 169 corr_data = np.array(corr_data) 170 171 rolled_array = np.rollaxis(corr_data, 0, 5) 172 173 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 174 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]): 175 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx]) 176 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx]) 177 matrix[si, sj, ci, cj] = CObs(real, imag) 178 179 return Npr_matrix(matrix, mom_in=mom)
Read hadrons ExternalLeg hdf5 file and output an array of CObs
Parameters
- path (str): path to the files to read
- filestem (str): namestem of the files to read
- ens_id (str): name of the ensemble, required for internal bookkeeping
- idl (range): If specified only configurations in the given range are read in.
View Source
182def read_Bilinear_hd5(path, filestem, ens_id, idl=None): 183 """Read hadrons Bilinear hdf5 file and output an array of CObs 184 185 Parameters 186 ---------- 187 path : str 188 path to the files to read 189 filestem : str 190 namestem of the files to read 191 ens_id : str 192 name of the ensemble, required for internal bookkeeping 193 idl : range 194 If specified only configurations in the given range are read in. 195 """ 196 197 files, idx = _get_files(path, filestem, idl) 198 199 mom_in = None 200 mom_out = None 201 202 corr_data = {} 203 for hd5_file in files: 204 file = h5py.File(path + '/' + hd5_file, "r") 205 for i in range(16): 206 name = file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['gamma'][0].decode('UTF-8') 207 if name not in corr_data: 208 corr_data[name] = [] 209 raw_data = file['Bilinear/Bilinear_' + str(i) + '/corr'][0][0].view('complex') 210 corr_data[name].append(raw_data) 211 if mom_in is None: 212 mom_in = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 213 if mom_out is None: 214 mom_out = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float) 215 216 file.close() 217 218 result_dict = {} 219 220 for key, data in corr_data.items(): 221 local_data = np.array(data) 222 223 rolled_array = np.rollaxis(local_data, 0, 5) 224 225 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 226 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]): 227 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx]) 228 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx]) 229 matrix[si, sj, ci, cj] = CObs(real, imag) 230 231 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out) 232 233 return result_dict
Read hadrons Bilinear hdf5 file and output an array of CObs
Parameters
- path (str): path to the files to read
- filestem (str): namestem of the files to read
- ens_id (str): name of the ensemble, required for internal bookkeeping
- idl (range): If specified only configurations in the given range are read in.
View Source
236def read_Fourquark_hd5(path, filestem, ens_id, idl=None, vertices=["VA", "AV"]): 237 """Read hadrons FourquarkFullyConnected hdf5 file and output an array of CObs 238 239 Parameters 240 ---------- 241 path : str 242 path to the files to read 243 filestem : str 244 namestem of the files to read 245 ens_id : str 246 name of the ensemble, required for internal bookkeeping 247 idl : range 248 If specified only configurations in the given range are read in. 249 vertices : list 250 Vertex functions to be extracted. 251 """ 252 253 files, idx = _get_files(path, filestem, idl) 254 255 mom_in = None 256 mom_out = None 257 258 vertex_names = [] 259 for vertex in vertices: 260 vertex_names += _get_lorentz_names(vertex) 261 262 corr_data = {} 263 264 tree = 'FourQuarkFullyConnected/FourQuarkFullyConnected_' 265 266 for hd5_file in files: 267 file = h5py.File(path + '/' + hd5_file, "r") 268 269 for i in range(32): 270 name = (file[tree + str(i) + '/info'].attrs['gammaA'][0].decode('UTF-8'), file[tree + str(i) + '/info'].attrs['gammaB'][0].decode('UTF-8')) 271 if name in vertex_names: 272 if name not in corr_data: 273 corr_data[name] = [] 274 raw_data = file[tree + str(i) + '/corr'][0][0].view('complex') 275 corr_data[name].append(raw_data) 276 if mom_in is None: 277 mom_in = np.array(str(file[tree + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 278 if mom_out is None: 279 mom_out = np.array(str(file[tree + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float) 280 281 file.close() 282 283 intermediate_dict = {} 284 285 for vertex in vertices: 286 lorentz_names = _get_lorentz_names(vertex) 287 for v_name in lorentz_names: 288 if v_name in [('SigmaXY', 'SigmaZT'), 289 ('SigmaXT', 'SigmaYZ'), 290 ('SigmaYZ', 'SigmaXT'), 291 ('SigmaZT', 'SigmaXY')]: 292 sign = -1 293 else: 294 sign = 1 295 if vertex not in intermediate_dict: 296 intermediate_dict[vertex] = sign * np.array(corr_data[v_name]) 297 else: 298 intermediate_dict[vertex] += sign * np.array(corr_data[v_name]) 299 300 result_dict = {} 301 302 for key, data in intermediate_dict.items(): 303 304 rolled_array = np.moveaxis(data, 0, 8) 305 306 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 307 for index in np.ndindex(rolled_array.shape[:-1]): 308 real = Obs([rolled_array[index].real], [ens_id], idl=[idx]) 309 imag = Obs([rolled_array[index].imag], [ens_id], idl=[idx]) 310 matrix[index] = CObs(real, imag) 311 312 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out) 313 314 return result_dict
Read hadrons FourquarkFullyConnected hdf5 file and output an array of CObs
Parameters
- path (str): path to the files to read
- filestem (str): namestem of the files to read
- ens_id (str): name of the ensemble, required for internal bookkeeping
- idl (range): If specified only configurations in the given range are read in.
- vertices (list): Vertex functions to be extracted.