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