wetb/hawc2/body_matrix_output_file.py

+# -*- coding: utf-8 -*-
+"""
+Read HAWC2 results files produced by the `new_htc_structure` commands.
+
+@author: ricriv
+"""
+
+# %% Import.
+
+from glob import glob
+import numpy as np
+
+
+# %% Functions.
+
+
+def read_body_matrix_output(file_name):
+    """
+    Read output of `new_htc_structure` / `body_matrix_output`, which is the mass, damping and stiffness matrices for each body.
+
+    Parameters
+    ----------
+    file_name : str
+        File name entered in the
+        `new_htc_structure` / `body_matrix_output`
+        command. A glob patter will be used to find all the files.
+
+    Returns
+    -------
+    bodies : dict
+        Dictionary containing the structural matrices for each body.
+        Its keys are:
+            mass : ndarray of shape (ndof, ndof)
+                Mass matrix.
+            damping : ndarray of shape (ndof, ndof)
+                Damping matrix.
+            stiffness : ndarray of shape (ndof, ndof)
+                Stiffness matrix.
+    """
+    # Get all files with the mass matrix.
+    files = glob(f"{file_name}*_m.bin")
+
+    # We get the body names by taking the part between
+    # the file_name and _m.bin.
+    names = [s[len(file_name) : -6] for s in files]
+
+    # Dict that will contain all the results.
+    bodies = dict.fromkeys(names)
+
+    # Loop over the bodies.
+    for name in names:
+        bodies[name] = {}
+
+        # Read mass matrix.
+        with open(f"{file_name}{name}_m.bin", "rb") as fid:
+            ndof = np.fromfile(fid, dtype=np.int32, count=2)[0]
+            bodies[name]["mass"] = np.fromfile(
+                fid, dtype=np.float64, count=ndof * ndof
+            ).reshape(ndof, ndof, order="F")
+
+        # Read damping matrix.
+        with open(f"{file_name}{name}_c.bin", "rb") as fid:
+            ndof = np.fromfile(fid, dtype=np.int32, count=2)[0]
+            bodies[name]["damping"] = np.fromfile(
+                fid, dtype=np.float64, count=ndof * ndof
+            ).reshape(ndof, ndof, order="F")
+
+        # Read stiffness matrix.
+        with open(f"{file_name}{name}_k.bin", "rb") as fid:
+            ndof = np.fromfile(fid, dtype=np.int32, count=2)[0]
+            bodies[name]["stiffness"] = np.fromfile(
+                fid, dtype=np.float64, count=ndof * ndof
+            ).reshape(ndof, ndof, order="F")
+
+    return bodies

wetb/hawc2/cmp_test_cases.py

-from __future__ import division
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import absolute_import
-from io import open
-from builtins import range
-from builtins import str
-from builtins import zip
-from future import standard_library
-standard_library.install_aliases()
 import string
 import shutil
 import numpy as np

wetb/hawc2/hawc2_input_writer.py

+import itertools
+import os
+from pathlib import Path
+import warnings
+
+from tqdm import tqdm
+import jinja2
+import pandas as pd
+from pandas.core.base import PandasObject
+import numpy as np
+
+from wetb.hawc2.htc_file import HTCFile
+
+
+class HAWC2InputWriter(object):
+    """Load a dlc and write the corresponding HAWC2 input files.
+
+    Parameters
+    -----------
+    base_htc_file : str or Pathlib.path
+        Path to base htc file from which the other htc files are created.
+
+    Notes
+    ------
+    You can load the dlc keywords and parameters using the `from_` methods presented
+    below. Any methods of the form `set_<tag>`(htc, **kwargs), where `<tag>` is one
+    of one of the keywords in the dlc, will be called when writing htc files.
+
+    You can write the htc files using the `write` or `write_all` methods presented below.
+    The `Name` tag must be defined when writing htc files.
+
+    Necessary/protected tags
+    -------------------------
+    Name : str
+        The name of the htc file (without extension) for each simulation. This will
+        be used in a `set_Name` call to the htc file, which will overwrite the
+        `simulation.logfile` and `output.filename` lines.
+    """
+
+    def __init__(self, base_htc_file, **kwargs):
+        self.base_htc_file = base_htc_file
+        self.contents = None
+        for k, v in kwargs.items():
+            setattr(self, k, v)
+
+    def __call__(self, path, **kwargs):
+        """Write single htc file if writer is called directly"""
+        return self.write(path, **kwargs)
+
+    def from_pandas(self, dataFrame):
+        """Load dlc contents from a PandasObject.
+
+        Parameters
+        ----------
+        dataFrame : Pandas.DataFrame or PandasObject
+            Dataframe of contents.
+
+        Returns
+        -------
+        self : HAWC2InputWriter
+            Instance of self to enable chaining.
+        """
+        if not isinstance(dataFrame, PandasObject) and hasattr(dataFrame, 'to_pandas'):
+            dataFrame = dataFrame.to_pandas()
+        self.contents = dataFrame
+        return self
+
+    def from_excel(self, excel_file):
+        """Load dlc contents from an excel file.
+
+        Parameters
+        ----------
+        excel_file : str or Pathlib.Path
+            Path to excel file.
+
+        Returns
+        -------
+        self : HAWC2InputWriter
+            Instance of self to enable chaining.
+        """
+
+        self.contents = pd.read_excel(excel_file, keep_default_na=False)
+        if 'Unnamed: 0' in self.contents:
+            self.contents = self.contents.drop(columns=['Unnamed: 0'])
+        return self
+
+    def from_CVF(self, constants, variables, functionals):
+        """Load dlc contents from dictionaries of constants, variables and functionals.
+
+        Parameters
+        ----------
+        constants : dict
+            (content, value) pairs for which the value remains constant for all
+            generated htc files. E.g.: {'time_stop': 700}.
+        variables : dict of lists
+            (content, list) pairs for which all combinations of the listed values are
+            generated in the htc files. E.g., {'wsp': [4, 6, 8]}.
+        functionals : dict of functions
+            (content, function) pairs for which the content is dependent on the
+            constants and variables contents. E.g., {'Name': lambda x: 'wsp' + str(x['wsp'])}.
+
+        Returns
+        -------
+        self : HAWC2InputWriter
+            Instance of self to enable chaining.
+        """
+        attributes = list(constants.keys()) + list(variables.keys()) + list(functionals.keys())
+        contents = []
+
+        for var in itertools.product(*list(variables.values())):
+            this_dict = dict(constants)
+            var_dict = dict(zip(variables.keys(), var))
+            this_dict.update(var_dict)
+
+            for key, func in functionals.items():
+                this_dict[key] = func(this_dict)
+
+            contents.append(list(this_dict.values()))
+
+        self.contents = pd.DataFrame(contents, columns=attributes)
+
+        return self
+
+    def set_Name(self, htc, **kwargs):
+        """Update the filename in simulation and output blocks of an htc file.
+
+        Notes
+        -----
+        This function is called during the `write` and `write_all` methods, thereby
+        overwriting the logfile and output filename. If you would like to specify a
+        subfolder within the `htc` and `res` directories where your simulations should
+        be placed, pass in the 'Folder' keyword into `kwargs`.
+
+        Parameters
+        ----------
+        htc : wetb.hawc2.htc_file.HTCFile
+            The htc file object to be modified.
+        name : str
+            The filename to be used in the logfile and simulation file (without
+            extension).
+        **kwargs
+            Keyword arguments. Must include keyword 'Name'. May optionally include
+            'Folder' to define a subfolder in your `htc` and `res` directories where
+            files will be stored.
+        """
+        htc.set_name(kwargs['Name'], subfolder=kwargs.get('Folder', ''))
+
+    def write(self, path, **kwargs):
+        """Write a single htc file.
+
+        Notes
+        -----
+        This function works both with the tagless and jinja file-writing systems. Any
+        methods of the form `set_<tag>` are called during file writing.
+
+        Parameters
+        ----------
+        path : str or pathlib.Path
+            The path to save the htc file to.
+        **kwargs : pd.DataFrame, pd.Series or dict
+            Keyword arguments. The tags to update/replace in the file.
+        """
+        htc = HTCFile(self.base_htc_file, jinja_tags=kwargs)
+        for k, v in kwargs.items():
+            k = k.replace('/', '.')
+            if '.' in k:  # directly replace tags like "wind.wsp"
+                command = k.split(".")[-1]
+                subsections = k.split(".")[:-1]
+                section = htc
+                for subsection in subsections:
+                    section = section.add_section(subsection)
+                if command in section:
+                    section[command].values = np.atleast_1d(v).tolist()
+                else:
+                    section[command] = v
+            else:  # otherwise, use the "set_" attribute
+                if hasattr(self, 'set_%s' % k):
+                    getattr(self, 'set_%s' % k)(htc, **kwargs)
+        htc.save(path)
+
+    def write_all(self, out_dir):
+        """Write all htc files for the dlc.
+
+        Notes
+        -----
+        This function works both with the tagless and jinja file-writing systems. Any
+        methods of the form `set_<tag>` are called during each htc write. DLC contents
+        must contain a 'Name' column.
+
+        Parameters
+        ----------
+        out_dir : str or pathlib.Path
+            The directory in which the htc files should be written.
+        """
+        try:
+            self.contents['Name']
+        except KeyError:
+            raise KeyError('"Name" not given in dlc contents! Cannot write files.')
+
+        out_dir = Path(out_dir)
+        out_dir.mkdir(parents=True, exist_ok=True)
+
+        N = len(self.contents)
+
+        print(f'Generating {N} htc files in directory: {out_dir}')
+
+        with tqdm(self.contents.iterrows(), total=N) as bar:
+            for _, row in bar:
+                if 'Folder' in row and row.Folder:
+                    path = out_dir / row.Folder
+                else:
+                    path = out_dir
+                self.write(path / (row.Name + '.htc'), **row.to_dict())
+
+
+class JinjaWriter(HAWC2InputWriter):
+    def __init__(self, base_htc_file, **kwargs):
+        """DEPRECATED!!!
+        Subclass of the HAWC2InputWriter object. Generates htc files using contents compatible with jinja2.
+        """
+        warnings.warn('The JinjaWriter is deprecated! Please switch to the HAWC2InputWriter',
+                      DeprecationWarning)
+        HAWC2InputWriter.__init__(self, base_htc_file, **kwargs)
+
+    def write(self, out_fn, **kwargs):
+        params = pd.Series(kwargs)
+        with open(self.base_htc_file) as f:
+            template = jinja2.Template(f.read())
+
+            out = template.render(params)
+
+            with open(out_fn, 'w') as f:
+                f.write(out)
+
+
+def main():
+    if __name__ == '__main__':
+        from wetb.hawc2.tests import test_files
+
+        # example of subclassing with custom set_<tag> method
+        path = os.path.dirname(test_files.__file__) + '/simulation_setup/DTU10MWRef6.0/'
+        base_htc_file = path + 'htc/DTU_10MW_RWT.htc'
+
+        class MyWriter(HAWC2InputWriter):
+            def set_time(self, htc, **kwargs):
+                htc.set_time(self.time_start, self.time_start + kwargs['time'])
+
+        myWriter = MyWriter(base_htc_file, time_start=100)
+        myWriter('tmp/t1.htc', Name='t1', **{'time': 600, 'wind.wsp': 10})
+
+        # example with constants, variables, and functionals
+        Constants = {
+            'time_stop': 100,
+        }
+
+        Variables = {
+            'wsp': [4, 6, 8, 10, 12, 14, 16, 18, 20, 22, 24],
+            'tint': [0.05, 0.10, 0.15],
+
+        }
+
+        Functionals = {
+            'Name': lambda x: 'dtu10mw_wsp' + str(x['wsp']) + '_ti' + str(x['tint']),
+        }
+
+        htc_template_fn = os.path.dirname(test_files.__file__) + \
+            '/simulation_setup/DTU10MWRef6.0/htc/DTU_10MW_RWT.htc.j2'
+        writer = HAWC2InputWriter(htc_template_fn)
+        writer.from_CVF(Constants, Variables, Functionals)
+        print(writer.contents)
+
+        # write all htc files to folder 'tmp'
+        writer.write_all('tmp')
+
+
+main()

wetb/hawc2/hawc2_pbs_file.py

+from wetb.utils.cluster_tools.pbsfile import PBSFile, Template
+import os
+from wetb.utils.cluster_tools.os_path import pjoin, relpath, abspath,\
+    cluster_path, repl
+
+template = Template("""
+
+[copy_hawc2]
+
+#===============================================================================
+echo copy input
+#===============================================================================
+
+cd "[modelpath]"
+(flock -x 200
+[copy_input_to_scratch]
+) 200>/scratch/$USER/$PBS_JOBID/[modelname]/lock_file_model
+cd "/scratch/$USER/$PBS_JOBID/[modelname]"
+[copy_input_to_exe_dir]
+
+
+#===============================================================================
+echo Run HAWC2
+#===============================================================================
+cd "/scratch/$USER/$PBS_JOBID/[modelname]/run_[jobname]/[rel_exe_dir]"
+[hawc2_cmd] [htc_file]
+
+#===============================================================================
+echo Copy output
+#===============================================================================
+cd "/scratch/$USER/$PBS_JOBID/[modelname]/run_[jobname]"
+[copy_output]
+
+rm -r "/scratch/$USER/$PBS_JOBID/[modelname]/run_[jobname]"
+
+echo Done
+""")
+
+
+def wine_cmd(platform='win32', hawc2='hawc2mb.exe', cluster='jess'):
+    wine_folder = {'jess': {'win32': '.wine32', 'win64': '.wine'}}[cluster][platform]
+    wine_prefix = "WINEARCH=%s WINEPREFIX=~/%s " % (platform, wine_folder)
+    if cluster == 'jess':
+        s = wine_prefix + "winefix\n"
+    else:
+        s = ""
+    return s + wine_prefix + "wine %s" % hawc2
+
+
+JESS_WINE32_HAWC2MB = wine_cmd()
+
+
+class HAWC2PBSFile(PBSFile):
+    def __init__(self, hawc2_path, hawc2_cmd, htc_file, exe_dir, input_files,
+                 output_files, queue='workq', walltime='00:10:00'):
+        self.hawc2_path = hawc2_path
+        self.hawc2_cmd = hawc2_cmd
+        self.htc_file = htc_file
+        self.exe_dir = exe_dir
+        self.queue = queue
+        self.walltime = walltime
+
+        if not os.path.isabs(htc_file):
+            htc_file = pjoin(exe_dir, htc_file)
+        else:
+            htc_file = repl(htc_file)
+
+        if htc_file not in input_files:
+            input_files.append(htc_file)
+        self.input_files = [abspath((pjoin(exe_dir, f), abspath(f))[os.path.isabs(f)])
+                            for f in input_files]
+
+        self.htc_file = relpath(htc_file, exe_dir)
+
+        self.output_files = [abspath((pjoin(exe_dir, f), abspath(f))[os.path.isabs(f)])
+                             for f in output_files]
+
+        drive = os.path.splitdrive(exe_dir)[0]
+#         print()
+#         print("\n".join(self.output_files))
+        common_prefix = os.path.commonprefix(self.input_files + self.output_files).rpartition("/")[0]
+        p = abspath(pjoin(exe_dir,
+                          relpath(common_prefix, exe_dir)))
+        self.modelpath = os.path.join(drive, os.path.splitdrive(p)[1])
+        self.modelname = os.path.basename(abspath(self.modelpath))
+        self.jobname = os.path.splitext(os.path.basename(htc_file))[0]
+
+        PBSFile.__init__(self, self.modelpath, self.jobname, self.commands, queue, walltime=walltime)
+
+    def commands(self):
+        rel_exe_dir = relpath(abspath(self.exe_dir), abspath(self.modelpath))
+        copy_input_to_scratch, copy_input_to_exe_dir = self.copy_input()
+        return template(copy_hawc2=self.copy_hawc2(),
+                        exe_dir=cluster_path(self.exe_dir),
+                        copy_input_to_scratch=copy_input_to_scratch,
+                        copy_input_to_exe_dir=copy_input_to_exe_dir,
+                        rel_exe_dir=rel_exe_dir,
+                        hawc2_cmd=self.hawc2_cmd,
+                        htc_file=self.htc_file,
+                        jobname=self.jobname,
+                        copy_output=self.copy_output(),
+                        modelpath=cluster_path(self.modelpath),
+                        modelname=self.modelname)
+
+    def copy_hawc2(self):
+        copy_hawc2 = Template('''#===============================================================================
+echo copy hawc2 to scratch
+#===============================================================================
+(flock -x 200
+mkdir -p "/scratch/$USER/$PBS_JOBID/hawc2/"
+unzip -u -o -q "[hawc2_path]/"*.zip -d "/scratch/$USER/$PBS_JOBID/hawc2/"
+find "[hawc2_path]/"* ! -name *.zip -exec cp -u -t "/scratch/$USER/$PBS_JOBID/hawc2/" {} +
+) 200>"/scratch/$USER/$PBS_JOBID/lock_file_hawc2"
+mkdir -p "/scratch/$USER/$PBS_JOBID/[modelname]/run_[jobname]/[rel_exe_dir]"
+cp "/scratch/$USER/$PBS_JOBID/hawc2/"* "/scratch/$USER/$PBS_JOBID/[modelname]/run_[jobname]/[rel_exe_dir]"''')
+        if self.hawc2_path is None:
+            return ""
+        else:
+            return copy_hawc2(hawc2_path=os.path.dirname(cluster_path(self.hawc2_path)))
+
+    def copy_input(self):
+        rel_input_files = [relpath(f, self.modelpath) for f in self.input_files]
+
+        copy_input = "\n".join(['mkdir -p "[TARGET]/%s" && cp -u -r "%s" "[TARGET]/%s"' % (os.path.dirname(f), f, os.path.dirname(f))
+                                for f in rel_input_files])
+        return (copy_input.replace("[TARGET]", "/scratch/$USER/$PBS_JOBID/[modelname]"),
+                copy_input.replace("[TARGET]", "/scratch/$USER/$PBS_JOBID/[modelname]/run_[jobname]"))
+
+    def copy_output(self):
+        rel_output_files = [relpath(f, self.modelpath) for f in self.output_files]
+        return "\n".join(['mkdir -p "[modelpath]/%s" && cp -u -r "%s" "[modelpath]/%s"' % (os.path.dirname(f), f, os.path.dirname(f))
+                          for f in rel_output_files])

wetb/hawc2/htc_contents.py

-'''
-Created on 20/01/2014
-
-@author: MMPE
-
-See documentation of HTCFile below
-
-'''
-from __future__ import division
-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import absolute_import
-from builtins import zip
-from builtins import int
-from builtins import str
-from future import standard_library
-import os
-standard_library.install_aliases()
-from collections import OrderedDict
-import collections
-
-
-class OrderedDict(collections.OrderedDict):
-    pass
-
-    def __str__(self):
-        return "\n".join(["%-30s\t %s" % ((str(k) + ":"), str(v)) for k, v in self.items()])
-
-
-def parse_next_line(lines):
-    _3to2list = list(lines.pop(0).split(";"))
-    line, comments, = _3to2list[:1] + [_3to2list[1:]]
-    comments = ";".join(comments).rstrip()
-    while lines and lines[0].lstrip().startswith(";"):
-        comments += "\n%s" % lines.pop(0).rstrip()
-    return line.strip(), comments
-
-def fmt_value(v):
-    try:
-        if int(float(v)) == float(v):
-            return int(float(v))
-        return float(v)
-    except ValueError:
-        return v
-
-c = 0
-class HTCContents(object):
-    lines = []
-    contents = None
-    name_ = ""
-    parent=None
-
-    def __setitem__(self, key, value):
-        if isinstance(key, str):
-            self.contents[key] = value
-        elif isinstance(key, int):
-            self.values[key] = value
-        else:
-            raise NotImplementedError
-        value.parent=self
-            
-
-    def __getitem__(self, key):
-        if isinstance(key, str):
-            key = key.replace(".", "/")
-            if "/" in key:
-                keys = key.split('/')
-                val = self.contents[keys[0]]
-                for k in keys[1:]:
-                    val = val[k]
-                return val
-            return self.contents[key]
-        else:
-            return self.values[key]
-
-    def __getattr__(self, *args, **kwargs):
-        if args[0] in ['__members__','__methods__']:
-            # fix python2 related issue. In py2, dir(self) calls __getattr__(('__members__',)), and this call must fail unhandled to work
-            return object.__getattribute__(self, *args, **kwargs)
-        try:
-            return object.__getattribute__(self, *args, **kwargs)
-        except:
-            return self.contents[args[0]]
-
-    def __setattr__(self, *args, **kwargs):
-        _3to2list1 = list(args)
-        k, v, = _3to2list1[:1] + _3to2list1[1:]
-        if k in dir(self):  # in ['section', 'filename', 'lines']:
-            return object.__setattr__(self, *args, **kwargs)
-        if isinstance(v, str):
-            v = [fmt_value(v) for v in v.split()]
-        if isinstance(v,HTCContents):
-            self.contents[k] = v
-        else:
-            if not isinstance(v, (list, tuple)):
-                v = [v]
-            self.contents[k] = HTCLine(k, v, "")
-
-    def __delattr__(self, *args, **kwargs):
-        k, = args
-        if k in self:
-            del self.contents[k]
-
-    def __iter__(self):
-        return iter(self.contents.values())
-
-
-    def __contains__(self, key):
-        if self.contents is None:
-            return False
-        return key in self.contents
-
-    def get(self, section, default=None):
-        try:
-            return self[section]
-        except KeyError:
-            return default
-
-    def keys(self):
-        return list(self.contents.keys())
-
-    def _add_contents(self, contents):
-        if contents.name_ not in self:
-            self[contents.name_] = contents
-        else:
-            ending = "__2"
-            while contents.name_ + ending in self:
-                ending = "__%d" % (1 + float("0%s" % ending.replace("__", "")))
-            self[contents.name_ + ending] = contents
-
-    def add_section(self, name, allow_duplicate=False):
-        if name in self and allow_duplicate is False:
-            return self[name]
-        section = HTCSection(name)
-        self._add_contents(section)
-        return section
-
-    def add_line(self, name, values, comments=""):
-        line = HTCLine(name, values, comments)
-        self._add_contents(line)
-        return line
-
-    def location(self):
-        if self.parent is None:
-            return os.path.basename(self.filename)
-        else:
-            return self.parent.location() + "/" + self.name_
-
-
-class HTCSection(HTCContents):
-    end_comments = ""
-    begin_comments = ""
-    def __init__(self, name, begin_comments="", end_comments=""):
-        self.name_ = name
-        self.begin_comments = begin_comments.strip(" \t")
-        self.end_comments = end_comments.strip(" \t")
-        self.contents = OrderedDict()
-
-    @staticmethod
-    def from_lines(lines):
-        line, begin_comments = parse_next_line(lines)
-        name = line[6:].lower()
-        if name == "output":
-            section = HTCOutputSection(name, begin_comments)
-        elif name.startswith("output_at_time"):
-            section = HTCOutputAtTimeSection(name, begin_comments)
-        else:
-            section = HTCSection(name, begin_comments)
-        while lines:
-            if lines[0].lower().startswith("begin"):
-                section._add_contents(HTCSection.from_lines(lines))
-            elif lines[0].lower().startswith("end"):
-                line, section.end_comments = parse_next_line(lines)
-                break
-            else:
-                section._add_contents(section.line_from_line(lines))
-        return section
-
-    def line_from_line(self, lines):
-        return HTCLine.from_lines(lines)
-
-    def __str__(self, level=0):
-        s = "%sbegin %s;%s\n" % ("  "*level, self.name_, (("", "\t" + self.begin_comments)[bool(self.begin_comments.strip())]).replace("\t\n","\n"))
-        s += "".join([c.__str__(level + 1) for c in self])
-        s += "%send %s;%s\n" % ("  "*level, self.name_, (("", "\t" + self.end_comments)[self.end_comments.strip() != ""]).replace("\t\n","\n"))
-        return s
-    
-    def get_subsection_by_name(self, name, field='name'):
-        lst = [s for s in self if field in s and s[field][0]==name]
-        if len(lst)==1:
-            return lst[0]
-        else:
-            if len(lst)==0:
-                raise ValueError("subsection '%s' not found"%name)
-            else:
-                raise NotImplementedError()
-
-class HTCLine(HTCContents):
-    values = None
-    comments = ""
-    def __init__(self, name, values, comments):
-        if "__" in name:
-            name = name[:name.index("__")]
-        self.name_ = name
-        self.values = list(values)
-        self.comments = comments.strip(" \t")
-
-    def __repr__(self):
-        return str(self)
-
-    def __str__(self, level=0):
-        if self.name_ == "":
-            return ""
-        return "%s%s%s;%s\n" % ("  "*(level), self.name_,
-                                  ("", "\t" + self.str_values())[bool(self.values)],
-                                  ("", "\t" + self.comments)[bool(self.comments.strip())])
-    def str_values(self):
-        return " ".join([str(v).lower() for v in self.values])
-
-    def __getitem__(self, key):
-        try:
-            return self.values[key]
-        except:
-            raise IndexError("Parameter %s does not exists for %s"%(key+1,self.location()))
-
-    @staticmethod
-    def from_lines(lines):
-        line, end_comments = parse_next_line(lines)
-        if len(line.split()) > 0:
-            _3to2list3 = list(line.split())
-            name, values, = _3to2list3[:1] + [_3to2list3[1:]]
-        else:
-            name = line
-            values = []
-
-        values = [fmt_value(v) for v in values]
-        return HTCLine(name, values, end_comments)
-
-
-    def remove(self):
-        self.name_ = ""
-        self.values = []
-        self.comments = ""
-
-
-class HTCOutputSection(HTCSection):
-    sensors = None
-    def __init__(self, name, begin_comments="", end_comments=""):
-        HTCSection.__init__(self, name, begin_comments=begin_comments, end_comments=end_comments)
-        self.sensors = []
-
-    def add_sensor(self, type, sensor, values=[], comment="", nr=None):
-        self._add_sensor(HTCSensor(type, sensor, values, comment), nr)
-
-
-    def _add_sensor(self, htcSensor, nr=None):
-        if nr is None:
-            nr = len(self.sensors)
-        self.sensors.insert(nr, htcSensor)
-
-
-    def line_from_line(self, lines):
-        while len(lines) and lines[0].strip() == "":
-            lines.pop(0)
-        name = lines[0].split()[0].strip()
-        if name in ['filename', 'data_format', 'buffer', 'time']:
-            return HTCLine.from_lines(lines)
-        else:
-            return HTCSensor.from_lines(lines)
-
-    def _add_contents(self, contents):
-        if isinstance(contents, HTCSensor):
-            self._add_sensor(contents)
-        else:
-            return HTCSection._add_contents(self, contents)
-
-
-    def __str__(self, level=0):
-        s = "%sbegin %s;%s\n" % ("  "*level, self.name_, ("", "\t" + self.begin_comments)[len(self.begin_comments.strip()) > 0])
-        s += "".join([c.__str__(level + 1) for c in self])
-        s += "".join([s.__str__(level + 1) for s in self.sensors])
-        s += "%send %s;%s\n" % ("  "*level, self.name_, ("", "\t" + self.end_comments)[self.end_comments.strip() != ""])
-        return s
-
-class HTCOutputAtTimeSection(HTCOutputSection):
-    type = None
-    time = None
-    def __init__(self, name, begin_comments="", end_comments=""):
-        if len(name.split()) < 3:
-            raise ValueError('"keyword" and "time" arguments required for output_at_time command:\n%s' % name)
-        name, self.type, time = name.split()
-        self.time = float(time)
-        HTCOutputSection.__init__(self, name, begin_comments=begin_comments, end_comments=end_comments)
-
-    def __str__(self, level=0):
-        s = "%sbegin %s %s %s;%s\n" % ("  "*level, self.name_, self.type, self.time, ("", "\t" + self.begin_comments)[len(self.begin_comments.strip())])
-        s += "".join([c.__str__(level + 1) for c in self])
-        s += "".join([s.__str__(level + 1) for s in self.sensors])
-        s += "%send %s;%s\n" % ("  "*level, self.name_, ("", "\t" + self.end_comments)[self.end_comments.strip() != ""])
-        return s
-
-class HTCSensor(HTCLine):
-    type = ""
-    sensor = ""
-    values = []
-
-    def __init__(self, type, sensor, values, comments):
-        self.type = type
-        self.sensor = sensor
-        self.values = values
-        self.comments = comments.strip(" \t")
-
-    @staticmethod
-    def from_lines(lines):
-        line, comments = parse_next_line(lines)
-        if len(line.split()) > 2:
-            _3to2list5 = list(line.split())
-            type, sensor, values, = _3to2list5[:2] + [_3to2list5[2:]]
-        elif len(line.split()) == 2:
-            type, sensor = line.split()
-            values = []
-        else:
-            type, sensor, values = "", "", []
-        def fmt(v):
-            try:
-                if int(float(v)) == float(v):
-                    return int(float(v))
-                return float(v)
-            except ValueError:
-                return v
-        values = [fmt(v) for v in values]
-        return HTCSensor(type, sensor, values, comments)
-
-    def __str__(self, level=0):
-        return "%s%s %s%s;%s\n" % ("  "*(level),
-                                self.type,
-                                self.sensor,
-                                ("", "\t" + self.str_values())[bool(self.values)],
-                                ("", "\t" + self.comments)[bool(self.comments.strip())])
-
-
+'''
+Created on 20/01/2014
+
+@author: MMPE
+
+See documentation of HTCFile below
+
+'''
+import os
+from collections import OrderedDict
+import collections
+
+
+class OrderedDict(collections.OrderedDict):
+    pass
+
+    def __str__(self):
+        return "\n".join(["%-30s\t %s" % ((str(k) + ":"), str(v)) for k, v in self.items()])
+
+
+def parse_next_line(lines):
+    _3to2list = list(lines.pop(0).split(";"))
+    line, comments, = _3to2list[:1] + [_3to2list[1:]]
+    comments = ";".join(comments).rstrip()
+    while lines and lines[0].lstrip().startswith(";"):
+        comments += "\n%s" % lines.pop(0).rstrip()
+    return line.strip(), comments
+
+
+def fmt_value(v):
+    try:
+        if int(float(v)) == float(v):
+            return int(float(v))
+        return float(v)
+    except ValueError:
+        return v.replace("\\", "/")
+
+
+class HTCContents(object):
+    lines = []
+    contents = None
+    name_ = ""
+    parent = None
+
+    def __getitem__(self, key):
+        if isinstance(key, str):
+            key = key.replace(".", "/")
+            if "/" in key:
+                keys = key.split('/')
+                val = self.contents[keys[0]]
+                for k in keys[1:]:
+                    val = val[k]
+                return val
+            return self.contents[key]
+        else:
+            return self.values[key]
+
+    def __getattr__(self, *args, **kwargs):
+        if args[0] in ['__members__', '__methods__']:
+            # fix python2 related issue. In py2, dir(self) calls
+            # __getattr__(('__members__',)), and this call must fail unhandled to work
+            return object.__getattribute__(self, *args, **kwargs)
+        try:
+            return object.__getattribute__(self, *args, **kwargs)
+        except Exception:
+            k = args[0]
+            if k.endswith("__1"):
+                k = k[:-3]
+            return self.contents[k]
+
+    def __setattr__(self, *args, **kwargs):
+        try:
+            object.__getattribute__(self, args[0])
+            return object.__setattr__(self, *args, **kwargs)
+        except AttributeError:
+            pass
+        _3to2list1 = list(args)
+        k, v, = _3to2list1[:1] + _3to2list1[1:]
+
+        if isinstance(v, str):
+            v = [fmt_value(v) for v in v.split()]
+        if not isinstance(v, HTCContents):
+            if not isinstance(v, (list, tuple)):
+                v = [v]
+            if k in self.contents:
+                self.contents[k].values = list(v)
+                return
+            v = HTCLine(k, v, "")
+        self.contents[k] = v
+        v.parent = self
+
+    def __delattr__(self, *args, **kwargs):
+        k, = args
+        if k in self:
+            del self.contents[k]
+
+    def __iter__(self):
+        # mainbodies must preceed constraints
+        values = ([v for v in self.contents.values() if v.name_ not in ['orientation', 'constraint']] +
+                  [v for v in self.contents.values() if v.name_ == 'orientation'] +
+                  [v for v in self.contents.values() if v.name_ == 'constraint'])
+        return iter(values)
+
+    def __contains__(self, key):
+        if self.contents is None:
+            return False
+        return key in self.contents
+
+    def get(self, section, default=None):
+        try:
+            return self[section]
+        except KeyError:
+            return default
+
+    def __call__(self, **kwargs):
+        """Allow accesing one of multiple subsections with same name, e.g. the main body where name=='shaft'
+        > htc.new_htc_structure.main_body(name='shaft')
+
+        or one of multiple lines with same name, e.g. the section in c2_def where value[0]==3
+        > htc.new_htc_structure.main_body.c2_def.sec(v0=3)
+        """
+        lst = [s for s in self.parent if s.name_ == self.name_ and (
+            all([k in s and s[k][0] == v for k, v in kwargs.items()]) or
+            (all([k[0] == 'v' for k in kwargs]) and all([s[int(k[1:])] == v for k, v in kwargs.items()]))
+        )]
+        assert len(lst) == 1
+        return lst[0]
+
+    def keys(self):
+        return list(self.contents.keys())
+
+    def _add_contents(self, contents):
+        if contents.name_ not in self:
+            self[contents.name_] = contents
+        else:
+            ending = "__2"
+            while contents.name_ + ending in self:
+                ending = "__%d" % (1 + float("0%s" % ending.replace("__", "")))
+            self[contents.name_ + ending] = contents
+        contents.parent = self
+
+    def add_section(self, section_name, members={}, section=None, allow_duplicate=False, **kwargs):
+        if isinstance(section_name, HTCSection):
+            section = section_name
+            section_name = section.name_
+        if section_name in self and allow_duplicate is False:
+            return self[section_name]
+        if section_name == "output":
+            section = HTCOutputSection(section_name)
+        elif section_name.startswith("output_at_time"):
+            section = HTCOutputAtTimeSection(section_name)
+        elif section is None:
+            section = HTCSection(section_name)
+        self._add_contents(section)
+        kwargs.update(members)
+        for k, v in kwargs.items():
+            section[k] = v
+        return section
+
+    def delete(self):
+        keys = [k for (k, v) in self.parent.contents.items() if v == self]
+        for k in keys:
+            del self.parent.contents[k]
+
+    def location(self):
+        if self.parent is None:
+            return os.path.basename(self.filename)
+        else:
+            name = [k for k in self.parent.keys() if self.parent[k] == self][0]
+            return self.parent.location() + "/" + name
+
+    def compare(self, other, compare_order=False):
+        my_keys = self.keys()
+        other_keys = other.keys()
+        s = ""
+        while my_keys or other_keys:
+            if my_keys:
+                if (my_keys[0] in other_keys):
+                    if compare_order:
+                        other_i = 0
+                    else:
+                        other_i = other_keys.index(my_keys[0])
+                    while other_keys[other_i] != my_keys[0]:
+                        s += "\n".join(["+ %s" % l for l in str(other[other_keys.pop(other_i)]
+                                                                ).strip().split("\n")]) + "\n\n"
+
+                    s += self[my_keys.pop(0)].compare(other[other_keys.pop(other_i)])
+
+                else:
+                    s += "\n".join(["- %s" % l for l in str(self[my_keys.pop(0)]).strip().split("\n")]) + "\n\n"
+            else:
+                s += "\n".join(["+ %s" % l for l in str(other[other_keys.pop(0)]).strip().split("\n")]) + "\n\n"
+        return s
+
+
+class HTCSection(HTCContents):
+    end_comments = ""
+    begin_comments = ""
+
+    def __init__(self, name, begin_comments="", end_comments=""):
+        self.name_ = name
+        self.begin_comments = begin_comments.strip(" \t")
+        self.end_comments = end_comments.strip(" \t")
+        self.contents = OrderedDict()
+        self.parent = None
+
+    @property
+    def section_name(self):
+        return self.name_
+
+    @section_name.setter
+    def section_name(self, value):
+        self.name_ = value
+
+    def add_line(self, name, values, comments=""):
+        line = HTCLine(name, values, comments)
+        self._add_contents(line)
+        return line
+
+    def __setitem__(self, key, value):
+        if isinstance(value, HTCContents):
+            self.contents[key] = value
+            value.parent = self
+        elif isinstance(value, (str, int, float)):
+            self.add_line(key, [value])
+        else:
+            self.add_line(key, value)
+
+    @staticmethod
+    def from_lines(lines):
+        line, begin_comments = parse_next_line(lines)
+        name = line[6:].lower().strip()
+        if name == "output":
+            section = HTCOutputSection(name, begin_comments)
+        elif name.startswith("output_at_time"):
+            section = HTCOutputAtTimeSection(name, begin_comments)
+        else:
+            section = HTCSection(name, begin_comments)
+        while lines:
+            if lines[0].strip() == "":
+                lines.pop(0)
+            if len(lines) == 0:
+                print()
+            if lines[0].lower().startswith("begin"):
+                section._add_contents(HTCSection.from_lines(lines))
+            elif lines[0].lower().startswith("end"):
+                line, section.end_comments = parse_next_line(lines)
+                break
+            elif lines:
+                section._add_contents(section.line_from_line(lines))
+        else:
+            raise Exception("Section '%s' has not end" % section.name_)
+        return section
+
+    def line_from_line(self, lines):
+        return HTCLine.from_lines(lines)
+
+    def __str__(self, level=0):
+        s = "%sbegin %s;%s\n" % ("  " * level, self.name_, (("", "\t" + self.begin_comments)
+                                                            [bool(self.begin_comments.strip())]).replace("\t\n", "\n"))
+        s += "".join([c.__str__(level + 1) for c in self])
+        s += "%send %s;%s\n" % ("  " * level, self.name_, (("", "\t" + self.end_comments)
+                                                           [self.end_comments.strip() != ""]).replace("\t\n", "\n"))
+        return s
+
+    def get_subsection_by_name(self, name, field='name'):
+        return self.get_section(name, field)
+
+    def get_section(self, name, field='name'):
+        lst = [s for s in self if field in s and s[field][0] == name]
+        if len(lst) == 1:
+            return lst[0]
+        elif len(lst) == 0:
+            raise ValueError("subsection with %s='%s' not found" % (field, name))
+        else:
+            raise ValueError("Multiple subsection with %s='%s' not found" % (field, name))
+
+    def get_element(self, key, value):
+        """Return subsection where subsection.<key>==value or line where line.values[key]==value"""
+        if isinstance(key, int):
+            lst = [s for s in self if s.values[key] == value]
+        elif isinstance(key, str):
+            lst = [s for s in self if key in s and s[key][0] == name]
+        else:
+            raise ValueError("Key argument must be int or str")
+        if len(lst) == 1:
+            return lst[0]
+        elif len(lst) == 0:
+            raise ValueError("contents with '%s=%s' not found" % (key, value))
+        else:
+            raise ValueError("Multiple contents with '%s=%s' not found" % (key, value))
+
+    def copy(self):
+        copy = self.__class__(name=self.name_, begin_comments=self.begin_comments, end_comments=self.end_comments)
+        for k, v in self.contents.items():
+            copy.contents[k] = v.copy()
+        return copy
+
+
+class HTCLine(HTCContents):
+    values = None
+    comments = ""
+
+    def __init__(self, name, values, comments):
+        if "__" in name:
+            name = name[:name.index("__")]
+        self.name_ = name
+        self.values = list(values)
+        self.comments = comments.strip(" \t")
+        self.parent = None
+
+    def __repr__(self):
+        return str(self)
+
+    def __str__(self, level=0):
+        if self.name_ == "":
+            return ""
+        return "%s%s%s;%s\n" % ("  " * (level), self.name_,
+                                ("", "\t" + self.str_values())[bool(self.values)],
+                                ("", "\t" + self.comments)[bool(self.comments.strip())])
+
+    def str_values(self):
+        return " ".join([str(v) for v in self.values])
+
+    def __getitem__(self, key):
+        try:
+            return self.values[key]
+        except Exception:
+            raise IndexError("Parameter %s does not exists for %s" % (key + 1, self.location()))
+
+    def __setitem__(self, key, value):
+        if isinstance(key, int):
+            self.values[key] = value
+        else:
+            raise NotImplementedError
+
+    @staticmethod
+    def from_lines(lines):
+        line, end_comments = parse_next_line(lines)
+        if len(line.split()) > 0:
+            _3to2list3 = list(line.split())
+            name, values, = _3to2list3[:1] + [_3to2list3[1:]]
+        else:
+            name = line
+            values = []
+
+        values = [fmt_value(v) for v in values]
+        return HTCLine(name, values, end_comments)
+
+    def compare(self, other):
+        s = ""
+        if self.values != other.values:
+            s += "\n".join(["+ %s" % l for l in str(self).strip().split("\n")]) + "\n"
+            s += "\n".join(["- %s" % l for l in str(other).strip().split("\n")]) + "\n"
+            s += "\n"
+        return s
+
+    def copy(self):
+        return HTCLine(name=self.name_, values=self.values, comments=self.comments)
+
+
+class HTCOutputSection(HTCSection):
+    sensors = None
+
+    def __init__(self, name, begin_comments="", end_comments=""):
+        HTCSection.__init__(self, name, begin_comments=begin_comments, end_comments=end_comments)
+        self.sensors = []
+
+    def add_sensor(self, type, sensor, values=[], comment="", nr=None):
+        self._add_sensor(HTCSensor(type, sensor, values, comment), nr)
+
+    def _add_sensor(self, htcSensor, nr=None):
+        if nr is None:
+            nr = len(self.sensors)
+        self.sensors.insert(nr, htcSensor)
+        htcSensor.parent = self
+
+    def line_from_line(self, lines):
+        name = lines[0].split()[0].strip()
+
+        if name in ['filename', 'data_format', 'buffer', 'time']:
+            return HTCLine.from_lines(lines)
+        else:
+            return HTCSensor.from_lines(lines)
+
+    def _add_contents(self, contents):
+        if isinstance(contents, HTCSensor):
+            self._add_sensor(contents)
+        else:
+            return HTCSection._add_contents(self, contents)
+
+    def __str__(self, level=0):
+        s = "%sbegin %s;%s\n" % ("  " * level, self.name_, ("", "\t" + self.begin_comments)
+                                 [len(self.begin_comments.strip()) > 0])
+        s += "".join([c.__str__(level + 1) for c in self])
+        s += "".join([s.__str__(level + 1) for s in self.sensors])
+        s += "%send %s;%s\n" % ("  " * level, self.name_, ("", "\t" + self.end_comments)
+                                [self.end_comments.strip() != ""])
+        return s
+
+    def compare(self, other):
+        s = HTCContents.compare(self, other)
+        for s1, s2 in zip(self.sensors, other.sensors):
+            s += s1.compare(s2)
+        for s1 in self.sensors[len(other.sensors):]:
+            s += "\n".join(["- %s" % l for l in str(s1).strip().split("\n")]) + "\n"
+        for s2 in self.sensors[len(self.sensors):]:
+            s += "\n".join(["- %s" % l for l in str(s2).strip().split("\n")]) + "\n"
+
+        return s
+
+    def copy(self):
+        copy = HTCSection.copy(self)
+        copy.sensors = list.copy(self.sensors)
+        return copy
+
+
+class HTCOutputAtTimeSection(HTCOutputSection):
+    type = None
+    time = None
+
+    def __init__(self, name, begin_comments="", end_comments=""):
+        if len(name.split()) < 3:
+            raise ValueError('"keyword" and "time" arguments required for output_at_time command:\n%s' % name)
+        name, self.type, time = name.split()
+        self.time = float(time)
+        HTCOutputSection.__init__(self, name, begin_comments=begin_comments, end_comments=end_comments)
+
+    def __str__(self, level=0):
+        s = "%sbegin %s %s %s;%s\n" % ("  " * level, self.name_, self.type, self.time,
+                                       ("", "\t" + self.begin_comments)[len(self.begin_comments.strip())])
+        s += "".join([c.__str__(level + 1) for c in self])
+        s += "".join([s.__str__(level + 1) for s in self.sensors])
+        s += "%send %s;%s\n" % ("  " * level, self.name_, ("", "\t" + self.end_comments)
+                                [self.end_comments.strip() != ""])
+        return s
+
+    def copy(self):
+        copy = HTCOutputAtTimeSection(name=f"{self.name_} {self.type} {self.time}",
+                                           begin_comments=self.begin_comments, end_comments=self.end_comments)
+        for k, v in self.contents.items():
+            copy.contents[k] = v.copy()
+        copy.sensors = list.copy(self.sensors)
+        return copy
+
+
+class HTCSensor(HTCLine):
+    type = ""
+    sensor = ""
+    values = []
+
+    def __init__(self, type, sensor, values, comments):
+        self.type = type
+        self.sensor = sensor
+        self.values = list(values)
+        self.comments = comments.strip(" \t")
+
+    @staticmethod
+    def from_lines(lines):
+        line, comments = parse_next_line(lines)
+        if len(line.split()) > 2:
+            _3to2list5 = list(line.split())
+            type, sensor, values, = _3to2list5[:2] + [_3to2list5[2:]]
+        elif len(line.split()) == 2:
+            type, sensor = line.split()
+            values = []
+        else:
+            type, sensor, values = "", "", []
+
+        def fmt(v):
+            try:
+                if int(float(v)) == float(v):
+                    return int(float(v))
+                return float(v)
+            except ValueError:
+                return v
+        values = [fmt(v) for v in values]
+        return HTCSensor(type, sensor, values, comments)
+
+    def __str__(self, level=0):
+        return "%s%s %s%s;%s\n" % ("  " * (level),
+                                   self.type,
+                                   self.sensor,
+                                   ("", "\t" + self.str_values())[bool(self.values)],
+                                   ("", "\t" + self.comments)[bool(self.comments.strip())])
+
+    def delete(self):
+        self.parent.sensors.remove(self)
+
+    def compare(self, other):
+        s = ""
+        if self.sensor != other.sensor or self.values != other.values:
+            s += "\n".join(["+ %s" % l for l in str(self).strip().split("\n")]) + "\n"
+            s += "\n".join(["- %s" % l for l in str(other).strip().split("\n")]) + "\n"
+            s += "\n"
+        return s

wetb/hawc2/htc_extensions.py

@@ -6,45 +6,20 @@ Created on 20/01/2014
 See documentation of HTCFile below
 
 '''
-from __future__ import division
-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import absolute_import
-from builtins import zip
-from builtins import int
-from builtins import str
-from future import standard_library
 import os
 
 
-standard_library.install_aliases()
-
-
-
 class HTCDefaults(object):
 
-
     empty_htc = """begin simulation;
         time_stop 600;
-        solvertype    1;    (newmark)
-        on_no_convergence continue;
-        convergence_limits 1E3 1.0 1E-7; ; . to run again, changed 07/11
+        solvertype    2;    (newmark)
         begin newmark;
           deltat    0.02;
         end newmark;
+        log_deltat 1;
     end simulation;
     ;
-    ;----------------------------------------------------------------------------------------------------------------------------------------------------------------
-    ;
-    begin new_htc_structure;
-      begin orientation;
-      end orientation;
-      begin constraint;
-      end constraint;
-    end new_htc_structure;
-    ;
-    ;----------------------------------------------------------------------------------------------------------------------------------------------------------------
-    ;
     begin wind ;
       density                 1.225 ;
       wsp                     10   ;
@@ -57,20 +32,14 @@ class HTCDefaults(object):
       tower_shadow_method     0     ;  0=none, 1=potential flow, 2=jet
     end wind;
     ;
-    ;----------------------------------------------------------------------------------------------------------------------------------------------------------------
-    ;
-    begin dll;
-    end dll;
-    ;
-    ;----------------------------------------------------------------------------------------------------------------------------------------------------------------
-    ;
     begin output;
+      filename    ./tmp;
+      data_format gtsdf;
       general time;
     end output;
     exit;"""
 
-
-    def add_mann_turbulence(self, L=29.4, ae23=1, Gamma=3.9, seed=1001, high_frq_compensation=True,
+    def add_mann_turbulence(self, L=33.6, ae23=1, Gamma=3.9, seed=1001, high_frq_compensation=True,
                             filenames=None,
                             no_grid_points=(16384, 32, 32), box_dimension=(6000, 100, 100),
                             dont_scale=False,
@@ -81,19 +50,22 @@ class HTCDefaults(object):
         if 'create_turb_parameters' in mann:
             mann.create_turb_parameters.values = [L, ae23, Gamma, seed, int(high_frq_compensation)]
         else:
-            mann.add_line('create_turb_parameters', [L, ae23, Gamma, seed, int(high_frq_compensation)], "L, alfaeps, gamma, seed, highfrq compensation")
+            mann.add_line('create_turb_parameters', [L, ae23, Gamma, seed, int(high_frq_compensation)],
+                          "L, alfaeps, gamma, seed, highfrq compensation")
         if filenames is None:
-            
+
             import numpy as np
-            dxyz = tuple(np.array(box_dimension) / no_grid_points)
+            dxyz = tuple(np.array(box_dimension) / np.array(no_grid_points))
             from wetb.wind.turbulence import mann_turbulence
-            filenames = ["./turb/" + mann_turbulence.name_format % ((L, ae23, Gamma, high_frq_compensation) + no_grid_points + dxyz + (seed, uvw)) for uvw in ['u', 'v', 'w']]
+            filenames = ["./turb/" + mann_turbulence.name_format %
+                         ((L, ae23, Gamma, high_frq_compensation) + no_grid_points + dxyz + (seed, uvw))
+                         for uvw in ['u', 'v', 'w']]
         if isinstance(filenames, str):
             filenames = ["./turb/%s_s%04d%s.bin" % (filenames, seed, c) for c in ['u', 'v', 'w']]
         for filename, c in zip(filenames, ['u', 'v', 'w']):
             setattr(mann, 'filename_%s' % c, filename)
         for c, n, dim in zip(['u', 'v', 'w'], no_grid_points, box_dimension):
-            setattr(mann, 'box_dim_%s' % c, "%d %.4f" % (n, dim / (n - 1)))
+            setattr(mann, 'box_dim_%s' % c, "%d %.4f" % (n, dim / n))
         if dont_scale:
             mann.dont_scale = 1
         else:
@@ -108,28 +80,24 @@ class HTCDefaults(object):
                 del mann.std_scaling
             except KeyError:
                 pass
-            
 
-
-    def add_turb_export(self, filename="export_%s.turb", samplefrq = None):
+    def add_turb_export(self, filename="export_%s.turb", samplefrq=None):
         exp = self.wind.add_section('turb_export', allow_duplicate=True)
         for uvw in 'uvw':
-            exp.add_line('filename_%s'%uvw, [filename%uvw])
-        sf = samplefrq or max(1,int( self.wind.mann.box_dim_u[1]/(self.wind.wsp[0] * self.deltat())))
+            exp.add_line('filename_%s' % uvw, [filename % uvw])
+        sf = samplefrq or max(1, int(self.wind.mann.box_dim_u[1] / (self.wind.wsp[0] * self.deltat())))
         exp.samplefrq = sf
         if "time" in self.output:
             exp.time_start = self.output.time[0]
         else:
             exp.time_start = 0
-        exp.nsteps = (self.simulation.time_stop[0]-exp.time_start[0]) / self.deltat()
+        exp.nsteps = (self.simulation.time_stop[0] - exp.time_start[0]) / self.deltat()
         for vw in 'vw':
-            exp.add_line('box_dim_%s'%vw, self.wind.mann['box_dim_%s'%vw].values)
-
-        
+            exp.add_line('box_dim_%s' % vw, self.wind.mann['box_dim_%s' % vw].values)
 
     def import_dtu_we_controller_input(self, filename):
         dtu_we_controller = [dll for dll in self.dll if dll.name[0] == 'dtu_we_controller'][0]
-        with open (filename) as fid:
+        with open(filename) as fid:
             lines = fid.readlines()
         K_r1 = float(lines[1].replace("K = ", '').replace("[Nm/(rad/s)^2]", ''))
         Kp_r2 = float(lines[4].replace("Kp = ", '').replace("[Nm/(rad/s)]", ''))
@@ -148,16 +116,23 @@ class HTCDefaults(object):
         cs.constant__21.values[1] = "%.6E" % KK1
         cs.constant__22.values[1] = "%.6E" % KK2
 
+    def add_hydro(self, mudlevel, mwl, gravity=9.81, rho=1027):
+        wp = self.add_section("hydro").add_section('water_properties')
+        wp.mudlevel = mudlevel
+        wp.mwl = mwl
+        wp.gravity = gravity
+        wp.rho = rho
+
 
 class HTCExtensions(object):
     def get_shear(self):
         shear_type, parameter = self.wind.shear_format.values
         z0 = -self.wind.center_pos0[2]
         wsp = self.wind.wsp[0]
-        if shear_type==1: #constant
-            return lambda z : wsp
-        elif shear_type==3:
+        if shear_type == 1:  # constant
+            return lambda z: wsp
+        elif shear_type == 3:
             from wetb.wind.shear import power_shear
             return power_shear(parameter, z0, wsp)
         else:
-            raise NotImplementedError
\ No newline at end of file
+            raise NotImplementedError

wetb/hawc2/htc_file.py

-'''
-Created on 20/01/2014
-
-@author: MMPE
-
-See documentation of HTCFile below
-
-'''
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import division
-from __future__ import absolute_import
-from io import open
-from builtins import str
-from future import standard_library
-from wetb.utils.process_exec import pexec
-from wetb.utils.cluster_tools.cluster_resource import unix_path_old
-standard_library.install_aliases()
-from collections import OrderedDict
-
-from wetb.hawc2.htc_contents import HTCContents, HTCSection, HTCLine
-from wetb.hawc2.htc_extensions import HTCDefaults, HTCExtensions
-import os
-from copy import copy
-
-
-def fmt_path(path):
-    return path.lower().replace("\\","/")
-
-class HTCFile(HTCContents, HTCDefaults, HTCExtensions):
-    """Wrapper for HTC files
-
-    Examples:
-    ---------
-    >>> htcfile = HTCFile('htc/test.htc')
-    >>> htcfile.wind.wsp = 10
-    >>> htcfile.save()
-    
-    #---------------------------------------------
-    >>> htc = HTCFile(filename=None, modelpath=None) # create minimal htcfile
-        
-    #Add section
-    >>> htc.add_section('hydro')
-    
-    #Add subsection
-    >>> htc.hydro.add_section("hydro_element")
-    
-    #Set values
-    >>> htc.hydro.hydro_element.wave_breaking = [2, 6.28, 1] # or
-    >>> htc.hydro.hydro_element.wave_breaking = 2, 6.28, 1
-    
-    #Set comments
-    >>> htc.hydro.hydro_element.wave_breaking.comments = "This is a comment"
-    
-    #Access section
-    >>> hydro_element = htc.hydro.hydro_element #or
-    >>> hydro_element = htc['hydro.hydro_element'] # or
-    >>> hydro_element = htc['hydro/hydro_element'] # or
-    >>> print (hydro_element.wave_breaking) #string represenation
-    wave_breaking    2 6.28 1;    This is a comment
-    >>> print (hydro_element.wave_breaking.name_) # command
-    wave_breaking
-    >>> print (hydro_element.wave_breaking.values) # values
-    [2, 6.28, 1
-    >>> print (hydro_element.wave_breaking.comments) # comments
-    This is a comment
-    >>> print (hydro_element.wave_breaking[0]) # first value
-    2
-    
-    #Delete element
-    del htc.simulation.logfile #Delete logfile line. Raise keyerror if not exists
-    """
-
-    filename = None
-    htc_inputfiles = []
-    level = 0
-    modelpath = "../"
-    initial_comments = None
-    _contents = None
-    def __init__(self, filename=None, modelpath=None):
-        """        
-        Parameters
-        ---------
-        filename : str
-            Absolute filename of htc file
-        modelpath : str
-            Model path relative to htc file 
-        """
-        
-        if filename is not None:
-            self.filename = filename
-        self.modelpath = modelpath or self.auto_detect_modelpath()
-        if filename and self.modelpath!="unknown" and not os.path.isabs(self.modelpath):
-            self.modelpath = os.path.realpath(os.path.join(os.path.dirname(self.filename), self.modelpath))
-        
-    
-                #assert 'simulation' in self.contents, "%s could not be loaded. 'simulation' section missing" % filename
-
-    def auto_detect_modelpath(self):
-        if self.filename is None:
-            return "../"
-        
-        #print (["../"*i for i in range(3)])
-        import numpy as np
-        input_files = HTCFile(self.filename, 'unknown').input_files()
-        rel_input_files = [f for f in input_files if not os.path.isabs(f)] 
-        found = ([np.sum([os.path.isfile(os.path.join(os.path.dirname(self.filename), "../"*i, f)) for f in rel_input_files]) for i in range(4)])
-        #for f in self.input_files():
-        #    print (os.path.isfile(os.path.join(os.path.dirname(self.filename), "../",f)), f)
-        if max(found)>0:
-            relpath = "../"* np.argmax(found)
-            return os.path.abspath(os.path.join(os.path.dirname(self.filename), relpath))
-        else:
-            raise ValueError("Modelpath cannot be autodetected for '%s'.\nInput files not found near htc file"%self.filename)
-        
-    def _load(self):
-        self.reset()
-        self.initial_comments = []
-        self.htc_inputfiles = []
-        self.contents = OrderedDict()
-        if self.filename is None:
-            lines = self.empty_htc.split("\n")
-        else:
-            lines = self.readlines(self.filename)
-
-        lines = [l.strip() for l in lines]
-
-        #lines = copy(self.lines)
-        while lines:
-            if lines[0].startswith(";"):
-                self.initial_comments.append(lines.pop(0).strip() + "\n")
-            elif lines[0].lower().startswith("begin"):
-                self._add_contents(HTCSection.from_lines(lines))
-            else:
-                line = HTCLine.from_lines(lines)
-                if line.name_ == "exit":
-                    break
-                self._add_contents(line)
-
-
-    def reset(self):
-        self._contents = None
-
-        
-    @property
-    def contents(self):
-        if self._contents is None:
-            self._load()
-        return self._contents
-    
-    @contents.setter
-    def contents(self, value):
-        self._contents = value
-        
-            
-    def readfilelines(self, filename):
-        with open(unix_path_old(filename), encoding='cp1252') as fid:
-            lines = list(fid.readlines())
-        if lines[0].encode().startswith(b'\xc3\xaf\xc2\xbb\xc2\xbf'):
-            lines[0] = lines[0][3:]
-        return lines
-
-    def readlines(self, filename):
-        self.htc_inputfiles.append(filename)
-        htc_lines = []
-        lines = self.readfilelines(filename)
-        for l in lines:
-            if l.lower().lstrip().startswith('continue_in_file'):
-                filename = l.lstrip().split(";")[0][len("continue_in_file"):].strip().lower()
-                
-                if self.modelpath=='unknown':
-                    self.htc_inputfiles.append(filename)
-                else:
-                    filename = os.path.join(self.modelpath, filename)
-                    for line in self.readlines(filename):
-                        if line.lstrip().lower().startswith('exit'):
-                            break
-                        htc_lines.append(line)
-            else:
-                htc_lines.append(l)
-        return htc_lines
-
-
-    def __setitem__(self, key, value):
-        self.contents[key] = value
-
-    def __str__(self):
-        self.contents #load
-        return "".join(self.initial_comments + [c.__str__(1) for c in self]+ ["exit;"])
-
-    def save(self, filename=None):
-        self.contents #load if not loaded
-        if filename is None:
-            filename = self.filename
-        else:
-            self.filename = filename
-        # exist_ok does not exist in Python27
-        if not os.path.exists(os.path.dirname(filename)):
-            os.makedirs(os.path.dirname(filename))  #, exist_ok=True)
-        with open(filename, 'w', encoding='cp1252') as fid:
-            fid.write(str(self))
-
-    def set_name(self, name, subfolder=''):
-        #if os.path.isabs(folder) is False and os.path.relpath(folder).startswith("htc" + os.path.sep):
-        self.contents #load if not loaded
-        fmt_folder = lambda folder, subfolder : "./" + os.path.relpath(os.path.join(folder, subfolder)).replace("\\", "/")
-        
-        self.filename = os.path.abspath(os.path.join(self.modelpath, fmt_folder('htc', subfolder), "%s.htc" % name)).replace("\\", "/")
-        if 'simulation' in self and 'logfile' in self.simulation:
-            self.simulation.logfile = os.path.join(fmt_folder('log', subfolder), "%s.log" % name).replace("\\", "/")
-            if 'animation' in self.simulation:
-                self.simulation.animation = os.path.join(fmt_folder('animation', subfolder), "%s.dat" % name).replace("\\", "/")
-            if 'visualization' in self.simulation:
-                self.simulation.visualization = os.path.join(fmt_folder('visualization', subfolder), "%s.hdf5" % name).replace("\\", "/")
-        elif 'test_structure' in self and 'logfile' in self.test_structure:  # hawc2aero
-            self.test_structure.logfile = os.path.join(fmt_folder('log', subfolder), "%s.log" % name).replace("\\", "/")
-        self.output.filename = os.path.join(fmt_folder('res', subfolder), "%s" % name).replace("\\", "/")
-
-    def set_time(self, start=None, stop=None, step=None):
-        self.contents # load if not loaded
-        if stop is not None:
-            self.simulation.time_stop = stop
-        else:
-            stop = self.simulation.time_stop[0]
-        if step is not None:
-            self.simulation.newmark.deltat = step
-        if start is not None:
-            self.output.time = start, stop
-            if "wind" in self:# and self.wind.turb_format[0] > 0:
-                self.wind.scale_time_start = start
-
-    def input_files(self):
-        self.contents # load if not loaded
-        if self.modelpath=="unknown":
-            files = [f.replace("\\","/") for f in self.htc_inputfiles]
-        else:
-            files = [os.path.abspath(f).replace("\\","/") for f in self.htc_inputfiles]
-        if 'new_htc_structure' in self:
-            for mb in [self.new_htc_structure[mb] for mb in self.new_htc_structure.keys() if mb.startswith('main_body')]:
-                if "timoschenko_input" in mb:
-                    files.append(mb.timoschenko_input.filename[0])
-                files.append(mb.get('external_bladedata_dll', [None, None, None])[2])
-        if 'aero' in self:
-            files.append(self.aero.ae_filename[0])
-            files.append(self.aero.pc_filename[0])
-            files.append(self.aero.get('external_bladedata_dll', [None, None, None])[2])
-            files.append(self.aero.get('output_profile_coef_filename', [None])[0])
-            if 'dynstall_ateflap' in self.aero:
-                files.append(self.aero.dynstall_ateflap.get('flap', [None] * 3)[2])
-            if 'bemwake_method' in self.aero:
-                files.append(self.aero.bemwake_method.get('a-ct-filename', [None] * 3)[0])
-        for dll in [self.dll[dll] for dll in self.get('dll', {}).keys() if 'filename' in self.dll[dll]]:
-            files.append(dll.filename[0])
-        if 'wind' in self:
-            files.append(self.wind.get('user_defined_shear', [None])[0])
-            files.append(self.wind.get('user_defined_shear_turbulence', [None])[0])
-        if 'wakes' in self:
-            files.append(self.wind.get('use_specific_deficit_file', [None])[0])
-            files.append(self.wind.get('write_ct_cq_file', [None])[0])
-            files.append(self.wind.get('write_final_deficits', [None])[0])
-        if 'hydro' in self:
-            if 'water_properties' in self.hydro:
-                files.append(self.hydro.water_properties.get('water_kinematics_dll', [None])[0])
-                files.append(self.hydro.water_properties.get('water_kinematics_dll', [None, None])[1])
-        if 'soil' in self:
-            if 'soil_element' in self.soil:
-                files.append(self.soil.soil_element.get('datafile', [None])[0])
-        try:
-            dtu_we_controller = self.dll.get_subsection_by_name('dtu_we_controller')
-            theta_min = dtu_we_controller.init.constant__5[1]
-            files.append(os.path.join(os.path.dirname(dtu_we_controller.filename[0]), "wpdata.%d"%theta_min).replace("\\","/"))
-        except:
-            pass
-            
-    
-        try:
-            files.append(self.force.dll.dll[0])
-        except:
-            pass
-
-        return [f for f in set(files) if f]
-
-    def output_files(self):
-        self.contents # load if not loaded
-        files = []
-        for k, index in [('simulation/logfile', 0),
-                         ('simulation/animation', 0),
-                         ('simulation/visualization', 0),
-                         ('new_htc_structure/beam_output_file_name', 0),
-                         ('new_htc_structure/body_output_file_name', 0),
-                         ('new_htc_structure/struct_inertia_output_file_name', 0),
-                         ('new_htc_structure/body_eigenanalysis_file_name', 0),
-                         ('new_htc_structure/constraint_output_file_name', 0),
-                         ('wind/turb_export/filename_u', 0),
-                         ('wind/turb_export/filename_v', 0),
-                         ('wind/turb_export/filename_w', 0)]:
-            line = self.get(k)
-            if line:
-                files.append(line[index])
-        if 'new_htc_structure' in self:
-            if 'system_eigenanalysis' in self.new_htc_structure:
-                f = self.new_htc_structure.system_eigenanalysis[0]
-                files.append(f)
-                files.append(os.path.join(os.path.dirname(f), 'mode*.dat').replace("\\", "/"))
-            if 'structure_eigenanalysis_file_name' in self.new_htc_structure:
-                f = self.new_htc_structure.structure_eigenanalysis_file_name[0]
-                files.append(f)
-                files.append(os.path.join(os.path.dirname(f), 'mode*.dat').replace("\\", "/"))
-        files.extend(self.res_file_lst())
-
-        for key in [k for k in self.contents.keys() if k.startswith("output_at_time")]:
-            files.append(self[key]['filename'][0] + ".dat")
-        return [f for f in files if f]
-
-    def turbulence_files(self):
-        self.contents # load if not loaded
-        if 'wind' not in self.contents.keys() or self.wind.turb_format[0] == 0:
-            return []
-        elif self.wind.turb_format[0] == 1:
-            files = [self.get('wind.mann.filename_%s' % comp, [None])[0] for comp in ['u', 'v', 'w']]
-        elif self.wind.turb_format[0] == 2:
-            files = [self.get('wind.flex.filename_%s' % comp, [None])[0] for comp in ['u', 'v', 'w']]
-        return [f for f in files if f]
-
-
-    def res_file_lst(self):
-        self.contents # load if not loaded
-        res = []
-        for output in [self[k] for k in self.keys() if self[k].name_=="output"]:
-            dataformat = output.get('data_format', 'hawc_ascii')
-            res_filename = output.filename[0]
-            if dataformat[0] == "gtsdf" or dataformat[0] == "gtsdf64":
-                res.append(res_filename + ".hdf5")
-            elif dataformat[0] == "flex_int":
-                res.extend([res_filename + ".int", os.path.join(os.path.dirname(res_filename), 'sensor')])
-            else:
-                res.extend([res_filename + ".sel", res_filename + ".dat"])
-        return res
-
-
-    def simulate(self, exe, skip_if_up_to_date=False):
-        self.contents # load if not loaded
-        if skip_if_up_to_date:
-            from os.path import isfile, getmtime, isabs
-            res_file = os.path.join(self.modelpath, self.res_file_lst()[0])
-            htc_file = os.path.join(self.modelpath, self.filename)
-            if isabs(exe):
-                exe_file = exe
-            else:
-                exe_file = os.path.join(self.modelpath, exe)
-            #print (from_unix(getmtime(res_file)), from_unix(getmtime(htc_file)))
-            if (isfile(htc_file) and isfile(res_file) and isfile(exe_file) and
-                str(HTCFile(htc_file))==str(self) and
-                getmtime(res_file) > getmtime(htc_file) and getmtime(res_file) > getmtime(exe_file)):
-                if "".join(self.readfilelines(htc_file)) == str(self):
-                        return
-
-        self.save()
-        htcfile = os.path.relpath(self.filename, self.modelpath)
-        hawc2exe = exe
-        errorcode, stdout, stderr, cmd = pexec([hawc2exe, htcfile], self.modelpath)
-
-        if "logfile" in self.simulation:
-            with open(os.path.join(self.modelpath, self.simulation.logfile[0])) as fid:
-                log = fid.read()
-        else:
-            log = stderr
-
-
-        if errorcode or 'Elapsed time' not in log:
-            raise Exception (str(stdout) + str(stderr))
-        
-    def deltat(self):
-        return self.simulation.newmark.deltat[0]
-    
-
-#     
-#     def get_body(self, name):
-#         lst = [b for b in self.new_htc_structure if b.name_=="main_body" and b.name[0]==name]
-#         if len(lst)==1:
-#             return lst[0]
-#         else:
-#             if len(lst)==0:
-#                 raise ValueError("Body '%s' not found"%name)
-#             else:
-#                 raise NotImplementedError()
-#         
-
-class H2aeroHTCFile(HTCFile):
-    def __init__(self, filename=None, modelpath=None):
-        HTCFile.__init__(self, filename=filename, modelpath=modelpath)
-
-    @property
-    def simulation(self):
-        return self.test_structure
-
-    def set_time(self, start=None, stop=None, step=None):
-        if stop is not None:
-            self.test_structure.time_stop = stop
-        else:
-            stop = self.simulation.time_stop[0]
-        if step is not None:
-            self.test_structure.deltat = step
-        if start is not None:
-            self.output.time = start, stop
-            if "wind" in self and self.wind.turb_format[0] > 0:
-                self.wind.scale_time_start = start
-
-
-if "__main__" == __name__:
-    f = HTCFile(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\htc\DTU_10MW_RWT_power_curve.htc", "../")
-    f.save(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\htc\DTU_10MW_RWT_power_curve.htc")
-
-    f = HTCFile(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\htc\DTU_10MW_RWT.htc", "../")
-    f.save(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\htc\DTU_10MW_RWT.htc")
-
+'''
+Created on 20/01/2014
+
+See documentation of HTCFile below
+
+'''
+from wetb.utils.process_exec import pexec
+from wetb.hawc2.hawc2_pbs_file import HAWC2PBSFile
+import jinja2
+from wetb.utils.cluster_tools.os_path import fixcase, abspath, pjoin
+from collections import OrderedDict
+from wetb.hawc2.htc_contents import HTCContents, HTCSection, HTCLine
+from wetb.hawc2.htc_extensions import HTCDefaults, HTCExtensions
+import os
+
+
+def fmt_path(path):
+    return path.lower().replace("\\", "/")
+
+
+class HTCFile(HTCContents, HTCDefaults, HTCExtensions):
+    """Wrapper for HTC files
+
+    Examples:
+    ---------
+    >>> htcfile = HTCFile('htc/test.htc')
+    >>> htcfile.wind.wsp = 10
+    >>> htcfile.save()
+
+    #---------------------------------------------
+    >>> htc = HTCFile(filename=None, modelpath=None) # create minimal htcfile
+
+    #Add section
+    >>> htc.add_section('hydro')
+
+    #Add subsection
+    >>> htc.hydro.add_section("hydro_element")
+
+    #Set values
+    >>> htc.hydro.hydro_element.wave_breaking = [2, 6.28, 1] # or
+    >>> htc.hydro.hydro_element.wave_breaking = 2, 6.28, 1
+
+    #Set comments
+    >>> htc.hydro.hydro_element.wave_breaking.comments = "This is a comment"
+
+    #Access section
+    >>> hydro_element = htc.hydro.hydro_element #or
+    >>> hydro_element = htc['hydro.hydro_element'] # or
+    >>> hydro_element = htc['hydro/hydro_element'] # or
+    >>> print (hydro_element.wave_breaking) #string represenation
+    wave_breaking    2 6.28 1;    This is a comment
+    >>> print (hydro_element.wave_breaking.name_) # command
+    wave_breaking
+    >>> print (hydro_element.wave_breaking.values) # values
+    [2, 6.28, 1
+    >>> print (hydro_element.wave_breaking.comments) # comments
+    This is a comment
+    >>> print (hydro_element.wave_breaking[0]) # first value
+    2
+
+    #Delete element
+    htc.simulation.logfile.delete()
+    #or
+    del htc.simulation.logfile #Delete logfile line. Raise keyerror if not exists
+
+    """
+
+    filename = None
+    jinja_tags = {}
+    htc_inputfiles = []
+    level = 0
+    modelpath = "../"
+    initial_comments = None
+
+    def __init__(self, filename=None, modelpath=None, jinja_tags={}):
+        """
+        Parameters
+        ---------
+        filename : str
+            Absolute filename of htc file
+        modelpath : str
+            Model path relative to htc file
+        """
+        if filename is not None:
+            try:
+                filename = fixcase(abspath(filename))
+                with self.open(str(filename)):
+                    pass
+            except Exception:
+                pass
+
+            self.filename = filename
+
+        self.jinja_tags = jinja_tags
+        self.modelpath = modelpath or self.auto_detect_modelpath()
+
+        if filename and self.modelpath != "unknown" and not os.path.isabs(self.modelpath):
+            drive, p = os.path.splitdrive(os.path.join(os.path.dirname(str(self.filename)), self.modelpath))
+            self.modelpath = os.path.join(drive, os.path.splitdrive(os.path.realpath(p))[1]).replace("\\", "/")
+        if self.modelpath != 'unknown' and self.modelpath[-1] != '/':
+            self.modelpath += "/"
+
+        self.load()
+
+    def auto_detect_modelpath(self):
+        if self.filename is None:
+            return "../"
+
+        #print (["../"*i for i in range(3)])
+        import numpy as np
+        input_files = HTCFile(self.filename, 'unknown').input_files()
+        if len(input_files) == 1:  # only input file is the htc file
+            return "../"
+        rel_input_files = [f for f in input_files if not os.path.isabs(f)]
+
+        def isfile_case_insensitive(f):
+            try:
+                f = fixcase(f)  # raises exception if not existing
+                return os.path.isfile(f)
+            except IOError:
+                return False
+        found = ([np.sum([isfile_case_insensitive(os.path.join(os.path.dirname(self.filename), "../" * i, f))
+                          for f in rel_input_files]) for i in range(4)])
+
+        if max(found) > 0:
+            relpath = "../" * np.argmax(found)
+            return abspath(pjoin(os.path.dirname(self.filename), relpath))
+        else:
+            raise ValueError(
+                "Modelpath cannot be autodetected for '%s'.\nInput files not found near htc file" % self.filename)
+
+    def load(self):
+        self.contents = OrderedDict()
+        self.initial_comments = []
+        self.htc_inputfiles = []
+        if self.filename is None:
+            lines = self.empty_htc.split("\n")
+        else:
+            lines = self.readlines(self.filename)
+
+        lines = [l.strip() for l in lines]
+
+        #lines = copy(self.lines)
+        while lines:
+            if lines[0].startswith(";"):
+                self.initial_comments.append(lines.pop(0).strip() + "\n")
+            elif lines[0].lower().startswith("begin"):
+                self._add_contents(HTCSection.from_lines(lines))
+            else:
+                line = HTCLine.from_lines(lines)
+                if line.name_ == "exit":
+                    break
+                self._add_contents(line)
+
+    def readfilelines(self, filename):
+        with self.open(self.unix_path(os.path.abspath(filename.replace('\\', '/'))), encoding='cp1252') as fid:
+            txt = fid.read()
+        if txt[:10].encode().startswith(b'\xc3\xaf\xc2\xbb\xc2\xbf'):
+            txt = txt[3:]
+        if self.jinja_tags:
+            template = jinja2.Template(txt)
+            txt = template.render(**self.jinja_tags)
+        return txt.replace("\r", "").split("\n")
+
+    def readlines(self, filename):
+        if filename != self.filename:  # self.filename may be changed by set_name/save. Added it when needed instead
+            self.htc_inputfiles.append(filename)
+        htc_lines = []
+        lines = self.readfilelines(filename)
+        for l in lines:
+            if l.lower().lstrip().startswith('continue_in_file'):
+                filename = l.lstrip().split(";")[0][len("continue_in_file"):].strip().lower()
+
+                if self.modelpath == 'unknown':
+                    p = os.path.dirname(self.filename)
+
+                    def isfile(f):
+                        try:
+                            return os.path.isfile(self.unix_path(f))
+                        except OSError:
+                            return False
+                    lu = [isfile(os.path.abspath(os.path.join(p, "../" * i, filename)))
+                          for i in range(4)].index(True)
+                    filename = os.path.join(p, "../" * lu, filename)
+                else:
+                    filename = os.path.join(self.modelpath, filename)
+                for line in self.readlines(filename):
+                    if line.lstrip().lower().startswith('exit'):
+                        break
+                    htc_lines.append(line)
+            else:
+                htc_lines.append(l)
+        return htc_lines
+
+    def __setitem__(self, key, value):
+        self.contents[key] = value
+
+    def __str__(self):
+        self.contents  # load
+        return "".join(self.initial_comments + [c.__str__(1) for c in self] + ["exit;"])
+
+    def save(self, filename=None):
+        """Saves the htc object to an htc file.
+
+        Args:
+            filename (str, optional): Specifies the filename of the htc file to be saved.
+            If the value is none, the filename attribute of the object will be used as the filename.
+            Defaults to None.
+        """
+        self.contents  # load if not loaded
+        if filename is None:
+            filename = self.filename
+        else:
+            self.filename = filename
+        # exist_ok does not exist in Python27
+        if not os.path.exists(os.path.dirname(filename)) and os.path.dirname(filename) != "":
+            os.makedirs(os.path.dirname(filename))  # , exist_ok=True)
+        with self.open(filename, 'w', encoding='cp1252') as fid:
+            fid.write(str(self))
+
+    def copy(self):
+        """
+        Copy this htc file.
+
+        Returns
+        -------
+        new_htc : HTCFile
+            Copy of this htc file.
+
+        """
+        new_htc = HTCFile()
+        for key in self.keys():
+            new_htc[key] = self[key].copy()
+        new_htc.modelpath = self.modelpath
+        new_htc.filename = self.filename[:-4] + "_copy.htc"
+        return new_htc
+
+    def set_name(self, name, subfolder=''):
+        """Sets the base filename of the simulation files.
+
+        Args:
+            name (str): Specifies name of the log file, dat file (for animation), hdf5 file (for visualization) and htc file.
+            subfolder (str, optional): Specifies the name of a subfolder to place the files in.
+                If the value is an empty string, no subfolders will be created.
+                Defaults to ''.
+
+        Returns:
+            None
+        """
+        # if os.path.isabs(folder) is False and os.path.relpath(folder).startswith("htc" + os.path.sep):
+        self.contents  # load if not loaded
+
+        def fmt_folder(folder, subfolder): return "./" + \
+            os.path.relpath(os.path.join(folder, subfolder)).replace("\\", "/")
+
+        self.filename = os.path.abspath(os.path.join(self.modelpath, fmt_folder(
+            'htc', subfolder), "%s.htc" % name)).replace("\\", "/")
+        if 'simulation' in self:
+            if 'logfile' in self.simulation:
+                self.simulation.logfile = os.path.join(fmt_folder('log', subfolder), "%s.log" % name).replace("\\", "/")
+            if 'animation' in self.simulation:
+                self.simulation.animation = os.path.join(fmt_folder(
+                    'animation', subfolder), "%s.dat" % name).replace("\\", "/")
+            if 'visualization' in self.simulation:
+                f = os.path.join(fmt_folder('visualization', subfolder), "%s.hdf5" % name).replace("\\", "/")
+                self.simulation.visualization[0] = f
+        elif 'test_structure' in self and 'logfile' in self.test_structure:  # hawc2aero
+            self.test_structure.logfile = os.path.join(fmt_folder('log', subfolder), "%s.log" % name).replace("\\", "/")
+        if 'output' in self:
+            self.output.filename = os.path.join(fmt_folder('res', subfolder), "%s" % name).replace("\\", "/")
+
+    def set_time(self, start=None, stop=None, step=None):
+        self.contents  # load if not loaded
+        if stop is not None:
+            self.simulation.time_stop = stop
+        else:
+            stop = self.simulation.time_stop[0]
+        if step is not None:
+            self.simulation.newmark.deltat = step
+        if start is not None:
+            self.output.time = start, stop
+            if "wind" in self:  # and self.wind.turb_format[0] > 0:
+                self.wind.scale_time_start = start
+
+    def expected_simulation_time(self):
+        return 600
+
+    def pbs_file(self, hawc2_path, hawc2_cmd, queue='workq', walltime=None,
+                 input_files=None, output_files=None, copy_turb=(True, True)):
+        walltime = walltime or self.expected_simulation_time() * 2
+        if len(copy_turb) == 1:
+            copy_turb_fwd, copy_turb_back = copy_turb, copy_turb
+        else:
+            copy_turb_fwd, copy_turb_back = copy_turb
+
+        input_files = input_files or self.input_files()
+        if copy_turb_fwd:
+            input_files += [f for f in self.turbulence_files() if os.path.isfile(f)]
+
+        output_files = output_files or self.output_files()
+        if copy_turb_back:
+            output_files += self.turbulence_files()
+
+        return HAWC2PBSFile(hawc2_path, hawc2_cmd, self.filename, self.modelpath,
+                            input_files, output_files,
+                            queue, walltime)
+
+    def input_files(self):
+        self.contents  # load if not loaded
+        if self.modelpath == "unknown":
+            files = [str(f).replace("\\", "/") for f in [self.filename] + self.htc_inputfiles]
+        else:
+            files = [os.path.abspath(str(f)).replace("\\", "/") for f in [self.filename] + self.htc_inputfiles]
+        if 'new_htc_structure' in self:
+            for mb in [self.new_htc_structure[mb]
+                       for mb in self.new_htc_structure.keys() if mb.startswith('main_body')]:
+                if "timoschenko_input" in mb:
+                    files.append(mb.timoschenko_input.filename[0])
+                files.append(mb.get('external_bladedata_dll', [None, None, None])[2])
+        if 'aero' in self:
+            files.append(self.aero.ae_filename[0])
+            files.append(self.aero.pc_filename[0])
+            files.append(self.aero.get('external_bladedata_dll', [None, None, None])[2])
+            files.append(self.aero.get('output_profile_coef_filename', [None])[0])
+            if 'dynstall_ateflap' in self.aero:
+                files.append(self.aero.dynstall_ateflap.get('flap', [None] * 3)[2])
+            if 'bemwake_method' in self.aero:
+                files.append(self.aero.bemwake_method.get('a-ct-filename', [None] * 3)[0])
+        for dll in [self.dll[dll] for dll in self.get('dll', {}).keys() if 'filename' in self.dll[dll]]:
+            files.append(dll.filename[0])
+            f, ext = os.path.splitext(dll.filename[0])
+            files.append(f + "_64" + ext)
+            files.append(f + ".so")
+        if 'wind' in self:
+            files.append(self.wind.get('user_defined_shear', [None])[0])
+            files.append(self.wind.get('user_defined_shear_turbulence', [None])[0])
+            if 'met_mast_wind' in self.wind:
+                if isinstance(self.wind.met_mast_wind, HTCLine):
+                    files.append(self.wind.met_mast_wind[0])
+                else:
+                    files.append(self.wind.met_mast_wind.get('filename', [None])[0])
+        if 'wakes' in self:
+            files.append(self.wind.get('use_specific_deficit_file', [None])[0])
+            files.append(self.wind.get('write_ct_cq_file', [None])[0])
+            files.append(self.wind.get('write_final_deficits', [None])[0])
+        if 'hydro' in self:
+            if 'water_properties' in self.hydro:
+                files.append(self.hydro.water_properties.get('water_kinematics_dll', [None])[0])
+                files.append(self.hydro.water_properties.get('water_kinematics_dll', [None, None])[1])
+        if 'soil' in self:
+            if 'soil_element' in self.soil:
+                files.append(self.soil.soil_element.get('datafile', [None])[0])
+        try:
+            dtu_we_controller = self.dll.get_subsection_by_name('dtu_we_controller')
+            theta_min = dtu_we_controller.init.constant__5[1]
+            if theta_min >= 90:
+                files.append(os.path.join(os.path.dirname(
+                    dtu_we_controller.filename[0]), "wpdata.%d" % theta_min).replace("\\", "/"))
+        except Exception:
+            pass
+
+        try:
+            files.append(self.force.dll.dll[0])
+        except Exception:
+            pass
+
+        def fix_path_case(f):
+            if os.path.isabs(f):
+                return self.unix_path(f)
+            elif self.modelpath != "unknown":
+                try:
+                    return "./" + os.path.relpath(self.unix_path(os.path.join(self.modelpath, f)),
+                                                  self.modelpath).replace("\\", "/")
+                except IOError:
+                    return f
+            else:
+                return f
+        return [fix_path_case(f) for f in set(files) if f]
+
+    def output_files(self):
+        self.contents  # load if not loaded
+        files = []
+        for k, index in [('simulation/logfile', 0),
+                         ('simulation/animation', 0),
+                         ('simulation/visualization', 0),
+                         ('new_htc_structure/beam_output_file_name', 0),
+                         ('new_htc_structure/body_output_file_name', 0),
+                         ('new_htc_structure/struct_inertia_output_file_name', 0),
+                         ('new_htc_structure/body_eigenanalysis_file_name', 0),
+                         ('new_htc_structure/constraint_output_file_name', 0),
+                         ('wind/turb_export/filename_u', 0),
+                         ('wind/turb_export/filename_v', 0),
+                         ('wind/turb_export/filename_w', 0)]:
+            line = self.get(k)
+            if line:
+                files.append(line[index])
+        if 'new_htc_structure' in self:
+            if 'system_eigenanalysis' in self.new_htc_structure:
+                f = self.new_htc_structure.system_eigenanalysis[0]
+                files.append(f)
+                files.append(os.path.join(os.path.dirname(f), 'mode*.dat').replace("\\", "/"))
+            if 'structure_eigenanalysis_file_name' in self.new_htc_structure:
+                f = self.new_htc_structure.structure_eigenanalysis_file_name[0]
+                files.append(f)
+                files.append(os.path.join(os.path.dirname(f), 'mode*.dat').replace("\\", "/"))
+        files.extend(self.res_file_lst())
+
+        for key in [k for k in self.contents.keys() if k.startswith("output_at_time")]:
+            files.append(self[key]['filename'][0] + ".dat")
+        return [f.lower() for f in files if f]
+
+    def turbulence_files(self):
+        self.contents  # load if not loaded
+        if 'wind' not in self.contents.keys() or self.wind.turb_format[0] == 0:
+            return []
+        elif self.wind.turb_format[0] == 1:
+            files = [self.get('wind.mann.filename_%s' % comp, [None])[0] for comp in ['u', 'v', 'w']]
+        elif self.wind.turb_format[0] == 2:
+            files = [self.get('wind.flex.filename_%s' % comp, [None])[0] for comp in ['u', 'v', 'w']]
+        return [f for f in files if f]
+
+    def res_file_lst(self):
+        self.contents  # load if not loaded
+        res = []
+        for output in [self[k] for k in self.keys()
+                       if self[k].name_.startswith("output") and not self[k].name_.startswith("output_at_time")]:
+            dataformat = output.get('data_format', 'hawc_ascii')
+            res_filename = output.filename[0]
+            if dataformat[0] == "gtsdf" or dataformat[0] == "gtsdf64":
+                res.append(res_filename + ".hdf5")
+            elif dataformat[0] == "flex_int":
+                res.extend([res_filename + ".int", os.path.join(os.path.dirname(res_filename), 'sensor')])
+            else:
+                res.extend([res_filename + ".sel", res_filename + ".dat"])
+        return res
+
+    def _simulate(self, exe, skip_if_up_to_date=False):
+        self.contents  # load if not loaded
+        if skip_if_up_to_date:
+            from os.path import isfile, getmtime, isabs
+            res_file = os.path.join(self.modelpath, self.res_file_lst()[0])
+            htc_file = os.path.join(self.modelpath, self.filename)
+            if isabs(exe):
+                exe_file = exe
+            else:
+                exe_file = os.path.join(self.modelpath, exe)
+            #print (from_unix(getmtime(res_file)), from_unix(getmtime(htc_file)))
+            if (isfile(htc_file) and isfile(res_file) and isfile(exe_file) and
+                str(HTCFile(htc_file)) == str(self) and
+                    getmtime(res_file) > getmtime(htc_file) and getmtime(res_file) > getmtime(exe_file)):
+                if "".join(self.readfilelines(htc_file)) == str(self):
+                    return
+
+        self.save()
+        htcfile = os.path.relpath(self.filename, self.modelpath)
+        assert any([os.path.isfile(os.path.join(f, exe)) for f in [''] + os.environ['PATH'].split(";")]), exe
+        return pexec([exe, htcfile], self.modelpath)
+
+    def simulate(self, exe, skip_if_up_to_date=False):
+        errorcode, stdout, stderr, cmd = self._simulate(exe, skip_if_up_to_date)
+        if ('simulation' in self.keys() and "logfile" in self.simulation and
+                os.path.isfile(os.path.join(self.modelpath, self.simulation.logfile[0]))):
+            with self.open(os.path.join(self.modelpath, self.simulation.logfile[0])) as fid:
+                log = fid.read()
+        else:
+            log = "%s\n%s" % (str(stdout), str(stderr))
+
+        if errorcode or 'Elapsed time' not in log:
+            log_lines = log.split("\n")
+            error_lines = [i for i, l in enumerate(log_lines) if 'error' in l.lower()]
+            if error_lines:
+                import numpy as np
+                line_i = np.r_[np.array([error_lines + i for i in np.arange(-3, 4)]).flatten(),
+                               np.arange(-5, 0) + len(log_lines)]
+                line_i = sorted(np.unique(np.maximum(np.minimum(line_i, len(log_lines) - 1), 0)))
+
+                lines = ["%04d %s" % (i, log_lines[i]) for i in line_i]
+                for jump in np.where(np.diff(line_i) > 1)[0]:
+                    lines.insert(jump, "...")
+
+                error_log = "\n".join(lines)
+            else:
+                error_log = log
+            raise Exception("\nError code: %s\nstdout:\n%s\n--------------\nstderr:\n%s\n--------------\nlog:\n%s\n--------------\ncmd:\n%s" %
+                            (errorcode, str(stdout), str(stderr), error_log, cmd))
+        return str(stdout) + str(stderr), log
+
+    def simulate_hawc2stab2(self, exe):
+        errorcode, stdout, stderr, cmd = self._simulate(exe, skip_if_up_to_date=False)
+
+        if errorcode:
+            raise Exception("\nstdout:\n%s\n--------------\nstderr:\n%s\n--------------\ncmd:\n%s" %
+                            (str(stdout), str(stderr), cmd))
+        return str(stdout) + str(stderr)
+
+    def deltat(self):
+        return self.simulation.newmark.deltat[0]
+
+    def compare(self, other):
+        if isinstance(other, str):
+            other = HTCFile(other)
+        return HTCContents.compare(self, other)
+
+    @property
+    def open(self):
+        return open
+
+    def unix_path(self, filename):
+        filename = os.path.realpath(str(filename)).replace("\\", "/")
+        ufn, rest = os.path.splitdrive(filename)
+        ufn += "/"
+        for f in rest[1:].split("/"):
+            f_lst = [f_ for f_ in os.listdir(ufn) if f_.lower() == f.lower()]
+            if len(f_lst) > 1:
+                # use the case sensitive match
+                f_lst = [f_ for f_ in f_lst if f_ == f]
+            if len(f_lst) == 0:
+                raise IOError("'%s' not found in '%s'" % (f, ufn))
+            else:  # one match found
+                ufn = os.path.join(ufn, f_lst[0])
+        return ufn.replace("\\", "/")
+
+
+#
+#     def get_body(self, name):
+#         lst = [b for b in self.new_htc_structure if b.name_=="main_body" and b.name[0]==name]
+#         if len(lst)==1:
+#             return lst[0]
+#         else:
+#             if len(lst)==0:
+#                 raise ValueError("Body '%s' not found"%name)
+#             else:
+#                 raise NotImplementedError()
+#
+
+class H2aeroHTCFile(HTCFile):
+    def __init__(self, filename=None, modelpath=None):
+        HTCFile.__init__(self, filename=filename, modelpath=modelpath)
+
+    @property
+    def simulation(self):
+        return self.test_structure
+
+    def set_time(self, start=None, stop=None, step=None):
+        if stop is not None:
+            self.test_structure.time_stop = stop
+        else:
+            stop = self.simulation.time_stop[0]
+        if step is not None:
+            self.test_structure.deltat = step
+        if start is not None:
+            self.output.time = start, stop
+            if "wind" in self and self.wind.turb_format[0] > 0:
+                self.wind.scale_time_start = start
+
+
+class SSH_HTCFile(HTCFile):
+    def __init__(self, ssh, filename=None, modelpath=None):
+        object.__setattr__(self, 'ssh', ssh)
+        HTCFile.__init__(self, filename=filename, modelpath=modelpath)
+
+    @property
+    def open(self):
+        return self.ssh.open
+
+    def unix_path(self, filename):
+        rel_filename = os.path.relpath(filename, self.modelpath).replace("\\", "/")
+        _, out, _ = self.ssh.execute("find -ipath ./%s" % rel_filename, cwd=self.modelpath)
+        out = out.strip()
+        if out == "":
+            raise IOError("'%s' not found in '%s'" % (rel_filename, self.modelpath))
+        elif "\n" in out:
+            raise IOError("Multiple '%s' found in '%s' (due to case senitivity)" % (rel_filename, self.modelpath))
+        else:
+            drive, path = os.path.splitdrive(os.path.join(self.modelpath, out))
+            path = os.path.realpath(path).replace("\\", "/")
+            return os.path.join(drive, os.path.splitdrive(path)[1])
+
+
+if "__main__" == __name__:
+    f = HTCFile(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\htc\DTU_10MW_RWT_power_curve.htc", "../")
+    print(f.input_files())
+#     f.save(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\htc\DTU_10MW_RWT_power_curve.htc")
+#
+#     f = HTCFile(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\htc\DTU_10MW_RWT.htc", "../")
+#     f.set_time = 0, 1, .1
+#     print(f.simulate(r"C:\mmpe\HAWC2\bin\HAWC2_12.8\hawc2mb.exe"))
+#
+#     f.save(r"C:\mmpe\HAWC2\models\DTU10MWRef6.0\htc\DTU_10MW_RWT.htc")

wetb/hawc2/htc_file_set.py

+import glob
+import os
+import copy
+from wetb.hawc2.hawc2_pbs_file import JESS_WINE32_HAWC2MB
+from wetb.hawc2.htc_file import HTCFile
+from wetb.utils.cluster_tools.pbsfile import PBSMultiRunner
+
+
+class HTCFileSet():
+    def __init__(self, model_path, htc_lst="**/*.htc"):
+        self.model_path = model_path
+
+        if not isinstance(htc_lst, list):
+            htc_lst = [htc_lst]
+
+        self.htc_files = []
+        for htc_path in htc_lst:
+            if os.path.isfile(htc_path):
+                self.htc_files.append(htc_path)
+            else:
+                if not os.path.isabs(htc_path):
+                    htc_path = os.path.join(model_path, htc_path)
+                for filename in glob.iglob(htc_path, recursive=True):
+                    self.htc_files.append(filename)
+
+    def pbs_files(self, hawc2_path, hawc2_cmd, queue='workq', walltime=None,
+                  input_files=None, output_files=None, copy_turb=(True, True)):
+
+        return (HTCFile(htc).pbs_file(hawc2_path, hawc2_cmd, queue=queue, walltime=walltime,
+                                      input_files=copy.copy(input_files),
+                                      output_files=copy.copy(output_files),
+                                      copy_turb=copy_turb) for htc in self.htc_files)
+
+    def save_pbs_files(self, hawc2_path=None, hawc2_cmd=JESS_WINE32_HAWC2MB, queue='workq', walltime=None,
+                       input_files=None, output_files=None, copy_turb=(True, True)):
+        for pbs in self.pbs_files(hawc2_path, hawc2_cmd, queue=queue, walltime=walltime,
+                                  input_files=input_files, output_files=output_files,
+                                  copy_turb=copy_turb):
+            pbs.save(self.model_path)
+
+
+if __name__ == '__main__':
+    #model_path = r'R:\HAWC2_tests\v12.6_mmpe3\win32\simple1'
+    model_path = "w:/simple1"
+    pbs_files = HTCFileSet(model_path).pbs_files(
+        hawc2_path=r"R:\HAWC2_tests\v12.6_mmpe3\hawc2\win32", hawc2_cmd=JESS_WINE32_HAWC2MB, input_files=['data/*'])
+    import pandas as pd
+    time_overview = pd.read_excel(
+        r'C:\mmpe\programming\Fortran\HAWC2_git\HAWC2\pytest_hawc2\release_tests\Time_overview.xlsx')
+    for pbs in pbs_files:
+        f = pbs.filename
+
+        pbs.walltime = time_overview.loc[f[:-3].replace("pbs_in/", 'simple1/')]['mean'] * 24 * 3600
+        pbs.save(model_path)
+    PBSMultiRunner(model_path, nodes=1, ppn=10).save()

wetb/hawc2/log_file.py

@@ -3,14 +3,6 @@ Created on 18/11/2015
 
 @author: MMPE
 '''
-from __future__ import print_function
-from __future__ import division
-from __future__ import unicode_literals
-from __future__ import absolute_import
-from io import open
-from builtins import int
-from future import standard_library
-standard_library.install_aliases()
 import os
 from wetb.hawc2.htc_file import HTCFile
 from collections import OrderedDict
@@ -23,6 +15,7 @@ INITIALIZATION = 'Initializing simulation'
 SIMULATING = "Simulating"
 DONE = "Simulation succeded"
 
+
 class LogInterpreter(object):
     def __init__(self, time_stop):
         self.time_stop = time_stop
@@ -60,7 +53,7 @@ class LogInterpreter(object):
         try:
             return float(time_line[time_line.index('=') + 1:time_line.index('Iter')])
         except:
-            print ("Cannot extract time from #" + time_line + "#")
+            print("Cannot extract time from #" + time_line + "#")
             pass
 
     def update_status(self, new_lines=""):
@@ -97,14 +90,16 @@ class LogInterpreter(object):
                     _3to2list1 = list(txt.split('Elapsed time'))
                     simulation_txt, rest, = _3to2list1[:1] + [_3to2list1[1:]]
                     if "*** ERROR ***" in simulation_txt:
-                        self.errors.extend([l.strip() for l in simulation_txt.strip().split("\n") if "error" in l.lower()])
+                        self.errors.extend([l.strip()
+                                            for l in simulation_txt.strip().split("\n") if "error" in l.lower()])
                     i1 = simulation_txt.rfind("Global time")
                     if i1 > -1:
                         self.current_time = self.extract_time(simulation_txt[i1:])
                     if self.current_time is not None and self.time_stop > 0:
                         self.pct = int(100 * self.current_time // self.time_stop)
                     try:
-                        self.remaining_time = (time.time() - self.start_time[1]) / (self.current_time - self.start_time[0]) * (self.time_stop - self.current_time)
+                        self.remaining_time = (
+                            time.time() - self.start_time[1]) / (self.current_time - self.start_time[0]) * (self.time_stop - self.current_time)
                     except:
                         pass
                     if rest:
@@ -118,7 +113,6 @@ class LogInterpreter(object):
             error_dict[error] = error_dict.get(error, 0) + 1
         return "\n".join([("%d x %s" % (v, k), k)[v == 1] for k, v in error_dict.items()])
 
-
     def remaining_time_str(self):
         if self.remaining_time:
             if self.remaining_time < 3600:
@@ -138,7 +132,6 @@ class LogFile(LogInterpreter):
         self.filename = log_filename
         LogInterpreter.__init__(self, time_stop)
 
-
     @staticmethod
     def from_htcfile(htcfile, modelpath=None):
         logfilename = htcfile.simulation.logfile[0]
@@ -150,13 +143,14 @@ class LogFile(LogInterpreter):
 
     def clear(self):
         # exist_ok does not exist in Python27
-        if not os.path.exists(os.path.dirname(self.filename)):
-            os.makedirs(os.path.dirname(self.filename))  #, exist_ok=True)
-        try:
-            with open(self.filename, 'w', encoding='utf-8'):
-                pass
-        except PermissionError as e:
-            raise PermissionError(str(e) + "\nLog file cannot be cleared. Check if it is open in another program")
+        # if not os.path.exists(os.path.dirname(self.filename)):
+        #    os.makedirs(os.path.dirname(self.filename))  #, exist_ok=True)
+        if os.path.isfile(self.filename):
+            try:
+                with open(self.filename, 'w', encoding='utf-8'):
+                    pass
+            except PermissionError as e:
+                raise PermissionError(str(e) + "\nLog file cannot be cleared. Check if it is open in another program")
         LogInterpreter.clear(self)
 
     def update_status(self):
@@ -174,6 +168,7 @@ class LogFile(LogInterpreter):
             if txt != "":
                 LogInterpreter.update_status(self, txt)
 
+
 class LogInfo(LogFile):
     def __init__(self, status, pct, remaining_time, lastline):
         self.status = status

wetb/hawc2/pc_file.py

@@ -3,17 +3,8 @@ Created on 24/04/2014
 
 @author: MMPE
 '''
-from __future__ import division
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import absolute_import
-from io import open
-from builtins import range
-from builtins import int
-from future import standard_library
-standard_library.install_aliases()
-
-from wetb.hawc2.ae_file import AEFile
+
+import os
 import numpy as np
 
 class PCFile(object):
@@ -29,11 +20,19 @@ class PCFile(object):
     >>> pcfile.CM(21,10) # CM for thickness 21% and AOA=10deg
     -0.1103
     """
-    def __init__(self, filename):
-        with open (filename) as fid:
-            lines = fid.readlines()
-        nsets = int(lines[0].split()[0])
+    def __init__(self, filename=None):
         self.pc_sets = {}
+        if filename is not None:
+            with open (filename) as fid:
+                lines = fid.readlines()
+            self._parse_lines(lines)
+            self.filename = filename
+        self.fmt = ' 19.015e'
+
+    def _parse_lines(self, lines):
+        """Read HAWC2 PC file (profile coefficient file).
+        """
+        nsets = int(lines[0].split()[0])
         lptr = 1
         for nset in range(1, nsets + 1):
             nprofiles = int(lines[lptr].split()[0])
@@ -62,7 +61,7 @@ class PCFile(object):
         Cx1 = np.interp(alpha, Cx1[:, 0], Cx1[:, column])
         th0, th1 = thicknesses[index - 1:index + 1]
         return Cx0 + (Cx1 - Cx0) * (thickness - th0) / (th1 - th0)
-    
+
     def _CxxxH2(self, thickness, alpha, column, pc_set_nr=1):
         thicknesses, profiles = self.pc_sets[pc_set_nr]
         index = np.searchsorted(thicknesses, thickness)
@@ -70,7 +69,7 @@ class PCFile(object):
             index = 1
 
         Cx0, Cx1 = profiles[index - 1:index + 1]
-        
+
         Cx0 = np.interp(np.arange(360), Cx0[:,0]+180, Cx0[:,column])
         Cx1 = np.interp(np.arange(360), Cx1[:,0]+180, Cx1[:,column])
         #Cx0 = np.interp(alpha, Cx0[:, 0], Cx0[:, column])
@@ -78,8 +77,6 @@ class PCFile(object):
         th0, th1 = thicknesses[index - 1:index + 1]
         cx = Cx0 + (Cx1 - Cx0) * (thickness - th0) / (th1 - th0)
         return np.interp(alpha+180, np.arange(360), cx)
-    
-        
 
     def CL(self, thickness, alpha, pc_set_nr=1):
         """Lift coefficient
@@ -99,10 +96,9 @@ class PCFile(object):
         """
         return self._Cxxx(thickness, alpha, 1, pc_set_nr)
 
-
     def CL_H2(self, thickness, alpha, pc_set_nr=1):
         return self._CxxxH2(thickness, alpha, 1, pc_set_nr)
-    
+
     def CD(self, thickness, alpha, pc_set_nr=1):
         """Drag coefficient
 
@@ -124,11 +120,41 @@ class PCFile(object):
     def CM(self, thickness, alpha, pc_set_nr=1):
         return self._Cxxx(thickness, alpha, 3, pc_set_nr)
 
+    def __str__(self, comments=None):
+        """This method will create a string that is formatted like a pc file
+        with the data in this class.
+        """
+
+        if comments is None:
+            comments = {}
+
+        cols = ['Angle of Attac', 'cl', 'cd', 'cm']
+        linefmt = ' '.join(['{%i:%s}' % (i, self.fmt) for i in range(len(cols))])
+
+        n_sets = len(self.pc_sets)
+        retval = str(n_sets) + '\n'
+        for idx_pc, (set_tcs, set_pcs) in self.pc_sets.items():
+            retval += str(len(set_tcs)) + '\n'
+            for i, (tc, pc) in enumerate(zip(set_tcs, set_pcs)):
+                nr = pc.shape[0]
+                retval += '%i %i %1.08f\n' % (i+1, nr, tc)
+                for line in pc:
+                    retval += linefmt.format(*line) + '\n'
+        return retval
+
+    def save(self, filename):
+        if not os.path.isdir(os.path.dirname(filename)):
+            # fails if dirname is empty string
+            if len(os.path.dirname(filename)) > 0:
+                os.makedirs(os.path.dirname(filename))
+        with open(filename, 'w') as fid:
+            fid.write(str(self))
+        self.filename = filename
+
+
 if __name__ == "__main__":
     pcfile = PCFile("tests/test_files/NREL_5MW_pc.txt")
-    aefile = AEFile("tests/test_files/NREL_5MW_ae.txt")
-    print (aefile.thickness(36))
-    
+
     print (pcfile.CL(21,10)) # CL for thickness 21% and AOA=10deg
     #1.358
     print (pcfile.CD(21,10)) # CD for thickness 21% and AOA=10deg

wetb/hawc2/sel_file.py

-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import division
-from __future__ import absolute_import
-from io import open
-from builtins import str
-from builtins import int
-from future import standard_library
-standard_library.install_aliases()
-from builtins import object
 from datetime import datetime
 import os
 import numpy as np

wetb/hawc2/sensor_names.py

+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Mon Mar 10 11:11:48 2025
+
+@author: dave
+"""
+
+import re
+import copy
+import pandas as pd
+
+from wetb.prepost import misc
+
+
+def unified_channel_names(ChInfo):
+    """Create consistant and unique channel names for a HAWC2 result file.
+
+    Parameters
+    ----------
+    ChInfo : List of list
+        The list of list as given by wetb.hawc2.ReadHawc2.ChInfo.
+
+    Make certain channels independent from their index.
+
+    The unified channel dictionary ch_dict holds consequently named
+    channels as the key, and the all information is stored in the value
+    as another dictionary.
+
+    The ch_dict key/values pairs are structured differently for different
+    type of channels. Currently supported channels are:
+
+    For forcevec, momentvec, state commands:
+        node numbers start with 0 at the root
+        element numbers start with 1 at the root
+        key:
+            coord-bodyname-pos-sensortype-component
+            global-tower-node-002-forcevec-z
+            local-blade1-node-005-momentvec-z
+            hub1-blade1-elem-011-zrel-1.00-state pos-z
+        value:
+            ch_dict[tag]['coord']
+            ch_dict[tag]['bodyname']
+            ch_dict[tag]['pos']
+            ch_dict[tag]['sensortype']
+            ch_dict[tag]['component']
+            ch_dict[tag]['chi']
+            ch_dict[tag]['sensortag']
+            ch_dict[tag]['units']
+
+    For the DLL's this is:
+        key:
+            DLL-dll_name-io-io_nr
+            DLL-yaw_control-outvec-3
+            DLL-yaw_control-inpvec-1
+        value:
+            ch_dict[tag]['dll_name']
+            ch_dict[tag]['io']
+            ch_dict[tag]['io_nr']
+            ch_dict[tag]['chi']
+            ch_dict[tag]['sensortag']
+            ch_dict[tag]['units']
+
+    For the bearings this is:
+        key:
+            bearing-bearing_name-output_type-units
+            bearing-shaft_nacelle-angle_speed-rpm
+        value:
+            ch_dict[tag]['bearing_name']
+            ch_dict[tag]['output_type']
+            ch_dict[tag]['chi']
+            ch_dict[tag]['units']
+
+    For many of the aero sensors:
+        'Cl', 'Cd', 'Alfa', 'Vrel'
+        key:
+            sensortype-blade_nr-pos
+            Cl-1-0.01
+        value:
+            ch_dict[tag]['sensortype']
+            ch_dict[tag]['blade_nr']
+            ch_dict[tag]['pos']
+            ch_dict[tag]['chi']
+            ch_dict[tag]['units']
+
+    Returns
+    -------
+    None.
+
+    """
+
+    # save them in a dictionary, use the new coherent naming structure
+    # as the key, and as value again a dict that hols all the different
+    # classifications: (chi, channel nr), (coord, coord), ...
+    ch_dict = dict()
+
+    # all columns for the output
+    cols = ['bearing_name', 'sensortag', 'bodyname', 'chi', 'component',
+            'pos', 'coord', 'sensortype', 'radius', 'blade_nr', 'units',
+            'output_type', 'io_nr', 'io', 'dll', 'azimuth', 'flap_nr',
+            'direction', 'wake_source_nr', 'center', 's', 'srel',
+            'radius_actual']
+    # so we can look up missing columns later
+    # don't use the set for building df_dict since the order will become random
+    colsset = set(cols)
+
+    # some channel ID's are unique, use them
+    ch_unique = set(['Omega', 'Ae rot. torque', 'Ae rot. power',
+                     'Ae rot. thrust', 'Time', 'Azi  1'])
+    ch_aero = set(['Cl', 'Cd', 'Cm', 'Alfa', 'Vrel', 'Tors_e', 'Alfa',
+                   'Lift', 'Drag'])
+    ch_aerogrid = set(['a_grid', 'am_grid', 'CT', 'CQ'])
+
+    # also safe as df
+    df_dict = {col: [] for col in cols}
+    df_dict['unique_ch_name'] = []
+
+    # -----------------------------------------------------------------
+    # REGEXes
+    # -----------------------------------------------------------------
+    # ESYS output: ESYS line3 SENSOR           66
+    re_esys = re.compile(r'ESYS (\w+) SENSOR\s*(\d*)')
+    # FORCE fext_damp   1
+    re_force = re.compile(r'FORCE (\w+) \s*(\d*)')
+
+    # scan through all channels and see which can be converted
+    # to sensible unified name
+    for ich in range(len(ChInfo[0])):
+
+        ch_id = ChInfo[0][ich]
+        ch_unit = ChInfo[1][ich]
+        ch_desc = ChInfo[2][ich]
+
+        # if len(ch_id) < 1 or len(ch_desc) < 1:
+        #     continue
+
+        items_ch0 = ch_id.split()
+        items_ch2 = ch_desc.split()
+
+        dll = False
+
+        # be carefull, identify only on the starting characters, because
+        # the signal tag can hold random text that in some cases might
+        # trigger a false positive
+
+        # -----------------------------------------------------------------
+        # check for all the unique channel descriptions
+        if ch_id.strip() in ch_unique:
+            tag = ch_id.strip()
+            channelinfo = {}
+            channelinfo['units'] = ch_unit
+            channelinfo['sensortag'] = ch_desc
+            channelinfo['chi'] = ich
+
+        # -----------------------------------------------------------------
+        # or in the long description:
+        #    0          1        2      3  4    5     6 and up
+        # MomentMz Mbdy:blade nodenr:   5 coo: blade  TAG TEXT
+        elif ch_desc.startswith('MomentM'):
+            coord = items_ch2[5]
+            bodyname = items_ch2[1].replace('Mbdy:', '')
+            # set nodenr to sortable way, include leading zeros
+            # node numbers start with 0 at the root
+            nodenr = '%03i' % int(items_ch2[3])
+            # skip the attached the component
+            # sensortype = items[0][:-2]
+            # or give the sensor type the same name as in HAWC2
+            sensortype = 'momentvec'
+            component = items_ch2[0][-1:len(items_ch2[0])]
+            # the tag only exists if defined
+            if len(items_ch2) > 6:
+                sensortag = ' '.join(items_ch2[6:])
+            else:
+                sensortag = ''
+
+            # and tag it
+            pos = 'node-%s' % nodenr
+            tagitems = (coord, bodyname, pos, sensortype, component)
+            tag = '%s-%s-%s-%s-%s' % tagitems
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['coord'] = coord
+            channelinfo['bodyname'] = bodyname
+            channelinfo['pos'] = pos
+            channelinfo['sensortype'] = sensortype
+            channelinfo['component'] = component
+            channelinfo['chi'] = ich
+            channelinfo['sensortag'] = sensortag
+            channelinfo['units'] = ch_unit
+
+        # -----------------------------------------------------------------
+        #   0    1      2        3       4  5     6     7 and up
+        # Force  Fx Mbdy:blade nodenr:   2 coo: blade  TAG TEXT
+        elif ch_desc.startswith('Force  F'):
+            coord = items_ch2[6]
+            bodyname = items_ch2[2].replace('Mbdy:', '')
+            nodenr = '%03i' % int(items_ch2[4])
+            # skipe the attached the component
+            # sensortype = items[0]
+            # or give the sensor type the same name as in HAWC2
+            sensortype = 'forcevec'
+            component = items_ch2[1][1]
+            if len(items_ch2) > 7:
+                sensortag = ' '.join(items_ch2[7:])
+            else:
+                sensortag = ''
+
+            # and tag it
+            pos = 'node-%s' % nodenr
+            tagitems = (coord, bodyname, pos, sensortype, component)
+            tag = '%s-%s-%s-%s-%s' % tagitems
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['coord'] = coord
+            channelinfo['bodyname'] = bodyname
+            channelinfo['pos'] = pos
+            channelinfo['sensortype'] = sensortype
+            channelinfo['component'] = component
+            channelinfo['chi'] = ich
+            channelinfo['sensortag'] = sensortag
+            channelinfo['units'] = ch_unit
+
+        # -----------------------------------------------------------------
+        #        0    1      2        3       4    5     6     7    8           9 and up
+        # Force_intp  Fz Mbdy:blade1 s=  11.87[m] s/S=   0.95 coo: local_aero center:default
+        # Moment_intp  Mx Mbdy:blade1 s=  11.87[m] s/S=   0.95 coo: local_aero center:default
+        elif items_ch2[0].endswith('_intp'):
+
+            sensortype = 'forcemomentvec_interp'
+
+            coord = items_ch2[8]
+            bodyname = items_ch2[2].replace('Mbdy:', '')
+            s = items_ch2[4].replace('[m]', '')
+            srel = items_ch2[6]
+            center = items_ch2[9].split(':')[1]
+            component = items_ch2[1]
+
+            if len(items_ch2) > 9:
+                sensortag = ' '.join(items_ch2[10:])
+            else:
+                sensortag = ''
+
+            # and tag it
+            pos = 's-%s' % (s)
+            tag = f'{sensortype}-{bodyname}-{center}-{coord}-{s}-{component}'
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['coord'] = coord
+            channelinfo['bodyname'] = bodyname
+            channelinfo['s'] = float(s)
+            channelinfo['srel'] = float(srel)
+            channelinfo['sensortype'] = sensortype
+            # channelinfo['output_type'] = output_type
+            channelinfo['component'] = component
+            channelinfo['center'] = center
+            channelinfo['chi'] = ich
+            channelinfo['sensortag'] = sensortag
+            channelinfo['units'] = ch_unit
+
+        # -----------------------------------------------------------------
+        # ELEMENT STATES: pos, vel, acc, rot, ang
+        #   0    1  2      3       4      5   6         7    8
+        # State pos x  Mbdy:blade E-nr:   1 Z-rel:0.00 coo: blade
+        #   0           1     2    3        4    5   6         7     8     9+
+        # State_rot proj_ang tx Mbdy:bname E-nr: 1 Z-rel:0.00 coo: cname  label
+        # State_rot omegadot tz Mbdy:bname E-nr: 1 Z-rel:1.00 coo: cname  label
+        elif ch_desc.startswith('State'):
+#                 or ch_details[ich,0].startswith('euler') \
+#                 or ch_details[ich,0].startswith('ax') \
+#                 or ch_details[ich,0].startswith('omega') \
+#                 or ch_details[ich,0].startswith('proj'):
+            coord = items_ch2[8]
+            bodyname = items_ch2[3].replace('Mbdy:', '')
+            # element numbers start with 1 at the root
+            elementnr = '%03i' % int(items_ch2[5])
+            zrel = '%04.2f' % float(items_ch2[6].replace('Z-rel:', ''))
+            # skip the attached the component
+            #sensortype = ''.join(items[0:2])
+            # or give the sensor type the same name as in HAWC2
+            tmp = ch_id.split(' ')
+            sensortype = tmp[0]
+            if sensortype.startswith('State'):
+                sensortype += ' ' + tmp[1]
+            component = items_ch2[2]
+            if len(items_ch2) > 8:
+                sensortag = ' '.join(items_ch2[9:])
+            else:
+                sensortag = ''
+
+            # and tag it
+            pos = 'elem-%s-zrel-%s' % (elementnr, zrel)
+            tagitems = (coord, bodyname, pos, sensortype, component)
+            tag = '%s-%s-%s-%s-%s' % tagitems
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['coord'] = coord
+            channelinfo['bodyname'] = bodyname
+            channelinfo['pos'] = pos
+            channelinfo['sensortype'] = sensortype
+            channelinfo['component'] = component
+            channelinfo['chi'] = ich
+            channelinfo['sensortag'] = sensortag
+            channelinfo['units'] = ch_unit
+
+        # -----------------------------------------------------------------
+        # statevec_new
+        #   0     1          2      3  4           5     6       7
+        # elastic Deflection blade1 Dx Mbdy:blade1 s=   0.00[m] s/S=
+        # 8     9     10     11
+        # 0.00 coo: blade1 center:c2def
+        # note that: 2 and 10 are the same
+        elif items_ch2[0] == 'elastic' or items_ch2[0] == 'absolute':
+            output_type = ' '.join(items_ch2[0:2])
+            bodyname = items_ch2[4].replace('Mbdy:', '')
+            s = '%06.02f' % float(items_ch2[6].replace('[m]', ''))
+            srel = '%04.02f' % float(items_ch2[8])
+            coord = items_ch2[10]
+            center = items_ch2[11].split(':')[1]
+            sensortype = 'statevec_new'
+
+            component = items_ch0[0]
+
+            if len(items_ch2) > 12:
+                sensortag = ' '.join(items_ch2[12:])
+            else:
+                sensortag = ''
+
+            # and tag it, allmost the same as in htc file here
+            tagitems = (sensortype, bodyname, center, coord, items_ch2[0],
+                        s, component)
+            tag = '-'.join(['%s']*7) % tagitems
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['coord'] = coord
+            channelinfo['bodyname'] = bodyname
+            channelinfo['s'] = float(s)
+            channelinfo['srel'] = float(srel)
+            channelinfo['sensortype'] = sensortype
+            channelinfo['output_type'] = output_type
+            channelinfo['component'] = component
+            channelinfo['center'] = center
+            channelinfo['chi'] = ich
+            channelinfo['sensortag'] = sensortag
+            channelinfo['units'] = ch_unit
+
+        # -----------------------------------------------------------------
+        # DLL CONTROL I/O
+        # there are two scenario's on how the channel description is formed
+        # the channel id is always the same though
+        # id for all three cases:
+        #          DLL out  1:  3
+        #          DLL inp  2:  3
+        # description case 1 ("dll type2_dll b2h2 inpvec 30" in htc output)
+        #               0         1    2   3     4+
+        #          yaw_control outvec  3  yaw_c input reference angle
+        # description case 2 ("dll inpvec 2 1" in htc output):
+        #           0  1 2     3  4  5  6+
+        #          DLL : 2 inpvec :  4  mgen hss
+        # description case 3
+        #           0         1     2       4
+        #          hawc_dll :echo outvec :  1
+        elif ch_id.startswith('DLL'):
+            # case 3
+            if items_ch2[0] == 'hawc_dll':
+                # hawc_dll named case (case 3) is polluted with colons
+                dll = items_ch2[1].replace(':', '')
+                io = items_ch2[2]
+                io_nr = items_ch2[4]
+                tag = 'DLL-%s-%s-%s' % (dll, io, io_nr)
+                sensortag = ''
+            # case 2: no reference to dll name
+            elif ch_desc.startswith('DLL'):
+                dll = items_ch2[2]
+                io = items_ch2[3]
+                io_nr = items_ch2[5]
+                sensortag = ' '.join(items_ch2[6:])
+                # and tag it
+                tag = 'DLL-%s-%s-%s' % (dll,io,io_nr)
+            # case 1: type2 dll name is given
+            else:
+                dll = items_ch2[0]
+                io = items_ch2[1]
+                io_nr = items_ch2[2]
+                sensortag = ' '.join(items_ch2[3:])
+                tag = 'DLL-%s-%s-%s' % (dll, io, io_nr)
+
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['dll'] = dll
+            channelinfo['io'] = io
+            channelinfo['io_nr'] = io_nr
+            channelinfo['chi'] = ich
+            channelinfo['sensortag'] = sensortag
+            channelinfo['units'] = ch_unit
+            channelinfo['sensortype'] = 'dll-io'
+
+        # -----------------------------------------------------------------
+        # BEARING OUTPUS
+        # bea1 angle_speed       rpm      shaft_nacelle angle speed
+        elif ch_id.startswith('bea'):
+            output_type = ch_id.split(' ')[1]
+            bearing_name = items_ch2[0]
+            units = ch_unit
+            # there is no label option for the bearing output
+
+            # and tag it
+            tag = 'bearing-%s-%s-%s' % (bearing_name, output_type, units)
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['bearing_name'] = bearing_name
+            channelinfo['output_type'] = output_type
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+
+        # -----------------------------------------------------------------
+        # AS DEFINED IN: ch_aero
+        # AERO CL, CD, CM, VREL, ALFA, LIFT, DRAG, etc
+        # Cl, R=  0.5     deg      Cl of blade  1 at radius   0.49
+        # Azi  1          deg      Azimuth of blade  1
+        #
+        # ch_desc:
+        # Angle of attack of blade   1 at radius   8.59 FOLLOWD BY USER LABEL
+        #
+        # NOTE THAT RADIUS FROM ch_id REFERS TO THE RADIUS
+        # YOU ASKED FOR, AND ch_desc IS WHAT YOU GET, which is
+        # still based on a mean radius (deflections change the game)
+        elif ch_id.split(',')[0] in ch_aero:
+            sensortype = ch_id.split(',')[0]
+
+            # sometimes the units for aero sensors are wrong!
+            units = ch_unit
+            # there is no label option
+
+            # Blade number is identified as the first integer in the string
+            # blade_nr = re.search(r'\d+', ch_desc).group()
+            # blade_nr = int(blade_nr)
+
+            # actual radius
+            rq = r'\.*of blade\s*(\d) at radius\s*([-+]?\d*\.\d+|\d+)'
+            s = ch_desc
+            blade_nr, radius_actual = re.findall(rq, s)[0]
+            blade_nr = int(blade_nr)
+
+            # radius what you asked for, identified as the last float in the string
+            s = ch_id
+            radius = float(re.findall(r"[-+]?\d*\.\d+|\d+", s)[-1])
+
+            # and tag it
+            tag = '%s-%s-%s' % (sensortype, blade_nr, radius)
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['sensortype'] = sensortype
+            channelinfo['radius'] = float(radius)
+            channelinfo['radius_actual'] = float(radius_actual)
+            channelinfo['blade_nr'] = blade_nr
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+
+        # -----------------------------------------------------------------
+        # for the induction grid over the rotor
+        # a_grid, azi    0.00 r   1.74
+        elif ch_id.split(',')[0] in ch_aerogrid:
+            items_ = ch_id.split(',')
+            sensortype = items_[0]
+            items2 = items_[1].split(' ')
+            items2 = misc.remove_items(items2, '')
+            azi = items2[1]
+            # radius what you asked for
+            radius = items2[3]
+            units = ch_unit
+            # and tag it
+            tag = '%s-azi-%s-r-%s' % (sensortype,azi,radius)
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['sensortype'] = sensortype
+            channelinfo['radius'] = float(radius)
+            channelinfo['azimuth'] = float(azi)
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+
+        # -----------------------------------------------------------------
+        # INDUCTION AT THE BLADE
+        # 0: Induc. Vz, rpco, R=  1.4
+        # 1: m/s
+        # 2: Induced wsp Vz of blade  1 at radius   1.37, RP. coo.
+        # Induc. Vx, locco, R=  1.4
+        #    Induced wsp Vx of blade  1 at radius   1.37, local ae coo.
+        # Induc. Vy, blco, R=  1.4
+        #    Induced wsp Vy of blade  1 at radius   1.37, local bl coo.
+        # Induc. Vz, glco, R=  1.4
+        #    Induced wsp Vz of blade  1 at radius   1.37, global coo.
+        # Induc. Vx, rpco, R=  8.4
+        #    Induced wsp Vx of blade  1 at radius   8.43, RP. coo.
+        elif ch_id.strip()[:5] == 'Induc':
+
+            coord = ch_desc.split(', ')[1].strip()
+            blade_nr = int(items_ch2[5])
+
+            # radius what you get
+            #  radius = float(items[8].replace(',', ''))
+            # radius what you asked for, identified as the last float in the string
+#                s = ch_desc
+#                radius = float(re.findall(r"[-+]?\d*\.\d+|\d+", s)[-1])
+            radius = float(items_ch2[8][:-1])
+
+            component = items_ch2[2]
+            units = ch_unit
+
+            # and tag it
+            rpl = (coord, blade_nr, component, radius)
+            tag = 'induc-%s-blade-%1i-%s-r-%03.01f' % rpl
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['blade_nr'] = blade_nr
+            channelinfo['sensortype'] = 'induction'
+            channelinfo['radius'] = radius
+            channelinfo['coord'] = coord
+            channelinfo['component'] = component
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+
+        # -----------------------------------------------------------------
+        # MORE AERO SENSORS
+        # Ae intfrc Fx, rpco, R=  0.0
+        #     Aero int. force Fx of blade  1 at radius   0.00, RP coo.
+        # Ae secfrc Fy, R= 25.0
+        #     Aero force  Fy of blade  1 at radius  24.11
+        # Ae pos x, glco, R= 88.2
+        #     Aero position x of blade  1 at radius  88.17, global coo.
+        elif ch_id.strip()[:2] == 'Ae':
+            units = ch_unit
+            # Blade number is identified as the first integer in the string
+            blade_nr = re.search(r'\d+', ch_desc).group()
+            blade_nr = int(blade_nr)
+            # radius what you get
+            tmp = ch_desc.split('radius ')[1].strip()
+            tmp = tmp.split(',')
+            # radius = float(tmp[0])
+            # radius what you asked for, identified as the last float in the string
+            s = ch_desc
+            radius = float(re.findall(r"[-+]?\d*\.\d+|\d+", s)[-1])
+
+            if len(tmp) > 1:
+                coord = tmp[1].strip()
+            else:
+                coord = 'aero'
+
+            sensortype = items_ch0[1]
+            component = items_ch0[2].replace(',', '')
+
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['blade_nr'] = blade_nr
+            channelinfo['sensortype'] = sensortype
+            channelinfo['radius'] = radius
+            channelinfo['coord'] = coord
+            channelinfo['component'] = component
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+
+            rpl = (coord, blade_nr, sensortype, component, radius)
+            tag = 'aero-%s-blade-%1i-%s-%s-r-%03.01f' % rpl
+
+        # TODO: wind speed
+        # some spaces have been trimmed here
+        # WSP gl. coo.,Vy          m/s
+        # // Free wind speed Vy, gl. coo, of gl. pos   0.00,  0.00,  -2.31
+        # WSP gl. coo.,Vdir_hor          deg
+        # Free wind speed Vdir_hor, gl. coo, of gl. pos  0.00,  0.00, -2.31
+
+        # -----------------------------------------------------------------
+        # Water surface gl. coo, at gl. coo, x,y=   0.00,   0.00
+        # Water vel., x-dir, at gl. coo, x,y,z=   0.00,   0.00,   1.00
+        elif ch_desc.startswith('Water'):
+            units = ch_unit
+            channelinfo = {}
+
+            # surface, vel or acc?
+            if items_ch2[1]=='surface':
+                # but remove the comma
+                x = items_ch2[-2][:-1]
+                y = items_ch2[-1]
+                # and tag it
+                tag = 'watersurface-global-%s-%s' % (x, y)
+                channelinfo['pos'] = (float(x), float(y))
+            else:
+                # but remove the comma
+                x = items_ch2[-3][:-1]
+                y = items_ch2[-2][:-1]
+                z = items_ch2[-1]
+                tag = f'{items_ch0[0]}-{x}-{y}-{z}'
+
+            # save all info in the dict
+            channelinfo['coord'] = 'global'
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+
+        # -----------------------------------------------------------------
+        # WIND SPEED
+        elif ch_desc.startswith('Free wind speed'):
+            units = ch_unit
+            direction = ch_id.split(',')[1]
+            # WSP gl. coo.,Vx
+            # Free wind speed Vx, gl. coo, of gl. pos    0.00,   0.00,  -6.00  LABEL
+            if ch_desc.startswith('Free '):
+                tmp = ch_desc.split('pos')[1]
+                x, y, z = tmp.split(',')
+                x, y, z = x.strip(), y.strip(), z.strip()
+                tmp = z.split('  ')
+                sensortag = ''
+                if len(tmp) == 2:
+                    z, sensortag = tmp
+                elif len(tmp) == 1:
+                    z = tmp[0]
+                pos = (float(x), float(y), float(z))
+                posstr = '%s-%s-%s' % (x, y, z)
+                coord = 'global'
+            else:
+                pos = items_ch2[6]
+                posstr = pos
+                coord = items_ch2[0].lower()
+                if len(items_ch2) > 6:
+                    sensortag = ' '.join(items_ch2[7:])
+
+            # and tag it
+            tag = 'windspeed-%s-%s-%s' % (coord, direction, posstr)
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['coord'] = 'global'
+            channelinfo['pos'] = pos
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+            channelinfo['sensortag'] = sensortag
+            # FIXME: direction is the same as component, right?
+            channelinfo['direction'] = direction
+            channelinfo['sensortype'] = 'wsp-global'
+
+        # WIND SPEED AT BLADE
+        # 0: WSP Vx, glco, R= 61.5
+        # 2: Wind speed Vx of blade  1 at radius  61.52, global coo.
+        elif ch_id.startswith('WSP V'):
+            units = ch_unit.strip()
+            tmp = ch_id.split(' ')[1].strip()
+            direction = tmp.replace(',', '')
+            coord = ch_desc.split(',')[1].split()[0]
+            # Blade number is identified as the first integer in the string
+            blade_nr = re.search(r'\d+', ch_desc).group()
+            blade_nr = int(blade_nr)
+
+            # radius what you get
+            # radius = ch_desc.split('radius')[1].split(',')[0]
+            # radius = radius.strip()
+            # radius what you asked for, identified as the last float in the string
+            s=ch_desc
+            radius=float(re.findall(r"[-+]?\d*\.\d+|\d+", s)[-1])
+
+            # and tag it
+            rpl = (direction, blade_nr, radius, coord)
+            tag = 'wsp-blade-%s-%s-%s-%s' % rpl
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['coord'] = coord
+            # FIXME: direction is the same as component, right?
+            channelinfo['direction'] = direction
+            channelinfo['blade_nr'] = blade_nr
+            channelinfo['radius'] = float(radius)
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+            channelinfo['sensortype'] = 'wsp-blade'
+
+        # FLAP ANGLE
+        # 2: Flap angle for blade  3 flap number  1
+        elif ch_id[:7] == 'setbeta':
+            units = ch_unit.strip()
+            # Blade number is identified as the first integer in the string
+            blade_nr = re.search(r'\d+', ch_desc).group()
+            blade_nr = int(blade_nr)
+            flap_nr = ch_desc.split(' ')[-1].strip()
+
+            # and tag it
+            tag = 'setbeta-bladenr-%s-flapnr-%s' % (blade_nr, flap_nr)
+            # save all info in the dict
+            channelinfo = {}
+            channelinfo['flap_nr'] = int(flap_nr)
+            channelinfo['blade_nr'] = blade_nr
+            channelinfo['units'] = units
+            channelinfo['chi'] = ich
+
+        # harmonic channel output
+        # Harmonic
+        # Harmonic sinus function
+        elif ch_id[:7] == 'Harmoni':
+
+            func_name = ' '.join(ch_unit.split(' ')[1:])
+
+            channelinfo = {}
+            channelinfo['output_type'] = func_name
+            channelinfo['sensortype'] = 'harmonic'
+            channelinfo['chi'] = ich
+
+            base = ch_desc.strip().lower().replace(' ', '_')
+            tag = base + '_0'
+            if tag in ch_dict:
+                tag_nr = int(tag.split('_')[-1]) + 1
+                tag = base + '_%i' % tag_nr
+
+        elif ch_id[:6] == 'a_norm':
+            channelinfo = {}
+            channelinfo['chi'] = ich
+            channelinfo['units'] = ch_unit.strip()
+            channelinfo['sensortype'] = 'aero'
+            tag = 'aero-induc_a_norm'
+
+        # wake   1 gl. pos pos_z  //  Wake pos_z of source   1, Met. coo.
+        elif ch_id[:4] == 'wake':
+            wake_nr = re.search(r'\d+', ch_id).group()
+            comp = re.search(r'pos_([xyz])', ch_id).group(1)
+
+            channelinfo = {}
+            channelinfo['output_type'] = 'wake_pos'
+            channelinfo['sensortype'] = 'wind_wake'
+            channelinfo['component'] = comp
+            channelinfo['units'] = ch_unit.strip()
+            channelinfo['chi'] = ich
+            channelinfo['wake_source_nr'] = int(wake_nr)
+            channelinfo['coord'] = 'met'
+
+            tag = 'wind_wake-wake_pos_%s_%s' % (comp, wake_nr)
+
+        # ESYS line1 SENSOR            1
+        elif ch_desc[:4] == 'ESYS':
+            # body = re.findall(regex, ch_desc)
+            body, outnr = re_esys.match(ch_desc).groups()
+
+            channelinfo = {}
+            channelinfo['output_type'] = 'esys'
+            channelinfo['sensortype'] = 'esys'
+            channelinfo['io_nr'] = int(outnr)
+            channelinfo['units'] = ch_unit.strip()
+            channelinfo['chi'] = ich
+
+            tag = 'esys-%s-%s' % (body, outnr)
+
+        elif ch_desc[:4] == 'FORC':
+            # body = re.findall(regex, ch_desc)
+            dllname, outnr = re_force.match(ch_desc).groups()
+
+            channelinfo = {}
+            channelinfo['output_type'] = 'force-dll'
+            channelinfo['sensortype'] = 'force-dll'
+            channelinfo['io_nr'] = int(outnr)
+            channelinfo['units'] = ch_unit.strip()
+            channelinfo['chi'] = ich
+
+            tag = 'force-%s-%s' % (dllname, outnr)
+
+        # -----------------------------------------------------------------
+        # If all this fails, just combine channel id and its tag if present
+        else:
+            tag = f'{ch_id}'
+            channelinfo = {}
+            channelinfo['chi'] = ich
+            channelinfo['units'] = ch_unit.strip()
+
+        # -----------------------------------------------------------------
+        # add a v_XXX tag in case the channel already exists
+        if tag in ch_dict:
+            jj = 1
+            while True:
+                tag_new = tag + '_v%i' % jj
+                if tag_new in ch_dict:
+                    jj += 1
+                else:
+                    tag = tag_new
+                    break
+
+        ch_dict[tag] = copy.copy(channelinfo)
+
+        # -----------------------------------------------------------------
+        # save in for DataFrame format
+        cols_ch = set(channelinfo.keys())
+        for col in cols_ch:
+            df_dict[col].append(channelinfo[col])
+        # the remainder columns we have not had yet. Fill in blank
+        for col in (colsset - cols_ch):
+            df_dict[col].append('')
+        df_dict['unique_ch_name'].append(tag)
+
+    ch_df = pd.DataFrame(df_dict)
+    ch_df.set_index('chi', inplace=True)
+
+    return ch_dict, ch_df
+
+
+def htc_channel_names(ChInfo):
+    """Incomplete prototype to give every channel a unique channel name which
+    is (nearly) identical to the channel names as defined in the htc output
+    section. Instead of spaces, use colon (;) to seperate the different commands.
+
+    THIS IS STILL A WIP
+
+    see also issue #11:
+    https://gitlab.windenergy.dtu.dk/toolbox/WindEnergyToolbox/issues/11
+    """
+
+    index = {}
+
+    names = {'htc_name':[], 'chi':[], 'label':[], 'unit':[], 'index':[],
+             'name':[], 'description':[]}
+    constraint_fmts = {'bea1':'constraint;bearing1',
+                       'bea2':'constraint;bearing2',
+                       'bea3':'constraint;bearing3',
+                       'bea4':'constraint;bearing4'}
+    # mbdy momentvec tower  1 1 global
+    force_fmts = {'F':'mbdy;forcevec;{body};{nodenr:03i};{coord};{comp}',
+                  'M':'mbdy;momentvec;{body};{nodenr:03i};{coord};{comp}'}
+    state_fmt = 'mbdy;{state};{typ};{body};{elnr:03i};{zrel:01.02f};{coord}'
+
+    wind_coord_map = {'Vx':'1', 'Vy':'2', 'Vz':'3'}
+    wind_fmt = 'wind;{typ};{coord};{x};{y};{z};{comp}'
+
+    for ch in range(len(ChInfo[0])):
+
+        name = ChInfo[0][ch]
+        name_items = misc.remove_items(name.split(' '), '')
+
+        description = ChInfo[2][ch]
+        descr_items = misc.remove_items(description.split(' '), '')
+
+        unit = ChInfo[1][ch]
+
+        # default names
+        htc_name = ' '.join(name_items+descr_items)
+        label = ''
+        coord = ''
+        typ = ''
+        elnr = ''
+        nodenr = ''
+        zrel = ''
+        state = ''
+
+        # CONSTRAINTS: BEARINGS
+        if name_items[0] in constraint_fmts:
+            htc_name = constraint_fmts[name_items[0]] + ';'
+            htc_name += (descr_items[0] + ';')
+            htc_name += unit
+
+        # MBDY FORCES/MOMENTS
+        elif name_items[0][0] in force_fmts:
+            comp = name_items[0]
+            if comp[0] == 'F':
+                i0 = 1
+            else:
+                i0 = 0
+            label = description.split('coo: ')[1].split('  ')[1]
+            coord = descr_items[i0+5]
+            body = descr_items[i0+1][5:]#.replace('Mbdy:', '')
+            nodenr = int(descr_items[i0+3])
+            htc_name = force_fmts[comp[0]].format(body=body, coord=coord,
+                                                  nodenr=nodenr, comp=comp)
+
+        # STATE: POS, VEL, ACC, STATE_ROT
+        elif descr_items[0][:5] == 'State':
+            if name_items[0] == 'State':
+                i0 = 1
+                state = 'state'
+            else:
+                i0 = 0
+                state = 'state_rot'
+            typ = name_items[i0+0]
+            comp = name_items[i0+1]
+            coord = name_items[i0+3]
+
+            body = descr_items[3][5:]#.replace('Mbdy:', '')
+            elnr = int(descr_items[5])
+            zrel = float(descr_items[6][6:])#.replace('Z-rel:', ''))
+            if len(descr_items) > 8:
+                label = ' '.join(descr_items[9:])
+            htc_name = state_fmt.format(typ=typ, body=body, elnr=elnr,
+                                        zrel=zrel, coord=coord,
+                                        state=state)
+
+        # WINDSPEED
+        elif description[:9] == 'Free wind':
+            if descr_items[4] == 'gl.':
+                coord = '1' # global
+            else:
+                coord = '2' # non-rotating rotor coordinates
+
+            try:
+                comp = wind_coord_map[descr_items[3][:-1]]
+                typ = 'free_wind'
+            except KeyError:
+                comp = descr_items[3]
+                typ = 'free_wind_hor'
+
+            tmp = description.split('pos')[1]
+            x, y, z = tmp.split(',')
+            # z might hold a label....
+            z_items  = z.split('  ')
+            if len(z_items) > 1:
+                label = '  '.join(z_items[1:])
+                z = z_items[0]
+            x, y, z = x.strip(), y.strip(), z.strip()
+
+            htc_name = wind_fmt.format(typ=typ, coord=coord, x=x, y=y, z=z,
+                                       comp=comp)
+
+
+        names['htc_name'].append(htc_name)
+        names['chi'].append(ch)
+        # this is the Channel column from the sel file, so the unique index
+        # which is dependent on the order of the channels
+        names['index'].append(ch+1)
+        names['unit'].append(unit)
+        names['name'].append(name)
+        names['description'].append(description)
+        names['label'].append(label)
+        names['state'].append(state)
+        names['type'].append(typ)
+        names['comp'].append(comp)
+        names['coord'].append(coord)
+        names['elnr'].append(coord)
+        names['nodenr'].append(coord)
+        names['zrel'].append(coord)
+        index[name] = ch
+
+    return names, index
+

wetb/hawc2/shear_file.py

@@ -3,15 +3,7 @@ Created on 24/04/2014
 
 @author: MMPE
 '''
-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import division
-from __future__ import absolute_import
-from builtins import range
-from io import open
-from future import standard_library
 from wetb.hawc2.htc_file import HTCFile
-standard_library.install_aliases()
 import numpy as np
 import os
 

wetb/hawc2/simulation.py

-from __future__ import absolute_import
-from __future__ import division
-from __future__ import print_function
-from __future__ import unicode_literals
-
-from builtins import str
-import glob
-from io import open
-import io
-import json
-import os
-import re
-import shutil
-import subprocess
-import sys
-from threading import Thread
-import time
-
-from future import standard_library
-from wetb.hawc2 import log_file
-from wetb.hawc2.htc_file import HTCFile, fmt_path
-from wetb.hawc2.log_file import LogFile
-
-
-standard_library.install_aliases()
-
-QUEUED = "queued"  #until start
-PREPARING = "Copy to host"  # during prepare simulation
-INITIALIZING = "Initializing"  #when starting
-SIMULATING = "Simulating"  # when logfile.status=simulating
-FINISHING = "Copy from host"  # during prepare simulation
-FINISH = "Simulation finish"  # when HAWC2 finish
-ERROR = "Error"  # when hawc2 returns error
-ABORTED = "Aborted"  # when stopped and logfile.status != Done
-CLEANED = "Cleaned"  # after copy back
-
-
-class Simulation(object):
-    """Class for doing hawc2 simulations
-
-
-
-
-    Examples
-    --------
-    >>> sim = Simulation("<path>/MyModel","htc/MyHtc")
-
-    Blocking inplace simulation
-    >>> sim.simulate()
-
-    Non-blocking distributed simulation(everything copied to temporary folder)\n
-    Starts two threads:
-    - non_blocking_simulation_thread:
-        - prepare_simulation() # copy to temporary folder
-        - simulate() # simulate
-        - finish_simulation # copy results back again
-    - updateSimStatusThread:
-        - update status every second
-    >>> sim.start()
-    >>> while sim.status!=CLEANED:
-    >>>     sim.show_status()
-
-
-    The default host is LocalSimulationHost. To simulate on pbs featured cluster
-    >>> sim.host = PBSClusterSimulationHost(sim, <hostname>, <username>, <password>)
-    >>> sim.start()
-    >>> while sim.status!=CLEANED:
-    >>>     sim.show_status()
-    """
-
-    is_simulating = False
-    is_done = False
-    status = QUEUED
-    def __init__(self, modelpath, htcfilename, hawc2exe="HAWC2MB.exe", copy_turbulence=True):
-        self.modelpath = os.path.abspath(modelpath) + "/"
-        if os.path.isabs(htcfilename):
-            htcfilename = os.path.relpath(htcfilename, modelpath)
-        if htcfilename.startswith("input/"):
-            htcfilename=htcfilename[6:]
-        exists = [os.path.isfile(os.path.join(modelpath, htcfilename)), 
-                  os.path.isfile(os.path.join(modelpath, "input/", htcfilename))]
-        if all(exists):
-            raise Exception("Both standard and input/output file structure available for %s in %s. Delete one of the options"%(htcfilename, modelpath) )
-        if not any(exists):
-            raise Exception ("%s not found in %s"%(htcfilename, modelpath))
-        else:
-            self.ios = exists[1] #input/output file structure
-            
-        if self.ios:
-            self.exepath = self.modelpath + "input/"
-        else:
-            self.exepath = self.modelpath 
-        htcfilename = fmt_path(htcfilename)
-        
-        self.tmp_modelpath = self.exepath
-        self.folder = os.path.dirname(htcfilename)
-        
-        self.filename = os.path.basename(htcfilename)
-        self.htcFile = HTCFile(os.path.join(self.exepath, htcfilename), self.exepath)
-        self.time_stop = self.htcFile.simulation.time_stop[0]
-        self.hawc2exe = hawc2exe
-        self.copy_turbulence = copy_turbulence
-        self.simulation_id = (htcfilename.replace("\\","/").replace("/", "_")[:50]+ "_%d" % id(self))
-        if self.simulation_id.startswith("input_"):
-            self.simulation_id = self.simulation_id[6:]
-        self.stdout_filename = fmt_path(os.path.join(os.path.relpath(self.exepath, self.modelpath), 
-                                                     (os.path.splitext(htcfilename)[0] + ".out").replace('htc', 'stdout', 1)))
-        if self.ios:
-            assert self.stdout_filename.startswith("input/")
-            self.stdout_filename = self.stdout_filename.replace("input/", "../output/")
-        #self.stdout_filename = "stdout/%s.out" % self.simulation_id
-        if 'logfile' in self.htcFile.simulation:
-            self.log_filename = self.htcFile.simulation.logfile[0]
-        else:
-            self.log_filename = self.stdout_filename
-        if os.path.isabs(self.log_filename):
-            self.log_filename = os.path.relpath(self.log_filename, self.modelpath)
-        else:
-            self.log_filename = os.path.relpath(self.log_filename)
-        self.log_filename = fmt_path(self.log_filename)
-        self.logFile = LogFile(os.path.join(self.exepath, self.log_filename), self.time_stop)
-        self.logFile.clear()
-        self.last_status = self.status
-        self.errors = []
-        self.non_blocking_simulation_thread = Thread(target=self.simulate_distributed)
-        self.updateSimStatusThread = UpdateSimStatusThread(self)
-        from wetb.hawc2.simulation_resources import LocalSimulationHost
-        self.host = LocalSimulationHost(self)
-        self.stdout = ""
-        self.returncode = 0
-
-    def start(self, update_interval=1):
-        """Start non blocking distributed simulation"""
-        self.is_simulating = True
-        if update_interval > 0:
-            self.updateSimStatusThread.interval = update_interval
-            self.updateSimStatusThread.start()
-        self.non_blocking_simulation_thread.start()
-
-    def wait(self):
-        self.non_blocking_simulation_thread.join()
-        self.update_status()
-
-    def abort(self, update_status=True):
-        self.status = ABORTED
-        self.is_simulating = False
-        self.is_done = True
-        self.host.stop()
-        if update_status:
-            self.update_status()
-
-#         if self.status != QUEUED:
-#             self.host.stop()
-#             for _ in range(50):
-#                 if self.is_simulating is False:
-#                     break
-#                 time.sleep(0.1)
-#         if self.logFile.status not in [log_file.DONE]:
-#             self.status = ABORTED
-#         self.is_simulating = False
-#         self.is_done = True
-#         if update_status:
-#             self.update_status()
-
-    def show_status(self):
-        #print ("log status:", self.logFile.status)
-        if self.logFile.status == log_file.SIMULATING:
-            if self.last_status != log_file.SIMULATING:
-                print ("|" + ("-"*50) + "|" + ("-"*49) + "|")
-                sys.stdout.write("|")
-            sys.stdout.write("."*(self.logFile.pct - getattr(self, 'last_pct', 0)))
-            sys.stdout.flush()
-            self.last_pct = self.logFile.pct
-        elif self.last_status == log_file.SIMULATING:
-            sys.stdout.write("."*(100 - self.last_pct) + "|")
-            sys.stdout.flush()
-            print ("\n")
-        elif self.logFile.status == log_file.UNKNOWN:
-            print (self.status)
-        else:
-            print (self.logFile.status)
-        if self.logFile.status != log_file.SIMULATING:
-            if self.logFile.errors:
-                print (self.logFile.errors)
-        self.last_status = self.logFile.status
-
-
-
-
-    def prepare_simulation(self):
-        self.status = PREPARING
-        self.tmp_modelpath = os.path.join(".hawc2launcher/%s/" % self.simulation_id)
-        self.tmp_exepath = os.path.join(self.tmp_modelpath, os.path.relpath(self.exepath, self.modelpath) ) + "/"
-        self.set_id(self.simulation_id, str(self.host), self.tmp_modelpath)
-
-        def fmt(src):
-            if os.path.isabs(src):
-                src = os.path.relpath(os.path.abspath(src), self.exepath)
-            else:
-                src = os.path.relpath (src)
-            assert not src.startswith(".."), "%s referes to a file outside the model path\nAll input files be inside model path" % src
-            return src
-        if self.ios:
-            input_folder_files = []
-            for root, _, filenames in os.walk(os.path.join(self.modelpath, "input/")):
-                for filename in filenames:
-                    input_folder_files.append(os.path.join(root, filename))
-             
-            input_patterns = [fmt(src) for src in input_folder_files + ([], self.htcFile.turbulence_files())[self.copy_turbulence] + self.additional_files().get('input', [])]
-        else:
-            input_patterns = [fmt(src) for src in self.htcFile.input_files() + ([], self.htcFile.turbulence_files())[self.copy_turbulence] + self.additional_files().get('input', [])]
-        input_files = set([f for pattern in input_patterns for f in glob.glob(os.path.join(self.exepath, pattern)) if os.path.isfile(f) and ".hawc2launcher" not in f])
-        if not os.path.isdir(os.path.dirname(self.exepath + self.stdout_filename)):
-            os.makedirs(os.path.dirname(self.exepath + self.stdout_filename))
-        self.host._prepare_simulation(input_files)
-
-
-#        return [fmt(src) for src in self.htcFile.input_files() + self.htcFile.turbulence_files() + self.additional_files().get('input', [])]
-#
-#        for src in self._input_sources():
-#            for src_file in glob.glob(os.path.join(self.modelpath, src)):
-#
-#
-#        self.host._prepare_simulation()
-
-    def simulate(self):
-        #starts blocking simulation
-        
-        self.is_simulating = True
-        self.errors = []
-        self.status = INITIALIZING
-        self.logFile.clear()
-        self.host._simulate()
-        self.returncode, self.stdout = self.host.returncode, self.host.stdout
-        if self.host.returncode or 'error' in self.host.stdout.lower():
-            if self.status==ABORTED:
-                return
-            if "error" in self.host.stdout.lower():
-                self.errors = (list(set([l for l in self.host.stdout.split("\n") if 'error' in l.lower()])))
-            self.status = ERROR
-        if  'HAWC2MB version:' not in self.host.stdout:
-            self.errors.append(self.host.stdout)
-            self.status = ERROR
-
-        self.logFile.update_status()
-        self.errors.extend(list(set(self.logFile.errors)))
-        self.update_status()
-        self.is_simulating = False
-        if self.host.returncode or self.errors:
-            raise Exception("Simulation error:\nReturn code: %d\n%s" % (self.host.returncode, "\n".join(self.errors)))
-        elif self.logFile.status != log_file.DONE or self.logFile.errors:
-            raise Warning("Simulation succeded with errors:\nLog status:%s\nErrors:\n%s" % (self.logFile.status, "\n".join(self.logFile.errors)))
-        else:
-            self.status = FINISH
-
-
-    def finish_simulation(self):
-        if self.status == ABORTED:
-            return
-        if self.status != ERROR:
-            self.status = FINISHING
-
-        def fmt(dst):
-            if os.path.isabs(dst):
-                dst = os.path.relpath(os.path.abspath(dst), self.exepath)
-            else:
-                dst = os.path.relpath (dst)
-            dst = fmt_path(dst)
-            assert not os.path.relpath(os.path.join(self.exepath, dst), self.modelpath).startswith(".."), "%s referes to a file outside the model path\nAll input files be inside model path" % dst
-            return dst
-        turb_files = [f for f in self.htcFile.turbulence_files() if self.copy_turbulence and not os.path.isfile(os.path.join(self.exepath, f))]
-        if self.ios:
-            output_patterns = ["../output/*", "../output/"] + turb_files + [os.path.join(self.exepath, self.stdout_filename)]
-        else:
-            output_patterns = self.htcFile.output_files() + turb_files + [os.path.join(self.exepath, self.stdout_filename)]
-        output_files = set(self.host.glob([fmt_path(os.path.join(self.tmp_exepath,fmt(p))) for p in output_patterns], recursive=self.ios))
-            
-            
-        try:
-            self.host._finish_simulation(output_files)
-            if self.status != ERROR:
-                self.status = CLEANED
-        except Warning as e:
-            self.errors.append(str(e))
-            if self.status != ERROR:
-                self.status = CLEANED
-        except Exception as e:
-            self.errors.append(str(e))
-            raise
-            
-        finally:
-            self.set_id(self.filename)
-            self.logFile.reset()
-            self.htcFile.reset()
-
-
-
-
-    def update_status(self, *args, **kwargs):
-        self.host.update_logFile_status()
-        if self.status in [INITIALIZING, SIMULATING]:
-            if self.logFile.status == log_file.SIMULATING:
-                self.status = SIMULATING
-            if self.logFile.status == log_file.DONE and self.is_simulating is False:
-                self.status = FINISH
-
-
-
-    def __str__(self):
-        return "Simulation(%s)" % self.filename
-
-
-
-
-    def additional_files(self):
-        additional_files_file = os.path.join(self.modelpath, 'additional_files.txt')
-        additional_files = {}
-        if os.path.isfile(additional_files_file):
-            with open(additional_files_file, encoding='utf-8') as fid:
-                additional_files = json.loads(fid.read().replace("\\", "/"))
-        return additional_files
-
-    def add_additional_input_file(self, file):
-        additional_files = self.additional_files()
-        additional_files['input'] = list(set(additional_files.get('input', []) + [file]))
-        additional_files_file = os.path.join(self.modelpath, 'additional_files.txt')
-        with open(additional_files_file, 'w', encoding='utf-8') as fid:
-                json.dump(additional_files, fid)
-
-
-    def simulate_distributed(self):
-        try:
-            self.prepare_simulation()
-            try:
-                self.simulate()
-            except Warning as e:
-                print ("simulation failed", str(self))
-                print ("Trying to finish")
-                raise
-            finally:
-                try:
-                    if self.status!=ABORTED:
-                        self.finish_simulation()
-                except:
-                    print ("finish_simulation failed", str(self))
-                    raise
-        except Exception as e:
-            self.status = ERROR
-            self.errors.append(str(e))
-            raise e
-        finally:
-            self.is_done = True
-
-
-
-
-    def fix_errors(self):
-        def confirm_add_additional_file(folder, file):
-            if os.path.isfile(os.path.join(self.modelpath, folder, file)):
-                filename = fmt_path(os.path.join(folder, file))
-                if self.get_confirmation("File missing", "'%s' seems to be missing in the temporary working directory. \n\nDo you want to add it to additional_files.txt" % filename):
-                    self.add_additional_input_file(filename)
-                    self.show_message("'%s' is now added to additional_files.txt.\n\nPlease restart the simulation" % filename)
-        for error in self.errors:
-            for regex in [r".*\*\*\* ERROR \*\*\* File '(.*)' does not exist in the (.*) folder",
-                          r".*\*\*\* ERROR \*\*\* DLL (.*)()"]:
-                m = re.compile(regex).match(error.strip())
-                if m is not None:
-                    file, folder = m.groups()
-                    confirm_add_additional_file(folder, file)
-                    continue
-            m = re.compile(r".*\*\*\* ERROR \*\*\* File '(.*)' does not exist in the working directory").match(error.strip())
-            if m is not None:
-                file = m.groups()[0]
-                for root, folder, files in os.walk(self.modelpath):
-                    if "__Thread" not in root and file in files:
-                        folder = os.path.relpath(root, self.modelpath)
-                        confirm_add_additional_file(folder, file)
-                continue
-
-    def get_confirmation(self, title, msg):
-        """override in subclass"""
-        return True
-
-    def show_message(self, msg, title="Information"):
-        print (msg)
-
-    def set_id(self, *args, **kwargs):
-        pass
-
-    def progress_callback(self,*args, **kwargs):
-        pass
-
-class UpdateSimStatusThread(Thread):
-    def __init__(self, simulation, interval=1):
-        Thread.__init__(self)
-        self.simulation = simulation
-        self.interval = interval
-
-    def start(self):
-        Thread.start(self)
-
-    def run(self):
-        while self.simulation.is_done is False:
-            self.simulation.update_status()
-            time.sleep(0.5)
-            t = time.time()
-            while self.simulation.is_simulating and time.time() < t + self.interval:
-                time.sleep(1)
+
+import glob
+import json
+import os
+import re
+import sys
+from threading import Thread
+import time
+
+from wetb.hawc2 import log_file
+from wetb.hawc2.htc_file import HTCFile, fmt_path
+from wetb.hawc2.log_file import LogFile
+import tempfile
+import stat
+import shutil
+
+
+QUEUED = "queued"  # until start
+PREPARING = "Copy to host"  # during prepare simulation
+INITIALIZING = "Initializing"  # when starting
+SIMULATING = "Simulating"  # when logfile.status=simulating
+FINISHING = "Copy from host"  # during prepare simulation
+FINISH = "Simulation finish"  # when HAWC2 finish
+ERROR = "Error"  # when hawc2 returns error
+ABORTED = "Aborted"  # when stopped and logfile.status != Done
+CLEANED = "Cleaned"  # after copy back
+
+
+class Simulation(object):
+    """Class for doing hawc2 simulations
+
+
+
+
+    Examples
+    --------
+    >>> sim = Simulation("<path>/MyModel","htc/MyHtc")
+
+    Blocking inplace simulation
+    >>> sim.simulate()
+
+    Non-blocking distributed simulation(everything copied to temporary folder)\n
+    Starts two threads:
+    - non_blocking_simulation_thread:
+        - prepare_simulation() # copy to temporary folder
+        - simulate() # simulate
+        - finish_simulation # copy results back again
+    - updateSimStatusThread:
+        - update status every second
+    >>> sim.start()
+    >>> while sim.status!=CLEANED:
+    >>>     sim.show_status()
+
+
+    The default host is LocalSimulationHost. To simulate on pbs featured cluster
+    >>> sim.host = PBSClusterSimulationHost(sim, <hostname>, <username>, <password>)
+    >>> sim.start()
+    >>> while sim.status!=CLEANED:
+    >>>     sim.show_status()
+    """
+
+    is_simulating = False
+    is_done = False
+    status = QUEUED
+
+    def __init__(self, modelpath, htcfile, hawc2exe="HAWC2MB.exe", copy_turbulence=True):
+        if isinstance(htcfile, HTCFile):
+            htcfilename = htcfile.filename
+        else:
+            htcfilename = htcfile
+
+        self.modelpath = os.path.abspath(modelpath).replace("\\", "/")
+        if self.modelpath[-1] != '/':
+            self.modelpath += "/"
+        if os.path.isabs(htcfilename):
+            htcfilename = os.path.relpath(htcfilename, modelpath)
+        if htcfilename.startswith("input/"):
+            htcfilename = htcfilename[6:]
+        exists = [os.path.isfile(os.path.join(modelpath, htcfilename)),
+                  os.path.isfile(os.path.join(modelpath, "input/", htcfilename))]
+        if all(exists):
+            raise Exception(
+                "Both standard and input/output file structure available for %s in %s. Delete one of the options" % (htcfilename, modelpath))
+        if not any(exists):
+            raise Exception("%s not found in %s" % (htcfilename, modelpath))
+        else:
+            self.ios = exists[1]  # input/output file structure
+
+        if self.ios:
+            self.exepath = self.modelpath + "input/"
+        else:
+            self.exepath = self.modelpath
+        # model_path: top level path containing all resources
+        # exepath: parent path for relative paths
+
+        htcfilename = fmt_path(htcfilename)
+
+        self.tmp_modelpath = self.exepath
+        self.folder = os.path.dirname(htcfilename)
+
+        self.filename = os.path.basename(htcfilename)
+        if isinstance(htcfile, HTCFile):
+            self.htcFile = htcfile
+        else:
+            self.htcFile = HTCFile(os.path.join(self.exepath, htcfilename), self.exepath)
+        self.time_stop = self.htcFile.simulation.time_stop[0]
+        self.hawc2exe = hawc2exe
+        self.copy_turbulence = copy_turbulence
+        self.simulation_id = (htcfilename.replace("\\", "/").replace("/", "_")[:50] + "_%d" % id(self))
+        if self.simulation_id.startswith("input_"):
+            self.simulation_id = self.simulation_id[6:]
+        self.stdout_filename = fmt_path(os.path.join(os.path.relpath(self.exepath, self.modelpath),
+                                                     (os.path.splitext(htcfilename)[0] + ".out").replace('htc', 'stdout', 1)))
+        if self.ios:
+            assert self.stdout_filename.startswith("input/")
+            self.stdout_filename = self.stdout_filename.replace("input/", "../output/")
+        # self.stdout_filename = "stdout/%s.out" % self.simulation_id
+        if 'logfile' in self.htcFile.simulation:
+            self.log_filename = self.htcFile.simulation.logfile[0]
+        else:
+            self.log_filename = self.stdout_filename
+        if os.path.isabs(self.log_filename):
+            self.log_filename = os.path.relpath(self.log_filename, self.modelpath)
+        else:
+            self.log_filename = os.path.relpath(self.log_filename)
+        self.log_filename = fmt_path(self.log_filename)
+        self.logFile = LogFile(os.path.join(self.exepath, self.log_filename), self.time_stop)
+        self.logFile.clear()
+        self.last_status = self.status
+        self.errors = []
+        from wetb.hawc2.simulation_resources import LocalSimulationHost
+        self.host = LocalSimulationHost(self)
+        self.stdout = ""
+        self.returncode = 0
+
+    def start(self, update_interval=1):
+        """Start non blocking distributed simulation"""
+        self.is_simulating = True
+        if update_interval > 0:
+            self.updateSimStatusThread = UpdateSimStatusThread(self)
+            self.updateSimStatusThread.interval = update_interval
+            self.updateSimStatusThread.start()
+        self.non_blocking_simulation_thread = Thread(
+            target=lambda: self.simulate_distributed(raise_simulation_errors=False))
+        self.non_blocking_simulation_thread.start()
+
+    def wait(self):
+        self.non_blocking_simulation_thread.join()
+        self.update_status()
+
+    def abort(self, update_status=True):
+        self.status = ABORTED
+        self.is_simulating = False
+        self.is_done = True
+        self.host.stop()
+        if update_status:
+            self.update_status()
+
+#         if self.status != QUEUED:
+#             self.host.stop()
+#             for _ in range(50):
+#                 if self.is_simulating is False:
+#                     break
+#                 time.sleep(0.1)
+#         if self.logFile.status not in [log_file.DONE]:
+#             self.status = ABORTED
+#         self.is_simulating = False
+#         self.is_done = True
+#         if update_status:
+#             self.update_status()
+
+    def show_status(self):
+        # print ("log status:", self.logFile.status)
+        if self.logFile.status == log_file.SIMULATING:
+            if self.last_status != log_file.SIMULATING:
+                print("|" + ("-" * 50) + "|" + ("-" * 49) + "|")
+                sys.stdout.write("|")
+            sys.stdout.write("." * (self.logFile.pct - getattr(self, 'last_pct', 0)))
+            sys.stdout.flush()
+            self.last_pct = self.logFile.pct
+        elif self.last_status == log_file.SIMULATING:
+            sys.stdout.write("." * (100 - self.last_pct) + "|")
+            sys.stdout.flush()
+            print("\n")
+        elif self.logFile.status == log_file.UNKNOWN:
+            print(self.status)
+        else:
+            print(self.logFile.status)
+        if self.logFile.status != log_file.SIMULATING:
+            if self.logFile.errors:
+                print(self.logFile.errors)
+        self.last_status = self.logFile.status
+
+    def prepare_simulation(self):
+        self.status = PREPARING
+        # self.tmp_modelpath = os.path.join(".hawc2launcher/%s/" % self.simulation_id)
+        # self.tmp_exepath = os.path.join(self.tmp_modelpath, os.path.relpath(self.exepath, self.modelpath) ) + "/"
+        self.set_id(self.simulation_id, str(self.host), self.tmp_modelpath)
+
+        def fmt(src):
+            if os.path.isabs(src):
+                src = os.path.relpath(os.path.abspath(src), self.exepath)
+            else:
+                src = os.path.relpath(src)
+            assert not src.startswith(
+                ".."), "%s referes to a file outside the model path\nAll input files be inside model path" % src
+            return src
+        if self.ios:
+            input_folder_files = []
+            for root, _, filenames in os.walk(os.path.join(self.modelpath, "input/")):
+                for filename in filenames:
+                    input_folder_files.append(os.path.join(root, filename))
+
+            input_patterns = [fmt(src) for src in input_folder_files + ([], self.htcFile.turbulence_files())
+                              [self.copy_turbulence] + self.additional_files().get('input', [])]
+        else:
+            input_patterns = [fmt(src) for src in self.htcFile.input_files(
+            ) + ([], self.htcFile.turbulence_files())[self.copy_turbulence] + self.additional_files().get('input', [])]
+        input_files = set([f for pattern in input_patterns for f in glob.glob(
+            os.path.join(self.exepath, pattern)) if os.path.isfile(f) and ".hawc2launcher" not in f])
+
+        self.host._prepare_simulation(input_files)
+
+
+#        return [fmt(src) for src in self.htcFile.input_files() + self.htcFile.turbulence_files() + self.additional_files().get('input', [])]
+#
+#        for src in self._input_sources():
+#            for src_file in glob.glob(os.path.join(self.modelpath, src)):
+#
+#
+#        self.host._prepare_simulation()
+
+    def simulate(self, cancel_event=None):
+        # starts blocking simulation
+
+        self.is_simulating = True
+        self.errors = []
+        self.status = INITIALIZING
+        self.logFile.clear()
+        self.host._simulate(cancel_event=cancel_event)
+        self.returncode, self.stdout = self.host.returncode, self.host.stdout
+        if self.host.returncode or ('elapsed time :' not in self.host.stdout.lower()):
+            if self.status == ABORTED:
+                return
+            if "error" in self.host.stdout.lower():
+                self.errors = (list(set([l for l in self.host.stdout.split(
+                    "\n") if '*** error ***' in l.lower() and not 'rms error' in l])))
+            self.status = ERROR
+        if 'HAWC2MB version:' not in self.host.stdout and 'Build information for HAWC2MB' not in self.host.stdout:
+            self.errors.append(self.host.stdout)
+            self.status = ERROR
+
+        self.logFile.update_status()
+        self.errors.extend(list(set(self.logFile.errors)))
+        self.update_status()
+        self.is_simulating = False
+        if self.host.returncode or self.errors:
+            raise Exception("Simulation error:\nReturn code: %d\n%s" % (self.host.returncode, "\n".join(self.errors)))
+        elif self.logFile.status != log_file.DONE or self.logFile.errors:
+            raise Warning("Simulation succeded with errors:\nLog status:%s\nErrors:\n%s" %
+                          (self.logFile.status, "\n".join(self.logFile.errors)))
+        else:
+            self.status = FINISH
+
+    def finish_simulation(self):
+        if self.status == ABORTED:
+            return
+        if self.status != ERROR:
+            self.status = FINISHING
+
+        def fmt(dst):
+            if os.path.isabs(dst):
+                dst = os.path.relpath(os.path.abspath(dst), self.exepath)
+            else:
+                dst = os.path.relpath(dst)
+            dst = fmt_path(dst)
+            assert not os.path.relpath(os.path.join(self.exepath, dst), self.modelpath).startswith(
+                ".."), "%s referes to a file outside the model path\nAll input files be inside model path" % dst
+            return dst
+        turb_files = [f for f in self.htcFile.turbulence_files()
+                      if self.copy_turbulence and not os.path.isfile(os.path.join(self.exepath, f))]
+        if self.ios:
+            output_patterns = ["../output/*", "../output/"] + turb_files + \
+                [os.path.join(self.exepath, self.stdout_filename)]
+        else:
+            output_patterns = self.htcFile.output_files() + turb_files + \
+                [os.path.join(self.exepath, self.stdout_filename)]
+        output_files = set(self.host.glob([fmt_path(os.path.join(self.tmp_exepath, fmt(p)))
+                                           for p in output_patterns], recursive=self.ios))
+        if not os.path.isdir(os.path.dirname(self.exepath + self.stdout_filename)):
+            os.makedirs(os.path.dirname(self.exepath + self.stdout_filename))
+
+        try:
+            self.host._finish_simulation(output_files)
+            if self.status != ERROR:
+                self.status = CLEANED
+        except Warning as e:
+            self.errors.append(str(e))
+            if self.status != ERROR:
+                self.status = CLEANED
+        except Exception as e:
+            self.errors.append(str(e))
+            raise
+
+        finally:
+            self.set_id(self.filename)
+            self.logFile.reset()
+            self.htcFile.reset()
+
+    def update_status(self, *args, **kwargs):
+        self.host.update_logFile_status()
+        if self.status in [INITIALIZING, SIMULATING]:
+            if self.logFile.status == log_file.SIMULATING:
+                self.status = SIMULATING
+            if self.logFile.status == log_file.DONE and self.is_simulating is False:
+                self.status = FINISH
+
+    def __str__(self):
+        return "Simulation(%s)" % self.filename
+
+    def additional_files(self):
+        additional_files_file = os.path.join(self.modelpath, 'additional_files.txt')
+        additional_files = {}
+        if os.path.isfile(additional_files_file):
+            with open(additional_files_file, encoding='utf-8') as fid:
+                additional_files = json.loads(fid.read().replace("\\", "/"))
+        return additional_files
+
+    def add_additional_input_file(self, file):
+        additional_files = self.additional_files()
+        additional_files['input'] = list(set(additional_files.get('input', []) + [file]))
+        additional_files_file = os.path.join(self.modelpath, 'additional_files.txt')
+        with open(additional_files_file, 'w', encoding='utf-8') as fid:
+            json.dump(additional_files, fid)
+
+    def simulate_distributed(self, raise_simulation_errors=True):
+        try:
+            with tempfile.TemporaryDirectory(prefix="h2launcher_%s_" % os.path.basename(self.filename)) as tmpdirname:
+                self.tmp_modelpath = tmpdirname + "/"
+                self.tmp_exepath = os.path.join(self.tmp_modelpath, os.path.relpath(self.exepath, self.modelpath)) + "/"
+                self.prepare_simulation()
+                try:
+                    self.simulate()
+                except Warning as e:
+                    print("simulation failed", str(self))
+                    print("Trying to finish")
+                    raise
+                finally:
+                    try:
+                        if self.status != ABORTED:
+                            self.finish_simulation()
+                    except Exception:
+                        print("finish_simulation failed", str(self))
+                        raise
+#                     finally:
+#                         def remove_readonly(fn, path, excinfo):
+#                             try:
+#                                 os.chmod(path, stat.S_IWRITE)
+#                                 fn(path)
+#                             except Exception as exc:
+#                                 print("Skipped:", path, "because:\n", exc)
+#
+#                         shutil.rmtree(tmpdirname, onerror=remove_readonly)
+
+        except Exception as e:
+            self.status = ERROR
+            self.errors.append(str(e))
+            if raise_simulation_errors:
+                raise e
+        finally:
+            self.is_done = True
+
+    def fix_errors(self):
+        def confirm_add_additional_file(folder, file):
+            if os.path.isfile(os.path.join(self.modelpath, folder, file)):
+                filename = fmt_path(os.path.join(folder, file))
+                if self.get_confirmation(
+                        "File missing", "'%s' seems to be missing in the temporary working directory. \n\nDo you want to add it to additional_files.txt" % filename):
+                    self.add_additional_input_file(filename)
+                    self.show_message(
+                        "'%s' is now added to additional_files.txt.\n\nPlease restart the simulation" % filename)
+        for error in self.errors:
+            for regex in [r".*\*\*\* ERROR \*\*\* File '(.*)' does not exist in the (.*) folder",
+                          r".*\*\*\* ERROR \*\*\* DLL (.*)()"]:
+                m = re.compile(regex).match(error.strip())
+                if m is not None:
+                    file, folder = m.groups()
+                    confirm_add_additional_file(folder, file)
+                    continue
+            m = re.compile(r".*\*\*\* ERROR \*\*\* File '(.*)' does not exist in the working directory").match(error.strip())
+            if m is not None:
+                file = m.groups()[0]
+                for root, folder, files in os.walk(self.modelpath):
+                    if "__Thread" not in root and file in files:
+                        folder = os.path.relpath(root, self.modelpath)
+                        confirm_add_additional_file(folder, file)
+                continue
+
+    def get_confirmation(self, title, msg):
+        """override in subclass"""
+        return True
+
+    def show_message(self, msg, title="Information"):
+        print(msg)
+
+    def set_id(self, *args, **kwargs):
+        pass
+
+    def progress_callback(self, *args, **kwargs):
+        pass
+
+
+class UpdateSimStatusThread(Thread):
+    def __init__(self, simulation, interval=1):
+        Thread.__init__(self)
+        self.simulation = simulation
+        self.interval = interval
+
+    def start(self):
+        Thread.start(self)
+
+    def run(self):
+        while self.simulation.is_done is False:
+            self.simulation.update_status()
+            time.sleep(0.5)
+            t = time.time()
+            while self.simulation.is_simulating and time.time() < t + self.interval:
+                time.sleep(1)

wetb/hawc2/simulation_resources.py

@@ -12,6 +12,7 @@ from subprocess import STDOUT
 import subprocess
 from threading import Thread
 import time
+from platform import system
 
 from wetb.hawc2 import log_file
 from wetb.hawc2.log_file import LogInfo, LogFile
@@ -24,36 +25,37 @@ from wetb.utils.cluster_tools.pbsjob import SSHPBSJob, NOT_SUBMITTED, DONE
 from wetb.hawc2.htc_file import fmt_path
 import numpy as np
 
+
 class SimulationHost(object):
     def __init__(self, simulation):
         self.sim = simulation
-    logFile = property(lambda self : self.sim.logFile, lambda self, v: setattr(self.sim, "logFile", v))
-    errors = property(lambda self : self.sim.errors)
-    modelpath = property(lambda self : self.sim.modelpath)
-    exepath = property(lambda self : self.sim.exepath)
-    tmp_modelpath = property(lambda self : self.sim.tmp_modelpath, lambda self, v: setattr(self.sim, "tmp_modelpath", v))
-    tmp_exepath = property(lambda self : self.sim.tmp_exepath, lambda self, v: setattr(self.sim, "tmp_exepath", v))
-    simulation_id = property(lambda self : self.sim.simulation_id)
-    stdout_filename = property(lambda self : self.sim.stdout_filename)
-    htcFile = property(lambda self : self.sim.htcFile)
-    additional_files = property(lambda self : self.sim.additional_files)
-    _input_sources = property(lambda self : self.sim._input_sources)
-    _output_sources = property(lambda self : self.sim._output_sources)
-    log_filename = property(lambda self : self.sim.log_filename)
-
-    status = property(lambda self : self.sim.status, lambda self, v: setattr(self.sim, "status", v))
-    is_simulating = property(lambda self : self.sim.is_simulating, lambda self, v: setattr(self.sim, "is_simulating", v))
+    logFile = property(lambda self: self.sim.logFile, lambda self, v: setattr(self.sim, "logFile", v))
+    errors = property(lambda self: self.sim.errors)
+    modelpath = property(lambda self: self.sim.modelpath)
+    exepath = property(lambda self: self.sim.exepath)
+    tmp_modelpath = property(lambda self: self.sim.tmp_modelpath, lambda self, v: setattr(self.sim, "tmp_modelpath", v))
+    tmp_exepath = property(lambda self: self.sim.tmp_exepath, lambda self, v: setattr(self.sim, "tmp_exepath", v))
+    simulation_id = property(lambda self: self.sim.simulation_id)
+    stdout_filename = property(lambda self: self.sim.stdout_filename)
+    htcFile = property(lambda self: self.sim.htcFile)
+    additional_files = property(lambda self: self.sim.additional_files)
+    _input_sources = property(lambda self: self.sim._input_sources)
+    _output_sources = property(lambda self: self.sim._output_sources)
+    log_filename = property(lambda self: self.sim.log_filename)
+
+    status = property(lambda self: self.sim.status, lambda self, v: setattr(self.sim, "status", v))
+    is_simulating = property(lambda self: self.sim.is_simulating, lambda self, v: setattr(self.sim, "is_simulating", v))
 
     def __str__(self):
         return self.resource.host
-    
+
+
 class LocalSimulationHost(SimulationHost):
     def __init__(self, simulation, resource=None):
         SimulationHost.__init__(self, simulation)
         if resource is None:
             resource = LocalResource(1)
         self.resource = resource
-        self.simulationThread = SimulationThread(self.sim)
 
     def get_datetime(self):
         return datetime.now()
@@ -68,36 +70,41 @@ class LocalSimulationHost(SimulationHost):
 
     def _prepare_simulation(self, input_files):
         # must be called through simulation object
-        self.tmp_modelpath = os.path.join(self.modelpath, self.tmp_modelpath)
-        self.tmp_exepath = os.path.join(self.tmp_modelpath, os.path.relpath(self.sim.exepath, self.sim.modelpath) ) + "/"
+        # self.tmp_modelpath = os.path.join(self.modelpath, self.tmp_modelpath)
+        # self.tmp_exepath = os.path.join(self.tmp_modelpath, os.path.relpath(self.sim.exepath, self.sim.modelpath) ) + "/"
         self.sim.set_id(self.simulation_id, 'Localhost', self.tmp_modelpath)
         for src_file in input_files:
             dst = os.path.join(self.tmp_modelpath, os.path.relpath(src_file, self.modelpath))
             # exist_ok does not exist in Python27
             if not os.path.exists(os.path.dirname(dst)):
-                os.makedirs(os.path.dirname(dst))  #, exist_ok=True)
+                os.makedirs(os.path.dirname(dst))  # , exist_ok=True)
             shutil.copy(src_file, dst)
             if not os.path.isfile(dst) or os.stat(dst).st_size != os.stat(src_file).st_size:
-                print ("error copy ", dst)
+                print("error copy ", dst)
 
         stdout_folder = os.path.join(self.tmp_exepath, os.path.dirname(self.sim.stdout_filename))
         if not os.path.exists(stdout_folder):
-            os.makedirs(stdout_folder)  #, exist_ok=True)
+            os.makedirs(stdout_folder)  # , exist_ok=True)
         self.logFile.filename = os.path.join(self.tmp_exepath, self.log_filename)
-        self.simulationThread.modelpath = self.tmp_modelpath
-        self.simulationThread.exepath = self.tmp_exepath
-
 
-    def _simulate(self):
-        #must be called through simulation object
+    def _simulate(self, cancel_event=None):
+        # must be called through simulation object
         self.returncode, self.stdout = 1, "Simulation failed"
+        self.simulationThread = SimulationThread(self.sim)
         self.simulationThread.start()
-        self.sim.set_id(self.sim.simulation_id, "Localhost(pid:%d)" % self.simulationThread.process.pid, self.tmp_modelpath)
+        self.sim.set_id(self.sim.simulation_id, "Localhost(pid:%d)" %
+                        self.simulationThread.process.pid, self.tmp_modelpath)
+        if cancel_event is not None:
+            while self.simulationThread.is_alive():
+                if cancel_event.isSet():
+                    self.simulationThread.stop()
+                    break
+                time.sleep(0.1)
         self.simulationThread.join()
         self.returncode, self.stdout = self.simulationThread.res
         self.logFile.update_status()
         self.errors.extend(list(set(self.logFile.errors)))
-
+        self.simulationThread = None
 
     def _finish_simulation(self, output_files):
         missing_result_files = []
@@ -106,39 +113,37 @@ class LocalSimulationHost(SimulationHost):
             # exist_ok does not exist in Python27
             try:
                 if not os.path.isdir(os.path.dirname(dst_file)):
-                    os.makedirs(os.path.dirname(dst_file))  #, exist_ok=True)
+                    os.makedirs(os.path.dirname(dst_file))  # , exist_ok=True)
                 if not os.path.isfile(dst_file) or os.path.getmtime(dst_file) != os.path.getmtime(src_file):
                     shutil.copy(src_file, dst_file)
-            except:
+            except Exception:
                 missing_result_files.append(dst_file)
 
         self.logFile.filename = os.path.join(self.sim.exepath, self.log_filename)
         if missing_result_files:
-            raise Warning("Failed to copy %s from %s"%(",".join(missing_result_files), self.resource.host))
-        try:
-            shutil.rmtree(self.tmp_modelpath, ignore_errors=False)
-        except (PermissionError, OSError) as e:
-            try:
-                #change permissions and delete
-                for root, folders, files in os.walk(self.tmp_modelpath):  
-                    for folder in folders:  
-                        os.chmod(os.path.join(root, folder), 0o666)
-                    for file in files:
-                        os.chmod(os.path.join(root, file), 0o666)
-                shutil.rmtree(self.tmp_modelpath)
-            except (PermissionError, OSError) as e:
-                raise Warning("Fail to remove temporary files and folders on %s\n%s"%(self.resource.host, str(e)))
-
+            raise Warning("Failed to copy %s from %s" % (",".join(missing_result_files), self.resource.host))
+#         try:
+#             shutil.rmtree(self.tmp_modelpath, ignore_errors=False)
+#         except (PermissionError, OSError) as e:
+#             try:
+#                 # change permissions and delete
+#                 for root, folders, files in os.walk(self.tmp_modelpath):
+#                     for folder in folders:
+#                         os.chmod(os.path.join(root, folder), 0o666)
+#                     for file in files:
+#                         os.chmod(os.path.join(root, file), 0o666)
+#                 shutil.rmtree(self.tmp_modelpath)
+#             except (PermissionError, OSError) as e:
+#                 raise Warning("Fail to remove temporary files and folders on %s\n%s" % (self.resource.host, str(e)))
 
     def update_logFile_status(self):
         self.logFile.update_status()
 
     def stop(self):
-            self.simulationThread.stop()
-            if self.simulationThread.is_alive():
-                self.simulationThread.join()
-            print ("simulatino_resources.stop joined")
-
+        self.simulationThread.stop()
+        if self.simulationThread.is_alive():
+            self.simulationThread.join()
+        print("simulatino_resources.stop joined")
 
 
 class SimulationThread(Thread):
@@ -151,41 +156,47 @@ class SimulationThread(Thread):
         self.res = [0, "", ""]
         self.low_priority = low_priority
 
-
     def start(self):
         CREATE_NO_WINDOW = 0x08000000
-        exepath = self.exepath #overwritten in _prepare_simulation
-        modelpath = self.modelpath #overwritten in _prepare_simulation
+        exepath = self.exepath  # overwritten in _prepare_simulation
+        modelpath = self.modelpath  # overwritten in _prepare_simulation
         htcfile = os.path.relpath(self.sim.htcFile.filename, self.sim.exepath)
-        
+
         hawc2exe = self.sim.hawc2exe
         stdout = self.sim.stdout_filename
         if not os.path.isdir(os.path.dirname(exepath + self.sim.stdout_filename)):
             os.makedirs(os.path.dirname(exepath + self.sim.stdout_filename))
-        
-        with open (os.path.join(exepath, stdout), 'wb') as stdout:
+
+        with open(os.path.join(exepath, stdout), 'wb') as stdout:
             if isinstance(hawc2exe, tuple):
                 wine, hawc2exe = hawc2exe
-                self.process = subprocess.Popen(" ".join([wine, hawc2exe, htcfile]), stdout=stdout, stderr=STDOUT, shell=True, cwd=exepath) #shell must be True to inwoke wine
+                # shell must be True to inwoke wine
+                self.process = subprocess.Popen(" ".join([wine, hawc2exe, htcfile]),
+                                                stdout=stdout, stderr=STDOUT, shell=True, cwd=exepath)
             else:
-                self.process = subprocess.Popen([hawc2exe, htcfile], stdout=stdout, stderr=STDOUT, shell=False, cwd=exepath, creationflags=CREATE_NO_WINDOW)
-            #self.process.communicate()
+                if system().lower() == "linux":
+                    self.process = subprocess.Popen([hawc2exe, htcfile], stdout=stdout,
+                                                    stderr=STDOUT, shell=False, cwd=exepath)
+                else:
+                    self.process = subprocess.Popen([hawc2exe, htcfile], stdout=stdout,
+                                                    stderr=STDOUT, shell=False, cwd=exepath, creationflags=CREATE_NO_WINDOW)
+                    
+            # self.process.communicate()
 
         import psutil
         try:
             self.sim.host.resource.process_name = psutil.Process(self.process.pid).name()
-        except:
+        except Exception:
             pass
         Thread.start(self)
 
-
     def run(self):
         import psutil
         p = psutil.Process(os.getpid())
         try:
             if self.low_priority:
                 p.set_nice(psutil.BELOW_NORMAL_PRIORITY_CLASS)
-        except:
+        except Exception:
             pass
         self.process.communicate()
         errorcode = self.process.returncode
@@ -198,9 +209,9 @@ class SimulationThread(Thread):
         if hasattr(self, 'process'):
             try:
                 subprocess.Popen("TASKKILL /F /PID {pid} /T".format(pid=self.process.pid))
-            except:
+            except Exception:
                 pass
-        
+
 
 class PBSClusterSimulationResource(SSHPBSClusterResource):
     def __init__(self, sshclient, min_cpu, min_free, init_cmd, wine_cmd, python_cmd, queue="workq"):
@@ -209,19 +220,18 @@ class PBSClusterSimulationResource(SSHPBSClusterResource):
         self.wine_cmd = wine_cmd
         self.python_cmd = python_cmd
         self.queue = queue
-        
-    
+
     def is_clean(self):
         return self.execute("find .hawc2launcher/ -type f | wc -l")[1] > 0
 
     def clean(self):
         try:
             self.execute('rm .hawc2launcher -r -f')
-        except:
+        except Exception:
             pass
         try:
             self.shared_ssh.close()
-        except:
+        except Exception:
             pass
 
     def update_resource_status(self):
@@ -229,31 +239,35 @@ class PBSClusterSimulationResource(SSHPBSClusterResource):
             _, out, _ = self.ssh.execute("find .hawc2launcher/ -name '*.out'")
             self.finished = set([f.split("/")[1] for f in out.split("\n") if "/" in f])
         except Exception as e:
-            print ("resource_manager.update_status, out", str(e))
+            print("resource_manager.update_status, out", str(e))
             pass
- 
+
         try:
             _, out, _ = self.ssh.execute("find .hawc2launcher -name 'status*' -exec cat {} \;")
-            self.loglines = {l.split(";")[0] : l.split(";")[1:] for l in out.split("\n") if ";" in l}
+            self.loglines = {l.split(";")[0]: l.split(";")[1:] for l in out.split("\n") if ";" in l}
         except Exception as e:
-            print ("resource_manager.update_status, status file", str(e))
+            print("resource_manager.update_status, status file", str(e))
             pass
         try:
             _, out, _ = self.ssh.execute("qstat -u %s" % self.username)
             self.is_executing = set([j.split(".")[0] for j in out.split("\n")[5:] if "." in j])
         except Exception as e:
-            print ("resource_manager.update_status, qstat", str(e))
+            print("resource_manager.update_status, qstat", str(e))
             pass
 
+
 class GormSimulationResource(PBSClusterSimulationResource):
     init_cmd = """export PATH=/home/python/miniconda3/bin:$PATH
 source activate wetb_py3"""
     queue = "workq"
     host = "gorm.risoe.dk"
+
     def __init__(self, username, password, wine_cmd="WINEARCH=win32 WINEPREFIX=~/.wine32 wine"):
 
         from wetb.utils.cluster_tools.ssh_client import SSHClient
-        PBSClusterSimulationResource.__init__(self, SSHClient(self.host, username, password, 22), 25, 100, self.init_cmd, wine_cmd, "python", self.queue)
+        PBSClusterSimulationResource.__init__(self, SSHClient(
+            self.host, username, password, 22), 25, 100, self.init_cmd, wine_cmd, "python", self.queue)
+
 
 class JessSimulationResource(PBSClusterSimulationResource):
     host = 'jess.dtu.dk'
@@ -261,12 +275,12 @@ class JessSimulationResource(PBSClusterSimulationResource):
 source activate wetb_py3
 WINEARCH=win32 WINEPREFIX=~/.wine32 winefix"""
     queue = "windq"
-    def __init__(self, username, password, wine_cmd="WINEARCH=win32 WINEPREFIX=~/.wine32 wine"):
-    
-        from wetb.utils.cluster_tools.ssh_client import SSHClient
-        PBSClusterSimulationResource.__init__(self, SSHClient(self.host, username, password, 22), 25, 600, self.init_cmd, wine_cmd, "python", self.queue)
 
+    def __init__(self, username, password, wine_cmd="WINEARCH=win32 WINEPREFIX=~/.wine32 wine"):
 
+        from wetb.utils.cluster_tools.ssh_client import SSHClient
+        PBSClusterSimulationResource.__init__(self, SSHClient(
+            self.host, username, password, 22), 25, 600, self.init_cmd, wine_cmd, "python", self.queue)
 
 
 class PBSClusterSimulationHost(SimulationHost):
@@ -276,22 +290,24 @@ class PBSClusterSimulationHost(SimulationHost):
         self.pbsjob = SSHPBSJob(resource.new_ssh_connection())
         self.resource = resource
 
-    hawc2exe = property(lambda self : os.path.basename(self.sim.hawc2exe))
+    hawc2exe = property(lambda self: os.path.basename(self.sim.hawc2exe))
 
     def glob(self, filepattern, cwd="", recursive=False):
         return self.ssh.glob(filepattern, cwd, recursive)
-    
+
     def get_datetime(self):
         v, out, err = self.ssh.execute('date "+%Y,%m,%d,%H,%M,%S"')
         if v == 0:
             return datetime.strptime(out.strip(), "%Y,%m,%d,%H,%M,%S")
 
-    #@print_time
     def _prepare_simulation(self, input_files):
         with self.ssh:
             self.ssh.execute(["mkdir -p .hawc2launcher/%s" % self.simulation_id], verbose=False)
+            self.simulationThread.modelpath = self.tmp_modelpath
+            self.simulationThread.exepath = self.tmp_exepath
             self.ssh.execute("mkdir -p %s%s" % (self.tmp_exepath, os.path.dirname(self.log_filename)))
-            self.ssh.upload_files(self.modelpath, self.tmp_modelpath, file_lst = [os.path.relpath(f, self.modelpath) for f in input_files], callback=self.sim.progress_callback("Copy to host"))
+            self.ssh.upload_files(self.modelpath, self.tmp_modelpath, file_lst=[os.path.relpath(
+                f, self.modelpath) for f in input_files], callback=self.sim.progress_callback("Copy to host"))
 #             for src_file in input_files:
 #                     dst = unix_path(self.tmp_modelpath + os.path.relpath(src_file, self.modelpath))
 #                     self.execute("mkdir -p %s" % os.path.dirname(dst), verbose=False)
@@ -303,39 +319,38 @@ class PBSClusterSimulationHost(SimulationHost):
             self.ssh.execute("mkdir -p %s%s" % (self.tmp_exepath, os.path.dirname(self.stdout_filename)))
             remote_log_filename = "%s%s" % (self.tmp_exepath, self.log_filename)
             self.ssh.execute("rm -f %s" % remote_log_filename)
-            
-    #@print_time
+
     def _finish_simulation(self, output_files):
         with self.ssh:
             download_failed = []
             try:
-                self.ssh.download_files(self.tmp_modelpath, self.modelpath, file_lst = [os.path.relpath(f, self.tmp_modelpath) for f in output_files], callback=self.sim.progress_callback("Copy from host") )
-            except:
-#                 
-#             for src_file in output_files:
-#                 try:
-#                     dst_file = os.path.join(self.modelpath, os.path.relpath(src_file, self.tmp_modelpath))
-#                     os.makedirs(os.path.dirname(dst_file), exist_ok=True)
-#                     self.download(src_file, dst_file, retry=10)
-#                 except Exception as e:
-#                     download_failed.append(dst_file)
-#             if download_failed:
-                raise Warning("Failed to download %s from %s"%(",".join(download_failed), self.ssh.host))
+                self.ssh.download_files(self.tmp_modelpath, self.modelpath, file_lst=[os.path.relpath(
+                    f, self.tmp_modelpath) for f in output_files], callback=self.sim.progress_callback("Copy from host"))
+            except Exception:
+                #
+                #             for src_file in output_files:
+                #                 try:
+                #                     dst_file = os.path.join(self.modelpath, os.path.relpath(src_file, self.tmp_modelpath))
+                #                     os.makedirs(os.path.dirname(dst_file), exist_ok=True)
+                #                     self.download(src_file, dst_file, retry=10)
+                #                 except Exception as e:
+                #                     download_failed.append(dst_file)
+                #             if download_failed:
+                raise Warning("Failed to download %s from %s" % (",".join(download_failed), self.ssh.host))
             else:
                 try:
                     self.ssh.execute('rm -r .hawc2launcher/%s' % self.simulation_id)
                 finally:
                     try:
                         self.ssh.execute('rm .hawc2launcher/status_%s' % self.simulation_id)
-                    except:
-                        raise Warning("Fail to remove temporary files and folders on %s"%self.ssh.host)
-
+                    except Exception:
+                        raise Warning("Fail to remove temporary files and folders on %s" % self.ssh.host)
 
     def _simulate(self):
         """starts blocking simulation"""
         self.sim.logFile = LogInfo(log_file.MISSING, 0, "None", "")
 
-        self.pbsjob.submit("%s.in" % self.simulation_id, self.tmp_exepath , self.sim.stdout_filename)
+        self.pbsjob.submit("%s.in" % self.simulation_id, self.tmp_exepath, self.sim.stdout_filename)
         sleeptime = 1
         while self.is_simulating:
             time.sleep(sleeptime)
@@ -353,11 +368,9 @@ class PBSClusterSimulationHost(SimulationHost):
             try:
                 self.ssh.download(self.tmp_exepath + self.log_filename, self.exepath + self.log_filename)
             except Exception:
-                raise Warning ("Logfile not found", self.tmp_modelpath + self.log_filename)
+                raise Warning("Logfile not found", self.tmp_modelpath + self.log_filename)
         self.sim.logFile = LogFile.from_htcfile(self.htcFile, self.exepath)
 
-
-
     def update_logFile_status(self):
         def pbsjob_status():
             if self.pbsjob._status in [NOT_SUBMITTED, DONE]:
@@ -374,7 +387,7 @@ class PBSClusterSimulationHost(SimulationHost):
                 self.logline = self.resource.loglines[self.simulation_id]
                 self.status = self.logline[0]
                 self.logFile = LogInfo(*self.logline[1:])
-                
+
         status = pbsjob_status()
         if status == pbsjob.NOT_SUBMITTED:
             pass
@@ -385,19 +398,16 @@ class PBSClusterSimulationHost(SimulationHost):
         else:
             set_status()
 
-
     def start(self):
         """Start non blocking distributed simulation"""
         self.non_blocking_simulation_thread.start()
 
-
-
     def abort(self):
         self.pbsjob.stop()
         self.stop()
         try:
             self.finish_simulation()
-        except:
+        except Exception:
             pass
         self.is_simulating = False
         self.is_done = True
@@ -407,22 +417,22 @@ class PBSClusterSimulationHost(SimulationHost):
     def stop(self):
         self.is_simulating = False
         self.pbsjob.stop()
-        
-
-
 
     def pbsjobfile(self, ios=False):
         cp_back = ""
-        for folder in set([fmt_path(os.path.relpath(os.path.dirname(f))) for f in self.htcFile.output_files() + self.htcFile.turbulence_files()]):
+        for folder in set([fmt_path(os.path.dirname(os.path.relpath(f)))
+                           for f in self.htcFile.output_files() + self.htcFile.turbulence_files()]):
+            if folder == '':
+                continue
             cp_back += "mkdir -p $PBS_O_WORKDIR/%s/. \n" % folder
             cp_back += "cp -R -f %s/. $PBS_O_WORKDIR/%s/.\n" % (folder, folder)
         rel_htcfilename = fmt_path(os.path.relpath(self.htcFile.filename, self.exepath))
         try:
             steps = self.htcFile.simulation.time_stop[0] / self.htcFile.simulation.newmark.deltat[0]
-            walltime = "%02d:00:00"%np.ceil(steps/500/60)
-        except:
+            walltime = "%02d:00:00" % np.ceil(steps / 500 / 60)
+        except Exception:
             walltime = "04:00:00"
-        init="""
+        init = """
 ### Standard Output
 #PBS -N h2l_%s
 ### merge stderr into stdout
@@ -436,21 +446,21 @@ class PBSClusterSimulationHost(SimulationHost):
 #PBS -q %s
 ### Create scratch directory and copy data to it
 cd $PBS_O_WORKDIR
-pwd"""% (self.simulation_id, self.stdout_filename, walltime, self.resource.queue)
-        copy_to="""
+pwd""" % (self.simulation_id, self.stdout_filename, walltime, self.resource.queue)
+        copy_to = """
 cp -R %s /scratch/$USER/$PBS_JOBID
 ### Execute commands on scratch nodes
 cd /scratch/$USER/$PBS_JOBID%s
-pwd"""%((".","../")[ios], ("", "/input")[ios])
-        run='''
+pwd""" % ((".", "../")[ios], ("", "/input")[ios])
+        run = '''
 %s
 ### modelpath: %s
-### htc: %s 
+### htc: %s
 echo "---------------------"
 %s -c "from wetb.hawc2.cluster_simulation import ClusterSimulation;ClusterSimulation('.','%s', ('%s','%s'))"
 echo "---------------------"
 echo $?
-echo "---------------------"'''% (self.resource.init_cmd, self.modelpath, self.htcFile.filename, self.resource.python_cmd, rel_htcfilename, self.resource.wine_cmd, self.hawc2exe)
+echo "---------------------"''' % (self.resource.init_cmd, self.modelpath, self.htcFile.filename, self.resource.python_cmd, rel_htcfilename, self.resource.wine_cmd, self.hawc2exe)
         copy_back = """
 ### Copy back from scratch directory
 cd /scratch/$USER/$PBS_JOBID%s
@@ -459,7 +469,4 @@ echo $PBS_JOBID
 cd /scratch/
 ### rm -r $PBS_JOBID
 exit""" % (("", "/input")[ios], cp_back)
-        return init+copy_to+run+copy_back
-
-
-
+        return init + copy_to + run + copy_back

wetb/hawc2/st_file.py

@@ -3,42 +3,70 @@ Created on 24/04/2014
 
 @author: MMPE
 '''
-from __future__ import print_function
-from __future__ import unicode_literals
-from __future__ import division
-from __future__ import absolute_import
-from io import open
-from builtins import range
-from builtins import int
-from future import standard_library
-standard_library.install_aliases()
-
+import types
+import os
 import numpy as np
 
+
+stc = "r m x_cg y_cg ri_x ri_y x_sh y_sh E G I_x I_y K k_x k_y A pitch x_e y_e"
+fpm = 'r m x_cg y_cg ri_x ri_y pitch x_e y_e K_11 K_12 K_13 K_14 K_15 K_16 K_22'
+fpm += ' K_23 K_24 K_25 K_26 K_33 K_34 K_35 K_36 K_44 K_45 K_46 K_55 K_56 K_66'
+
 class StFile(object):
     """Read HAWC2 St (beam element structural data) file
 
-    Methods are autogenerated for:
-
-    - r : curved length distance from main_body node 1 [m]
-    - m : mass per unit length [kg/m]
-    - x_cg : xc2-coordinate from C1/2 to mass center [m]
-    - y_cg : yc2-coordinate from C1/2 to mass center [m]
-    - ri_x : radius of gyration related to elastic center. Corresponds to rotation about principal bending xe axis [m]
-    - ri_y : radius of gyration related to elastic center. Corresponds to rotation about principal bending ye axis [m]
-    - xs : xc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion.
-    - ys : yc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion.
-    - E : modulus of elasticity [N/m2]
-    - G : shear modulus of elasticity [N/m2]
-    - Ix : area moment of inertia with respect to principal bending xe axis [m4]. This is the principal bending axis most parallel to the xc2 axis
-    - Iy : area moment of inertia with respect to principal bending ye axis [m4]
-    - K : torsional stiffness constant with respect to ze axis at the shear center [m4/rad]. For a circular section only this is identical to the polar moment of inertia.
-    - kx : shear factor for force in principal bending xe direction [-]
-    - ky : shear factor for force in principal bending ye direction [-]
-    - A : cross sectional area [m2]
-    - pitch : structural pitch about z_c2 axis. This is the angle between the xc2 -axis defined with the c2_def command and the main principal bending axis xe.
-    - xe : xc2-coordinate from C1/2 to center of elasticity [m]. The elastic center is the point where radial force (in the z-direction) does not contribute to bending around the x or y directions.
-    - ye : yc2-coordinate from C1/2 to center of elasticity [m]. The elastic center is the point where radial force (in the
+    Methods are autogenerated for the classical Timoshenko format:
+
+    - r: curved length distance from main_body node 1 [m]
+    - m: mass per unit length [kg/m]
+    - x_cg: xc2-coordinate from C1/2 to mass center [m]
+    - y_cg: yc2-coordinate from C1/2 to mass center [m]
+    - ri_x: radius of gyration related to elastic center. Corresponds to rotation about principal bending xe axis [m]
+    - ri_y: radius of gyration related to elastic center. Corresponds to rotation about principal bending ye axis [m]
+    - x_sh: xc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion.
+    - y_sh: yc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion.
+    - E: modulus of elasticity [N/m^2]
+    - G: shear modulus of elasticity [N/m^2]
+    - I_x: area moment of inertia with respect to principal bending xe axis [m^4]. This is the principal bending axis most parallel to the xc2 axis
+    - I_y: area moment of inertia with respect to principal bending ye axis [m^4]
+    - K: torsional stiffness constant with respect to ze axis at the shear center [m^4/rad]. For a circular section only this is identical to the polar moment of inertia.
+    - k_x: shear factor for force in principal bending xe direction [-]
+    - k_y: shear factor for force in principal bending ye direction [-]
+    - A: cross sectional area [m2]
+    - pitch: structural pitch about z_c2 axis. This is the angle between the xc2 -axis defined with the c2_def command and the main principal bending axis xe.
+    - x_e: xc2-coordinate from C1/2 to center of elasticity [m]. The elastic center is the point where radial force (in the z-direction) does not contribute to bending around the x or y directions.
+    - y_e: yc2-coordinate from C1/2 to center of elasticity [m]. The elastic center is the point where radial force (in the z-direction) does not contribute to bending around the x or y directions.
+
+    as well as the Fully-Populated Matrix format:
+
+    - r: curved length distance from main_body node 1 [m].
+    - m: mass per unit length [kg/m].
+    - x_cg: xc2-coordinate from C1/2 to mass center [m].
+    - y_cg: yc2-coordinate from C1/2 to mass center [m].
+    - ri_x: radius of gyration related to elastic center. Corresponds to rotation about principal bending xe axis [m].
+    - ri_y: radius of gyration related to elastic center. Corresponds to rotation about principal bending ye axis [m].
+    - pitch: structural pitch about z_c2 axis. This is the angle between the xc2 -axis defined with the c2_def command and the main principal bending axis xe.
+    - K_11: element of the cross sectional stiffness matrix [N].
+    - K_12: element of the cross sectional stiffness matrix [N].
+    - K_13: element of the cross sectional stiffness matrix [N].
+    - K_14: element of the cross sectional stiffness matrix [N*m].
+    - K_15: element of the cross sectional stiffness matrix [N*m].
+    - K_16: element of the cross sectional stiffness matrix [N*m].
+    - K_22: element of the cross sectional stiffness matrix [N].
+    - K_23: element of the cross sectional stiffness matrix [N].
+    - K_24: element of the cross sectional stiffness matrix [N*m].
+    - K_25: element of the cross sectional stiffness matrix [N*m].
+    - K_26: element of the cross sectional stiffness matrix [N*m].
+    - K_33: element of the cross sectional stiffness matrix [N].
+    - K_34: element of the cross sectional stiffness matrix [N*m].
+    - K_35: element of the cross sectional stiffness matrix [N*m].
+    - K_36: element of the cross sectional stiffness matrix [N*m].
+    - K_44: element of the cross sectional stiffness matrix [N*m^2].
+    - K_45: element of the cross sectional stiffness matrix [N*m^2].
+    - K_46: element of the cross sectional stiffness matrix [N*m^2].
+    - K_55: element of the cross sectional stiffness matrix [N*m^2].
+    - K_56: element of the cross sectional stiffness matrix [N*m^2].
+    - K_66: element of the cross sectional stiffness matrix [N*m^2].
 
     The autogenerated methods have the following structure
 
@@ -58,18 +86,38 @@ class StFile(object):
     Examples
     --------
     >>> stfile = StFile(r"tests/test_files/DTU_10MW_RWT_Blade_st.dat")
-    >>> print (stfile.m()) # Interpolated mass at radius 36
+    >>> print (stfile.m()) # Mass at nodes
     [ 1189.51054664  1191.64291781  1202.76694262  ... 15.42438683]
     >>> print (st.E(radius=36, mset=1, set=1))  # Elasticity interpolated to radius 36m
     8722924514.652649
     >>> print (st.E(radius=36, mset=1, set=2))  # Same for stiff blade set
     8.722924514652648e+17
     """
-    def __init__(self, filename):
-        with open (filename) as fid:
+
+    def __init__(self, filename=None):
+        """
+        Instantiate a new `StFile`.
+
+        Parameters
+        ----------
+        filename : str, optional
+            Path to the st file. If provided, this file will be immediately read. The default is None.
+
+        Returns
+        -------
+        None.
+
+        """
+
+        # in case the user wants to create a new non-existing st file
+        if filename is None:
+            self.main_data_sets = {}
+            return
+
+        with open(filename) as fid:
             txt = fid.read()
 #         Some files starts with first set ("#1...") with out specifying number of sets
-#         no_maindata_sets = int(txt.strip()[0]) 
+#         no_maindata_sets = int(txt.strip()[0])
 #         assert no_maindata_sets == txt.count("#")
         self.main_data_sets = {}
         for mset in txt.split("#")[1:]:
@@ -80,16 +128,32 @@ class StFile(object):
                 set_lines = set_txt.split("\n")
                 set_nr, no_rows = map(int, set_lines[0].split()[:2])
                 assert set_nr not in set_data_dict
-                set_data_dict[set_nr] = np.array([set_lines[i].split() for i in range(1, no_rows + 1)], dtype=np.float)
+                try:
+                    # HAWC2 will ignore everything after the 19th element,
+                    # some users have placed comments here after a ;
+                    linelst = [set_lines[i].split(';')[0].split() for i in range(1, no_rows + 1)]
+                except Exception as e:
+                    print('it went wrong at (set/subset):', mset_nr, set_nr,
+                          'with', no_rows, 'rows')
+                    raise e
+                set_data_dict[set_nr] = np.array(linelst, dtype=float)
             self.main_data_sets[mset_nr] = set_data_dict
 
-        for i, name in enumerate("r m x_cg y_cg ri_x ri_y x_sh y_sh E G I_x I_y I_p k_x k_y A pitch x_e y_e".split()):
-            setattr(self, name, lambda radius=None, mset=1, set=1, column=i: self._value(radius, column, mset, set))
+        if len(linelst[0])==len(stc.split()):
+            self.cols = stc.split()
+        elif len(linelst[0])==len(fpm.split()):
+            self.cols = fpm.split()
+        else:
+            raise TypeError('wrong number of columns in st file')
 
-    def _value(self, radius, column, mset_nr=1, set_nr=1):
-        st_data = self.main_data_sets[mset_nr][set_nr]
+        for i, name in enumerate(self.cols):
+            setattr(self, name, lambda radius=None, mset=1, set=1,
+                    column=i: self._value(radius, column, mset, set))
+
+    def _value(self, radius, column, mset=1, set=1):
+        st_data = self.main_data_sets[mset][set]
         if radius is None:
-            radius = self.radius_st(None, mset_nr, set_nr)
+            radius = self.radius_st(None, mset, set)
         return np.interp(radius, st_data[:, 0], st_data[:, column])
 
     def radius_st(self, radius=None, mset=1, set=1):
@@ -98,31 +162,210 @@ class StFile(object):
             return r
         return r[np.argmin(np.abs(r - radius))]
 
-    def to_str(self, mset=1, set=1):
+    def to_str(self, mset=1, set=1, precision='%12.5e '):
         d = self.main_data_sets[mset][set]
-        return "\n".join([("%12.5e "*d.shape[1]) % tuple(row) for row in d])
+        return '\n'.join([(precision * d.shape[1]) % tuple(row) for row in d])
+
+    def set_value(self, mset, set, **kwargs):
+        """
+        Set the value of one column in the st file.
+
+        Example: `st.set_value(1, 1, m=mass_array, x_sh=shear_offset)`.
+
+        Parameters
+        ----------
+        mset : int
+            Main set number.
+        set : int
+            Sub set number.
+        **kwargs : dict
+            An arbitrary number of columns, in the form `col_name = value`.
+            The column name must be as written in the class documentation.
+            The `value` argument must be a 1D array with size equal to the number of radial stations.
+
+        Returns
+        -------
+        None.
+
+        """
+        for k, v in kwargs.items():
+            column = self.cols.index(k)
+            self.main_data_sets[mset][set][:, column] = v
+
+    def save(self, filename, precision='%15.07e', encoding='utf-8'):
+        """Save all data defined in main_data_sets to st file.
+        """
+
+        # when creating empty st object, cols is not set yet
+        if not hasattr(self, 'cols'):
+            if self.main_data_sets[1][1].shape[1]==19:
+                self.cols = stc.split()
+            elif self.main_data_sets[1][1].shape[1]==30:
+                self.cols = fpm.split()
+            else:
+                c = self.main_data_sets[1][1].shape[1]
+                raise TypeError(f'St input needs 19 (iso) or 30 (aniso/fpm) cols, not {c}')
+
+        colwidth = len(precision % 1)
+        sep = '=' * colwidth * len(self.cols) + '\n'
+        colhead = ''.join([k.center(colwidth) for k in self.cols]) + '\n'
+
+        nsets = len(self.main_data_sets)
+        # fails if dirname is empty string
+        if len(os.path.dirname(filename)) > 0:
+            os.makedirs(os.path.dirname(filename), exist_ok=True)
+        with open(filename, 'w', encoding=encoding) as fid:
+            fid.write('%i ; number of sets, Nset\n' % nsets)
+            for mset, set_data_dict in self.main_data_sets.items():
+                fid.write('#%i ; set number\n' % mset)
+                for set, set_array in set_data_dict.items():
+                    dstr = self.to_str(mset=mset, set=set, precision=precision)
+                    npoints = self.main_data_sets[mset][set].shape[0]
+                    fid.write(sep + colhead + sep)
+                    fid.write('$%i %i\n' % (set, npoints))
+                    fid.write(dstr + '\n')
+
+    def element_stiffnessmatrix(self, radius, mset, set, length):
+        """Compute the element stiffness matrix
+
+        Parameters
+        ----------
+        radius : float
+            radius of element (used of obtain element properties).
+        mset : int
+            Main set number.
+        set : int
+            Sub set number.
+        length : float
+            Element length
+        """
+
+        # not supported for FPM format
+        if len(self.cols)==30:
+            return
+
+        K = np.zeros((13, 13))
+        "r m x_cg y_cg ri_x ri_y x_sh y_sh E G I_x I_y K k_x k_y A pitch x_e y_e"
+        ES1, ES2, EMOD, GMOD, IX, IY, IZ, KX, KY, A = [getattr(self, n)(radius, mset, set)
+                                                       for n in "x_sh,y_sh,E,G,I_x,I_y,K,k_x,k_y,A".split(",")]
+        ELLGTH = length
+
+        ETAX = EMOD * IX / (KY * GMOD * A * ELLGTH**2)
+        ETAY = EMOD * IY / (KX * GMOD * A * ELLGTH**2)
+        ROX = 1 / (1 + 12 * ETAX)
+        ROY = 1 / (1 + 12 * ETAY)
+        ROY = 1 / (1 + 12 * ETAX)
+        K[1, 1] = 12 * EMOD * IY * ROY / ELLGTH**3
+        K[1, 5] = 6 * EMOD * IY * ROY / ELLGTH**2
+        K[1, 6] = -K[1, 1] * ES2
+        K[1, 7] = -K[1, 1]
+        K[1, 11] = K[1, 5]
+        K[1, 12] = -K[1, 6]
+        K[2, 2] = 12 * EMOD * IX * ROX / ELLGTH**3
+        K[2, 4] = -6 * EMOD * IX * ROX / ELLGTH**2
+        K[2, 6] = K[2, 2] * ES1
+        K[2, 8] = -K[2, 2]
+        K[2, 10] = K[2, 4]
+        K[2, 12] = -K[2, 6]
+        K[3, 3] = A * EMOD / ELLGTH
+        K[3, 9] = -K[3, 3]
+        K[4, 4] = 4 * EMOD * IX * (1 + 3 * ETAX) * ROX / ELLGTH
+        K[4, 6] = K[2, 4] * ES1
+        K[4, 8] = -K[2, 4]
+        K[4, 10] = 2 * EMOD * IX * (1 - 6 * ETAX) * ROX / ELLGTH
+        K[4, 12] = -K[4, 6]
+        K[5, 5] = 4 * EMOD * IY * (1 + 3 * ETAY) * ROY / ELLGTH
+        K[5, 6] = -K[1, 5] * ES2
+        K[5, 7] = -K[1, 5]
+        K[5, 11] = 2 * EMOD * IY * (1 - 6 * ETAY) * ROY / ELLGTH
+        K[5, 12] = -K[5, 6]
+        K[6, 6] = GMOD * IZ / ELLGTH + 12 * EMOD * (IX * ES1**2 * ROX + IY * ES2**2 * ROY) / ELLGTH**3
+        K[6, 7] = K[1, 12]
+        K[6, 8] = K[2, 12]
+        K[6, 10] = -K[4, 12]
+        K[6, 11] = -K[5, 12]
+        K[6, 12] = -K[6, 6]
+        K[7, 7] = K[1, 1]
+        K[7, 11] = -K[1, 5]
+        K[7, 12] = K[1, 6]
+        K[8, 8] = K[2, 2]
+        K[8, 10] = -K[2, 4]
+        K[8, 12] = K[2, 6]
+        K[9, 9] = K[3, 3]
+        K[10, 10] = K[4, 4]
+        K[10, 12] = -K[4, 6]
+        K[11, 11] = K[5, 5]
+        K[11, 12] = -K[5, 6]
+        K[12, 12] = K[6, 6]
+        K = K[1:, 1:]
+        K = K + K.T - np.eye(12) * K
+        return K
+
+    def shape_function_ori(self, radius, mset, set, length, z):
+
+        # not supported for FPM format
+        if len(self.cols)==30:
+            return
+
+        XSC, YSC, EMOD, GMOD, IX, IY, IZ, KX, KY, AREA = [getattr(self, n)(radius, mset, set)
+                                                          for n in "x_sh,y_sh,E,G,I_x,I_y,K,k_x,k_y,A".split(",")]
+
+        etax = EMOD * IX / KY / GMOD / AREA / (length**2)
+        etay = EMOD * IY / KX / GMOD / AREA / (length**2)
+        rhox = 1 / (1 + 12 * etax)
+        rhoy = 1 / (1 + 12 * etay)
+
+        f1 = z / length
+        f2 = 1 - f1
+        f3x = f1 * (3 * f1 - 2 * f1**2 + 12 * etax) * rhox
+        f4x = f2 * (3 * f2 - 2 * f2**2 + 12 * etax) * rhox
+        f5x = length * (f1**2 - (1 - 6 * etax) * f1 - 6 * etax) * f1 * rhox
+        f6x = length * (f2**2 - (1 - 6 * etax) * f2 - 6 * etax) * f2 * rhox
+        f7x = 6 / length * f1 * f2 * rhox
+        f8x = f1 * (3 * f1 - 2 * (1 - 6 * etax)) * rhox
+        f9x = f2 * (3 * f2 - 2 * (1 - 6 * etax)) * rhox
+
+        f3y = f1 * (3 * f1 - 2 * f1**2 + 12 * etay) * rhoy
+        f4y = f2 * (3 * f2 - 2 * f2**2 + 12 * etay) * rhoy
+        f5y = length * (f1**2 - (1 - 6 * etay) * f1 - 6 * etay) * f1 * rhoy
+        f6y = length * (f2**2 - (1 - 6 * etay) * f2 - 6 * etay) * f2 * rhoy
+        f7y = 6 / length * f1 * f2 * rhoy
+        f8y = f1 * (3 * f1 - 2 * (1 - 6 * etay)) * rhoy
+        f9y = f2 * (3 * f2 - 2 * (1 - 6 * etay)) * rhoy
+
+        return np.array([[f4y, 0, 0, 0, -f7y, 0],
+                         [0,  f4x, 0,  f7x, 0, 0],
+                         [0, 0,   f2, 0, 0, 0],
+                         [0, f6x, 0, f9x, 0, 0],
+                         [-f6y, 0, 0, 0, f9y, 0],
+                         [(f2 - f4y) * YSC, -(f2 - f4x) * XSC, 0, f7x * XSC, f7y * YSC, f2],
+                         [f3y, 0, 0, 0, f7y, 0],
+                         [0, f3x, 0, -f7x, 0, 0],
+                         [0, 0,    f1, 0, 0, 0],
+                         [0, -f5x, 0, f8x, 0, 0],
+                         [f5y, 0, 0, 0, f8y, 0],
+                         [(f1 - f3y) * YSC, -(f1 - f3x) * XSC, 0, -f7x * XSC, -f7y * YSC, f1]]).T
 
 
 if __name__ == "__main__":
     import os
-    st = StFile(os.path.dirname(__file__) + r"/tests/test_files/DTU_10MW_RWT_Blade_st.dat")
-    print (st.m())
-    print (st.E(radius=36, mset=1, set=1))  # Elastic blade
-    print (st.E(radius=36, mset=1, set=2))  # stiff blade
+    cwd = os.path.dirname(__file__)
+    st = StFile(os.path.join(cwd, r'tests/test_files/DTU_10MW_RWT_Blade_st.dat'))
+    print(st.m())
+    print(st.E(radius=36, mset=1, set=1))  # Elastic blade
+    print(st.E(radius=36, mset=1, set=2))  # stiff blade
     #print (st.radius())
     xyz = np.array([st.x_e(), st.y_e(), st.r()]).T[:40]
     n = 2
     xyz = np.array([st.x_e(None, 1, n), st.y_e(None, 1, n), st.r(None, 1, n)]).T[:40]
     #print (xyz)
-    print (np.sqrt(np.sum((xyz[1:] - xyz[:-1]) ** 2, 1)).sum())
-    print (xyz[-1, 2])
-    print (np.sqrt(np.sum((xyz[1:] - xyz[:-1]) ** 2, 1)).sum() - xyz[-1, 2])
-    print (st.x_e(67.8883), st.y_e(67.8883))
+    print(np.sqrt(np.sum((xyz[1:] - xyz[:-1]) ** 2, 1)).sum())
+    print(xyz[-1, 2])
+    print(np.sqrt(np.sum((xyz[1:] - xyz[:-1]) ** 2, 1)).sum() - xyz[-1, 2])
+    print(st.x_e(67.8883), st.y_e(67.8883))
     #print (np.sqrt(np.sum(np.diff(xyz, 0) ** 2, 1)))
-    print (st.pitch(67.8883 - 0.01687))
-    print (st.pitch(23.2446))
-
-
+    print(st.pitch(67.8883 - 0.01687))
+    print(st.pitch(23.2446))
 
-    #print (st.)
-    #print (st.)
+    # print (st.)
+    # print (st.)

wetb/hawc2/tests/test_AtTimeFile.py

@@ -3,12 +3,6 @@ Created on 17/07/2014
 
 @author: MMPE
 '''
-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import division
-from __future__ import absolute_import
-from future import standard_library
-standard_library.install_aliases()
 import unittest
 from wetb.hawc2.at_time_file import AtTimeFile
 import numpy as np
@@ -62,7 +56,6 @@ class TestAtTimeFile(unittest.TestCase):
 
     def test_rotor_name_null(self):
         atfile = AtTimeFile(self.testfilepath + "at_time/test_rotor_name_null.dat", bladetip_radius=20.501)  # load file
-        print (atfile.radius_s())
 
 
 if __name__ == "__main__":

wetb/hawc2/tests/test_ae_file.py

@@ -3,20 +3,22 @@ Created on 05/11/2015
 
 @author: MMPE
 '''
-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import division
-from __future__ import absolute_import
-from future import standard_library
-standard_library.install_aliases()
+
+import os
 import unittest
+
+import mock
+
+from wetb.hawc2 import ae_file
 from wetb.hawc2.ae_file import AEFile
-import os
+import numpy as np
+
 
 testfilepath = os.path.join(os.path.dirname(__file__), 'test_files/')  # test file path
-class TestAEFile(unittest.TestCase):
 
 
+class TestAEFile(unittest.TestCase):
+
     def test_aefile(self):
         ae = AEFile(testfilepath + "NREL_5MW_ae.txt")
         self.assertEqual(ae.thickness(38.950), 21)
@@ -29,7 +31,56 @@ class TestAEFile(unittest.TestCase):
         self.assertEqual(ae.chord(32), 3.673)
         self.assertEqual(ae.pc_set_nr(32), 1)
 
+    def test_ae_file_main(self):
+        def no_print(s):
+            pass
+        with mock.patch.object(ae_file, "__name__", "__main__"):
+            with mock.patch.object(ae_file, "print", no_print):
+                getattr(ae_file, 'main')()
+
+    def test_add_set(self):
+        ae = AEFile(testfilepath + "NREL_5MW_ae.txt")
+
+        ae.add_set(radius=[0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
+                   chord=[1.1, 1.0, 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, 0.1],
+                   thickness=[100.0, 100.0, 90.0, 80.0, 70.0, 60.0, 50.0, 40.0, 30.0, 20.0, 10.0],
+                   pc_set_id=[1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0])
+
+        ae.add_set(radius=[0.0, 0.1],
+                   chord=[1.1, 1.0],
+                   thickness=[100.0, 100.0],
+                   pc_set_id=[1.0, 1.0],
+                   set_id=4)
+
+        self.assertEqual(ae.thickness(38.950), 21)
+        self.assertEqual(ae.chord(38.950), 3.256)
+        self.assertEqual(ae.pc_set_nr(38.950), 1)
+        np.testing.assert_array_equal(ae.chord(None, 2), [1.1, 1.0, 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, 0.1])
+        np.testing.assert_array_equal(ae.ae_sets[2][:2], ae.ae_sets[4])
+
+    def test_str(self):
+        ae = AEFile(testfilepath + "NREL_5MW_ae.txt")
+        ref = """1 r[m]           Chord[m]    T/C[%]  Set no.
+1 19
+  0.00000000000000000e+00   3.54199999999999982e+00   1.00000000000000000e+02     1"""
+
+        self.assertEqual(str(ae)[:len(ref)], ref)
+
+    def test_save(self):
+        ae = AEFile()
+
+        ae.add_set(radius=[0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0],
+                   chord=[1.1, 1.0, 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, 0.1],
+                   thickness=[100.0, 100.0, 90.0, 80.0, 70.0, 60.0, 50.0, 40.0, 30.0, 20.0, 10.0],
+                   pc_set_id=[1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0])
+        fn = testfilepath + "tmp/ae_file.txt"
+        ae.save(fn)
+        ae2 = AEFile(fn)
+        assert str(ae) == str(ae2)
 
+    def test_multiple_sets(self):
+        ae = AEFile(testfilepath + 'ae_files/HAWC2_ae.dat')
+        self.assertEqual(len(ae.ae_sets), 2)
 
 
 if __name__ == "__main__":

wetb/hawc2/tests/test_body_matrix_output_file.py

+# -*- coding: utf-8 -*-
+"""
+Test body_matrix_output_file module.
+
+@author: ricriv
+"""
+
+# %% Import.
+
+import os
+import unittest
+import numpy.testing as npt
+from wetb.hawc2.body_matrix_output_file import read_body_matrix_output
+
+
+# %% Test.
+
+# Test file path.
+tfp = os.path.join(os.path.dirname(__file__), "test_files/")
+
+
+class TestBodyOutput(unittest.TestCase):
+    def test_all(self):
+        # Read files.
+        bodies = read_body_matrix_output(f"{tfp}/body_matrix_output/body_")
+
+        # Have we read all files?
+        self.assertEqual(set(bodies.keys()), {"blade1_1", "tower"})
+        self.assertEqual(
+            set(bodies["tower"].keys()), {"mass", "damping", "stiffness"}
+        )
+        self.assertEqual(
+            set(bodies["blade1_1"].keys()), {"mass", "damping", "stiffness"}
+        )
+
+        # Check matrices shape.
+        npt.assert_array_equal(bodies["tower"]["mass"].shape, (66, 66))
+        npt.assert_array_equal(bodies["blade1_1"]["mass"].shape, (18, 18))
+
+        # Check 1 value.
+        npt.assert_almost_equal(
+            bodies["tower"]["mass"][0, 0], 627396.222002882
+        )
+
+
+if __name__ == "__main__":
+    unittest.main()

wetb/hawc2/tests/test_compare_test_cases.py

@@ -3,12 +3,6 @@ Created on 17/07/2014
 
 @author: MMPE
 '''
-from __future__ import unicode_literals
-from __future__ import print_function
-from __future__ import division
-from __future__ import absolute_import
-from future import standard_library
-standard_library.install_aliases()
 import unittest
 from wetb.hawc2.cmp_test_cases import CompareTestCases
 import numpy as np