From 6c41810e8182e5d62e5037871d379d5587abf0cf Mon Sep 17 00:00:00 2001 From: Justus Date: Mon, 15 Nov 2021 15:55:26 +0100 Subject: [PATCH] added Qtop extraction for oQCD1.2 --- pyerrors/input/openQCD.py | 128 ++++++++++++++++++++++++++++++++++++++ pyerrors/input/sfcf.py | 63 ------------------- 2 files changed, 128 insertions(+), 63 deletions(-) diff --git a/pyerrors/input/openQCD.py b/pyerrors/input/openQCD.py index 5e1c8d49..2483baa9 100644 --- a/pyerrors/input/openQCD.py +++ b/pyerrors/input/openQCD.py @@ -343,3 +343,131 @@ def _read_array_openQCD2(fp): arr = _parse_array_openQCD2(d, n, size, tmp, quadrupel=True) return {'d': d, 'n': n, 'size': size, 'arr': arr} + + +def read_qtop(path, prefix, version = "1.2",**kwargs): + """Read qtop format from given folder structure. + + Parameters + ---------- + target -- specifies the topological sector to be reweighted to (default 0) + full -- if true read the charge instead of the reweighting factor. + """ + dtr_cnfg = 4 + L = 20 + c = 0.35 + target = 0 + full = False + + if 'target' in kwargs: + target = kwargs.get('target') + + + if kwargs.get('full'): + full = True + + if "r_start" in kwargs: + r_start = kwargs.get("r_start") + if "r_stop" in kwargs: + r_stop = kwargs.get("r_stop") + #if one wants to read specific files with this method... + if "files" in kwargs: + files = kwargs.get("files") + else: + #find files in path + found = [] + files = [] + for (dirpath, dirnames, filenames) in os.walk(path+"/"): + #print(filenames) + found.extend(filenames) + break + for f in found: + if fnmatch.fnmatch(f, prefix+"*"+".ms.dat"): + files.append(f) + print(files) + #now that we found our files, we dechiffer them... + rep_names = [] + + deltas = [] + for rep,file in enumerate(files): + + with open(path+"/"+file, "rb") as fp: + #this, for now, is only for version 1.2 + #header + t = fp.read(12) + header = struct.unpack('iii', t) + dn = header[0] + nn = header[1] + tmax = header[2] + print('dn:', dn) + print('nn:', nn) + print('tmax:', tmax) + t = fp.read(8) + eps = struct.unpack('d', t)[0] + print('eps:', eps) + + Q = [] + i = 1 + while 0 < 1: + t = fp.read(4) + if(len(t) < 4): + break + nc = struct.unpack('i',t)[0] + if(nc != i): + print("WARNING: possible missing config:" +str(i)) + #raise Exception('Config missing?') + else: + t = fp.read(8 * tmax * (nn + 1)) + t = fp.read(8 * tmax * (nn + 1)) + t = fp.read(8 * tmax * (nn + 1)) + tmpd = struct.unpack('d' * tmax * (nn + 1), t) + Q.append(tmpd) + i += 1 + #print(tmp) + + print('max_t:', dn * (nn) * eps) + + t_aim = (c * L) ** 2 / 8 + + print('t_aim:', t_aim) + index_aim = round(t_aim / eps / dn) + print('index_aim:', index_aim) + + + Q_sum = [] + for i, item in enumerate(Q): + Q_sum.append([sum(item[current:current + tmax]) for current in range(0, len(item), tmax)]) + Q_round = [] + for i in range(len(Q) // dtr_cnfg): + Q_round.append(round(Q_sum[dtr_cnfg * i][index_aim])) + + replica = len(files) + + tmp = [] + for q in Q_round: + #floats = list(map(float, line.split())) + if full: + tmp.append(q) #round(Q_sum[dtr_cnfg * i][index_aim]) + else: + if int(q) == target: #round(Q_sum[dtr_cnfg * i][index_aim]) + tmp.append(1.0) + else: + tmp.append(0.0) + + truncated_file = file[:-7] #as seen in previous examples, this could lead to some weird behaviour... maybe -7 fixes this. + print(truncated_file) + idx = truncated_file.index('r') + #print(truncated_file[idx:]) + # this might be a quite fishy way to find out which replicum we are actually talking about... + if "r_start" in kwargs: + tmp = tmp[r_start[int(truncated_file[idx+1:])-1]:] + if "r_stop" in kwargs: + tmp = tmp[:r_stop[int(truncated_file[idx+1:])-1]] + + rep_names.append(truncated_file[:idx] + '|' + truncated_file[idx:]) + + deltas.append(np.array(tmp)) + + + result = Obs(deltas, rep_names) + return result \ No newline at end of file diff --git a/pyerrors/input/sfcf.py b/pyerrors/input/sfcf.py index e48bdd16..706e26a9 100644 --- a/pyerrors/input/sfcf.py +++ b/pyerrors/input/sfcf.py @@ -229,66 +229,3 @@ def read_sfcf_c(path, prefix, name, quarks='.*', noffset=0, wf=0, wf2=0, **kwarg result.append(Obs(deltas[t], new_names)) return result - -def read_qtop(path, prefix, **kwargs): - """Read qtop format from given folder structure. - - Parameters - ---------- - target -- specifies the topological sector to be reweighted to (default 0) - full -- if true read the charge instead of the reweighting factor. - """ - - if 'target' in kwargs: - target = kwargs.get('target') - else: - target = 0 - - if kwargs.get('full'): - full = 1 - else: - full = 0 - - ls = [] - for (dirpath, dirnames, filenames) in os.walk(path): - ls.extend(filenames) - break - - if not ls: - raise Exception('Error, directory not found') - - # Exclude files with different names - for exc in ls: - if not fnmatch.fnmatch(exc, prefix + '*'): - ls = list(set(ls) - set([exc])) - if len(ls) > 1: - ls.sort(key=lambda x: int(re.findall(r'\d+', x[len(prefix):])[0])) # New version, to cope with ids, etc. - replica = len(ls) - print('Read Q_top from', prefix[:-1], ',', replica, 'replica') - - deltas = [] - - for rep in range(replica): - tmp = [] - with open(path + '/' + ls[rep]) as fp: - for k, line in enumerate(fp): - floats = list(map(float, line.split())) - if full == 1: - tmp.append(floats[1]) - else: - if int(floats[1]) == target: - tmp.append(1.0) - else: - tmp.append(0.0) - - deltas.append(np.array(tmp)) - - rep_names = [] - for entry in ls: - truncated_entry = entry.split('.')[0] - idx = truncated_entry.index('r') - rep_names.append(truncated_entry[:idx] + '|' + truncated_entry[idx:]) - - result = Obs(deltas, rep_names) - - return result