noe/pipeline
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

Time average, other improvments

Browse Source
This commit is contained in:
Noe Brucy
2019-05-20 10:56:31 +02:00
parent 2f9bdf4bf9
commit aa8ba76162
2 changed files with 99 additions and 59 deletions
Download Patch File Download Diff File Expand all files Collapse all files
disk_postprocess.py
+77 -47
View File
@@ -152,8 +152,8 @@ def make_image_aux(
directory = path directory = path
# Checking for existing files # Checking for existing files
name = directory + "/coldens_z" + "_" + tag + "_" + format(num, "05") + out_ext name_save = directory + "/maps_disk" + "_" + tag + "_" + format(num, "05") + ".save"
if len(glob.glob(name)) == 1 and not force: if len(glob.glob(name_save)) == 1 and not force:
return return
# Prepare saving data # Prepare saving data
@@ -453,9 +453,6 @@ def make_image_aux(
P.close() P.close()
if save_data: if save_data:
name_save = (
directory + "/maps_disk" + "_" + tag + "_" + format(num, "05") + ".save"
)
f = open(name_save, "w") f = open(name_save, "w")
pickle.dump(maps_disk, f) pickle.dump(maps_disk, f)
f.close() f.close()
@@ -961,7 +958,7 @@ def disk_pdf(
P.yscale("log") P.yscale("log")
P.ylim([0.5 / nb_cells, 1.0]) P.ylim([0.5 / nb_cells, 1.0])
P.xlabel(r"$\log(N / \bar{N})$") P.xlabel(r"$\log(N / \bar{N})$")
P.ylabel("Probability density") P.ylabel(r"\mathcal{P}_N")
if put_title: if put_title:
P.title(title) P.title(title)
values, edges, _ = P.hist( values, edges, _ = P.hist(
@@ -974,19 +971,20 @@ def disk_pdf(
# Compute the slope of the right part of the histogramm # Compute the slope of the right part of the histogramm
mask_fit = (centers > 0.25) & (centers < 1.25) & (values > 0) mask_fit = (centers > 0.25) & (centers < 1.25) & (values > 0)
(a, b, rho, _, _) = linregress(centers[mask_fit], np.log10(values[mask_fit])) if np.sum(mask_fit > 0):
print("a=%e, b=%e, rho=%e" % (a, b, rho)) (a, b, rho, _, _) = linregress(centers[mask_fit], np.log10(values[mask_fit]))
fit = { P.plot(centers, 10 ** (a * centers + b), "--")
"beta": int(tag.split("_")[1][4:]), print("a=%e, b=%e, rho=%e" % (a, b, rho))
"slope": a, fit = {
"origin": b, "beta": int(tag.split("_")[1][4:]),
"correlation": rho, "slope": a,
} "origin": b,
P.plot(centers, 10 ** (a * centers + b), "--") "correlation": rho,
f = open(name_prop, "w") }
prop_disk["fit"] = fit f = open(name_prop, "w")
pickle.dump(prop_disk, f) prop_disk["fit"] = fit
f.close() pickle.dump(prop_disk, f)
f.close()
if interactive: if interactive:
pass pass
@@ -1020,6 +1018,7 @@ def compare(
path, path,
runs, runs,
nums, nums,
path_out=None,
force=False, force=False,
interactive=False, interactive=False,
Q_in_name=True, Q_in_name=True,
@@ -1039,7 +1038,7 @@ def compare(
if type(nums) == int: if type(nums) == int:
nums = [nums] nums = [nums]
nums_name = "_".join(str(num).zfill(5) for num in nums) nums_name = "_".join(str(num).zfill(5) for num in [nums[0], nums[-1]])
# Initialize arrays # Initialize arrays
alpha_rey = np.zeros(len(runs)) alpha_rey = np.zeros(len(runs))
@@ -1053,25 +1052,49 @@ def compare(
for i, run in enumerate(runs): for i, run in enumerate(runs):
path_run = path + "/" + run path_run = path + "/" + run
nb_outputs = 0
for j, num in enumerate(nums): for num in nums:
# Load property file try:
name_save = path_run + "/prop_disk_" + str(num).zfill(5) + ".save" # Load property file
name_save = path_run + "/prop_disk_" + str(num).zfill(5) + ".save"
# Check if the properties file exists # Check if the properties file exists
if len(glob.glob(name_save)) == 0: if len(glob.glob(name_save)) == 0:
raise ("no pickle file for disk properties. Run disk_prop() first") raise IOError(
f = open(name_save, "r") "no pickle file for disk properties {}. Run disk_prop() first".format(
prop_disk = pickle.load(f) name_save
f.close() )
)
f = open(name_save, "r")
prop_disk = pickle.load(f)
f.close()
alpha_rey[i] = alpha_rey[i] + prop_disk["alpha_rey_mean"] / len(nums) alpha_rey[i] = alpha_rey[i] + abs(prop_disk["alpha_rey_mean"])
alpha_grav[i] = alpha_grav[i] + prop_disk["alpha_grav_mean"] / len(nums) alpha_grav[i] = alpha_grav[i] + abs(prop_disk["alpha_grav_mean"])
Q[i] = Q[i] + prop_disk["Q_min"] / len(nums) Q[i] = Q[i] + prop_disk["Q_min"]
if pdf:
fit = prop_disk["fit"]
beta[i] = fit["beta"]
slope[i] = slope[i] + fit["slope"]
nb_outputs = nb_outputs + 1
print(run, num, nb_outputs)
except (IOError, KeyError) as e:
print(run, num, e, nb_outputs)
if nb_outputs > 0:
alpha_rey[i] = alpha_rey[i] / nb_outputs
alpha_grav[i] = alpha_grav[i] / nb_outputs
Q[i] = Q[i] / nb_outputs
if pdf: if pdf:
fit = prop_disk["fit"] slope[i] = slope[i] / nb_outputs
beta[i] = fit["beta"] else:
slope[i] = slope[i] + fit["slope"] / len(nums) for array in [alpha_rey, alpha_grav, Q]:
array[i] = np.nan
if pdf:
slope[i] = np.nan
if Q_in_name: if Q_in_name:
Q0[i] = float(run.split("_")[2][1:]) Q0[i] = float(run.split("_")[2][1:])
@@ -1096,7 +1119,7 @@ def compare(
if interactive: if interactive:
P.figure() P.figure()
else: else:
P.savefig(path + "/alphaQ_" + nums_name + out_ext) P.savefig(path_out + "/alphaQ_" + nums_name + out_ext)
P.close() P.close()
if Q_in_name: if Q_in_name:
@@ -1115,7 +1138,7 @@ def compare(
if interactive: if interactive:
P.figure() P.figure()
else: else:
P.savefig(path + "/alphaQ0_" + nums_name + out_ext) P.savefig(path_out + "/alphaQ0_" + nums_name + out_ext)
P.close() P.close()
if pdf: if pdf:
@@ -1130,7 +1153,7 @@ def compare(
if interactive: if interactive:
P.figure() P.figure()
else: else:
P.savefig(path + "/dcolslopebeta_" + nums_name + out_ext) P.savefig(path_out + "/dcolslopebeta_" + nums_name + out_ext)
P.close() P.close()
# alpha = f(beta) # alpha = f(beta)
@@ -1152,7 +1175,7 @@ def compare(
if interactive: if interactive:
P.figure() P.figure()
else: else:
P.savefig(path + "/alphabeta_" + nums_name + out_ext) P.savefig(path_out + "/alphabeta_" + nums_name + out_ext)
P.close() P.close()
@@ -1190,14 +1213,20 @@ def evolution(path, nums, force=False, interactive=False, pdf=False):
f.close() f.close()
time[i] = prop_disk["time"] time[i] = prop_disk["time"]
alpha_rey[i] = prop_disk["alpha_rey_mean"]
alpha_grav[i] = prop_disk["alpha_grav_mean"] try:
Qmin[i] = prop_disk["Q_min"] alpha_rey[i] = prop_disk["alpha_rey_mean"]
Qmean[i] = prop_disk["Q_mean"] alpha_grav[i] = prop_disk["alpha_grav_mean"]
mass_disk[i] = prop_disk["mass_disk"] Qmin[i] = prop_disk["Q_min"]
mass_box[i] = prop_disk["mass_box"] Qmean[i] = prop_disk["Q_mean"]
mass_disk[i] = prop_disk["mass_disk"]
mass_box[i] = prop_disk["mass_box"]
except:
for array in [alpha_rey, alpha_grav, Qmin, Qmean, mass_disk, mass_box]:
array[i] = np.nan
if pdf: if pdf:
slope = prop_disk["fit"]["slope"] slope[i] = prop_disk["fit"]["slope"]
# Check if the output file exists, and exit if it is the case # Check if the output file exists, and exit if it is the case
name_save = path + "/alpha_time" + out_ext name_save = path + "/alpha_time" + out_ext
@@ -1255,9 +1284,10 @@ def evolution(path, nums, force=False, interactive=False, pdf=False):
# PDF # PDF
if pdf: if pdf:
P.grid() P.grid()
print(time, slope)
P.plot(time, slope, "o-.") P.plot(time, slope, "o-.")
P.legend() P.legend()
P.ylabel(r"$d\log\% / d\logN$") P.ylabel(r"$d\log\mathcal{P}_{N} / d\logN$")
P.xlabel(r"time (code)") P.xlabel(r"time (code)")
if interactive: if interactive:
P.figure() P.figure()
pipeline_disk.py
+22 -12
View File
@@ -238,14 +238,15 @@ if args.compare:
path = storage_out + path_suffix path = storage_out + path_suffix
if args.time_avg: if args.time_avg:
# Select output availables for all runs # Select output availables for at least on run
output_ok = reduce(np.intersect1d, [run_succeded[run] for run in runs]) output_ok = reduce(np.union1d, [run_succeded[run] for run in runs]).astype(int)
print(output_ok) print(output_ok)
if len(output_ok) >= 1: if len(output_ok) >= 1:
dp.compare( dp.compare(
path, path,
runs, runs,
output_ok, output_ok,
path_out=path + "/comp",
force=args.force_redo, force=args.force_redo,
interactive=args.interactive, interactive=args.interactive,
Q_in_name=(not args.pdf), Q_in_name=(not args.pdf),
@@ -255,13 +256,22 @@ if args.compare:
for i in range(first, last + 1, step): for i in range(first, last + 1, step):
# Select usable runs # Select usable runs
runs_ok = [run for run in runs if i in run_succeded[run]] runs_ok = [run for run in runs if i in run_succeded[run]]
if len(runs_ok) > 1:
dp.compare( try:
path, if len(runs_ok) > 1:
runs_ok, dp.compare(
i, path,
force=args.force_redo, runs_ok,
interactive=args.interactive, i,
Q_in_name=(not args.pdf), path_out=path + "/comp",
pdf=args.pdf, force=args.force_redo,
) interactive=args.interactive,
Q_in_name=(not args.pdf),
pdf=args.pdf,
)
except (KeyError, IOError) as e:
print(e)
if args.skip:
pass
else:
raise