pyerrors.input.hadrons
1import os 2import warnings 3from collections import Counter 4import h5py 5import numpy as np 6from ..obs import Obs, CObs 7from ..correlators import Corr 8from ..dirac import epsilon_tensor_rank4 9 10 11def _get_files(path, filestem, idl): 12 ls = os.listdir(path) 13 14 # Clean up file list 15 files = list(filter(lambda x: x.startswith(filestem + "."), ls)) 16 17 if not files: 18 raise Exception('No files starting with', filestem, 'in folder', path) 19 20 def get_cnfg_number(n): 21 return int(n.replace(".h5", "")[len(filestem) + 1:]) # From python 3.9 onward the safer 'removesuffix' method can be used. 22 23 # Sort according to configuration number 24 files.sort(key=get_cnfg_number) 25 26 cnfg_numbers = [] 27 filtered_files = [] 28 for line in files: 29 no = get_cnfg_number(line) 30 if idl: 31 if no in list(idl): 32 filtered_files.append(line) 33 cnfg_numbers.append(no) 34 else: 35 filtered_files.append(line) 36 cnfg_numbers.append(no) 37 38 if idl: 39 if Counter(list(idl)) != Counter(cnfg_numbers): 40 raise Exception("Not all configurations specified in idl found, configurations " + str(list(Counter(list(idl)) - Counter(cnfg_numbers))) + " are missing.") 41 42 # Check that configurations are evenly spaced 43 dc = np.unique(np.diff(cnfg_numbers)) 44 if np.any(dc < 0): 45 raise Exception("Unsorted files") 46 if len(dc) == 1: 47 idx = range(cnfg_numbers[0], cnfg_numbers[-1] + dc[0], dc[0]) 48 elif idl: 49 idx = idl 50 warnings.warn("Configurations are not evenly spaced.", RuntimeWarning) 51 else: 52 raise Exception("Configurations are not evenly spaced.") 53 54 return filtered_files, idx 55 56 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 = h5file[tree + '/' + meson + '/corr'] 85 real_data = raw_data[:]["re"].astype(np.double) 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 100 101 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) 140 141 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) 180 181 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 234 235 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 315 316 317def _get_lorentz_names(name): 318 lorentz_index = ['X', 'Y', 'Z', 'T'] 319 320 res = [] 321 322 if name == "TT": 323 for i in range(4): 324 for j in range(i + 1, 4): 325 res.append(("Sigma" + lorentz_index[i] + lorentz_index[j], "Sigma" + lorentz_index[i] + lorentz_index[j])) 326 return res 327 328 if name == "TTtilde": 329 for i in range(4): 330 for j in range(i + 1, 4): 331 for k in range(4): 332 for o in range(k + 1, 4): 333 fac = epsilon_tensor_rank4(i, j, k, o) 334 if not np.isclose(fac, 0.0): 335 res.append(("Sigma" + lorentz_index[i] + lorentz_index[j], "Sigma" + lorentz_index[k] + lorentz_index[o])) 336 return res 337 338 assert len(name) == 2 339 340 if 'S' in name or 'P' in name: 341 if not set(name) <= set(['S', 'P']): 342 raise Exception("'" + name + "' is not a Lorentz scalar") 343 344 g_names = {'S': 'Identity', 345 'P': 'Gamma5'} 346 347 res.append((g_names[name[0]], g_names[name[1]])) 348 349 else: 350 if not set(name) <= set(['V', 'A']): 351 raise Exception("'" + name + "' is not a Lorentz scalar") 352 353 for ind in lorentz_index: 354 res.append(('Gamma' + ind + (name[0] == 'A') * 'Gamma5', 355 'Gamma' + ind + (name[1] == 'A') * 'Gamma5')) 356 357 return res
58def read_meson_hd5(path, filestem, ens_id, meson='meson_0', idl=None): 59 """Read hadrons meson hdf5 file and extract the meson labeled 'meson' 60 61 Parameters 62 ----------------- 63 path : str 64 path to the files to read 65 filestem : str 66 namestem of the files to read 67 ens_id : str 68 name of the ensemble, required for internal bookkeeping 69 meson : str 70 label of the meson to be extracted, standard value meson_0 which 71 corresponds to the pseudoscalar pseudoscalar two-point function. 72 idl : range 73 If specified only configurations in the given range are read in. 74 """ 75 76 files, idx = _get_files(path, filestem, idl) 77 78 tree = meson.rsplit('_')[0] 79 corr_data = [] 80 infos = [] 81 for hd5_file in files: 82 h5file = h5py.File(path + '/' + hd5_file, "r") 83 if not tree + '/' + meson in h5file: 84 raise Exception("Entry '" + meson + "' not contained in the files.") 85 raw_data = h5file[tree + '/' + meson + '/corr'] 86 real_data = raw_data[:]["re"].astype(np.double) 87 corr_data.append(real_data) 88 if not infos: 89 for k, i in h5file[tree + '/' + meson].attrs.items(): 90 infos.append(k + ': ' + i[0].decode()) 91 h5file.close() 92 corr_data = np.array(corr_data) 93 94 l_obs = [] 95 for c in corr_data.T: 96 l_obs.append(Obs([c], [ens_id], idl=[idx])) 97 98 corr = Corr(l_obs) 99 corr.tag = r", ".join(infos) 100 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.
103class Npr_matrix(np.ndarray): 104 105 def __new__(cls, input_array, mom_in=None, mom_out=None): 106 obj = np.asarray(input_array).view(cls) 107 obj.mom_in = mom_in 108 obj.mom_out = mom_out 109 return obj 110 111 @property 112 def g5H(self): 113 """Gamma_5 hermitean conjugate 114 115 Uses the fact that the propagator is gamma5 hermitean, so just the 116 in and out momenta of the propagator are exchanged. 117 """ 118 return Npr_matrix(self, 119 mom_in=self.mom_out, 120 mom_out=self.mom_in) 121 122 def _propagate_mom(self, other, name): 123 s_mom = getattr(self, name, None) 124 o_mom = getattr(other, name, None) 125 if s_mom is not None and o_mom is not None: 126 if not np.allclose(s_mom, o_mom): 127 raise Exception(name + ' does not match.') 128 return o_mom if o_mom is not None else s_mom 129 130 def __matmul__(self, other): 131 return self.__new__(Npr_matrix, 132 super().__matmul__(other), 133 self._propagate_mom(other, 'mom_in'), 134 self._propagate_mom(other, 'mom_out')) 135 136 def __array_finalize__(self, obj): 137 if obj is None: 138 return 139 self.mom_in = getattr(obj, 'mom_in', None) 140 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
143def read_ExternalLeg_hd5(path, filestem, ens_id, idl=None): 144 """Read hadrons ExternalLeg hdf5 file and output an array of CObs 145 146 Parameters 147 ---------- 148 path : str 149 path to the files to read 150 filestem : str 151 namestem of the files to read 152 ens_id : str 153 name of the ensemble, required for internal bookkeeping 154 idl : range 155 If specified only configurations in the given range are read in. 156 """ 157 158 files, idx = _get_files(path, filestem, idl) 159 160 mom = None 161 162 corr_data = [] 163 for hd5_file in files: 164 file = h5py.File(path + '/' + hd5_file, "r") 165 raw_data = file['ExternalLeg/corr'][0][0].view('complex') 166 corr_data.append(raw_data) 167 if mom is None: 168 mom = np.array(str(file['ExternalLeg/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 169 file.close() 170 corr_data = np.array(corr_data) 171 172 rolled_array = np.rollaxis(corr_data, 0, 5) 173 174 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 175 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]): 176 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx]) 177 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx]) 178 matrix[si, sj, ci, cj] = CObs(real, imag) 179 180 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.
183def read_Bilinear_hd5(path, filestem, ens_id, idl=None): 184 """Read hadrons Bilinear hdf5 file and output an array of CObs 185 186 Parameters 187 ---------- 188 path : str 189 path to the files to read 190 filestem : str 191 namestem of the files to read 192 ens_id : str 193 name of the ensemble, required for internal bookkeeping 194 idl : range 195 If specified only configurations in the given range are read in. 196 """ 197 198 files, idx = _get_files(path, filestem, idl) 199 200 mom_in = None 201 mom_out = None 202 203 corr_data = {} 204 for hd5_file in files: 205 file = h5py.File(path + '/' + hd5_file, "r") 206 for i in range(16): 207 name = file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['gamma'][0].decode('UTF-8') 208 if name not in corr_data: 209 corr_data[name] = [] 210 raw_data = file['Bilinear/Bilinear_' + str(i) + '/corr'][0][0].view('complex') 211 corr_data[name].append(raw_data) 212 if mom_in is None: 213 mom_in = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 214 if mom_out is None: 215 mom_out = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float) 216 217 file.close() 218 219 result_dict = {} 220 221 for key, data in corr_data.items(): 222 local_data = np.array(data) 223 224 rolled_array = np.rollaxis(local_data, 0, 5) 225 226 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 227 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]): 228 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx]) 229 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx]) 230 matrix[si, sj, ci, cj] = CObs(real, imag) 231 232 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out) 233 234 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.
237def read_Fourquark_hd5(path, filestem, ens_id, idl=None, vertices=["VA", "AV"]): 238 """Read hadrons FourquarkFullyConnected hdf5 file and output an array of CObs 239 240 Parameters 241 ---------- 242 path : str 243 path to the files to read 244 filestem : str 245 namestem of the files to read 246 ens_id : str 247 name of the ensemble, required for internal bookkeeping 248 idl : range 249 If specified only configurations in the given range are read in. 250 vertices : list 251 Vertex functions to be extracted. 252 """ 253 254 files, idx = _get_files(path, filestem, idl) 255 256 mom_in = None 257 mom_out = None 258 259 vertex_names = [] 260 for vertex in vertices: 261 vertex_names += _get_lorentz_names(vertex) 262 263 corr_data = {} 264 265 tree = 'FourQuarkFullyConnected/FourQuarkFullyConnected_' 266 267 for hd5_file in files: 268 file = h5py.File(path + '/' + hd5_file, "r") 269 270 for i in range(32): 271 name = (file[tree + str(i) + '/info'].attrs['gammaA'][0].decode('UTF-8'), file[tree + str(i) + '/info'].attrs['gammaB'][0].decode('UTF-8')) 272 if name in vertex_names: 273 if name not in corr_data: 274 corr_data[name] = [] 275 raw_data = file[tree + str(i) + '/corr'][0][0].view('complex') 276 corr_data[name].append(raw_data) 277 if mom_in is None: 278 mom_in = np.array(str(file[tree + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float) 279 if mom_out is None: 280 mom_out = np.array(str(file[tree + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float) 281 282 file.close() 283 284 intermediate_dict = {} 285 286 for vertex in vertices: 287 lorentz_names = _get_lorentz_names(vertex) 288 for v_name in lorentz_names: 289 if v_name in [('SigmaXY', 'SigmaZT'), 290 ('SigmaXT', 'SigmaYZ'), 291 ('SigmaYZ', 'SigmaXT'), 292 ('SigmaZT', 'SigmaXY')]: 293 sign = -1 294 else: 295 sign = 1 296 if vertex not in intermediate_dict: 297 intermediate_dict[vertex] = sign * np.array(corr_data[v_name]) 298 else: 299 intermediate_dict[vertex] += sign * np.array(corr_data[v_name]) 300 301 result_dict = {} 302 303 for key, data in intermediate_dict.items(): 304 305 rolled_array = np.moveaxis(data, 0, 8) 306 307 matrix = np.empty((rolled_array.shape[:-1]), dtype=object) 308 for index in np.ndindex(rolled_array.shape[:-1]): 309 real = Obs([rolled_array[index].real], [ens_id], idl=[idx]) 310 imag = Obs([rolled_array[index].imag], [ens_id], idl=[idx]) 311 matrix[index] = CObs(real, imag) 312 313 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out) 314 315 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.