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
def read_meson_hd5(path, filestem, ens_id, meson='meson_0', idl=None)
 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.
class Npr_matrix(numpy.ndarray):
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 (See ndarray.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 type int16 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). The base 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__:

  1. If buffer is None, then only shape, dtype, and order are used.
  2. 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])
Npr_matrix()
g5H

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
def read_ExternalLeg_hd5(path, filestem, ens_id, idl=None)
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.
def read_Bilinear_hd5(path, filestem, ens_id, idl=None)
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.
def read_Fourquark_hd5(path, filestem, ens_id, idl=None, vertices=['VA', 'AV'])
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.