Add extractor for turb rms

This commit is contained in:
Noe Brucy
2019-12-10 17:03:00 +01:00
parent 10600f53df
commit 8e5813f6e2
4 changed files with 233 additions and 173 deletions
+15 -3
View File
@@ -109,6 +109,13 @@ class Plotter(Aggregator, BaseProcessor):
else: else:
super(Plotter, self)._not_self_dep(name, dep, dep_arg, overwrite, **kwargs) super(Plotter, self)._not_self_dep(name, dep, dep_arg, overwrite, **kwargs)
def _needs_computation(self, overwrite, plot_filename):
return (
self.pp_params.out.interactive
or overwrite
or not os.path.exists(plot_filename)
)
def _process_rule(self, name, rule, arg, overwrite=False, **kwargs): def _process_rule(self, name, rule, arg, overwrite=False, **kwargs):
if not arg is None: if not arg is None:
name_full = name + "_" + str(arg) name_full = name + "_" + str(arg)
@@ -158,7 +165,7 @@ class Plotter(Aggregator, BaseProcessor):
def _plot_rule( def _plot_rule(
self, rule, save, arg, plot_filename, overwrite, open_figure=True, **kwargs self, rule, save, arg, plot_filename, overwrite, open_figure=True, **kwargs
): ):
if overwrite or not os.path.exists(plot_filename): if self._needs_computation(overwrite, plot_filename):
if open_figure: if open_figure:
P.figure() P.figure()
rule.plot(save, arg, **kwargs) rule.plot(save, arg, **kwargs)
@@ -217,12 +224,17 @@ class Plotter(Aggregator, BaseProcessor):
prop_value = self.comp.get_nml(nml_key, run) prop_value = self.comp.get_nml(nml_key, run)
if prop_name in self.value_convert: if prop_name in self.value_convert:
prop_value_str = self.value_convert[prop_name](prop_value) prop_value_str = self.value_convert[prop_name](prop_value)
else: elif type(prop_value) in [int, float]:
prop_value_str = "${:.6g}$".format(prop_value) prop_value_str = "${:.6g}$".format(prop_value)
else:
prop_value_str = str(prop_value)
return r"{} = {}".format(prop_label, prop_value_str) return r"{} = {}".format(prop_label, prop_value_str)
if nml_key is None and label is None: if nml_key is None and label is None:
label_run = r"{}".format(self.save.root._v_attrs.attrs[node.name].label) if "attrs" in self.save.root._v_attrs:
label_run = r"{}".format(self.save.root._v_attrs.attrs[run].label)
else:
label_run = run
elif not nml_key is None: elif not nml_key is None:
if not type(nml_key) == list: if not type(nml_key) == list:
nml_key = [nml_key] nml_key = [nml_key]
+51 -40
View File
@@ -20,7 +20,6 @@ from pymses.analysis import ScalarOperator, FractionOperator, MaxLevelOperator
from mypool import MyPool as Pool from mypool import MyPool as Pool
from functools import partial from functools import partial
from abc import ABCMeta, abstractmethod from abc import ABCMeta, abstractmethod
import contextlib
import bunch import bunch
from run_selector import * from run_selector import *
@@ -167,7 +166,7 @@ class BaseProcessor:
self._log("Dependency {} for {} is unknown".format(dep, name), "ERROR") self._log("Dependency {} for {} is unknown".format(dep, name), "ERROR")
def _needs_computation(self, overwrite, name_full): def _needs_computation(self, overwrite, name_full):
return overwrite or not (name_full in self.save) return overwrite
def _process_rule(self, name, rule, arg, overwrite=False, **kwargs): def _process_rule(self, name, rule, arg, overwrite=False, **kwargs):
if not arg is None: if not arg is None:
@@ -265,6 +264,9 @@ class HDF5Container(BaseProcessor):
self.save.close() self.save.close()
self.opened = False self.opened = False
def _needs_computation(self, overwrite, name_full):
return overwrite or not (name_full in self.save)
def _process_rule(self, name, rule, arg, overwrite, **kwargs): def _process_rule(self, name, rule, arg, overwrite, **kwargs):
self.open() self.open()
try: try:
@@ -349,15 +351,7 @@ class PostProcessor(HDF5Container):
cells_loaded = False cells_loaded = False
def __init__( def __init__(self, path=None, num=None, path_out=None, pp_params=None, tag=None):
self,
path=None,
num=None,
path_out=None,
pp_params=None,
tag=None,
variables=["rho", "vel", "Br", "Bl", "P", "g", "phi"],
):
""" """
Creates the basic structures needed for the outputs Creates the basic structures needed for the outputs
""" """
@@ -381,9 +375,10 @@ class PostProcessor(HDF5Container):
self.path = path self.path = path
self.run = os.path.basename(path) self.run = os.path.basename(path)
self.num = num self.num = num
self._ro = pymses.RamsesOutput(path, num, order=pp_params.pymses.order) self._ro = pymses.RamsesOutput(
self.variables = variables path, num, order=pp_params.pymses.order, verbose=pp_params.pymses.verbose
self._amr = self._ro.amr_source(self.variables) )
self._amr = self._ro.amr_source(self.pp_params.pymses.variables)
self.info = self._ro.info.copy() self.info = self._ro.info.copy()
# Density operator # Density operator
@@ -442,47 +437,55 @@ class PostProcessor(HDF5Container):
map_max_size=pp_params.out.map_size, map_max_size=pp_params.out.map_size,
) )
self._add_metadata()
self.close() self.close()
self.log_id = "[{}, {}] ".format(self.run, self.num) self.log_id = "[{}, {}] ".format(self.run, self.num)
self.def_rules() self.def_rules()
def _add_metadata(self):
"""
Add additional metadata to the file
"""
# Label of the run in the label.txt file
label_filename = self.path + "/" + self.pp_params.input.label_filename
if os.path.exists(label_filename):
label_file = open(label_filename, "r")
self.label = label_file.readline()
label_file.close()
else:
self.label = self.run
self.save.root._v_attrs.label = self.label
# def open_pymlog(self):
# if self.pp_params.pymses.verbose:
# return sys.stdout
# else:
# return open(os.devnull, "w")
def load_cells(self): def load_cells(self):
"""
Load all cells from the source file in the memory.
Cells will be accessible trough self.cells
(/!\ Long and memory heavy)
"""
if not self.cells_loaded: if not self.cells_loaded:
# with self.open_pymlog() as f, contextlib.redirect_stdout(f):
cell_source = CellsToPoints(self._amr) cell_source = CellsToPoints(self._amr)
self.cells = cell_source.flatten() self.cells = cell_source.flatten()
self.cells_loaded = True self.cells_loaded = True
def unload_cells(self): def unload_cells(self):
"""
Free space in the memory by telling the garbage collectors that
self.cells is not needed
"""
if self.cells_loaded: if self.cells_loaded:
del self.cells del self.cells
self.cells_loaded = False self.cells_loaded = False
def _slice(self, getter, ax_los="z", z=0, unit=cst.none): def _slice(self, getter, ax_los="z", z=0, unit=cst.none):
"""
Slice process function.
Return a slice of the source box.
Parameters
----------
getter : callable
A callable that extract the wanted data from a pymses dataset
ax_los : string
The axis perpendicular to the slice plane
z : float
Coordinate of the slice on the ax_los axis
unit : cst.Unit
Unit of the resulting dataset
Returns
-------
A numpy array containing the slice
"""
op = ScalarOperator(getter, unit) op = ScalarOperator(getter, unit)
datamap = slicing.SliceMap(self._amr, self._cam[ax_los], op, z=z) datamap = slicing.SliceMap(self._amr, self._cam[ax_los], op, z=z)
return datamap.map.T return datamap.map.T
@@ -490,6 +493,9 @@ class PostProcessor(HDF5Container):
def _ax_avg( def _ax_avg(
self, getter, ax_los, unit=cst.none, mass_weighted=True, surf_qty=False self, getter, ax_los, unit=cst.none, mass_weighted=True, surf_qty=False
): ):
"""
Map of the average of a quantity (given by getter) along an axis (ax_los)
"""
if mass_weighted: if mass_weighted:
def num(cells): def num(cells):
@@ -559,7 +565,7 @@ class PostProcessor(HDF5Container):
return dmap_P / dmap_rho return dmap_P / dmap_rho
def _levels(self, ax_los): def _levels(self, ax_los):
self._amr.set_read_levelmax(20) self._amr.set_read_levelmax(self.pp_params.pymses.levelmax)
level_op = MaxLevelOperator() level_op = MaxLevelOperator()
rt_level = raytracing.RayTracer(self._amr, self._ro.info, level_op) rt_level = raytracing.RayTracer(self._amr, self._ro.info, level_op)
datamap = rt_level.process(self._cam[ax_los], surf_qty=True) datamap = rt_level.process(self._cam[ax_los], surf_qty=True)
@@ -632,6 +638,9 @@ class PostProcessor(HDF5Container):
return map_Q return map_Q
def _radial_bins(self, _): def _radial_bins(self, _):
"""
Computes radial bins (for disk)
"""
pos_star = self.pp_params.disk.pos_star pos_star = self.pp_params.disk.pos_star
im_extent = self.save.root.maps._v_attrs.im_extent im_extent = self.save.root.maps._v_attrs.im_extent
@@ -660,6 +669,9 @@ class PostProcessor(HDF5Container):
return rad_bins return rad_bins
def _rr(self, _): def _rr(self, _):
"""
Computes the radius from the center
"""
im_extent = self.save.root.maps._v_attrs.im_extent im_extent = self.save.root.maps._v_attrs.im_extent
map_size = self.pp_params.out.map_size map_size = self.pp_params.out.map_size
pos_star = self.pp_params.disk.pos_star pos_star = self.pp_params.disk.pos_star
@@ -1010,7 +1022,6 @@ class Comparator(Aggregator, HDF5Container):
# Get postprocesor objets for each run # Get postprocesor objets for each run
self.pp_runs = {} self.pp_runs = {}
attrs = {}
for run in self.runs: for run in self.runs:
path_run = path + "/" + run path_run = path + "/" + run
@@ -1041,7 +1052,7 @@ class Comparator(Aggregator, HDF5Container):
else: else:
saved_nums = self.save.get_node(name_full)._v_attrs.nums saved_nums = self.save.get_node(name_full)._v_attrs.nums
missing_runs = len([run for run in self.nums if not run in saved_nums]) > 0 missing_runs = len([run for run in self.nums if not run in saved_nums]) > 0
missing_nums = missing_runs and all( missing_nums = missing_runs or all(
[ [
len([num for num in self.nums[run] if not num in saved_nums[run]]) len([num for num in self.nums[run] if not num in saved_nums[run]])
> 0 > 0
+52
View File
@@ -0,0 +1,52 @@
plot : # Plot parameters
out_ext : '.jpeg' # extension for plots
put_title : False # Add a title to plot
ntick : 6 # Number of ticks for maps
vel_red : 40 # Take point each vel_red for velocities
disk: # Disk speficic parameters
enable : False # Enable specific disk analysis
pos_star : [1., 1., 1.] # Position of the central star
binning : "log" # Kind of binning (lin : linear, log : logarithmic)
nb_bin : 100 # Number of bins for averaged quantities
bin_in : 1e-3 # Outer radius of the inner bin
bin_out : 0.25 # Inner radius of the outer bin
rmin_pdf : 0.075 # Inner radius for PDF computation
rmax_pdf : 0.3 # Outer radius for PDF computation
pdf: # parameters for probability density functions
nb_bin : 50 # Number of bins for the PDF
xmin_fit : 0. # Lower boundary of the fit
xmax_fit : 1.25 # Upper boundary of the fit
pymses: # Parameters for Pymses reader
order : '<' # In which order the output are read
fft : False # Quick and dirty rendering using FFT
verbose : False # Let pymses write on standart output
variables : ["rho","vel","Br","Bl","P", "g", "phi"] # Read these variables
levelmax : 20 # Maximal AMR level visited when looking levels
input: # Parameters on how to look for input files (: output from Ramses)
log_prefix : "run.log" # Prefix of the log file
label_filename : "label.txt" # Name of the label file
nml_filename : "run.nml" # name of the namelist file
out: # Parameters for post processing
center : [0.5, 0.5, 0.5] # Center of the image
zoom : 1. # Zoom of the image
map_size : 1024 # Size of the computed maps in pixel
tag : "" # Tag for the image
interactive : False # Interactive mode (do not save the plots on the disk)
process:
verbose : True
num_process : 1
+44 -59
View File
@@ -7,33 +7,24 @@ from pp_params import *
import f90nml import f90nml
class RunSelector: def __init__(self, path_in,
def __init__( in_runs=None, in_nums="all", pp_params=default_params(),
self, name_run='*', namelist_cond = dict(),
path_in, sort_run_by=None, time_min=None, time_max=None):
in_runs=None,
in_nums="all",
pp_params=default_params(),
name_run="*",
namelist_cond={},
sort_run_by=None,
time_min=None,
time_max=None,
):
self.path_in = path_in self.path_in = path_in
self.pp_params = pp_params self.pp_params = pp_params
self.namelist = {} self.namelist = dict()
self.runs = self.get_runs(in_runs, name_run, namelist_cond, sort_run_by) self.runs = self.get_runs(in_runs, name_run, namelist_cond, sort_run_by)
self.info = {} self.info = dict()
for run in self.runs: for run in self.runs:
self.info[run] = {} self.info[run] = dict()
self.nums = {} self.nums = dict()
if not type(in_nums) == dict: if not type(in_nums) == dict :
nums_temp = in_nums nums_temp = in_nums
in_nums = {} in_nums = dict()
for run in self.runs: for run in self.runs:
in_nums[run] = nums_temp in_nums[run] = nums_temp
@@ -42,50 +33,53 @@ class RunSelector:
def load_namelist(self, run): def load_namelist(self, run):
path_run = self.path_in + "/" + run path_run = self.path_in + "/" + run
path_nml = path_run + "/" + self.pp_params.input.nml_filename path_nml = path_run + '/' + self.pp_params.input.nml_filename
return f90nml.read(path_nml) return f90nml.read(path_nml)
def get_nml_value(self, nml_key, run): def get_nml_value(self, nml_key, run):
res = self.namelist[run] res = self.namelist[run]
for key in nml_key.split("/"): for key in nml_key.split('/'):
if key in res:
res = res[key] res = res[key]
elif key == nml_key.split('/')[-1]:
res = None
else:
raise KeyError(key)
return res return res
def get_runs(self, in_runs=None, name_run="*", namelist_cond={}, sort_run_by=None): def get_runs(self, in_runs=None, name_run='*', namelist_cond=dict(), sort_run_by=None):
def try_load_nml(run): def try_load_nml(run):
try: try :
self.namelist[run] = self.load_namelist(run) self.namelist[run] = self.load_namelist(run)
success = True success = True
except IOError: except IOError:
success = False success = False
return success return success
runs = map( runs = map(os.path.basename, filter(os.path.isdir, glob.glob(self.path_in + "/" + name_run)))
os.path.basename, if not in_runs is None :
filter(os.path.isdir, glob.glob(self.path_in + "/" + name_run)), runs = filter(lambda n : n in runs, in_runs)
)
if not in_runs is None:
runs = filter(lambda n: n in runs, in_runs)
runs = filter(try_load_nml, runs) runs = filter(try_load_nml, runs)
if type(namelist_cond) == tuple: if type(namelist_cond) == tuple:
namelist_cond = [namelist_cond] namelist_cond = [namelist_cond]
for (nml_key, operator, operand) in namelist_cond: for (nml_key, operator, operand) in namelist_cond:
value = {} value = dict()
for run in runs: for run in runs:
value[run] = self.get_nml_value(nml_key, run) value[run] = self.get_nml_value(nml_key, run)
if operator == "=": if operator == '=':
runs = filter(lambda r: value[r] == operand, runs) runs = filter(lambda r: value[r] == operand, runs)
if operator == "!=": if operator == '!=':
runs = filter(lambda r: not value[r] == operand, runs) runs = filter(lambda r: not value[r] == operand, runs)
elif operator == ">": elif operator == '>':
runs = filter(lambda r: value[r] > operand, runs) runs = filter(lambda r: value[r] > operand, runs)
elif operator == "<": elif operator == '<':
runs = filter(lambda r: value[r] < operand, runs) runs = filter(lambda r: value[r] < operand, runs)
elif operator == "in": elif operator == 'in':
runs = filter(lambda r: value[r] in operand, runs) runs = filter(lambda r: value[r] in operand, runs)
# Sort by the value in the namelist of sort_run_by # Sort by the value in the namelist of sort_run_by
if not sort_run_by is None: if not sort_run_by is None:
if type(sort_run_by) == str: if type(sort_run_by) == str:
@@ -96,20 +90,10 @@ class RunSelector:
return runs return runs
def load_info(self, run, num): def load_info(self, run, num):
info_file = open( info_file = open(self.path_in + "/" + run + "/" +
self.path_in "output_" + str(num).zfill(5) + "/" +
+ "/" "info_" + str(num).zfill(5) + ".txt", "r")
+ run info = dict()
+ "/"
+ "output_"
+ str(num).zfill(5)
+ "/"
+ "info_"
+ str(num).zfill(5)
+ ".txt",
"r",
)
info = {}
for line in info_file.readlines(): for line in info_file.readlines():
parsed = yaml.safe_load(line.replace("=", ":")) parsed = yaml.safe_load(line.replace("=", ":"))
if type(parsed) == dict: if type(parsed) == dict:
@@ -119,41 +103,42 @@ class RunSelector:
def get_nums(self, run, in_nums=None, time_min=None, time_max=None): def get_nums(self, run, in_nums=None, time_min=None, time_max=None):
def try_load_info(num): def try_load_info(num):
try: try :
self.info[run][num] = self.load_info(run, num) self.info[run][num] = self.load_info(run, num)
success = True success = True
except IOError: except IOError:
success = False success = False
return success return success
names = glob.glob( names = glob.glob(self.path_in + "/" + run + "/output_[0-9][0-9][0-9][0-9][0-9]")
self.path_in + "/" + run + "/output_[0-9][0-9][0-9][0-9][0-9]" nums = map(lambda n : int(n.split('/')[-1].split('_')[1]), names)
)
nums = map(lambda n: int(n.split("/")[-1].split("_")[1]), names)
if type(in_nums) == int: if type(in_nums) == int:
in_nums = [in_nums] in_nums = [in_nums]
if type(in_nums) == list: if type(in_nums) == list:
nums = filter(lambda n: n in nums, in_nums) nums = filter(lambda n : n in nums, in_nums)
nums = np.sort(nums) nums = np.sort(nums)
if in_nums == "first": if in_nums == "first":
i = 0 i = 0
while i < len(nums) and not try_load_info(nums[i]): while i < len(nums) - 1 and not try_load_info(nums[i]):
i = i + 1 i = i + 1
nums = [nums[i]] nums = [nums[i]]
elif in_nums == "last": elif in_nums == "last":
i = len(nums) - 1 i = len(nums) - 1
while i >= 0 and not try_load_info(nums[i]): while i > 0 and not try_load_info(nums[i]):
i = i - 1 i = i - 1
nums = [nums[i]] nums = [nums[i]]
else: else:
nums = filter(try_load_info, nums) nums = filter(try_load_info, nums)
if not time_min is None: if not time_min is None:
nums = filter(lambda n: self.info[run][n]["time"] >= time_min, nums) nums = filter(lambda n: self.info[run][n]['time'] >= time_min, nums)
if not time_max is None: if not time_max is None:
nums = filter(lambda n: self.info[run][n]["time"] <= time_max, nums) nums = filter(lambda n: self.info[run][n]['time'] <= time_max, nums)
return nums return nums