From 8e5813f6e28153de8b5701bd015059ba51457b2c Mon Sep 17 00:00:00 2001 From: Noe Brucy Date: Tue, 10 Dec 2019 17:03:00 +0100 Subject: [PATCH] Add extractor for turb rms --- plotter.py | 18 +++- postprocessor.py | 91 ++++++++++-------- pp_params.yml | 52 ++++++++++ run_selector.py | 245 ++++++++++++++++++++++------------------------- 4 files changed, 233 insertions(+), 173 deletions(-) create mode 100644 pp_params.yml diff --git a/plotter.py b/plotter.py index cb63b46..d4d6dd6 100644 --- a/plotter.py +++ b/plotter.py @@ -109,6 +109,13 @@ class Plotter(Aggregator, BaseProcessor): else: 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): if not arg is None: name_full = name + "_" + str(arg) @@ -158,7 +165,7 @@ class Plotter(Aggregator, BaseProcessor): def _plot_rule( 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: P.figure() rule.plot(save, arg, **kwargs) @@ -217,12 +224,17 @@ class Plotter(Aggregator, BaseProcessor): prop_value = self.comp.get_nml(nml_key, run) if prop_name in self.value_convert: 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) + else: + prop_value_str = str(prop_value) return r"{} = {}".format(prop_label, prop_value_str) 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: if not type(nml_key) == list: nml_key = [nml_key] diff --git a/postprocessor.py b/postprocessor.py index 9fdd0b2..6c159a5 100644 --- a/postprocessor.py +++ b/postprocessor.py @@ -20,7 +20,6 @@ from pymses.analysis import ScalarOperator, FractionOperator, MaxLevelOperator from mypool import MyPool as Pool from functools import partial from abc import ABCMeta, abstractmethod -import contextlib import bunch from run_selector import * @@ -167,7 +166,7 @@ class BaseProcessor: self._log("Dependency {} for {} is unknown".format(dep, name), "ERROR") 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): if not arg is None: @@ -265,6 +264,9 @@ class HDF5Container(BaseProcessor): self.save.close() 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): self.open() try: @@ -349,15 +351,7 @@ class PostProcessor(HDF5Container): cells_loaded = False - def __init__( - self, - path=None, - num=None, - path_out=None, - pp_params=None, - tag=None, - variables=["rho", "vel", "Br", "Bl", "P", "g", "phi"], - ): + def __init__(self, path=None, num=None, path_out=None, pp_params=None, tag=None): """ Creates the basic structures needed for the outputs """ @@ -381,9 +375,10 @@ class PostProcessor(HDF5Container): self.path = path self.run = os.path.basename(path) self.num = num - self._ro = pymses.RamsesOutput(path, num, order=pp_params.pymses.order) - self.variables = variables - self._amr = self._ro.amr_source(self.variables) + self._ro = pymses.RamsesOutput( + path, num, order=pp_params.pymses.order, verbose=pp_params.pymses.verbose + ) + self._amr = self._ro.amr_source(self.pp_params.pymses.variables) self.info = self._ro.info.copy() # Density operator @@ -442,47 +437,55 @@ class PostProcessor(HDF5Container): map_max_size=pp_params.out.map_size, ) - self._add_metadata() self.close() self.log_id = "[{}, {}] ".format(self.run, self.num) 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): + """ + 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: - # with self.open_pymlog() as f, contextlib.redirect_stdout(f): cell_source = CellsToPoints(self._amr) self.cells = cell_source.flatten() self.cells_loaded = True def unload_cells(self): + """ + Free space in the memory by telling the garbage collectors that + self.cells is not needed + """ if self.cells_loaded: del self.cells self.cells_loaded = False 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) datamap = slicing.SliceMap(self._amr, self._cam[ax_los], op, z=z) return datamap.map.T @@ -490,6 +493,9 @@ class PostProcessor(HDF5Container): def _ax_avg( 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: def num(cells): @@ -559,7 +565,7 @@ class PostProcessor(HDF5Container): return dmap_P / dmap_rho def _levels(self, ax_los): - self._amr.set_read_levelmax(20) + self._amr.set_read_levelmax(self.pp_params.pymses.levelmax) level_op = MaxLevelOperator() rt_level = raytracing.RayTracer(self._amr, self._ro.info, level_op) datamap = rt_level.process(self._cam[ax_los], surf_qty=True) @@ -632,6 +638,9 @@ class PostProcessor(HDF5Container): return map_Q def _radial_bins(self, _): + """ + Computes radial bins (for disk) + """ pos_star = self.pp_params.disk.pos_star im_extent = self.save.root.maps._v_attrs.im_extent @@ -660,6 +669,9 @@ class PostProcessor(HDF5Container): return rad_bins def _rr(self, _): + """ + Computes the radius from the center + """ im_extent = self.save.root.maps._v_attrs.im_extent map_size = self.pp_params.out.map_size pos_star = self.pp_params.disk.pos_star @@ -1010,7 +1022,6 @@ class Comparator(Aggregator, HDF5Container): # Get postprocesor objets for each run self.pp_runs = {} - attrs = {} for run in self.runs: path_run = path + "/" + run @@ -1041,7 +1052,7 @@ class Comparator(Aggregator, HDF5Container): else: 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_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]]) > 0 diff --git a/pp_params.yml b/pp_params.yml new file mode 100644 index 0000000..ef3c5ae --- /dev/null +++ b/pp_params.yml @@ -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 diff --git a/run_selector.py b/run_selector.py index ab1a693..3e5f5f4 100644 --- a/run_selector.py +++ b/run_selector.py @@ -7,153 +7,138 @@ from pp_params import * import f90nml -class RunSelector: - def __init__( - self, - path_in, - 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.pp_params = pp_params + def __init__(self, path_in, + in_runs=None, in_nums="all", pp_params=default_params(), + name_run='*', namelist_cond = dict(), + sort_run_by=None, time_min=None, time_max=None): + self.path_in = path_in + self.pp_params = pp_params - self.namelist = {} - self.runs = self.get_runs(in_runs, name_run, namelist_cond, sort_run_by) + self.namelist = dict() + self.runs = self.get_runs(in_runs, name_run, namelist_cond, sort_run_by) - self.info = {} - for run in self.runs: - self.info[run] = {} - self.nums = {} + self.info = dict() + for run in self.runs: + self.info[run] = dict() + self.nums = dict() - if not type(in_nums) == dict: - nums_temp = in_nums - in_nums = {} - for run in self.runs: - in_nums[run] = nums_temp + if not type(in_nums) == dict : + nums_temp = in_nums + in_nums = dict() + for run in self.runs: + in_nums[run] = nums_temp - for run in self.runs: - self.nums[run] = self.get_nums(run, in_nums[run], time_min, time_max) + for run in self.runs: + self.nums[run] = self.get_nums(run, in_nums[run], time_min, time_max) - def load_namelist(self, run): - path_run = self.path_in + "/" + run - path_nml = path_run + "/" + self.pp_params.input.nml_filename - return f90nml.read(path_nml) + def load_namelist(self, run): + path_run = self.path_in + "/" + run + path_nml = path_run + '/' + self.pp_params.input.nml_filename + return f90nml.read(path_nml) - def get_nml_value(self, nml_key, run): - res = self.namelist[run] - for key in nml_key.split("/"): - res = res[key] - return res + def get_nml_value(self, nml_key, run): + res = self.namelist[run] + for key in nml_key.split('/'): + if key in res: + res = res[key] + elif key == nml_key.split('/')[-1]: + res = None + else: + raise KeyError(key) + return res - def get_runs(self, in_runs=None, name_run="*", namelist_cond={}, sort_run_by=None): - def try_load_nml(run): - try: - self.namelist[run] = self.load_namelist(run) - success = True - except IOError: - success = False - return success + def get_runs(self, in_runs=None, name_run='*', namelist_cond=dict(), sort_run_by=None): + def try_load_nml(run): + try : + self.namelist[run] = self.load_namelist(run) + success = True + except IOError: + success = False + return success - runs = map( - os.path.basename, - filter(os.path.isdir, glob.glob(self.path_in + "/" + name_run)), - ) - if not in_runs is None: - runs = filter(lambda n: n in runs, in_runs) - runs = filter(try_load_nml, runs) + runs = map(os.path.basename, filter(os.path.isdir, glob.glob(self.path_in + "/" + name_run))) + if not in_runs is None : + runs = filter(lambda n : n in runs, in_runs) + runs = filter(try_load_nml, runs) - if type(namelist_cond) == tuple: - namelist_cond = [namelist_cond] + if type(namelist_cond) == tuple: + namelist_cond = [namelist_cond] - for (nml_key, operator, operand) in namelist_cond: - value = {} - for run in runs: - value[run] = self.get_nml_value(nml_key, run) - if operator == "=": - runs = filter(lambda r: value[r] == operand, runs) - if operator == "!=": - runs = filter(lambda r: not value[r] == operand, runs) - elif operator == ">": - runs = filter(lambda r: value[r] > operand, runs) - elif operator == "<": - runs = filter(lambda r: value[r] < operand, runs) - elif operator == "in": - runs = filter(lambda r: value[r] in operand, runs) + for (nml_key, operator, operand) in namelist_cond: + value = dict() + for run in runs: + value[run] = self.get_nml_value(nml_key, run) + if operator == '=': + runs = filter(lambda r: value[r] == operand, runs) + if operator == '!=': + runs = filter(lambda r: not value[r] == operand, runs) + elif operator == '>': + runs = filter(lambda r: value[r] > operand, runs) + elif operator == '<': + runs = filter(lambda r: value[r] < operand, runs) + elif operator == 'in': + runs = filter(lambda r: value[r] in operand, runs) - # Sort by the value in the namelist of sort_run_by - if not sort_run_by is None: - if type(sort_run_by) == str: - sort_run_by = [sort_run_by] - for nml_key in reversed(sort_run_by): - runs.sort(key=partial(self.get_nml_value, nml_key)) - return runs + # Sort by the value in the namelist of sort_run_by + if not sort_run_by is None: + if type(sort_run_by) == str: + sort_run_by = [sort_run_by] + for nml_key in reversed(sort_run_by): + runs.sort(key=partial(self.get_nml_value, nml_key)) - def load_info(self, run, num): - info_file = open( - self.path_in - + "/" - + run - + "/" - + "output_" - + str(num).zfill(5) - + "/" - + "info_" - + str(num).zfill(5) - + ".txt", - "r", - ) - info = {} - for line in info_file.readlines(): - parsed = yaml.safe_load(line.replace("=", ":")) - if type(parsed) == dict: - info.update(parsed) - info_file.close() - return info + return runs - def get_nums(self, run, in_nums=None, time_min=None, time_max=None): - def try_load_info(num): - try: - self.info[run][num] = self.load_info(run, num) - success = True - except IOError: - success = False - return success + def load_info(self, run, num): + info_file = open(self.path_in + "/" + run + "/" + + "output_" + str(num).zfill(5) + "/" + + "info_" + str(num).zfill(5) + ".txt", "r") + info = dict() + for line in info_file.readlines(): + parsed = yaml.safe_load(line.replace("=", ":")) + if type(parsed) == dict: + info.update(parsed) + info_file.close() + return info - names = glob.glob( - 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) + def get_nums(self, run, in_nums=None, time_min=None, time_max=None): + def try_load_info(num): + try : + self.info[run][num] = self.load_info(run, num) + success = True + except IOError: + success = False + return success - if type(in_nums) == int: - in_nums = [in_nums] - if type(in_nums) == list: - nums = filter(lambda n: n in nums, in_nums) + names = glob.glob(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 = np.sort(nums) - if in_nums == "first": - i = 0 - while i < len(nums) and not try_load_info(nums[i]): - i = i + 1 - nums = [nums[i]] - elif in_nums == "last": - i = len(nums) - 1 - while i >= 0 and not try_load_info(nums[i]): - i = i - 1 - nums = [nums[i]] - else: - nums = filter(try_load_info, nums) + if type(in_nums) == int: + in_nums = [in_nums] + if type(in_nums) == list: + nums = filter(lambda n : n in nums, in_nums) - if not time_min is None: - nums = filter(lambda n: self.info[run][n]["time"] >= time_min, nums) - if not time_max is None: - nums = filter(lambda n: self.info[run][n]["time"] <= time_max, nums) + nums = np.sort(nums) - return nums + if in_nums == "first": + i = 0 + while i < len(nums) - 1 and not try_load_info(nums[i]): + i = i + 1 + nums = [nums[i]] + elif in_nums == "last": + i = len(nums) - 1 + while i > 0 and not try_load_info(nums[i]): + i = i - 1 + nums = [nums[i]] + else: + nums = filter(try_load_info, nums) + + if not time_min is None: + nums = filter(lambda n: self.info[run][n]['time'] >= time_min, nums) + if not time_max is None: + nums = filter(lambda n: self.info[run][n]['time'] <= time_max, nums) + + + + return nums