diff --git a/docs/pyerrors/input/hadrons.html b/docs/pyerrors/input/hadrons.html
index e6ff3b26..18060141 100644
--- a/docs/pyerrors/input/hadrons.html
+++ b/docs/pyerrors/input/hadrons.html
@@ -150,415 +150,414 @@
48 idx = range(cnfg_numbers[0], cnfg_numbers[-1] + dc[0], dc[0])
49 elif idl:
50 idx = idl
- 51 warnings.warn("Configurations are not evenly spaced.", RuntimeWarning)
- 52 else:
- 53 raise Exception("Configurations are not evenly spaced. Provide an idl if you want to proceed with this set of configurations.")
- 54
- 55 return filtered_files, idx
+ 51 else:
+ 52 raise Exception("Configurations are not evenly spaced. Provide an idl if you want to proceed with this set of configurations.")
+ 53
+ 54 return filtered_files, idx
+ 55
56
- 57
- 58def read_meson_hd5(path, filestem, ens_id, meson='meson_0', idl=None, gammas=None):
- 59 r'''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 gammas : tuple of strings
- 73 Instrad of a meson label one can also provide a tuple of two strings
- 74 indicating the gamma matrices at source and sink.
- 75 ("Gamma5", "Gamma5") corresponds to the pseudoscalar pseudoscalar
- 76 two-point function. The gammas argument dominateds over meson.
- 77 idl : range
- 78 If specified only configurations in the given range are read in.
- 79 '''
- 80
- 81 files, idx = _get_files(path, filestem, idl)
- 82
- 83 tree = meson.rsplit('_')[0]
- 84 if gammas is not None:
- 85 h5file = h5py.File(path + '/' + files[0], "r")
- 86 found_meson = None
- 87 for key in h5file[tree].keys():
- 88 if gammas[0] == h5file[tree][key].attrs["gamma_snk"][0].decode() and h5file[tree][key].attrs["gamma_src"][0].decode() == gammas[1]:
- 89 found_meson = key
- 90 break
- 91 h5file.close()
- 92 if found_meson:
- 93 meson = found_meson
- 94 else:
- 95 raise Exception("Source Sink combination " + str(gammas) + " not found.")
- 96
- 97 corr_data = []
- 98 infos = []
- 99 for hd5_file in files:
-100 h5file = h5py.File(path + '/' + hd5_file, "r")
-101 if not tree + '/' + meson in h5file:
-102 raise Exception("Entry '" + meson + "' not contained in the files.")
-103 raw_data = h5file[tree + '/' + meson + '/corr']
-104 real_data = raw_data[:]["re"].astype(np.double)
-105 corr_data.append(real_data)
-106 if not infos:
-107 for k, i in h5file[tree + '/' + meson].attrs.items():
-108 infos.append(k + ': ' + i[0].decode())
-109 h5file.close()
-110 corr_data = np.array(corr_data)
-111
-112 l_obs = []
-113 for c in corr_data.T:
-114 l_obs.append(Obs([c], [ens_id], idl=[idx]))
-115
-116 corr = Corr(l_obs)
-117 corr.tag = r", ".join(infos)
-118 return corr
+ 57def read_meson_hd5(path, filestem, ens_id, meson='meson_0', idl=None, gammas=None):
+ 58 r'''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 gammas : tuple of strings
+ 72 Instrad of a meson label one can also provide a tuple of two strings
+ 73 indicating the gamma matrices at source and sink.
+ 74 ("Gamma5", "Gamma5") corresponds to the pseudoscalar pseudoscalar
+ 75 two-point function. The gammas argument dominateds over meson.
+ 76 idl : range
+ 77 If specified only configurations in the given range are read in.
+ 78 '''
+ 79
+ 80 files, idx = _get_files(path, filestem, idl)
+ 81
+ 82 tree = meson.rsplit('_')[0]
+ 83 if gammas is not None:
+ 84 h5file = h5py.File(path + '/' + files[0], "r")
+ 85 found_meson = None
+ 86 for key in h5file[tree].keys():
+ 87 if gammas[0] == h5file[tree][key].attrs["gamma_snk"][0].decode() and h5file[tree][key].attrs["gamma_src"][0].decode() == gammas[1]:
+ 88 found_meson = key
+ 89 break
+ 90 h5file.close()
+ 91 if found_meson:
+ 92 meson = found_meson
+ 93 else:
+ 94 raise Exception("Source Sink combination " + str(gammas) + " not found.")
+ 95
+ 96 corr_data = []
+ 97 infos = []
+ 98 for hd5_file in files:
+ 99 h5file = h5py.File(path + '/' + hd5_file, "r")
+100 if not tree + '/' + meson in h5file:
+101 raise Exception("Entry '" + meson + "' not contained in the files.")
+102 raw_data = h5file[tree + '/' + meson + '/corr']
+103 real_data = raw_data[:]["re"].astype(np.double)
+104 corr_data.append(real_data)
+105 if not infos:
+106 for k, i in h5file[tree + '/' + meson].attrs.items():
+107 infos.append(k + ': ' + i[0].decode())
+108 h5file.close()
+109 corr_data = np.array(corr_data)
+110
+111 l_obs = []
+112 for c in corr_data.T:
+113 l_obs.append(Obs([c], [ens_id], idl=[idx]))
+114
+115 corr = Corr(l_obs)
+116 corr.tag = r", ".join(infos)
+117 return corr
+118
119
-120
-121def read_DistillationContraction_hd5(path, ens_id, diagrams=["direct"], idl=None):
-122 """Read hadrons DistillationContraction hdf5 files in given directory structure
-123
-124 Parameters
-125 -----------------
-126 path : str
-127 path to the directories to read
-128 ens_id : str
-129 name of the ensemble, required for internal bookkeeping
-130 diagrams : list
-131 List of strings of the diagrams to extract, e.g. ["direct", "box", "cross"].
-132 idl : range
-133 If specified only configurations in the given range are read in.
-134 """
-135
-136 res_dict = {}
-137
-138 directories, idx = _get_files(path, "data", idl)
-139
-140 explore_path = Path(path + "/" + directories[0])
-141
-142 for explore_file in explore_path.iterdir():
-143 if explore_file.is_file():
-144 stem = explore_file.with_suffix("").with_suffix("").as_posix().split("/")[-1]
-145 else:
-146 continue
-147
-148 file_list = []
-149 for dir in directories:
-150 tmp_path = Path(path + "/" + dir)
-151 file_list.append((tmp_path / stem).as_posix() + tmp_path.suffix + ".h5")
-152
-153 corr_data = {}
-154
-155 for diagram in diagrams:
-156 corr_data[diagram] = []
-157
-158 for n_file, (hd5_file, n_traj) in enumerate(zip(file_list, list(idx))):
-159 h5file = h5py.File(hd5_file)
-160
-161 if n_file == 0:
-162 if h5file["DistillationContraction/Metadata"].attrs.get("TimeSources")[0].decode() != "0...":
-163 raise Exception("Routine is only implemented for files containing inversions on all timeslices.")
-164
-165 Nt = h5file["DistillationContraction/Metadata"].attrs.get("Nt")[0]
-166
-167 identifier = []
-168 for in_file in range(len(h5file["DistillationContraction/Metadata/DmfInputFiles"].attrs.keys()) - 1):
-169 encoded_info = h5file["DistillationContraction/Metadata/DmfInputFiles"].attrs.get("DmfInputFiles_" + str(in_file))
-170 full_info = encoded_info[0].decode().split("/")[-1].replace(".h5", "").split("_")
-171 my_tuple = (full_info[0], full_info[1][1:], full_info[2], full_info[3])
-172 identifier.append(my_tuple)
-173 identifier = tuple(identifier)
-174 # "DistillationContraction/Metadata/DmfSuffix" contains info about different quarks, irrelevant in the SU(3) case.
-175
-176 for diagram in diagrams:
-177 real_data = np.zeros(Nt)
-178 for x0 in range(Nt):
-179 raw_data = h5file["DistillationContraction/Correlators/" + diagram + "/" + str(x0)][:]["re"].astype(np.double)
-180 real_data += np.roll(raw_data, -x0)
-181 real_data /= Nt
-182
-183 corr_data[diagram].append(real_data)
-184 h5file.close()
-185
-186 res_dict[str(identifier)] = {}
-187
-188 for diagram in diagrams:
-189
-190 tmp_data = np.array(corr_data[diagram])
-191
-192 l_obs = []
-193 for c in tmp_data.T:
-194 l_obs.append(Obs([c], [ens_id], idl=[idx]))
-195
-196 corr = Corr(l_obs)
-197 corr.tag = str(identifier)
-198
-199 res_dict[str(identifier)][diagram] = corr
-200
-201 return res_dict
+120def read_DistillationContraction_hd5(path, ens_id, diagrams=["direct"], idl=None):
+121 """Read hadrons DistillationContraction hdf5 files in given directory structure
+122
+123 Parameters
+124 -----------------
+125 path : str
+126 path to the directories to read
+127 ens_id : str
+128 name of the ensemble, required for internal bookkeeping
+129 diagrams : list
+130 List of strings of the diagrams to extract, e.g. ["direct", "box", "cross"].
+131 idl : range
+132 If specified only configurations in the given range are read in.
+133 """
+134
+135 res_dict = {}
+136
+137 directories, idx = _get_files(path, "data", idl)
+138
+139 explore_path = Path(path + "/" + directories[0])
+140
+141 for explore_file in explore_path.iterdir():
+142 if explore_file.is_file():
+143 stem = explore_file.with_suffix("").with_suffix("").as_posix().split("/")[-1]
+144 else:
+145 continue
+146
+147 file_list = []
+148 for dir in directories:
+149 tmp_path = Path(path + "/" + dir)
+150 file_list.append((tmp_path / stem).as_posix() + tmp_path.suffix + ".h5")
+151
+152 corr_data = {}
+153
+154 for diagram in diagrams:
+155 corr_data[diagram] = []
+156
+157 for n_file, (hd5_file, n_traj) in enumerate(zip(file_list, list(idx))):
+158 h5file = h5py.File(hd5_file)
+159
+160 if n_file == 0:
+161 if h5file["DistillationContraction/Metadata"].attrs.get("TimeSources")[0].decode() != "0...":
+162 raise Exception("Routine is only implemented for files containing inversions on all timeslices.")
+163
+164 Nt = h5file["DistillationContraction/Metadata"].attrs.get("Nt")[0]
+165
+166 identifier = []
+167 for in_file in range(len(h5file["DistillationContraction/Metadata/DmfInputFiles"].attrs.keys()) - 1):
+168 encoded_info = h5file["DistillationContraction/Metadata/DmfInputFiles"].attrs.get("DmfInputFiles_" + str(in_file))
+169 full_info = encoded_info[0].decode().split("/")[-1].replace(".h5", "").split("_")
+170 my_tuple = (full_info[0], full_info[1][1:], full_info[2], full_info[3])
+171 identifier.append(my_tuple)
+172 identifier = tuple(identifier)
+173 # "DistillationContraction/Metadata/DmfSuffix" contains info about different quarks, irrelevant in the SU(3) case.
+174
+175 for diagram in diagrams:
+176 real_data = np.zeros(Nt)
+177 for x0 in range(Nt):
+178 raw_data = h5file["DistillationContraction/Correlators/" + diagram + "/" + str(x0)][:]["re"].astype(np.double)
+179 real_data += np.roll(raw_data, -x0)
+180 real_data /= Nt
+181
+182 corr_data[diagram].append(real_data)
+183 h5file.close()
+184
+185 res_dict[str(identifier)] = {}
+186
+187 for diagram in diagrams:
+188
+189 tmp_data = np.array(corr_data[diagram])
+190
+191 l_obs = []
+192 for c in tmp_data.T:
+193 l_obs.append(Obs([c], [ens_id], idl=[idx]))
+194
+195 corr = Corr(l_obs)
+196 corr.tag = str(identifier)
+197
+198 res_dict[str(identifier)][diagram] = corr
+199
+200 return res_dict
+201
202
-203
-204class Npr_matrix(np.ndarray):
-205
-206 def __new__(cls, input_array, mom_in=None, mom_out=None):
-207 obj = np.asarray(input_array).view(cls)
-208 obj.mom_in = mom_in
-209 obj.mom_out = mom_out
-210 return obj
-211
-212 @property
-213 def g5H(self):
-214 """Gamma_5 hermitean conjugate
-215
-216 Uses the fact that the propagator is gamma5 hermitean, so just the
-217 in and out momenta of the propagator are exchanged.
-218 """
-219 return Npr_matrix(self,
-220 mom_in=self.mom_out,
-221 mom_out=self.mom_in)
-222
-223 def _propagate_mom(self, other, name):
-224 s_mom = getattr(self, name, None)
-225 o_mom = getattr(other, name, None)
-226 if s_mom is not None and o_mom is not None:
-227 if not np.allclose(s_mom, o_mom):
-228 raise Exception(name + ' does not match.')
-229 return o_mom if o_mom is not None else s_mom
-230
-231 def __matmul__(self, other):
-232 return self.__new__(Npr_matrix,
-233 super().__matmul__(other),
-234 self._propagate_mom(other, 'mom_in'),
-235 self._propagate_mom(other, 'mom_out'))
-236
-237 def __array_finalize__(self, obj):
-238 if obj is None:
-239 return
-240 self.mom_in = getattr(obj, 'mom_in', None)
-241 self.mom_out = getattr(obj, 'mom_out', None)
+203class Npr_matrix(np.ndarray):
+204
+205 def __new__(cls, input_array, mom_in=None, mom_out=None):
+206 obj = np.asarray(input_array).view(cls)
+207 obj.mom_in = mom_in
+208 obj.mom_out = mom_out
+209 return obj
+210
+211 @property
+212 def g5H(self):
+213 """Gamma_5 hermitean conjugate
+214
+215 Uses the fact that the propagator is gamma5 hermitean, so just the
+216 in and out momenta of the propagator are exchanged.
+217 """
+218 return Npr_matrix(self,
+219 mom_in=self.mom_out,
+220 mom_out=self.mom_in)
+221
+222 def _propagate_mom(self, other, name):
+223 s_mom = getattr(self, name, None)
+224 o_mom = getattr(other, name, None)
+225 if s_mom is not None and o_mom is not None:
+226 if not np.allclose(s_mom, o_mom):
+227 raise Exception(name + ' does not match.')
+228 return o_mom if o_mom is not None else s_mom
+229
+230 def __matmul__(self, other):
+231 return self.__new__(Npr_matrix,
+232 super().__matmul__(other),
+233 self._propagate_mom(other, 'mom_in'),
+234 self._propagate_mom(other, 'mom_out'))
+235
+236 def __array_finalize__(self, obj):
+237 if obj is None:
+238 return
+239 self.mom_in = getattr(obj, 'mom_in', None)
+240 self.mom_out = getattr(obj, 'mom_out', None)
+241
242
-243
-244def read_ExternalLeg_hd5(path, filestem, ens_id, idl=None):
-245 """Read hadrons ExternalLeg hdf5 file and output an array of CObs
-246
-247 Parameters
-248 ----------
-249 path : str
-250 path to the files to read
-251 filestem : str
-252 namestem of the files to read
-253 ens_id : str
-254 name of the ensemble, required for internal bookkeeping
-255 idl : range
-256 If specified only configurations in the given range are read in.
-257 """
-258
-259 files, idx = _get_files(path, filestem, idl)
-260
-261 mom = None
-262
-263 corr_data = []
-264 for hd5_file in files:
-265 file = h5py.File(path + '/' + hd5_file, "r")
-266 raw_data = file['ExternalLeg/corr'][0][0].view('complex')
-267 corr_data.append(raw_data)
-268 if mom is None:
-269 mom = np.array(str(file['ExternalLeg/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
-270 file.close()
-271 corr_data = np.array(corr_data)
-272
-273 rolled_array = np.rollaxis(corr_data, 0, 5)
-274
-275 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
-276 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]):
-277 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx])
-278 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx])
-279 matrix[si, sj, ci, cj] = CObs(real, imag)
-280
-281 return Npr_matrix(matrix, mom_in=mom)
+243def read_ExternalLeg_hd5(path, filestem, ens_id, idl=None):
+244 """Read hadrons ExternalLeg hdf5 file and output an array of CObs
+245
+246 Parameters
+247 ----------
+248 path : str
+249 path to the files to read
+250 filestem : str
+251 namestem of the files to read
+252 ens_id : str
+253 name of the ensemble, required for internal bookkeeping
+254 idl : range
+255 If specified only configurations in the given range are read in.
+256 """
+257
+258 files, idx = _get_files(path, filestem, idl)
+259
+260 mom = None
+261
+262 corr_data = []
+263 for hd5_file in files:
+264 file = h5py.File(path + '/' + hd5_file, "r")
+265 raw_data = file['ExternalLeg/corr'][0][0].view('complex')
+266 corr_data.append(raw_data)
+267 if mom is None:
+268 mom = np.array(str(file['ExternalLeg/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
+269 file.close()
+270 corr_data = np.array(corr_data)
+271
+272 rolled_array = np.rollaxis(corr_data, 0, 5)
+273
+274 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
+275 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]):
+276 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx])
+277 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx])
+278 matrix[si, sj, ci, cj] = CObs(real, imag)
+279
+280 return Npr_matrix(matrix, mom_in=mom)
+281
282
-283
-284def read_Bilinear_hd5(path, filestem, ens_id, idl=None):
-285 """Read hadrons Bilinear hdf5 file and output an array of CObs
-286
-287 Parameters
-288 ----------
-289 path : str
-290 path to the files to read
-291 filestem : str
-292 namestem of the files to read
-293 ens_id : str
-294 name of the ensemble, required for internal bookkeeping
-295 idl : range
-296 If specified only configurations in the given range are read in.
-297 """
-298
-299 files, idx = _get_files(path, filestem, idl)
-300
-301 mom_in = None
-302 mom_out = None
-303
-304 corr_data = {}
-305 for hd5_file in files:
-306 file = h5py.File(path + '/' + hd5_file, "r")
-307 for i in range(16):
-308 name = file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['gamma'][0].decode('UTF-8')
-309 if name not in corr_data:
-310 corr_data[name] = []
-311 raw_data = file['Bilinear/Bilinear_' + str(i) + '/corr'][0][0].view('complex')
-312 corr_data[name].append(raw_data)
-313 if mom_in is None:
-314 mom_in = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
-315 if mom_out is None:
-316 mom_out = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float)
-317
-318 file.close()
-319
-320 result_dict = {}
-321
-322 for key, data in corr_data.items():
-323 local_data = np.array(data)
-324
-325 rolled_array = np.rollaxis(local_data, 0, 5)
-326
-327 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
-328 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]):
-329 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx])
-330 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx])
-331 matrix[si, sj, ci, cj] = CObs(real, imag)
-332
-333 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out)
-334
-335 return result_dict
+283def read_Bilinear_hd5(path, filestem, ens_id, idl=None):
+284 """Read hadrons Bilinear hdf5 file and output an array of CObs
+285
+286 Parameters
+287 ----------
+288 path : str
+289 path to the files to read
+290 filestem : str
+291 namestem of the files to read
+292 ens_id : str
+293 name of the ensemble, required for internal bookkeeping
+294 idl : range
+295 If specified only configurations in the given range are read in.
+296 """
+297
+298 files, idx = _get_files(path, filestem, idl)
+299
+300 mom_in = None
+301 mom_out = None
+302
+303 corr_data = {}
+304 for hd5_file in files:
+305 file = h5py.File(path + '/' + hd5_file, "r")
+306 for i in range(16):
+307 name = file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['gamma'][0].decode('UTF-8')
+308 if name not in corr_data:
+309 corr_data[name] = []
+310 raw_data = file['Bilinear/Bilinear_' + str(i) + '/corr'][0][0].view('complex')
+311 corr_data[name].append(raw_data)
+312 if mom_in is None:
+313 mom_in = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
+314 if mom_out is None:
+315 mom_out = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float)
+316
+317 file.close()
+318
+319 result_dict = {}
+320
+321 for key, data in corr_data.items():
+322 local_data = np.array(data)
+323
+324 rolled_array = np.rollaxis(local_data, 0, 5)
+325
+326 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
+327 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]):
+328 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx])
+329 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx])
+330 matrix[si, sj, ci, cj] = CObs(real, imag)
+331
+332 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out)
+333
+334 return result_dict
+335
336
-337
-338def read_Fourquark_hd5(path, filestem, ens_id, idl=None, vertices=["VA", "AV"]):
-339 """Read hadrons FourquarkFullyConnected hdf5 file and output an array of CObs
-340
-341 Parameters
-342 ----------
-343 path : str
-344 path to the files to read
-345 filestem : str
-346 namestem of the files to read
-347 ens_id : str
-348 name of the ensemble, required for internal bookkeeping
-349 idl : range
-350 If specified only configurations in the given range are read in.
-351 vertices : list
-352 Vertex functions to be extracted.
-353 """
-354
-355 files, idx = _get_files(path, filestem, idl)
-356
-357 mom_in = None
-358 mom_out = None
-359
-360 vertex_names = []
-361 for vertex in vertices:
-362 vertex_names += _get_lorentz_names(vertex)
-363
-364 corr_data = {}
-365
-366 tree = 'FourQuarkFullyConnected/FourQuarkFullyConnected_'
-367
-368 for hd5_file in files:
-369 file = h5py.File(path + '/' + hd5_file, "r")
-370
-371 for i in range(32):
-372 name = (file[tree + str(i) + '/info'].attrs['gammaA'][0].decode('UTF-8'), file[tree + str(i) + '/info'].attrs['gammaB'][0].decode('UTF-8'))
-373 if name in vertex_names:
-374 if name not in corr_data:
-375 corr_data[name] = []
-376 raw_data = file[tree + str(i) + '/corr'][0][0].view('complex')
-377 corr_data[name].append(raw_data)
-378 if mom_in is None:
-379 mom_in = np.array(str(file[tree + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
-380 if mom_out is None:
-381 mom_out = np.array(str(file[tree + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float)
-382
-383 file.close()
-384
-385 intermediate_dict = {}
-386
-387 for vertex in vertices:
-388 lorentz_names = _get_lorentz_names(vertex)
-389 for v_name in lorentz_names:
-390 if v_name in [('SigmaXY', 'SigmaZT'),
-391 ('SigmaXT', 'SigmaYZ'),
-392 ('SigmaYZ', 'SigmaXT'),
-393 ('SigmaZT', 'SigmaXY')]:
-394 sign = -1
-395 else:
-396 sign = 1
-397 if vertex not in intermediate_dict:
-398 intermediate_dict[vertex] = sign * np.array(corr_data[v_name])
-399 else:
-400 intermediate_dict[vertex] += sign * np.array(corr_data[v_name])
-401
-402 result_dict = {}
-403
-404 for key, data in intermediate_dict.items():
-405
-406 rolled_array = np.moveaxis(data, 0, 8)
-407
-408 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
-409 for index in np.ndindex(rolled_array.shape[:-1]):
-410 real = Obs([rolled_array[index].real], [ens_id], idl=[idx])
-411 imag = Obs([rolled_array[index].imag], [ens_id], idl=[idx])
-412 matrix[index] = CObs(real, imag)
-413
-414 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out)
-415
-416 return result_dict
+337def read_Fourquark_hd5(path, filestem, ens_id, idl=None, vertices=["VA", "AV"]):
+338 """Read hadrons FourquarkFullyConnected hdf5 file and output an array of CObs
+339
+340 Parameters
+341 ----------
+342 path : str
+343 path to the files to read
+344 filestem : str
+345 namestem of the files to read
+346 ens_id : str
+347 name of the ensemble, required for internal bookkeeping
+348 idl : range
+349 If specified only configurations in the given range are read in.
+350 vertices : list
+351 Vertex functions to be extracted.
+352 """
+353
+354 files, idx = _get_files(path, filestem, idl)
+355
+356 mom_in = None
+357 mom_out = None
+358
+359 vertex_names = []
+360 for vertex in vertices:
+361 vertex_names += _get_lorentz_names(vertex)
+362
+363 corr_data = {}
+364
+365 tree = 'FourQuarkFullyConnected/FourQuarkFullyConnected_'
+366
+367 for hd5_file in files:
+368 file = h5py.File(path + '/' + hd5_file, "r")
+369
+370 for i in range(32):
+371 name = (file[tree + str(i) + '/info'].attrs['gammaA'][0].decode('UTF-8'), file[tree + str(i) + '/info'].attrs['gammaB'][0].decode('UTF-8'))
+372 if name in vertex_names:
+373 if name not in corr_data:
+374 corr_data[name] = []
+375 raw_data = file[tree + str(i) + '/corr'][0][0].view('complex')
+376 corr_data[name].append(raw_data)
+377 if mom_in is None:
+378 mom_in = np.array(str(file[tree + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
+379 if mom_out is None:
+380 mom_out = np.array(str(file[tree + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float)
+381
+382 file.close()
+383
+384 intermediate_dict = {}
+385
+386 for vertex in vertices:
+387 lorentz_names = _get_lorentz_names(vertex)
+388 for v_name in lorentz_names:
+389 if v_name in [('SigmaXY', 'SigmaZT'),
+390 ('SigmaXT', 'SigmaYZ'),
+391 ('SigmaYZ', 'SigmaXT'),
+392 ('SigmaZT', 'SigmaXY')]:
+393 sign = -1
+394 else:
+395 sign = 1
+396 if vertex not in intermediate_dict:
+397 intermediate_dict[vertex] = sign * np.array(corr_data[v_name])
+398 else:
+399 intermediate_dict[vertex] += sign * np.array(corr_data[v_name])
+400
+401 result_dict = {}
+402
+403 for key, data in intermediate_dict.items():
+404
+405 rolled_array = np.moveaxis(data, 0, 8)
+406
+407 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
+408 for index in np.ndindex(rolled_array.shape[:-1]):
+409 real = Obs([rolled_array[index].real], [ens_id], idl=[idx])
+410 imag = Obs([rolled_array[index].imag], [ens_id], idl=[idx])
+411 matrix[index] = CObs(real, imag)
+412
+413 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out)
+414
+415 return result_dict
+416
417
-418
-419def _get_lorentz_names(name):
-420 lorentz_index = ['X', 'Y', 'Z', 'T']
-421
-422 res = []
-423
-424 if name == "TT":
-425 for i in range(4):
-426 for j in range(i + 1, 4):
-427 res.append(("Sigma" + lorentz_index[i] + lorentz_index[j], "Sigma" + lorentz_index[i] + lorentz_index[j]))
-428 return res
-429
-430 if name == "TTtilde":
-431 for i in range(4):
-432 for j in range(i + 1, 4):
-433 for k in range(4):
-434 for o in range(k + 1, 4):
-435 fac = epsilon_tensor_rank4(i, j, k, o)
-436 if not np.isclose(fac, 0.0):
-437 res.append(("Sigma" + lorentz_index[i] + lorentz_index[j], "Sigma" + lorentz_index[k] + lorentz_index[o]))
-438 return res
-439
-440 assert len(name) == 2
-441
-442 if 'S' in name or 'P' in name:
-443 if not set(name) <= set(['S', 'P']):
-444 raise Exception("'" + name + "' is not a Lorentz scalar")
-445
-446 g_names = {'S': 'Identity',
-447 'P': 'Gamma5'}
-448
-449 res.append((g_names[name[0]], g_names[name[1]]))
-450
-451 else:
-452 if not set(name) <= set(['V', 'A']):
-453 raise Exception("'" + name + "' is not a Lorentz scalar")
-454
-455 for ind in lorentz_index:
-456 res.append(('Gamma' + ind + (name[0] == 'A') * 'Gamma5',
-457 'Gamma' + ind + (name[1] == 'A') * 'Gamma5'))
-458
-459 return res
+418def _get_lorentz_names(name):
+419 lorentz_index = ['X', 'Y', 'Z', 'T']
+420
+421 res = []
+422
+423 if name == "TT":
+424 for i in range(4):
+425 for j in range(i + 1, 4):
+426 res.append(("Sigma" + lorentz_index[i] + lorentz_index[j], "Sigma" + lorentz_index[i] + lorentz_index[j]))
+427 return res
+428
+429 if name == "TTtilde":
+430 for i in range(4):
+431 for j in range(i + 1, 4):
+432 for k in range(4):
+433 for o in range(k + 1, 4):
+434 fac = epsilon_tensor_rank4(i, j, k, o)
+435 if not np.isclose(fac, 0.0):
+436 res.append(("Sigma" + lorentz_index[i] + lorentz_index[j], "Sigma" + lorentz_index[k] + lorentz_index[o]))
+437 return res
+438
+439 assert len(name) == 2
+440
+441 if 'S' in name or 'P' in name:
+442 if not set(name) <= set(['S', 'P']):
+443 raise Exception("'" + name + "' is not a Lorentz scalar")
+444
+445 g_names = {'S': 'Identity',
+446 'P': 'Gamma5'}
+447
+448 res.append((g_names[name[0]], g_names[name[1]]))
+449
+450 else:
+451 if not set(name) <= set(['V', 'A']):
+452 raise Exception("'" + name + "' is not a Lorentz scalar")
+453
+454 for ind in lorentz_index:
+455 res.append(('Gamma' + ind + (name[0] == 'A') * 'Gamma5',
+456 'Gamma' + ind + (name[1] == 'A') * 'Gamma5'))
+457
+458 return res
@@ -574,67 +573,67 @@
- 59def read_meson_hd5(path, filestem, ens_id, meson='meson_0', idl=None, gammas=None):
- 60 r'''Read hadrons meson hdf5 file and extract the meson labeled 'meson'
- 61
- 62 Parameters
- 63 -----------------
- 64 path : str
- 65 path to the files to read
- 66 filestem : str
- 67 namestem of the files to read
- 68 ens_id : str
- 69 name of the ensemble, required for internal bookkeeping
- 70 meson : str
- 71 label of the meson to be extracted, standard value meson_0 which
- 72 corresponds to the pseudoscalar pseudoscalar two-point function.
- 73 gammas : tuple of strings
- 74 Instrad of a meson label one can also provide a tuple of two strings
- 75 indicating the gamma matrices at source and sink.
- 76 ("Gamma5", "Gamma5") corresponds to the pseudoscalar pseudoscalar
- 77 two-point function. The gammas argument dominateds over meson.
- 78 idl : range
- 79 If specified only configurations in the given range are read in.
- 80 '''
- 81
- 82 files, idx = _get_files(path, filestem, idl)
- 83
- 84 tree = meson.rsplit('_')[0]
- 85 if gammas is not None:
- 86 h5file = h5py.File(path + '/' + files[0], "r")
- 87 found_meson = None
- 88 for key in h5file[tree].keys():
- 89 if gammas[0] == h5file[tree][key].attrs["gamma_snk"][0].decode() and h5file[tree][key].attrs["gamma_src"][0].decode() == gammas[1]:
- 90 found_meson = key
- 91 break
- 92 h5file.close()
- 93 if found_meson:
- 94 meson = found_meson
- 95 else:
- 96 raise Exception("Source Sink combination " + str(gammas) + " not found.")
- 97
- 98 corr_data = []
- 99 infos = []
-100 for hd5_file in files:
-101 h5file = h5py.File(path + '/' + hd5_file, "r")
-102 if not tree + '/' + meson in h5file:
-103 raise Exception("Entry '" + meson + "' not contained in the files.")
-104 raw_data = h5file[tree + '/' + meson + '/corr']
-105 real_data = raw_data[:]["re"].astype(np.double)
-106 corr_data.append(real_data)
-107 if not infos:
-108 for k, i in h5file[tree + '/' + meson].attrs.items():
-109 infos.append(k + ': ' + i[0].decode())
-110 h5file.close()
-111 corr_data = np.array(corr_data)
-112
-113 l_obs = []
-114 for c in corr_data.T:
-115 l_obs.append(Obs([c], [ens_id], idl=[idx]))
-116
-117 corr = Corr(l_obs)
-118 corr.tag = r", ".join(infos)
-119 return corr
+ 58def read_meson_hd5(path, filestem, ens_id, meson='meson_0', idl=None, gammas=None):
+ 59 r'''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 gammas : tuple of strings
+ 73 Instrad of a meson label one can also provide a tuple of two strings
+ 74 indicating the gamma matrices at source and sink.
+ 75 ("Gamma5", "Gamma5") corresponds to the pseudoscalar pseudoscalar
+ 76 two-point function. The gammas argument dominateds over meson.
+ 77 idl : range
+ 78 If specified only configurations in the given range are read in.
+ 79 '''
+ 80
+ 81 files, idx = _get_files(path, filestem, idl)
+ 82
+ 83 tree = meson.rsplit('_')[0]
+ 84 if gammas is not None:
+ 85 h5file = h5py.File(path + '/' + files[0], "r")
+ 86 found_meson = None
+ 87 for key in h5file[tree].keys():
+ 88 if gammas[0] == h5file[tree][key].attrs["gamma_snk"][0].decode() and h5file[tree][key].attrs["gamma_src"][0].decode() == gammas[1]:
+ 89 found_meson = key
+ 90 break
+ 91 h5file.close()
+ 92 if found_meson:
+ 93 meson = found_meson
+ 94 else:
+ 95 raise Exception("Source Sink combination " + str(gammas) + " not found.")
+ 96
+ 97 corr_data = []
+ 98 infos = []
+ 99 for hd5_file in files:
+100 h5file = h5py.File(path + '/' + hd5_file, "r")
+101 if not tree + '/' + meson in h5file:
+102 raise Exception("Entry '" + meson + "' not contained in the files.")
+103 raw_data = h5file[tree + '/' + meson + '/corr']
+104 real_data = raw_data[:]["re"].astype(np.double)
+105 corr_data.append(real_data)
+106 if not infos:
+107 for k, i in h5file[tree + '/' + meson].attrs.items():
+108 infos.append(k + ': ' + i[0].decode())
+109 h5file.close()
+110 corr_data = np.array(corr_data)
+111
+112 l_obs = []
+113 for c in corr_data.T:
+114 l_obs.append(Obs([c], [ens_id], idl=[idx]))
+115
+116 corr = Corr(l_obs)
+117 corr.tag = r", ".join(infos)
+118 return corr
122def read_DistillationContraction_hd5(path, ens_id, diagrams=["direct"], idl=None):
-123 """Read hadrons DistillationContraction hdf5 files in given directory structure
-124
-125 Parameters
-126 -----------------
-127 path : str
-128 path to the directories to read
-129 ens_id : str
-130 name of the ensemble, required for internal bookkeeping
-131 diagrams : list
-132 List of strings of the diagrams to extract, e.g. ["direct", "box", "cross"].
-133 idl : range
-134 If specified only configurations in the given range are read in.
-135 """
-136
-137 res_dict = {}
-138
-139 directories, idx = _get_files(path, "data", idl)
-140
-141 explore_path = Path(path + "/" + directories[0])
-142
-143 for explore_file in explore_path.iterdir():
-144 if explore_file.is_file():
-145 stem = explore_file.with_suffix("").with_suffix("").as_posix().split("/")[-1]
-146 else:
-147 continue
-148
-149 file_list = []
-150 for dir in directories:
-151 tmp_path = Path(path + "/" + dir)
-152 file_list.append((tmp_path / stem).as_posix() + tmp_path.suffix + ".h5")
-153
-154 corr_data = {}
-155
-156 for diagram in diagrams:
-157 corr_data[diagram] = []
-158
-159 for n_file, (hd5_file, n_traj) in enumerate(zip(file_list, list(idx))):
-160 h5file = h5py.File(hd5_file)
-161
-162 if n_file == 0:
-163 if h5file["DistillationContraction/Metadata"].attrs.get("TimeSources")[0].decode() != "0...":
-164 raise Exception("Routine is only implemented for files containing inversions on all timeslices.")
-165
-166 Nt = h5file["DistillationContraction/Metadata"].attrs.get("Nt")[0]
-167
-168 identifier = []
-169 for in_file in range(len(h5file["DistillationContraction/Metadata/DmfInputFiles"].attrs.keys()) - 1):
-170 encoded_info = h5file["DistillationContraction/Metadata/DmfInputFiles"].attrs.get("DmfInputFiles_" + str(in_file))
-171 full_info = encoded_info[0].decode().split("/")[-1].replace(".h5", "").split("_")
-172 my_tuple = (full_info[0], full_info[1][1:], full_info[2], full_info[3])
-173 identifier.append(my_tuple)
-174 identifier = tuple(identifier)
-175 # "DistillationContraction/Metadata/DmfSuffix" contains info about different quarks, irrelevant in the SU(3) case.
-176
-177 for diagram in diagrams:
-178 real_data = np.zeros(Nt)
-179 for x0 in range(Nt):
-180 raw_data = h5file["DistillationContraction/Correlators/" + diagram + "/" + str(x0)][:]["re"].astype(np.double)
-181 real_data += np.roll(raw_data, -x0)
-182 real_data /= Nt
-183
-184 corr_data[diagram].append(real_data)
-185 h5file.close()
-186
-187 res_dict[str(identifier)] = {}
-188
-189 for diagram in diagrams:
-190
-191 tmp_data = np.array(corr_data[diagram])
-192
-193 l_obs = []
-194 for c in tmp_data.T:
-195 l_obs.append(Obs([c], [ens_id], idl=[idx]))
-196
-197 corr = Corr(l_obs)
-198 corr.tag = str(identifier)
-199
-200 res_dict[str(identifier)][diagram] = corr
-201
-202 return res_dict
+ 121def read_DistillationContraction_hd5(path, ens_id, diagrams=["direct"], idl=None):
+122 """Read hadrons DistillationContraction hdf5 files in given directory structure
+123
+124 Parameters
+125 -----------------
+126 path : str
+127 path to the directories to read
+128 ens_id : str
+129 name of the ensemble, required for internal bookkeeping
+130 diagrams : list
+131 List of strings of the diagrams to extract, e.g. ["direct", "box", "cross"].
+132 idl : range
+133 If specified only configurations in the given range are read in.
+134 """
+135
+136 res_dict = {}
+137
+138 directories, idx = _get_files(path, "data", idl)
+139
+140 explore_path = Path(path + "/" + directories[0])
+141
+142 for explore_file in explore_path.iterdir():
+143 if explore_file.is_file():
+144 stem = explore_file.with_suffix("").with_suffix("").as_posix().split("/")[-1]
+145 else:
+146 continue
+147
+148 file_list = []
+149 for dir in directories:
+150 tmp_path = Path(path + "/" + dir)
+151 file_list.append((tmp_path / stem).as_posix() + tmp_path.suffix + ".h5")
+152
+153 corr_data = {}
+154
+155 for diagram in diagrams:
+156 corr_data[diagram] = []
+157
+158 for n_file, (hd5_file, n_traj) in enumerate(zip(file_list, list(idx))):
+159 h5file = h5py.File(hd5_file)
+160
+161 if n_file == 0:
+162 if h5file["DistillationContraction/Metadata"].attrs.get("TimeSources")[0].decode() != "0...":
+163 raise Exception("Routine is only implemented for files containing inversions on all timeslices.")
+164
+165 Nt = h5file["DistillationContraction/Metadata"].attrs.get("Nt")[0]
+166
+167 identifier = []
+168 for in_file in range(len(h5file["DistillationContraction/Metadata/DmfInputFiles"].attrs.keys()) - 1):
+169 encoded_info = h5file["DistillationContraction/Metadata/DmfInputFiles"].attrs.get("DmfInputFiles_" + str(in_file))
+170 full_info = encoded_info[0].decode().split("/")[-1].replace(".h5", "").split("_")
+171 my_tuple = (full_info[0], full_info[1][1:], full_info[2], full_info[3])
+172 identifier.append(my_tuple)
+173 identifier = tuple(identifier)
+174 # "DistillationContraction/Metadata/DmfSuffix" contains info about different quarks, irrelevant in the SU(3) case.
+175
+176 for diagram in diagrams:
+177 real_data = np.zeros(Nt)
+178 for x0 in range(Nt):
+179 raw_data = h5file["DistillationContraction/Correlators/" + diagram + "/" + str(x0)][:]["re"].astype(np.double)
+180 real_data += np.roll(raw_data, -x0)
+181 real_data /= Nt
+182
+183 corr_data[diagram].append(real_data)
+184 h5file.close()
+185
+186 res_dict[str(identifier)] = {}
+187
+188 for diagram in diagrams:
+189
+190 tmp_data = np.array(corr_data[diagram])
+191
+192 l_obs = []
+193 for c in tmp_data.T:
+194 l_obs.append(Obs([c], [ens_id], idl=[idx]))
+195
+196 corr = Corr(l_obs)
+197 corr.tag = str(identifier)
+198
+199 res_dict[str(identifier)][diagram] = corr
+200
+201 return res_dict
205class Npr_matrix(np.ndarray):
-206
-207 def __new__(cls, input_array, mom_in=None, mom_out=None):
-208 obj = np.asarray(input_array).view(cls)
-209 obj.mom_in = mom_in
-210 obj.mom_out = mom_out
-211 return obj
-212
-213 @property
-214 def g5H(self):
-215 """Gamma_5 hermitean conjugate
-216
-217 Uses the fact that the propagator is gamma5 hermitean, so just the
-218 in and out momenta of the propagator are exchanged.
-219 """
-220 return Npr_matrix(self,
-221 mom_in=self.mom_out,
-222 mom_out=self.mom_in)
-223
-224 def _propagate_mom(self, other, name):
-225 s_mom = getattr(self, name, None)
-226 o_mom = getattr(other, name, None)
-227 if s_mom is not None and o_mom is not None:
-228 if not np.allclose(s_mom, o_mom):
-229 raise Exception(name + ' does not match.')
-230 return o_mom if o_mom is not None else s_mom
-231
-232 def __matmul__(self, other):
-233 return self.__new__(Npr_matrix,
-234 super().__matmul__(other),
-235 self._propagate_mom(other, 'mom_in'),
-236 self._propagate_mom(other, 'mom_out'))
-237
-238 def __array_finalize__(self, obj):
-239 if obj is None:
-240 return
-241 self.mom_in = getattr(obj, 'mom_in', None)
-242 self.mom_out = getattr(obj, 'mom_out', None)
+ 204class Npr_matrix(np.ndarray):
+205
+206 def __new__(cls, input_array, mom_in=None, mom_out=None):
+207 obj = np.asarray(input_array).view(cls)
+208 obj.mom_in = mom_in
+209 obj.mom_out = mom_out
+210 return obj
+211
+212 @property
+213 def g5H(self):
+214 """Gamma_5 hermitean conjugate
+215
+216 Uses the fact that the propagator is gamma5 hermitean, so just the
+217 in and out momenta of the propagator are exchanged.
+218 """
+219 return Npr_matrix(self,
+220 mom_in=self.mom_out,
+221 mom_out=self.mom_in)
+222
+223 def _propagate_mom(self, other, name):
+224 s_mom = getattr(self, name, None)
+225 o_mom = getattr(other, name, None)
+226 if s_mom is not None and o_mom is not None:
+227 if not np.allclose(s_mom, o_mom):
+228 raise Exception(name + ' does not match.')
+229 return o_mom if o_mom is not None else s_mom
+230
+231 def __matmul__(self, other):
+232 return self.__new__(Npr_matrix,
+233 super().__matmul__(other),
+234 self._propagate_mom(other, 'mom_in'),
+235 self._propagate_mom(other, 'mom_out'))
+236
+237 def __array_finalize__(self, obj):
+238 if obj is None:
+239 return
+240 self.mom_in = getattr(obj, 'mom_in', None)
+241 self.mom_out = getattr(obj, 'mom_out', None)
245def read_ExternalLeg_hd5(path, filestem, ens_id, idl=None):
-246 """Read hadrons ExternalLeg hdf5 file and output an array of CObs
-247
-248 Parameters
-249 ----------
-250 path : str
-251 path to the files to read
-252 filestem : str
-253 namestem of the files to read
-254 ens_id : str
-255 name of the ensemble, required for internal bookkeeping
-256 idl : range
-257 If specified only configurations in the given range are read in.
-258 """
-259
-260 files, idx = _get_files(path, filestem, idl)
-261
-262 mom = None
-263
-264 corr_data = []
-265 for hd5_file in files:
-266 file = h5py.File(path + '/' + hd5_file, "r")
-267 raw_data = file['ExternalLeg/corr'][0][0].view('complex')
-268 corr_data.append(raw_data)
-269 if mom is None:
-270 mom = np.array(str(file['ExternalLeg/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
-271 file.close()
-272 corr_data = np.array(corr_data)
-273
-274 rolled_array = np.rollaxis(corr_data, 0, 5)
-275
-276 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
-277 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]):
-278 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx])
-279 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx])
-280 matrix[si, sj, ci, cj] = CObs(real, imag)
-281
-282 return Npr_matrix(matrix, mom_in=mom)
+ 244def read_ExternalLeg_hd5(path, filestem, ens_id, idl=None):
+245 """Read hadrons ExternalLeg hdf5 file and output an array of CObs
+246
+247 Parameters
+248 ----------
+249 path : str
+250 path to the files to read
+251 filestem : str
+252 namestem of the files to read
+253 ens_id : str
+254 name of the ensemble, required for internal bookkeeping
+255 idl : range
+256 If specified only configurations in the given range are read in.
+257 """
+258
+259 files, idx = _get_files(path, filestem, idl)
+260
+261 mom = None
+262
+263 corr_data = []
+264 for hd5_file in files:
+265 file = h5py.File(path + '/' + hd5_file, "r")
+266 raw_data = file['ExternalLeg/corr'][0][0].view('complex')
+267 corr_data.append(raw_data)
+268 if mom is None:
+269 mom = np.array(str(file['ExternalLeg/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
+270 file.close()
+271 corr_data = np.array(corr_data)
+272
+273 rolled_array = np.rollaxis(corr_data, 0, 5)
+274
+275 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
+276 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]):
+277 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx])
+278 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx])
+279 matrix[si, sj, ci, cj] = CObs(real, imag)
+280
+281 return Npr_matrix(matrix, mom_in=mom)
285def read_Bilinear_hd5(path, filestem, ens_id, idl=None):
-286 """Read hadrons Bilinear hdf5 file and output an array of CObs
-287
-288 Parameters
-289 ----------
-290 path : str
-291 path to the files to read
-292 filestem : str
-293 namestem of the files to read
-294 ens_id : str
-295 name of the ensemble, required for internal bookkeeping
-296 idl : range
-297 If specified only configurations in the given range are read in.
-298 """
-299
-300 files, idx = _get_files(path, filestem, idl)
-301
-302 mom_in = None
-303 mom_out = None
-304
-305 corr_data = {}
-306 for hd5_file in files:
-307 file = h5py.File(path + '/' + hd5_file, "r")
-308 for i in range(16):
-309 name = file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['gamma'][0].decode('UTF-8')
-310 if name not in corr_data:
-311 corr_data[name] = []
-312 raw_data = file['Bilinear/Bilinear_' + str(i) + '/corr'][0][0].view('complex')
-313 corr_data[name].append(raw_data)
-314 if mom_in is None:
-315 mom_in = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
-316 if mom_out is None:
-317 mom_out = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float)
-318
-319 file.close()
-320
-321 result_dict = {}
-322
-323 for key, data in corr_data.items():
-324 local_data = np.array(data)
-325
-326 rolled_array = np.rollaxis(local_data, 0, 5)
-327
-328 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
-329 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]):
-330 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx])
-331 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx])
-332 matrix[si, sj, ci, cj] = CObs(real, imag)
-333
-334 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out)
-335
-336 return result_dict
+ 284def read_Bilinear_hd5(path, filestem, ens_id, idl=None):
+285 """Read hadrons Bilinear hdf5 file and output an array of CObs
+286
+287 Parameters
+288 ----------
+289 path : str
+290 path to the files to read
+291 filestem : str
+292 namestem of the files to read
+293 ens_id : str
+294 name of the ensemble, required for internal bookkeeping
+295 idl : range
+296 If specified only configurations in the given range are read in.
+297 """
+298
+299 files, idx = _get_files(path, filestem, idl)
+300
+301 mom_in = None
+302 mom_out = None
+303
+304 corr_data = {}
+305 for hd5_file in files:
+306 file = h5py.File(path + '/' + hd5_file, "r")
+307 for i in range(16):
+308 name = file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['gamma'][0].decode('UTF-8')
+309 if name not in corr_data:
+310 corr_data[name] = []
+311 raw_data = file['Bilinear/Bilinear_' + str(i) + '/corr'][0][0].view('complex')
+312 corr_data[name].append(raw_data)
+313 if mom_in is None:
+314 mom_in = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
+315 if mom_out is None:
+316 mom_out = np.array(str(file['Bilinear/Bilinear_' + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float)
+317
+318 file.close()
+319
+320 result_dict = {}
+321
+322 for key, data in corr_data.items():
+323 local_data = np.array(data)
+324
+325 rolled_array = np.rollaxis(local_data, 0, 5)
+326
+327 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
+328 for si, sj, ci, cj in np.ndindex(rolled_array.shape[:-1]):
+329 real = Obs([rolled_array[si, sj, ci, cj].real], [ens_id], idl=[idx])
+330 imag = Obs([rolled_array[si, sj, ci, cj].imag], [ens_id], idl=[idx])
+331 matrix[si, sj, ci, cj] = CObs(real, imag)
+332
+333 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out)
+334
+335 return result_dict
339def read_Fourquark_hd5(path, filestem, ens_id, idl=None, vertices=["VA", "AV"]):
-340 """Read hadrons FourquarkFullyConnected hdf5 file and output an array of CObs
-341
-342 Parameters
-343 ----------
-344 path : str
-345 path to the files to read
-346 filestem : str
-347 namestem of the files to read
-348 ens_id : str
-349 name of the ensemble, required for internal bookkeeping
-350 idl : range
-351 If specified only configurations in the given range are read in.
-352 vertices : list
-353 Vertex functions to be extracted.
-354 """
-355
-356 files, idx = _get_files(path, filestem, idl)
-357
-358 mom_in = None
-359 mom_out = None
-360
-361 vertex_names = []
-362 for vertex in vertices:
-363 vertex_names += _get_lorentz_names(vertex)
-364
-365 corr_data = {}
-366
-367 tree = 'FourQuarkFullyConnected/FourQuarkFullyConnected_'
-368
-369 for hd5_file in files:
-370 file = h5py.File(path + '/' + hd5_file, "r")
-371
-372 for i in range(32):
-373 name = (file[tree + str(i) + '/info'].attrs['gammaA'][0].decode('UTF-8'), file[tree + str(i) + '/info'].attrs['gammaB'][0].decode('UTF-8'))
-374 if name in vertex_names:
-375 if name not in corr_data:
-376 corr_data[name] = []
-377 raw_data = file[tree + str(i) + '/corr'][0][0].view('complex')
-378 corr_data[name].append(raw_data)
-379 if mom_in is None:
-380 mom_in = np.array(str(file[tree + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
-381 if mom_out is None:
-382 mom_out = np.array(str(file[tree + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float)
-383
-384 file.close()
-385
-386 intermediate_dict = {}
-387
-388 for vertex in vertices:
-389 lorentz_names = _get_lorentz_names(vertex)
-390 for v_name in lorentz_names:
-391 if v_name in [('SigmaXY', 'SigmaZT'),
-392 ('SigmaXT', 'SigmaYZ'),
-393 ('SigmaYZ', 'SigmaXT'),
-394 ('SigmaZT', 'SigmaXY')]:
-395 sign = -1
-396 else:
-397 sign = 1
-398 if vertex not in intermediate_dict:
-399 intermediate_dict[vertex] = sign * np.array(corr_data[v_name])
-400 else:
-401 intermediate_dict[vertex] += sign * np.array(corr_data[v_name])
-402
-403 result_dict = {}
-404
-405 for key, data in intermediate_dict.items():
-406
-407 rolled_array = np.moveaxis(data, 0, 8)
-408
-409 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
-410 for index in np.ndindex(rolled_array.shape[:-1]):
-411 real = Obs([rolled_array[index].real], [ens_id], idl=[idx])
-412 imag = Obs([rolled_array[index].imag], [ens_id], idl=[idx])
-413 matrix[index] = CObs(real, imag)
-414
-415 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out)
-416
-417 return result_dict
+ 338def read_Fourquark_hd5(path, filestem, ens_id, idl=None, vertices=["VA", "AV"]):
+339 """Read hadrons FourquarkFullyConnected hdf5 file and output an array of CObs
+340
+341 Parameters
+342 ----------
+343 path : str
+344 path to the files to read
+345 filestem : str
+346 namestem of the files to read
+347 ens_id : str
+348 name of the ensemble, required for internal bookkeeping
+349 idl : range
+350 If specified only configurations in the given range are read in.
+351 vertices : list
+352 Vertex functions to be extracted.
+353 """
+354
+355 files, idx = _get_files(path, filestem, idl)
+356
+357 mom_in = None
+358 mom_out = None
+359
+360 vertex_names = []
+361 for vertex in vertices:
+362 vertex_names += _get_lorentz_names(vertex)
+363
+364 corr_data = {}
+365
+366 tree = 'FourQuarkFullyConnected/FourQuarkFullyConnected_'
+367
+368 for hd5_file in files:
+369 file = h5py.File(path + '/' + hd5_file, "r")
+370
+371 for i in range(32):
+372 name = (file[tree + str(i) + '/info'].attrs['gammaA'][0].decode('UTF-8'), file[tree + str(i) + '/info'].attrs['gammaB'][0].decode('UTF-8'))
+373 if name in vertex_names:
+374 if name not in corr_data:
+375 corr_data[name] = []
+376 raw_data = file[tree + str(i) + '/corr'][0][0].view('complex')
+377 corr_data[name].append(raw_data)
+378 if mom_in is None:
+379 mom_in = np.array(str(file[tree + str(i) + '/info'].attrs['pIn'])[3:-2].strip().split(), dtype=float)
+380 if mom_out is None:
+381 mom_out = np.array(str(file[tree + str(i) + '/info'].attrs['pOut'])[3:-2].strip().split(), dtype=float)
+382
+383 file.close()
+384
+385 intermediate_dict = {}
+386
+387 for vertex in vertices:
+388 lorentz_names = _get_lorentz_names(vertex)
+389 for v_name in lorentz_names:
+390 if v_name in [('SigmaXY', 'SigmaZT'),
+391 ('SigmaXT', 'SigmaYZ'),
+392 ('SigmaYZ', 'SigmaXT'),
+393 ('SigmaZT', 'SigmaXY')]:
+394 sign = -1
+395 else:
+396 sign = 1
+397 if vertex not in intermediate_dict:
+398 intermediate_dict[vertex] = sign * np.array(corr_data[v_name])
+399 else:
+400 intermediate_dict[vertex] += sign * np.array(corr_data[v_name])
+401
+402 result_dict = {}
+403
+404 for key, data in intermediate_dict.items():
+405
+406 rolled_array = np.moveaxis(data, 0, 8)
+407
+408 matrix = np.empty((rolled_array.shape[:-1]), dtype=object)
+409 for index in np.ndindex(rolled_array.shape[:-1]):
+410 real = Obs([rolled_array[index].real], [ens_id], idl=[idx])
+411 imag = Obs([rolled_array[index].imag], [ens_id], idl=[idx])
+412 matrix[index] = CObs(real, imag)
+413
+414 result_dict[key] = Npr_matrix(matrix, mom_in=mom_in, mom_out=mom_out)
+415
+416 return result_dict