Source code for pywasp.wasp.config

"""
WAsP Configuration module

This module contains the Config class that is used to set various parameters used in the
WAsP modelling. The config class contains two elements, terrain and climate to control
their respective models.
"""

import json
import pprint
import warnings
import textwrap

import pandas as pd
import tabulate
import numpy as np

from .._errors import PywaspError
from ..core import Rvea0287, Rvea0288
from .parameters_metadata import _DICT_TEMPLATE, CLIMATE_DESC, TERRAIN_DESC

_DEFAULT_PAR_SET = "WAsP_12.8"


def _dict2json(dict_in, filepath):
    """Exports nested dictionary to human-readable json.

    Parameters
    ----------
    dict_in : dict
        Nested dictionary containing metadata.
    filepath : str or Path
        Path to json file
    """
    with open(filepath, "w") as outfile:
        json.dump(dict_in, outfile, indent=4)


def _json2dict(filepath):
    """Loads json file to nested dictionary.

    Parameters
    ----------
    filepath : str or Path
        Path to json file

    Returns
    -------
    dict
        Nested dictionary containing metadata.
    """
    with open(filepath) as handle:
        dictdump = json.loads(handle.read())
    return dictdump


def _print_desc(par):  # pragma: no cover
    """
    Prints description of parameter

    Parameters
    ----------
    par : dict
        PyWAsP formated dictionary of parameter, see _DICT_TEMPLATE for details
    """

    short_fmt = "{0:15}: {1:16}\n"
    long_fmt = "\n{0}:\n{1}\n"

    _long = ["notes", "definition", "references"]
    _not_print = ["from_windows_gui", "investigation_required", "internal_parameter"]

    out = "================================================================================\n"
    out += "{}\n".format(par["short_name"])
    out += "================================================================================\n"

    for key in par.keys():
        if key in _long and key not in _not_print and par[key] != "":
            if len(par[key]) > 80:
                formated_str = pprint.pformat(par[key], indent=0)
                formated_str = (
                    formated_str.replace("('", "")
                    .replace("')", "")
                    .replace("'\n '", "\n")
                )
            else:
                formated_str = par[key]
            out += long_fmt.format(key, formated_str)
        elif key not in _not_print and par[key] != "":
            out += short_fmt.format(key, par[key])

    out += "================================================================================\n"

    print(out)


class _Parameters:
    """This is a base class for manipulating with WAsP configuration parameters.

    Parameters
    ----------
    Rvea : module
        module that links fotran deconf.90 to python via f2py
    vocab_dict : dict
        Dictionary of metadata that contains additional information for each parameter
    par_set : str
        Named parameter set with different defaults
    """

    def __init__(self, Rvea, vocab_dict):
        super().__init__()
        self.Rvea = Rvea
        self.param, par_desc, self.nodatavalue = self.Rvea.deconf()
        self._vocab_dict = vocab_dict

        # converting byte strings to Unicode (suffix 180 indicates max no of characters)
        if len(par_desc.shape) == 2:
            par_desc = par_desc.view(
                "S180"
            ).squeeze()  # Convert from S1 to S180 and squeeze character dimension
        par_desc = np.char.strip(par_desc).astype("U180")
        self.par_desc = dict((i + 1, desc) for i, desc in enumerate(par_desc))

        # temporary make short names until we build a vocab for config parameters
        # those should be part of fotran code
        par_names = []
        for i in range(0, len(self.param)):
            if (i + 1) in self._vocab_dict:
                new_par_name = self._vocab_dict[i + 1]["short_name"]
                if new_par_name == "":
                    new_par_name = "short_name_" + str(i + 1)
                par_names += [new_par_name]
                self.par_desc[i + 1] = self._vocab_dict[i + 1]["definition"]
            else:
                par_names += ["short_name_" + str(i + 1)]

        self.par_values = dict((i + 1, value) for i, value in enumerate(self.param))
        self.par_names = dict((i + 1, par_name) for i, par_name in enumerate(par_names))

        # switch keys and values
        self.alias = {y: x for x, y in self.par_names.items()}

    def __getitem__(self, key):
        """This is a modified dict method __getitem__
        which allows getting values from the dict
        either via parameter short name or its index number.
        """
        # checks if key exists as a parameter index
        if key in self.alias.values():
            return self.par_values[key]
        # checks if key exists as a parameter short-name
        elif key in self.alias:
            return self.par_values[self.alias[key]]
        # otherwise it raise KeyError
        else:
            raise KeyError(f"Key '{key}' does not exist in dictionary")

    def __setitem__(self, key, value):
        """This is a modified dict method __setitem__
        which allows setting values to the dict
        either via parameter short name or its index number.
        """
        if key in self.alias.values():
            self.par_values[key] = value
            self.param[key - 1] = value  # numpy array, index starts at 0
        elif key in self.alias:
            self.par_values[self.alias[key]] = value
            self.param[self.alias[key] - 1] = value  # numpy array, index starts at 0
        else:
            raise KeyError(f"Key '{key}' does not exist in dictionary")

    def __eq__(self, comp):
        if not isinstance(comp, type(self)):
            raise TypeError("Comparison classes don't match.")
        return all(np.isclose(self.param, comp.param))

    def keys(self, ret_type=int):
        """Return the keys of the Parameter in requested type

        Parameters
        ----------
        ret_type : type, optional
            Return type of the keys, either int or str, by default int

        Returns
        -------
        list
            list containing keys in requested format
        """
        if ret_type == int:
            return list(self.par_names.keys())
        elif ret_type == str:
            return list(self.par_names.values())
        else:
            raise ValueError("ret_type must be either int or str.")

    def getpar(self, indexbasedone):  # pragma: no cover
        """Get parameter with base 1 index

        Parameters
        ----------
        indexbasedone: int
            One based integer of the parameter of interest

        Returns
        -------
        param: float
            Value of the requested parameter
        """
        msg = (
            "getpar is deprecated, use square-brackets [] and list-like syntax instead."
        )
        warnings.warn(msg, FutureWarning)
        return self.__getitem__(indexbasedone)

    def putpar(self, indexbasedone, newval):  # pragma: no cover
        """Set parameter with index to new value

        Parameters
        ----------
        indexbasedone: int
            One based integer of the parameter of interest
        newval: float
            New value of the parameter

        """
        msg = (
            "putpar is deprecated, use square-brackets [] and list-like syntax instead."
        )
        warnings.warn(msg, FutureWarning)
        self.__setitem__(indexbasedone, newval)

    def get_desc(self, par):  # pragma: no cover
        """Get parameter description

        Parameters
        ----------
        par: str, int
            Short name or index of parameter

        Returns
        -------
        desc: str
            Description of the requested parameter
        """
        if type(par) == str:
            if par in self.alias.keys() and self.alias[par] in self._vocab_dict:
                par_desc = self._vocab_dict[self.alias[par]]
                return _print_desc(par_desc)
            else:
                print("Requested parameter does not exist")
        elif type(par) == int:
            if par in self.par_names.keys() and par in self._vocab_dict:
                par_desc = self._vocab_dict[par]
                return _print_desc(par_desc)
            else:
                print("Requested parameter does not exist")

    @classmethod
    def from_dict(cls, d):
        """Initialize class from dictionary

        This is most often used to create a parameter set saved from file.

        Parameters
        ----------
        d : dict
            Stacked dictionary of the configuration parameters

        Returns
        -------
        self : self
            Object filled with the requested fields
        """
        self = cls()
        keys = [int(i) for i in d.keys()]
        keys.sort()

        param = []
        par_desc = {}
        par_names = {}
        par_values = {}
        for key in keys:
            try:
                val = d[key]
            except KeyError:  # JSON encodes keys as strings so need to handle that
                val = d[str(key)]
            par_desc.update({key: val["definition"]})
            par_names.update({key: val["short_name"]})
            par_values.update({key: float(val["default_value"])})
            param += [float(val["default_value"])]

        param = np.asarray(param)

        self.param = param
        self.par_desc = par_desc
        self.par_names = par_names
        self.par_values = par_values
        self.alias = {y: x for x, y in self.par_names.items()}

        return self

    def to_dict(self):
        """Create stacked dictionary with all information needed to recreate object

        Returns
        -------
        dict
            Stacked dictionary that can recreate entire object
        """
        dict_out = {}
        for key in self.par_values.keys():
            dict_instance = _DICT_TEMPLATE.copy()
            dict_instance["index"] = key
            dict_instance["default_value"] = float(self.par_values[key])
            dict_instance["definition"] = self.par_desc[key]
            dict_instance["short_name"] = self.par_names[key]
            dict_out.update({key: dict_instance})
        return dict_out

    def compare(self, comp, exact=False):
        """Print values that differ between two parameter sets

        Parameters
        ----------
        comp : Climate or Terrain
            Same class as self to use for comparison
        exact : bool, optional
            use an equality comparison, by default use numpy.isclose
        """
        out = "Index  Parameter        Own           Comp          Description\n"
        out += "                        value         value\n"
        out += "=====  ===============  ============  ===========   ==================================================\n"
        fmt = "{0:5d}  {1:14}  {2:13f}  {3:12}  {4}\n"  # first digit indicated order of the input elements

        if not isinstance(comp, type(self)):
            raise TypeError("Comparison classes don't match.")
        for i, val in enumerate(self.param):
            if exact:
                is_eq = val == comp.param[i]
            else:
                is_eq = np.isclose(val, comp.param[i])
            if not is_eq:
                idx = i + 1
                out += fmt.format(
                    idx,
                    self.par_names[idx],
                    self.param[i],
                    comp.param[i],
                    self.par_desc[idx],
                )
        print(out)

    def __str__(self):
        # Set header and other information
        out = "Index  Parameter                   Current       Default       Description\n"
        out += "                                   value         value\n"
        out += "=====  ==========================  ============  ===========   ==================================================\n"
        fmt = "{0:5d}  {1:25}  {2:13f}  {3:12}  {4}\n"  # first digit indicated order of the input elements

        # Loop over all values and print the ones that are set
        for i in range(0, len(self.par_desc)):
            if self._vocab_dict.get(i + 1, None) is None:
                continue
            if (
                self.par_desc[i + 1] != ""
                and self._vocab_dict[i + 1]["internal_parameter"] == False
                and self._vocab_dict[i + 1]["investigation_required"] == False
            ):
                out += fmt.format(
                    i + 1,
                    self.par_names[i + 1],
                    self.param[i],
                    self._vocab_dict[i + 1]["default_value"],
                    self.par_desc[i + 1],
                )

        return out

    def __repr__(self):
        return self.__str__()


class Climate(_Parameters):
    """Parameters for WAsP's climate module

    NOTE: This class is usually accessed from inside of the Config class

    Parameters
    ----------
    par_set : str
        Named parameter set with different defaults
    """

    def __init__(self, par_set=_DEFAULT_PAR_SET):
        super().__init__(Rvea0288, CLIMATE_DESC)
        if par_set in ["WAsP_12.6", "WAsP_12.7"]:
            self.set_profile_model(1)
        elif par_set == "WAsP_11.4":
            self.set_profile_model(-1)
        elif par_set == "WAsP_12.8":
            self.set_profile_model(3)
        else:
            raise PywaspError("Unknown WAsP Parameter set.")

    def set_profile_model(self, model):
        """Sets vertical profile model Parameters

        .. warning::
            This function is experimental and its signature may change.

        Parameters
        ----------
        model: int
            Integer which sets corresponds to:

            * *-1* WAsP 11 profile model using c_g=1.65 and powerlaw estimation
              for reversal height

            * *0* WAsP 12 profile model with c_g solved by iteration and dA/dmu
              set such that it corresponds with the old profile model

            * *1* WAsP 12 profile model, as version 0 but including baroclinicity

            * *2* WAsP 12 profile model, including sectorwise stability as described in mesoclimate class

            * *3* WAsP 12 profile model, including sectorwise stability and baroclinicity as described in mesoclimate class
        """
        # Get defaults from deconf in case they have been changed by the user
        defaults, _, _ = self.Rvea.deconf()  # 0-indexed
        if model == -1:
            self[10] = 1.8
            self[11] = 4.5
            self[44] = 10
            self[48] = -1
            self[54] = 0.002
            self[55] = 1.65
            self[62] = 0.6
            self[64] = 16
            self[86] = 1e-4
            self[87] = 0.1653
            self[90] = 0.0
            self[91] = 0.0
            self[95] = 0.0
            self[96] = 0.1653
            self[100] = 0.0
            self[101] = 1.515152
            self[102] = -1
            self[104] = 360
        elif model == 0:
            self[10] = 1.8
            self[11] = 4.5
            self[44] = 10
            self[48] = 0
            self[54] = -1
            self[55] = defaults[55 - 1]
            self[62] = 0.6
            self[64] = 16
            self[86] = 1e-4
            self[87] = 0.1653
            self[90] = 0.0
            self[91] = 0.0
            self[95] = 0.0
            self[96] = 0.1653
            self[100] = 0.0
            self[101] = 2.5 / 1.65
            self[102] = -1
            self[104] = 360
        elif model == 2:
            # wasp sectorwise tstar stability
            self[10] = defaults[10 - 1]
            self[11] = defaults[11 - 1]
            self[44] = defaults[44 - 1]  # latitude below which coriolis is constant.
            self[48] = 2
            self[54] = defaults[54 - 1]
            self[55] = defaults[55 - 1]
            self[62] = defaults[62 - 1]
            self[64] = defaults[64 - 1]
            self[86] = defaults[
                86 - 1
            ]  # set coriolis varying with latitude above 45 degrees lat
            self[90] = 0.0
            self[91] = 0.0
            self[95] = self[11] / self[8]
            self[96] = defaults[96 - 1]
            self[100] = defaults[100 - 1]
            self[101] = defaults[101 - 1]
            self[102] = defaults[102 - 1]  # slope of psi function in stable conditions
            self[104] = defaults[104 - 1]
        elif model == 3:
            # wasp sectorwise tstar stability + baro
            self[10] = defaults[10 - 1]
            self[11] = defaults[11 - 1]
            self[44] = defaults[44 - 1]  # latitude below which coriolis is constant.
            self[48] = defaults[48 - 1]
            self[54] = defaults[54 - 1]
            self[55] = defaults[55 - 1]
            self[62] = defaults[62 - 1]
            self[64] = defaults[64 - 1]  # minimum geostrophic wind
            self[86] = defaults[86 - 1]
            self[87] = defaults[87 - 1]
            self[90] = defaults[
                90 - 1
            ]  # form factor for linear geostrophic shear with height
            self[91] = defaults[
                91 - 1
            ]  # form factor for linear geostrophic shear with height
            self[
                95
            ] = (
                -1
            )  # correction of wind profile depends on varying B based on baroclinicity
            self[96] = defaults[
                96 - 1
            ]  # boundary layer height from external source is used
            self[100] = defaults[100 - 1]
            self[101] = defaults[101 - 1]
            self[102] = defaults[102 - 1]
            self[104] = defaults[104 - 1]
        elif model == 1:
            # wasp 12.1 baroclinic, current default
            self[10] = 1.8
            self[11] = 4.5
            self[44] = 10
            self[48] = 1
            self[54] = -1
            self[55] = defaults[55 - 1]
            self[62] = 0.6
            self[64] = 16
            self[86] = 1e-4
            self[87] = 0.3
            self[90] = 0.5
            self[91] = 0.5
            self[95] = 0.0
            self[96] = 0.1653
            self[100] = 0.0
            self[101] = 2.5 / 1.65
            self[102] = -1
            self[104] = 360
        else:
            raise ValueError(f"Option '{model}' is not valid.")


class Terrain(_Parameters):
    """Parameters of for the orography and roughness models

    NOTE: This class is usually accessed from inside of the Config class

    Parameters
    ----------
    par_set : str
        Named parameter set with different defaults
    """

    def __init__(self, par_set=_DEFAULT_PAR_SET):
        super().__init__(Rvea0287, TERRAIN_DESC)
        if par_set in ["WAsP_11.4", "WAsP_12.6"]:
            self.set_terrain_analysis(0)
        elif par_set in ["WAsP_12.7", "WAsP_12.8"]:
            self.set_terrain_analysis(1)
        else:
            raise PywaspError("Unknown WAsP Parameter set.")

    def set_terrain_analysis(self, indexterr):
        """Set terrain analysis

        Set the terrain analysis to classic (0) or spider analysis (1).
        See WAsP Core for references and documentation.

        .. warning::
            This function is experimental and its signature may change.

        Parameters
        ----------
        indexterr: int
            Terrain analysis classic (0) or spider analysis (1)

        """
        if indexterr == 1:
            # the spider grid is supposed to deal with forest
            # so the BZ model was modified. See "introduction of canopies in wasp"
            # available at rofl
            self[82] = 1
            self[26] = -0.29
            self[30] = 2
        else:
            self[82] = 0
            self[26] = 0.27
            self[30] = 0.67

    def use_displacement_height(self, indexterr=True):
        """Use displacement height

        Enable (1) or disable (0) adding the displacement height to the orographic grid.

        .. warning::
            This function is experimental and its signature may change.

        Parameters
        ----------
        indexterr: boolean
            Enable (True) or disable (False) adding the displacement height to the
            orographic grid. Defaults to True.
        """

        if indexterr == True:
            self[69] = 1
        else:
            self[69] = 0


def _params_to_table(params, columns=None):
    """Convert params to .rst table"""
    df = pd.DataFrame(params.values())
    if columns is None:
        columns = df.columns
    return df[columns].to_markdown(tablefmt="rst", index=False)


[docs] class Config: """Configuration class for the WAsP model parameters Due to the structure of the original fortran testbench code, parameters from the WAsP model are set in a long array with integers in a Fortran data block. Here we initialize both the climate and terrain parameters. The obstacle model does not, yet, have a similar structure to that of the climate and terrain models. Parameter can be set with a parameter-suite name, like "WAsP_12.8" for the parameters corresponding to WAsP version 12.8. Individual parameters can also be set manually through the methods of the ``.terrain`` and ``.climate`` class attributes (see below). By default, "WAsP_12.8" is used. Parameters ---------- par_set : str Named parameter set with different defaults: * WAsP_11.4 -- Profile model -1; terrain analysis 0 * WAsP_12.6 -- Profile model 1; terrain analysis 0 * WAsP_12.7 -- Profile model 1; terrain analysis 1 * WAsP_12.8 -- Profile model 3; terrain analysis 1 The above numbers are the parameter when calling ``Config.climate.set_profile_model`` and ``Config.terrain.set_terrain_analysis`` respectively. Notes ----- The ``.terrain`` and ``.climate`` class attributes each respectively hold the parameters for the terrain and climate analysis in WAsP. The attributes are themselves classes with associated methods described below. To interact with the parameters, both ``terrain`` and ``.climate`` acts as lists, allowing to get and set parameters by indicies:: from pywasp.wasp import Config conf = Config() conf.climate[10] = 1.8 # Set the climate parameter 10 (A0) value to 1.8 decay_length = conf.terrain[31] # get the decay length parameter from the terrain model **Terrain model config** The terrain model allows users to set the terrain analysis parameter-set and whether toggle whether displacement heights are used. This is done through the two methods: 1. ``.set_terrain_analysis(index)``, which set the terrain analysis. It has two options: - 0: for classic analysis - 1: for "spider" analysis. 2. ``.use_displacement_height(flag)``, which sets the displacement heights. Options are True or False. The full list of parameters for the terrain model is: {terrain_params_table} **Climate model config** The climate model config allows users to set the profile model through the ``.set_profile_model(model)``. Currently, the options are: - -1: WAsP 11 profile model using c_g=1.65 and powerlaw estimation for reversal height - 0: WAsP 12 profile model with c_g solved by iteration and dA/dmu set such that it corresponds with the old profile model - 1: WAsP 12 profile model, as version 0 but including baroclinicity - 2: WAsP 12 profile model, including sectorwise stability as described in mesoclimate class - 3: WAsP 12 profile model, including sectorwise stability and baroclinicity as described in mesoclimate class The full list of parameters for the climate model is: {climate_params_table} Examples -------- Config objects can be instantiated with a standard parameter-suite, for example WAsP 12.8: >>> import pywasp as pw >>> conf = pw.wasp.Config("WAsP_12.8") >>> print(conf) Config >>> conf.climate[10] = 1.8 # Set the climate parameter 10 (A0) value to 1.8 >>> conf.set_profile_model(-1) # Use the WAsP 11 profile model >>> conf.terrain.set_terrain_analysis # Use spider analysis terrain model ... >>> # Use the Config in a WAsP model-call >>> wwc = pw.wasp.downscale(gwc, topo_map, output_locs, conf=conf) """ # Update the docstring with current terrain and climate parameter tables # in rst format (requires the tabulate package) __doc__ = __doc__.format( terrain_params_table=_params_to_table( TERRAIN_DESC, columns=[ "index", "short_name", "valid_min", "valid_max", "default_value", "definition", ], ).replace( "\n", "\n ", # needed to match docstring indentation ), climate_params_table=_params_to_table( CLIMATE_DESC, columns=[ "index", "short_name", "valid_min", "valid_max", "default_value", "definition", ], ).replace( "\n", "\n ", # needed to match docstring indentation ), )
[docs] def __init__(self, par_set=_DEFAULT_PAR_SET): self.climate = Climate(par_set) self.terrain = Terrain(par_set)
def __str__(self): string = "WAsP Config" string += "\n\nClimate parameters\n==================\n" string += self.climate.__str__() string += "\nTerrain parameters\n==================\n" string += self.terrain.__str__() return string def __eq__(self, comp): return (self.climate == comp.climate) and (self.terrain == comp.terrain)
[docs] def to_dict(self): """Export WAsP configuration parameters to dictionary. Returns ------- dict Dictionary of all parameters and their metadata """ dict_out = {} dict_out["climate"] = self.climate.to_dict() dict_out["terrain"] = self.terrain.to_dict() return dict_out
[docs] def to_file(self, filepath): """Export WAsP configuration parameters to file. .. warning:: This function is experimental and its signature may change. Parameters ---------- filepath : str or Path Path to file that will store WAsP configuration parameters """ warnings.warn( "Writing pywasp configuration to files may lead to inaccurate values when reading them back, use with caution" ) _dict2json(self.to_dict(), filepath)
[docs] @classmethod def from_dict(cls, d): """Initialize class from dictionary This is most often used to create a parameter set saved from file. Parameters ---------- d : dict Stacked dictionary of the configuration parameters Returns ------- self : self Object filled with the requested fields """ self = cls() self.climate = self.climate.from_dict(d["climate"]) self.terrain = self.terrain.from_dict(d["terrain"]) return self
[docs] @classmethod def read_file(cls, filepath): """Initialize class from exported JSON file .. warning:: This function is experimental and its signature may change. Parameters ---------- filepath : str or pathlib.Path Path to file Returns ------- self : self Object filled with the requested fields """ warnings.warn( "Reading pywasp configuration from files may lead to inaccurate values, use with caution" ) d = _json2dict(filepath) return cls.from_dict(d)