pyerrors.input.hadrons
1import os 2import warnings 3from collections import Counter 4import h5py 5from pathlib import Path 6import numpy as np 7from ..obs import Obs, CObs 8from ..correlators import Corr 9from ..dirac import epsilon_tensor_rank4 10 11 12def _get_files(path, filestem, idl): 13 ls = os.listdir(path) 14 15 # Clean up file list 16 files = list(filter(lambda x: x.startswith(filestem + "."), ls)) 17 18 if not files: 19 raise Exception('No files starting with', filestem, 'in folder', path) 20 21 def get_cnfg_number(n): 22 return int(n.replace(".h5", "")[len(filestem) + 1:]) # From python 3.9 onward the safer 'removesuffix' method can be used. 23 24 # Sort according to configuration number 25 files.sort(key=get_cnfg_number) 26 27 cnfg_numbers = [] 28 filtered_files = [] 29 for line in files: 30 no = get_cnfg_number(line) 31 if idl: 32 if no in list(idl): 33 filtered_files.append(line) 34 cnfg_numbers.append(no) 35 else: 36 filtered_files.append(line) 37 cnfg_numbers.append(no) 38 39 if idl: 40 if Counter(list(idl)) != Counter(cnfg_numbers): 41 raise Exception("Not all configurations specified in idl found, configurations " + str(list(Counter(list(idl)) - Counter(cnfg_numbers))) + " are missing.") 42 43 # Check that configurations are evenly spaced 44 dc = np.unique(np.diff(cnfg_numbers)) 45 if np.any(dc < 0): 46 raise Exception("Unsorted files") 47 if len(dc) == 1: 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 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 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 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) 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) 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 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 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
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
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.
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
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.
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)
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
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)
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.
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
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.
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
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.